fixing the vfp running with flow_ebos.

opm/autodiff/StandardWell.hpp

@@ -99,7 +99,6 @@ namespace Opm
 
         virtual void init(const PhaseUsage* phase_usage_arg,
                           const std::vector<bool>* active_arg,
-                          const VFPProperties* vfp_properties_arg,
                           const std::vector<double>& depth_arg,
                           const double gravity_arg,
                           const int num_cells);
@@ -161,6 +160,7 @@ namespace Opm
                                            std::vector<double>& well_potentials) const;
     protected:
 
+        // protected functions from the Base class
         using Base::phaseUsage;
         using Base::active;
         using Base::flowToEbosPvIdx;
@@ -172,7 +172,7 @@ namespace Opm
         using Base::wellHasTHPConstraints;
         using Base::mostStrictBhpFromBhpLimits;
 
-        // TODO: decide wether to use member function to refer to private member later
+        // protected member variables from the Base class
         using Base::name_;
         using Base::vfp_properties_;
         using Base::gravity_;
@@ -201,8 +201,6 @@ namespace Opm
         // pressure drop between different perforations
         std::vector<double> perf_pressure_diffs_;
 
-        // TODO: probably, they should be moved to the WellInterface, then
-        // we decide the template paramters.
         // two off-diagonal matrices
         OffDiagMatWell duneB_;
         OffDiagMatWell duneC_;
@@ -214,6 +212,7 @@ namespace Opm
         mutable BVectorWell invDrw_;
         mutable BVector scaleAddRes_;
 
+        // residuals of the well equations
         BVectorWell resWell_;
 
         std::vector<EvalWell> well_variables_;
@@ -235,7 +234,7 @@ namespace Opm
 
         EvalWell extendEval(const Eval& in) const;
 
-        // TODO: maybe this function can go to some helper file.
+        // TODO: maybe this type of function can go to some helper file.
         void localInvert(DiagMatWell& istlA) const;
 
         // xw = inv(D)*(rw - C*x)
@@ -250,7 +249,7 @@ namespace Opm
                                                          std::vector<double>& rvmax_perf,
                                                          std::vector<double>& surf_dens_perf) const;
 
-        // TODO: not total sure whether it is a good idea to put here
+        // TODO: not total sure whether it is a good idea to put this function here
         // the major reason to put here is to avoid the usage of Wells struct
         void computeConnectionDensities(const std::vector<double>& perfComponentRates,
                                         const std::vector<double>& b_perf,
@@ -270,7 +269,7 @@ namespace Opm
                                      const ModelParameters& param,
                                      WellState& well_state);
 
-        // TODO: maybe we should provide a light version of computeWellFlux, which does not include the
+        // TODO: maybe we should provide a light version of computePerfRate, which does not include the
         // calculation of the derivatives
         void computeWellRatesWithBhp(const Simulator& ebosSimulator,
                                      const EvalWell& bhp,

opm/autodiff/StandardWell_impl.hpp

@@ -45,14 +45,12 @@ namespace Opm
     StandardWell<TypeTag>::
     init(const PhaseUsage* phase_usage_arg,
          const std::vector<bool>* active_arg,
-         const VFPProperties* vfp_properties_arg,
          const std::vector<double>& depth_arg,
          const double gravity_arg,
          const int num_cells)
     {
         Base::init(phase_usage_arg, active_arg,
-                   vfp_properties_arg, depth_arg,
+                   depth_arg, gravity_arg, num_cells);
-                   gravity_arg, num_cells);
 
         perf_depth_.resize(number_of_perforations_, 0.);
         for (int perf = 0; perf < number_of_perforations_; ++perf) {
@@ -166,8 +164,7 @@ namespace Opm
 
             // pick the density in the top layer
             const double rho = perf_densities_[0];
-            // TODO: not sure whether it is always correct
+            const double well_ref_depth = ref_depth_;
-            const double well_ref_depth = perf_depth_[0];
             const double dp = wellhelpers::computeHydrostaticCorrection(well_ref_depth, vfp_ref_depth, rho, gravity_);
             bhp -= dp;
             return bhp;
@@ -989,7 +986,7 @@ namespace Opm
                     const WellType& well_type = well_type_;
                     // pick the density in the top layer
                     const double rho = perf_densities_[0];
-                    const double well_ref_depth = perf_depth_[0];
+                    const double well_ref_depth = ref_depth_;
 
                     if (well_type == INJECTOR) {
                         const double vfp_ref_depth = vfp_properties_->getInj()->getTable(vfp)->getDatumDepth();
@@ -1082,7 +1079,7 @@ namespace Opm
 
                     const WellType& well_type = well_type_;
                     const double rho = perf_densities_[0];
-                    const double well_ref_depth = perf_depth_[0];
+                    const double well_ref_depth = ref_depth_;
                     if (well_type == INJECTOR) {
                         const double vfp_ref_depth = vfp_properties_->getInj()->getTable(table_id)->getDatumDepth();
 
@@ -1184,7 +1181,7 @@ namespace Opm
 
             // pick the density in the top layer
             const double rho = perf_densities_[0];
-            const double well_ref_depth = perf_depth_[0];
+            const double well_ref_depth = ref_depth_;
 
             // TODO: make the following a function and we call it so many times.
             if (well_type_ == INJECTOR) {
@@ -2027,7 +2024,7 @@ namespace Opm
                     // Calculating the BHP value based on THP
                     // TODO: check whether it is always correct to do calculation based on the depth of the first perforation.
                     const double rho = perf_densities_[0]; // TODO: this item is the one keeping the function from WellInterface
-                    const double well_ref_depth = perf_depth_[0];
+                    const double well_ref_depth = ref_depth_;
                     if (well_type_ == INJECTOR) {
                         const double vfp_ref_depth = vfp_properties_->getInj()->getTable(vfp)->getDatumDepth();
                         const double dp = wellhelpers::computeHydrostaticCorrection(well_ref_depth, vfp_ref_depth, rho, gravity_);

opm/autodiff/StandardWellsDense.hpp

@@ -258,8 +258,6 @@ namespace Opm {
 
             PhaseUsage phase_usage_;
             std::vector<bool>  active_;
-            const VFPProperties* vfp_properties_;
-            double gravity_;
             const RateConverterType& rate_converter_;
 
             std::vector<double> pv_;

opm/autodiff/StandardWellsDense_impl.hpp

@@ -56,7 +56,6 @@ namespace Opm {
 
         phase_usage_ = phase_usage_arg;
         active_ = active_arg;
-        gravity_ = gravity_arg;
         pv_ = pv_arg;
 
         calculateEfficiencyFactors();
@@ -84,7 +83,7 @@ namespace Opm {
         // TODO: to see whether we can postpone of the intialization of the well containers to
         // optimize the usage of the following several member variables
         for (auto& well : well_container_) {
-            well->init(&phase_usage_, &active_, vfp_properties_, depth_arg, gravity_, nc);
+            well->init(&phase_usage_, &active_, depth_arg, gravity_arg, nc);
         }
     }
 
@@ -97,7 +96,9 @@ namespace Opm {
     StandardWellsDense<TypeTag>::
     setVFPProperties(const VFPProperties*  vfp_properties_arg)
     {
-        vfp_properties_ = vfp_properties_arg;
+        for (auto& well : well_container_) {
+            well->setVFPProperties(vfp_properties_arg);
+        }
     }
 
 
@@ -781,6 +782,8 @@ namespace Opm {
             WellControls* wc = well_container_[w]->wellControls();
             const int control = well_controls_get_current(wc);
             well_state.currentControls()[w] = control;
+            // TODO: for VFP control, the perf_densities are still zero here, investigate better
+            // way to handle it later.
             well_container_[w]->updateWellStateWithTarget(control, well_state);
 
             // The wells are not considered to be newly added

opm/autodiff/WellInterface.hpp

@@ -96,9 +96,10 @@ namespace Opm
         /// Well controls
         WellControls* wellControls() const;
 
+        void setVFPProperties(const VFPProperties* vfp_properties_arg);
+
         virtual void init(const PhaseUsage* phase_usage_arg,
                           const std::vector<bool>* active_arg,
-                          const VFPProperties* vfp_properties_arg,
                           const std::vector<double>& depth_arg,
                           const double gravity_arg,
                           const int num_cells);

opm/autodiff/WellInterface_impl.hpp

@@ -97,14 +97,12 @@ namespace Opm
     WellInterface<TypeTag>::
     init(const PhaseUsage* phase_usage_arg,
          const std::vector<bool>* active_arg,
-         const VFPProperties* vfp_properties_arg,
          const std::vector<double>& /* depth_arg */,
          const double gravity_arg,
          const int /* num_cells */)
     {
         phase_usage_ = phase_usage_arg;
         active_ = active_arg;
-        vfp_properties_ = vfp_properties_arg;
         gravity_ = gravity_arg;
     }
 
@@ -112,6 +110,18 @@ namespace Opm
 
 
 
+    template<typename TypeTag>
+    void
+    WellInterface<TypeTag>::
+    setVFPProperties(const VFPProperties* vfp_properties_arg)
+    {
+        vfp_properties_ = vfp_properties_arg;
+    }
+
+
+
+
+
     template<typename TypeTag>
     const std::string&
     WellInterface<TypeTag>::