OilfieldToolsNet/opm-simulators
4
0
Fork 0
mirror of https://github.com/OPM/opm-simulators.git synced 2025-02-25 18:55:30 -06:00
Code Issues Packages Projects Releases Wiki Activity

Start using face geometry computed the Ecl way

Browse Source 
Face centers are computed using the cell corners. With this
implementation the face center seen from a cell may be different from the
face center seen from its neighbour.

Face normals with area lenghts are calculated using the face corners
directly not using a triangulation point in the center of the face.

This gives transmissibility almost equal to eclipse.
This commit is contained in:
Tor Harald Sandve
2017-01-19 08:27:01 +01:00
parent 5955cb07ba
commit e4f345d124
1 changed files with 6 additions and 17 deletions
Download Patch File Download Diff File Expand all files Collapse all files

23
opm/autodiff/GeoProps.hpp
View File

@@ -615,27 +615,15 @@ namespace Opm
double sgn = 2.0 * (faceCells(faceIdx, 0) == cellIdx) - 1;
const int dim = Opm::UgGridHelpers::dimensions(grid);
const double* faceNormal = Opm::UgGridHelpers::faceNormal(grid, faceIdx);
#if HAVE_OPM_GRID
assert( dim <= 3 );
Dune::FieldVector< double, 3 > scaledFaceNormal( 0 );
for (int indx = 0; indx < dim; ++indx) {
scaledFaceNormal[ indx ] = faceNormal[ indx ];
}
// compute unit normal incase the normal is already scaled
scaledFaceNormal /= scaledFaceNormal.two_norm();
// compute proper normal scaled with face area
scaledFaceNormal *= Opm::UgGridHelpers::faceArea(grid, faceIdx);
#else
const double* scaledFaceNormal = faceNormal;
#endif
int cartesianCellIdx = AutoDiffGrid::globalCell(grid)[cellIdx];
auto cellCenter = eclGrid.getCellCenter(cartesianCellIdx);
const auto& faceCenter = Opm::UgGridHelpers::faceCenterEcl(grid, cellIdx, faceTag);
const auto& faceAreaNormalEcl = Opm::UgGridHelpers::faceAreaNormalEcl(grid, faceIdx);
for (int indx = 0; indx < dim; ++indx) {
const double Ci = Opm::UgGridHelpers::faceCentroid(grid, faceIdx)[indx] - cellCenter[indx];
const double Ci = faceCenter[indx] - cellCenter[indx];
dist += Ci*Ci;
cn += sgn * Ci * scaledFaceNormal[ indx ];
cn += sgn * Ci * faceAreaNormalEcl[ indx ];
}
if (cn < 0){
@@ -656,6 +644,7 @@ namespace Opm
cn = -cn;
}
hTrans[cellFaceIdx] = perm[cellIdx*dim*dim + d] * cn / dist;
}
}