From e7e5d946288850a26b3ed99ce3afbd6b54c40d1c Mon Sep 17 00:00:00 2001
From: =?UTF-8?q?B=C3=A5rd=20Skaflestad?= <Bard.Skaflestad@sintef.no>
Date: Tue, 20 Dec 2011 02:13:50 +0100
Subject: [PATCH] Compute perforation fluxes (no gravity adjustments).

Requires additional parameter 'pmobc' (phase mobilities per cell) in
public function cfs_tpfa_res_flux().
---
 opm/core/pressure/tpfa/cfs_tpfa_residual.c | 43 ++++++++++++++++++++--
 opm/core/pressure/tpfa/cfs_tpfa_residual.h |  1 +
 2 files changed, 40 insertions(+), 4 deletions(-)

diff --git a/opm/core/pressure/tpfa/cfs_tpfa_residual.c b/opm/core/pressure/tpfa/cfs_tpfa_residual.c
index 317732eed..f09066ebb 100644
--- a/opm/core/pressure/tpfa/cfs_tpfa_residual.c
+++ b/opm/core/pressure/tpfa/cfs_tpfa_residual.c
@@ -987,6 +987,38 @@ compute_flux(grid_t       *G,
 }
 
 
+static void
+compute_wflux(int                        np    ,
+              struct cfs_tpfa_res_wells *W     ,
+              const double              *pmobc ,
+              const double              *cpress,
+              const double              *wpress,
+              double                    *wflux )
+{
+    int           w, c, i, p;
+    double        pw, dp, t;
+    const double *pmob;
+
+    for (w = i = 0; w < W->conn->number_of_wells; w++) {
+        pw = wpress[w];
+
+        for (; i < W->conn->well_connpos[w + 1]; i++) {
+            c  = W->conn->well_cells[ i ];
+            dp = pw - cpress[c];
+
+            if (dp > 0) { pmob = W->data->phasemob + (i * np); } /* w->c */
+            else        { pmob = pmobc             + (c * np); } /* c->w */
+
+            for (p = 0, t = 0.0; p < np; p++) {
+                t += pmob[ p ];
+            }
+
+            wflux[ i ] = W->data->WI[i] * t * dp;
+        }
+    }
+}
+
+
 static size_t
 maxconn(grid_t *G)
 {
@@ -1141,6 +1173,7 @@ cfs_tpfa_res_flux(grid_t                     *G        ,
                   struct cfs_tpfa_res_forces *forces   ,
                   int                         np       ,
                   const double               *trans    ,
+                  const double               *pmobc    ,
                   const double               *pmobf    ,
                   const double               *gravcap_f,
                   const double               *cpress   ,
@@ -1149,11 +1182,13 @@ cfs_tpfa_res_flux(grid_t                     *G        ,
                   double                     *wflux    )
 /* ---------------------------------------------------------------------- */
 {
-    (void) forces;
-    (void) wpress;
-    (void) wflux;
-
     compute_flux(G, np, trans, pmobf, gravcap_f, cpress, fflux);
+
+    if ((forces != NULL) && (forces->W != NULL) &&
+        (wpress != NULL) && (wflux     != NULL)) {
+
+        compute_wflux(np, forces->W, pmobc, cpress, wpress, wflux);
+    }
 }
 
 
diff --git a/opm/core/pressure/tpfa/cfs_tpfa_residual.h b/opm/core/pressure/tpfa/cfs_tpfa_residual.h
index b2a60de5d..1401b5479 100644
--- a/opm/core/pressure/tpfa/cfs_tpfa_residual.h
+++ b/opm/core/pressure/tpfa/cfs_tpfa_residual.h
@@ -79,6 +79,7 @@ cfs_tpfa_res_flux(grid_t                     *G        ,
                   struct cfs_tpfa_res_forces *forces   ,
                   int                         np       ,
                   const double               *trans    ,
+                  const double               *pmobc    ,
                   const double               *pmobf    ,
                   const double               *gravcap_f,
                   const double               *cpress   ,