1 #include <opm/common/utility/numeric/blas_lapack.h>
3 #include <opm/grid/GridHelpers.hpp>
14 namespace UgGridHelpers
24 inline const double* multiplyFaceNormalWithArea(
const Dune::CpGrid& grid,
int face_index,
const double* in)
26 int d=Opm::UgGridHelpers::dimensions(grid);
27 double* out=
new double[d];
28 double area=Opm::UgGridHelpers::faceArea(grid, face_index);
35 inline void maybeFreeFaceNormal(
const Dune::CpGrid&,
const double* array)
40 inline const double* multiplyFaceNormalWithArea(
const UnstructuredGrid&,
int,
const double* in)
57 using namespace Opm::UgGridHelpers;
59 double s, dist, denom;
69 MAT_SIZE_T nrows, ncols, ldA, incx, incy;
74 nrows = ncols = ldA = d;
78 for (
int c =0, i = 0; c < numCells(*G); c++) {
79 K = perm + (c * d * d);
82 FaceRow faces = c2f[c];
84 for(
typename FaceRow::const_iterator f=faces.begin(), end=faces.end();
87 s = 2.0*(face_cells(*f, 0) == c) - 1.0;
88 n = faceNormal(*G, *f);
89 const double* nn=multiplyFaceNormalWithArea(*G, *f, n);
90 const double* fc = &(faceCentroid(*G, *f)[0]);
91 dgemv_(
"No Transpose", &nrows, &ncols,
92 &a1, K, &ldA, nn, &incx, &a2, &Kn[0], &incy);
93 maybeFreeFaceNormal(*G, nn);
95 htrans[i] = denom = 0.0;
96 for (j = 0; j < d; j++) {
97 dist = fc[j] - getCoordinate(cc, j);
99 htrans[i] += s * dist * Kn[j];
100 denom += dist * dist;
105 htrans[i] = std::abs(htrans[i]);
108 cc = increment(cc, 1, d);
119 using namespace Opm::UgGridHelpers;
121 for (
int f = 0; f < numFaces(*G); f++) {
127 for (
int c = 0, i = 0; c < numCells(*G); c++) {
129 FaceRow faces = c2f[c];
131 for(
typename FaceRow::const_iterator f=faces.begin(), end=faces.end();
134 trans[*f] += 1.0 / htrans[i];
138 for (
int f = 0; f < numFaces(*G); f++) {
139 trans[f] = 1.0 / trans[f];
148 const double *totmob,
149 const double *htrans,
153 using namespace Opm::UgGridHelpers;
155 for (
int f = 0; f < numFaces(*G); f++) {
161 for (
int c = 0, i = 0; c < numCells(*G); c++) {
163 FaceRow faces = c2f[c];
165 for(
typename FaceRow::const_iterator f=faces.begin(), end=faces.end();
168 trans[*f] += 1.0 / (totmob[c] * htrans[i]);
172 for (
int f = 0; f < numFaces(*G); f++) {
173 trans[f] = 1.0 / trans[f];
[ provides Dune::Grid ]
Definition: CpGrid.hpp:207
Copyright 2019 Equinor AS.
Definition: CartesianIndexMapper.hpp:10
T area(const Point< T, 2 > *c)
Computes the area of a 2-dimensional triangle.
Definition: Volumes.hpp:118
Holds the implementation of the CpGrid as a pimple.
Definition: CellQuadrature.hpp:29
Maps the grid type to the associated type of the cell to faces mapping.
Definition: GridHelpers.hpp:277
Traits of the cell centroids of a grid.
Definition: GridHelpers.hpp:126
Traits of the face to attached cell mappping of a grid.
Definition: GridHelpers.hpp:334
Data structure for an unstructured grid, unstructured meaning that any cell may have an arbitrary num...
Definition: UnstructuredGrid.h:99
Routines to assist in the calculation of two-point transmissibilities.
void tpfa_eff_trans_compute(struct UnstructuredGrid *G, const double *totmob, const double *htrans, double *trans)
Calculate effective two-point transmissibilities from one-sided, total mobility weighted,...
Definition: trans_tpfa.c:103
void tpfa_htrans_compute(struct UnstructuredGrid *G, const double *perm, double *htrans)
Calculate static, one-sided transmissibilities for use in the two-point flux approximation method.
Definition: trans_tpfa.c:19
void tpfa_trans_compute(struct UnstructuredGrid *G, const double *htrans, double *trans)
Compute two-point transmissibilities from one-sided transmissibilities.
Definition: trans_tpfa.c:74