From 2dffbb3a70577a006fdec1c13f5c3ce0ec14e0d6 Mon Sep 17 00:00:00 2001
From: =?UTF-8?q?B=C3=A5rd=20Skaflestad?= <Bard.Skaflestad@sintef.no>
Date: Fri, 24 May 2013 17:22:35 +0200
Subject: [PATCH] Add gather/scatter support for wells

---
 opm/autodiff/FullyImplicitBlackoilSolver.cpp | 29 ++++++++++++++++++++
 opm/autodiff/FullyImplicitBlackoilSolver.hpp |  7 +++++
 2 files changed, 36 insertions(+)

diff --git a/opm/autodiff/FullyImplicitBlackoilSolver.cpp b/opm/autodiff/FullyImplicitBlackoilSolver.cpp
index 55a40a8fa..f9a477b2b 100644
--- a/opm/autodiff/FullyImplicitBlackoilSolver.cpp
+++ b/opm/autodiff/FullyImplicitBlackoilSolver.cpp
@@ -147,6 +147,7 @@ namespace Opm {
         , canph_ (active2Canonical(fluid.phaseUsage()))
         , cells_ (buildAllCells(grid.number_of_cells))
         , ops_   (grid)
+        , wops_  (wells)
         , grav_  (gravityOperator(grid_, ops_, geo_))
         , rq_    (fluid.numPhases())
     {
@@ -215,6 +216,34 @@ namespace Opm {
     }
 
 
+    FullyImplicitBlackoilSolver::
+    WellOps::WellOps(const Wells& wells)
+        : w2p(wells.well_connpos[ wells.number_of_wells ],
+              wells.number_of_wells)
+        , p2w(wells.number_of_wells,
+              wells.well_connpos[ wells.number_of_wells ])
+    {
+        const int        nw   = wells.number_of_wells;
+        const int* const wpos = wells.well_connpos;
+
+        typedef Eigen::Triplet<double> Tri;
+
+        std::vector<Tri> scatter, gather;
+        scatter.reserve(wpos[nw]);
+        gather .reserve(wpos[nw]);
+
+        for (int w = 0, i = 0; w < nw; ++w) {
+            for (; i < wpos[ w + 1 ]; ++i) {
+                scatter.push_back(Tri(i, w, 1.0));
+                gather .push_back(Tri(w, i, 1.0));
+            }
+        }
+
+        w2p.setFromTriplets(scatter.begin(), scatter.end());
+        p2w.setFromTriplets(gather .begin(), gather .end());
+    }
+
+
     void
     FullyImplicitBlackoilSolver::allocateResidual()
     {
diff --git a/opm/autodiff/FullyImplicitBlackoilSolver.hpp b/opm/autodiff/FullyImplicitBlackoilSolver.hpp
index e1f539f0a..7773e50bf 100644
--- a/opm/autodiff/FullyImplicitBlackoilSolver.hpp
+++ b/opm/autodiff/FullyImplicitBlackoilSolver.hpp
@@ -80,6 +80,12 @@ namespace Opm {
             ADB              Rs;
         };
 
+        struct WellOps {
+            WellOps(const Wells& wells);
+            M w2p;              // well -> perf (scatter)
+            M p2w;              // perf -> well (gather)
+        };
+
         enum { Water = BlackoilPropsAdInterface::Water,
                Oil   = BlackoilPropsAdInterface::Oil  ,
                Gas   = BlackoilPropsAdInterface::Gas  };
@@ -95,6 +101,7 @@ namespace Opm {
         const std::vector<int>          canph_;
         const std::vector<int>          cells_;  // All grid cells
         HelperOps                       ops_;
+        const WellOps                   wops_;
         const M                         grav_;
 
         std::vector<ReservoirResidualQuant> rq_;