[satfunc] Relperm scaling

[berland]

User feature request: Combine endpoint scaling data from INIT and relperm-data from DATA files, and then provide an API to give out the scaled relperm data as a dataframe.

It is also possible to extract the relperm from the binary output files from Eclipse, that is not supported yet.

[berland]

Old proof-of-concept for extracting relperm from INIT files.

#!/usr/bin/env python2.7
# This code is licensed under GPLv3 [vkip]

import argparse
import sys

from numpy import array, logical_and
from matplotlib.mlab import find

from ecl.ecl.ecl_file import EclFile

## INTEHEAD INDICES
IPHS = 14
OWG = [ [True, False, False],
        [False, True, False],
        [True, True, False], 
        [False, False, True], 
        [True, False, True], 
        [False, True, True],
        [True, True, True] ]

## TABDIMS INDICES
RPWBASE = 20
RPWNN   = 21
RPWNT   = 22

RPGBASE = 23
RPGNN   = 24
RPGNT   = 25

RPOBASE = 26
RPONN   = 28
RPONT   = 29

def fix_tab(tabs, ix):
    t = tabs[:,ix,:]
    iok = find(logical_and(t[:,0] >=-1.e21, t[:,0] <= 1.0e21))
    return t[iok,:]

def print_tabs(tabs, headers):
    h = reduce(lambda a,b: a+b, [ "%-21s " % (x) for x in headers ])
    for i,t in enumerate(tabs):
        print "-- TABLE NUMBER %d " % (i+1)
        print "-- %s" % (h)
        for j in range(t.shape[0]):
            vstr = reduce(lambda a,b: a+b, [ "%.15e " % (x) for x in t[j,:] ])
            print "   %s" % (vstr)
        print "-- "
        print ""

def main(argv = sys.argv):
    parser = argparse.ArgumentParser(description="Extract relative permeabilities from init file")
    parser.add_argument("init_file", help="INIT file name (including suffix)")
    args = parser.parse_args(argv[1:])

    init = EclFile(args.init_file)
    ih = array(init['INTEHEAD'][0])
    lh = array(init['LOGIHEAD'][0])
    tabdims = array(init['TABDIMS'][0])
    tab = array(init['TAB'][0])

    oil, water, gas = OWG[ih[IPHS]-1]
    threephase = oil and water and gas

    otabs = wtabs = gtabs = None

    if water: 
        ncols = 3
        rpwbase = tabdims[RPWBASE]
        rpwnn = tabdims[RPWNN]
        rpwnt = tabdims[RPWNT]
        wtab_ = array(tab[(rpwbase-1):(rpwbase+ncols*rpwnn*rpwnt)-1]).reshape(rpwnn, rpwnt, ncols, order='F')
        wtabs = [ fix_tab(wtab_, i) for i in range(rpwnt) ]
        print "######################################################################"
        print "# WATER "
        print "######################################################################"
        print_tabs(wtabs, ["Sw", "Krw", "Pcow"])

    if gas: 
        ncols = 3
        rpgbase = tabdims[RPGBASE]
        rpgnn = tabdims[RPGNN]
        rpgnt = tabdims[RPGNT]
        gtab_ = array(tab[(rpgbase-1):(rpgbase+ncols*rpgnn*rpgnt)-1]).reshape(rpgnn, rpgnt, ncols, order='F')
        gtabs = [ fix_tab(gtab_, i) for i in range(rpgnt) ]
        print "######################################################################"
        print "# GAS "
        print "######################################################################"
        print_tabs(gtabs, ["Sg", "Krg", "Pcog"])

    if oil:
        ncols = 3 if threephase else 2
        rpobase = tabdims[RPOBASE]
        rponn = tabdims[RPONN]
        rpont = tabdims[RPONT]
        otab_ = array(tab[(rpobase-1):(rpobase+ncols*rponn*rpont)-1]).reshape(rponn, rpont, ncols, order='F')
        otabs = [ fix_tab(otab_, i) for i in range(rpont) ]
        headers = [ "So", "Krow", "Krog" ] if threephase else [ "So", "Kro" ]
        print "######################################################################"
        print "# OIL "
        print "######################################################################"
        print_tabs(otabs, headers)

if __name__ == '__main__':
    sys.exit(main())