Made IncompTpfa receive wells-output-arguments

examples/spu_2p.cpp

@@ -510,6 +510,7 @@ main(int argc, char** argv)
         } else {
             computeTotalMobility(*props, allcells, state.saturation(), totmob);
         }
+        std::vector<double> empty_vector_for_wells;
         pressure_timer.start();
         if (rock_comp->isActive()) {
             rc.resize(num_cells);
@@ -521,8 +522,10 @@ main(int argc, char** argv)
                     rc[cell] = rock_comp->rockComp(state.pressure()[cell]);
                 }
                 state.pressure() = initial_pressure;
+                
+
                 psolver.solve(totmob, omega, src, bcs.c_bcs(), porevol, rc, simtimer.currentStepLength(),
-                              state.pressure(), state.faceflux());
+                              state.pressure(), state.faceflux(), empty_vector_for_wells, empty_vector_for_wells);
                 double max_change = 0.0;
                 for (int cell = 0; cell < num_cells; ++cell) {
                     max_change = std::max(max_change, std::fabs(state.pressure()[cell] - prev_pressure[cell]));
@@ -534,7 +537,8 @@ main(int argc, char** argv)
             }
             computePorevolume(*grid->c_grid(), *props, *rock_comp, state.pressure(), porevol);
         } else {
-            psolver.solve(totmob, omega, src, bcs.c_bcs(), state.pressure(), state.faceflux());
+            psolver.solve(totmob, omega, src, bcs.c_bcs(), state.pressure(), state.faceflux(), empty_vector_for_wells,
+                          empty_vector_for_wells);
         }
         pressure_timer.stop();
         double pt = pressure_timer.secsSinceStart();

opm/core/ProductionSpecification.cpp

@@ -10,7 +10,8 @@ namespace Opm
       oil_max_rate_(-1.0),
       water_production_target_(-1.0),
       fluid_volume_max_rate_(-1.0),
-      BHP_limit_(-1.0)
+      BHP_limit_(-1.0),
+            guide_rate_(-1.0)
     {
     }
 

opm/core/ProductionSpecification.hpp

@@ -18,6 +18,11 @@ namespace Opm
         {
             NONE_P, RATE, WELL
         };
+        
+        enum GuideRateType
+        {
+            OIL, RAT
+        };
 
         ProductionSpecification();
  
@@ -30,6 +35,8 @@ namespace Opm
         double BHP_limit_;
         double gas_max_rate_;
         double liquid_max_rate_;
+        double guide_rate_;
+        GuideRateType guide_rate_type_;
 
     };
 }

opm/core/WellsGroup.cpp

@@ -97,9 +97,7 @@ namespace Opm
             }
         }
         
-        
+        return true;
-        
-        
     }
 
     WellsGroupInterface* WellNode::findGroup(std::string name_of_node)

opm/core/WellsManager.cpp

@@ -527,6 +527,11 @@ namespace Opm
             // We know that getLeafNodes() is ordered the same way as they're indexed in w_
             node->setWellsPointer(w_, i);
         }
+        
+        // Set the guide rates:
+        if(deck.hasField("WCONPROD")) {
+            WCONPROD wconprod = deck.getWCONPROD();
+        }
     }
 
 

opm/core/pressure/IncompTpfa.cpp

@@ -91,12 +91,18 @@ namespace Opm
     ///                        If null, noflow conditions are assumed.
     /// \param[out] pressure   Will contain N cell-pressure values.
     /// \param[out] faceflux   Will contain F signed face flux values.
+    /// \param[out] well_bhp   Will contain bhp values for each well passed
+    ///                        in the constructor.
+    /// \param[out] well_rate  Will contain rate values for each well passed
+    ///                        in the constructor.
     void IncompTpfa::solve(const std::vector<double>& totmob,
 			   const std::vector<double>& omega,
 			   const std::vector<double>& src,
 			   const FlowBoundaryConditions* bcs,
 			   std::vector<double>& pressure,
-			   std::vector<double>& faceflux)
+			   std::vector<double>& faceflux,
+                           std::vector<double>& well_bhp,
+                           std::vector<double>& well_rate)
     {
 	UnstructuredGrid* gg = const_cast<UnstructuredGrid*>(&grid_);
 	tpfa_eff_trans_compute(gg, &totmob[0], &htrans_[0], &trans_[0]);
@@ -150,6 +156,10 @@ namespace Opm
     /// \param[in]  dt         Timestep.
     /// \param[out] pressure   Will contain N cell-pressure values.
     /// \param[out] faceflux   Will contain F signed face flux values.
+    /// \param[out] well_bhp   Will contain bhp values for each well passed
+    ///                        in the constructor
+    /// \param[out] well_rate  Will contain rate values for each well passed
+    ///                        in the constructor
     void IncompTpfa::solve(const std::vector<double>& totmob,
                            const std::vector<double>& omega,
                            const std::vector<double>& src,
@@ -158,7 +168,9 @@ namespace Opm
                            const std::vector<double>& rock_comp,
                            const double dt,
                            std::vector<double>& pressure,
-                           std::vector<double>& faceflux)
+                           std::vector<double>& faceflux,
+                           std::vector<double>& well_bhp,
+                           std::vector<double>& well_rate)
     {
 	UnstructuredGrid* gg = const_cast<UnstructuredGrid*>(&grid_);
 	tpfa_eff_trans_compute(gg, &totmob[0], &htrans_[0], &trans_[0]);

opm/core/pressure/IncompTpfa.hpp

@@ -72,12 +72,18 @@ namespace Opm
 	///                        If null, noflow conditions are assumed.
 	/// \param[out] pressure   Will contain N cell-pressure values.
 	/// \param[out] faceflux   Will contain F signed face flux values.
+        /// \param[out] well_bhp   Will contain bhp values for each well passed
+        ///                        in the constructor
+        /// \param[out] well_rate  Will contain rate values for each well passed
+        ///                        in the constructor
 	void solve(const std::vector<double>& totmob,
 		   const std::vector<double>& omega,
 		   const std::vector<double>& src,
 		   const FlowBoundaryConditions* bcs,
 		   std::vector<double>& pressure,
-		   std::vector<double>& faceflux);
+		   std::vector<double>& faceflux, 
+                   std::vector<double>& well_bhp,
+                   std::vector<double>& well_rate);
 
         /// Assemble and solve pressure system with rock compressibility (assumed constant per cell).
         /// \param[in]  totmob     Must contain N total mobility values (one per cell).
@@ -97,6 +103,10 @@ namespace Opm
         /// \param[in]  dt         Timestep.
         /// \param[out] pressure   Will contain N cell-pressure values.
         /// \param[out] faceflux   Will contain F signed face flux values.
+        /// \param[out] well_bhp   Will contain bhp values for each well passed
+        ///                        in the constructor
+        /// \param[out] well_rate  Will contain rate values for each well passed
+        ///                        in the constructor
         void solve(const std::vector<double>& totmob,
                    const std::vector<double>& omega,
                    const std::vector<double>& src,
@@ -105,7 +115,9 @@ namespace Opm
                    const std::vector<double>& rock_comp,
                    const double dt,
                    std::vector<double>& pressure,
-                   std::vector<double>& faceflux);
+                   std::vector<double>& faceflux,
+                   std::vector<double>& well_bhp,
+                   std::vector<double>& well_rate);
 
         /// Expose read-only reference to internal half-transmissibility.
         const ::std::vector<double>& getHalfTrans() const { return htrans_; }