Merge pull request #1388 from tskille/pybind_eclipse_grid

additional python bindings for eclipse_grid

python/cxx/eclipse_grid.cpp

@@ -4,6 +4,8 @@
 #include <opm/parser/eclipse/EclipseState/Grid/Fault.hpp>
 #include <opm/parser/eclipse/EclipseState/Grid/FaceDir.hpp>
 
+#include <pybind11/stl.h>
+#include "converters.hpp"
 #include "export.hpp"
 
 
@@ -39,6 +41,73 @@ namespace {
       return grid.getCellVolume(i_idx, j_idx, k_idx);
     }
 
+    py::array cellVolumeAll( const EclipseGrid& grid)
+    {
+        std::vector<double> cellVol;
+        std::array<int, 3> dims = grid.getNXYZ();        
+        size_t nCells = dims[0]*dims[1]*dims[2];
+        cellVol.reserve(nCells);
+        
+        for (size_t n = 0; n < nCells; n++)
+            cellVol.push_back(grid.getCellVolume(n));
+        
+        return convert::numpy_array(cellVol);
+    }
+
+    py::array cellVolumeMask( const EclipseGrid& grid, std::vector<int>& mask)
+    {
+        std::array<int, 3> dims = grid.getNXYZ();        
+        size_t nCells = dims[0]*dims[1]*dims[2];
+    
+        if (nCells != mask.size()) 
+            throw std::logic_error("size of input mask doesn't match size of grid");
+        
+        std::vector<double> cellVol(nCells, 0.0);
+
+        for (size_t n = 0; n < nCells; n++)
+            if (mask[n]==1)
+                cellVol[n] = grid.getCellVolume(n);
+                
+        return convert::numpy_array(cellVol);
+    }
+    
+    double cellDepth1G( const EclipseGrid& grid, size_t glob_idx) {
+      return grid.getCellDepth(glob_idx);
+    }
+
+    double cellDepth3( const EclipseGrid& grid, size_t i_idx, size_t j_idx, size_t k_idx) {
+      return grid.getCellDepth(i_idx, j_idx, k_idx);
+    }
+    
+    py::array cellDepthAll( const EclipseGrid& grid)
+    {
+        std::vector<double> cellDepth;
+        std::array<int, 3> dims = grid.getNXYZ();        
+        size_t nCells = dims[0]*dims[1]*dims[2];
+        cellDepth.reserve(nCells);
+        
+        for (size_t n = 0; n < nCells; n++)
+            cellDepth.push_back(grid.getCellDepth(n));
+        
+        return convert::numpy_array(cellDepth);
+    }
+
+    py::array cellDepthMask( const EclipseGrid& grid, std::vector<int>& mask)
+    {
+        std::array<int, 3> dims = grid.getNXYZ();        
+        size_t nCells = dims[0]*dims[1]*dims[2];
+    
+        if (nCells != mask.size()) 
+            throw std::logic_error("size of input mask doesn't match size of grid");
+        
+        std::vector<double> cellDepth(nCells, 0.0);
+
+        for (size_t n = 0; n < nCells; n++)
+            if (mask[n]==1)
+                cellDepth[n] = grid.getCellDepth(n);
+                
+        return convert::numpy_array(cellDepth);
+    }
 }
 
 void python::common::export_EclipseGrid(py::module& module) {
@@ -54,6 +123,12 @@ void python::common::export_EclipseGrid(py::module& module) {
         .def( "getIJK",         &getIJK )
         .def( "getCellVolume",  &cellVolume1G)
         .def( "getCellVolume",   &cellVolume3)
+        .def( "getCellVolume",   &cellVolumeAll)
+        .def( "getCellVolume",   &cellVolumeMask)
+        .def( "getCellDepth",  &cellDepth1G)
+        .def( "getCellDepth",  &cellDepth3)
+        .def( "getCellDepth",  &cellDepthAll)
+        .def( "getCellDepth",  &cellDepthMask)
       ;
 
 }

python/tests/test_props.py

@@ -1,5 +1,6 @@
 import unittest
 import opm.io
+import numpy as np 
 
 from opm.io.parser import Parser
 from opm.io.ecl_state import EclipseState, test_field_props
@@ -61,17 +62,66 @@ class TestProps(unittest.TestCase):
                     self.assertClose(darcy, perm)
 
     def test_volume(self):
+        exp = 293.3 * 293.3 * 30  # cubicfeet = 73 078.6084 cubic meter
+        exp *= (12*0.0254)**3  # cubic feet to cubic meter
         grid  = self.spe3.grid()
         for i in range(grid.nx):
             for j in range(grid.ny):
                 for k in range(grid.nz):
                     g_idx = grid.globalIndex(i,j,k)
-                    exp = 293.3 * 293.3 * 30  # cubicfeet = 73 078.6084 cubic meter
-                    exp *= (12*0.0254)**3  # cubic feet to cubic meter
                     if k == 0:
                         self.assertClose(exp, grid.getCellVolume(g_idx))
                     self.assertEqual(grid.getCellVolume(g_idx), grid.getCellVolume(i, j, k))
 
+        celVol1 = grid.getCellVolume()
+        self.assertTrue(isinstance(celVol1, np.ndarray))
+        self.assertEqual(celVol1.dtype, "float64")
+        self.assertEqual(len(celVol1), grid.nx*grid.ny*grid.nz)
 
+        mask = [0] * (grid.nx*grid.ny*grid.nz)
+        for ind in [0,10,11,12,15]:
+            mask[ind]=1
+            
+        self.assertEqual(len(mask), grid.nx*grid.ny*grid.nz)
+        self.assertEqual(sum(mask), 5)
+
+        celVol2 = grid.getCellVolume(mask)
+        self.assertEqual(len(celVol2), grid.nx*grid.ny*grid.nz)
+        self.assertClose(exp, sum(celVol2)/5.0)
+
+
+    def test_depth(self):
+        
+        refDepth = np.array([7330.0, 7360.0, 7400.0, 7450.0], dtype="float64")
+        
+        grid  = self.spe3.grid()
+        for i in range(grid.nx):
+            for j in range(grid.ny):
+                for k in range(grid.nz):
+                    g_idx = grid.globalIndex(i,j,k)
+                    exp = refDepth[k] * 12*0.0254
+                    g_idx = grid.globalIndex(i,j,k)
+                    self.assertClose(exp, grid.getCellDepth(g_idx))
+                    self.assertClose(exp, grid.getCellDepth(i, j, k))
+        
+
+        depth1 = grid.getCellDepth()
+        self.assertTrue(isinstance(depth1, np.ndarray))
+        self.assertEqual(depth1.dtype, "float64")
+        self.assertEqual(len(depth1), grid.nx*grid.ny*grid.nz)
+
+        mask = [0] * (grid.nx*grid.ny*grid.nz)
+        for k in range(grid.nz):
+            mask[grid.globalIndex(0, 0, k)]=1
+
+        self.assertEqual(len(mask), grid.nx*grid.ny*grid.nz)
+        self.assertEqual(sum(mask), grid.nz)
+
+        depth2 = grid.getCellDepth(mask)
+        self.assertEqual(len(depth2), grid.nx*grid.ny*grid.nz)
+        self.assertClose(sum(refDepth)*12*0.0254, sum(depth2))
+        
+
+               
 if __name__ == "__main__":
     unittest.main()