[bugfix] Take normalized face normals into account when calculating tpfa.

At least for Cpgrid the face normal is normalized which has to be taken into
account in tpfa_htrans_compute.
We therefore multiply the facenormal with the face area in cases (read Cpgrid)
where this needed.

opm/core/grid/GridHelpers.cpp

@@ -58,6 +58,11 @@ const double* faceNormal(const UnstructuredGrid& grid, int face_index)
     return grid.face_normals+face_index*grid.dimensions;
 }
 
+double faceArea(const UnstructuredGrid& grid, int face_index)
+{
+    return grid.face_areas[face_index];
+}
+
 SparseTableView cell2Faces(const UnstructuredGrid& grid)
 {
     return SparseTableView(grid.cell_faces, grid.cell_facepos, numCells(grid));

opm/core/grid/GridHelpers.hpp

@@ -176,6 +176,10 @@ faceCentroid(const UnstructuredGrid& grid, int face_index);
 /// \param face_index The index of the face in the grid.
 const double* faceNormal(const UnstructuredGrid& grid, int face_index);
 
+/// \brief Get the area of a face
+/// \param grid The grid that the face is part of.
+/// \param face_index The index of the face in the grid.
+double faceArea(const UnstructuredGrid& grid, int face_index);
 
 /// \brief Maps the grid type to the associated type of the cell to faces mapping.
 ///

opm/core/pressure/tpfa/TransTpfa_impl.hpp

@@ -8,6 +8,48 @@
 #include <opm/core/pressure/tpfa/trans_tpfa.h>
 #include <opm/core/grid/GridHelpers.hpp>
 
+namespace Dune
+{
+class CpGrid;
+}
+
+namespace Opm
+{
+namespace UgGridHelpers
+{
+int dimensions(const Dune::CpGrid&);
+
+double faceArea(const Dune::CpGrid&, int);
+}
+}
+
+namespace
+{
+const double* multiplyFaceNormalWithArea(const Dune::CpGrid& grid, int face_index, const double* in)
+{
+    int d=Opm::UgGridHelpers::dimensions(grid);
+    double* out=new double[d];
+    double area=Opm::UgGridHelpers::faceArea(grid, face_index);
+    
+    for(int i=0;i<d;++i)
+        out[i]=in[i]*area;
+    return out;
+}
+
+inline const double* multiplyFaceNormalWithArea(const UnstructuredGrid& grid, int face_index, const double* in)
+{
+    return in;
+}
+
+inline void maybeFreeFaceNormal(const Dune::CpGrid&, const double* array)
+{
+    delete[] array;
+}
+
+inline void maybeFreeFaceNormal(const UnstructuredGrid&, const double*)
+{}
+}
+
 /* ---------------------------------------------------------------------- */
 /* htrans <- sum(C(:,i) .* K(cellNo,:) .* N(:,j), 2) ./ sum(C.*C, 2) */
 /* ---------------------------------------------------------------------- */
@@ -47,12 +89,12 @@ tpfa_htrans_compute(const Grid* G, const double *perm, double *htrans)
             f!=end; ++f, ++i)
         {
             s = 2.0*(face_cells(*f, 0) == c) - 1.0;
-
             n = faceNormal(*G, *f);
+            const double* nn=multiplyFaceNormalWithArea(*G, *f, n);
             const double* fc = &(faceCentroid(*G, *f)[0]);
-
             dgemv_("No Transpose", &nrows, &ncols,
-                   &a1, K, &ldA, n, &incx, &a2, &Kn[0], &incy);
+                   &a1, K, &ldA, nn, &incx, &a2, &Kn[0], &incy);
+            maybeFreeFaceNormal(*G, nn);
             
             htrans[i] = denom = 0.0;
             for (j = 0; j < d; j++) {