makd PolymerProperties handle PLYSHEAR keyword.

examples/sim_poly2p_comp_reorder.cpp

@@ -175,9 +175,15 @@ try
         ads_vals[2] = 0.0025;
         // ads_vals[1] = 0.0;
         // ads_vals[2] = 0.0;
+        std::vector<double> water_vel_vals(2, -1e100);
+        water_vel_vals[0] = 0.0;
+        water_vel_vals[1] = 10.0;
+        std::vector<double> shear_vrf_vals(2, -1e100);
+        shear_vrf_vals[0] = 1.0;
+        shear_vrf_vals[1] = 1.0;
         poly_props.set(c_max, mix_param, rock_density, dead_pore_vol, res_factor, c_max_ads,
                        static_cast<Opm::PolymerProperties::AdsorptionBehaviour>(ads_index),
-                       c_vals_visc,  visc_mult_vals, c_vals_ads, ads_vals);
+                       c_vals_visc,  visc_mult_vals, c_vals_ads, ads_vals, water_vel_vals, shear_vrf_vals);
     }
 
     bool use_gravity = (gravity[0] != 0.0 || gravity[1] != 0.0 || gravity[2] != 0.0);

examples/sim_poly2p_incomp_reorder.cpp

@@ -173,9 +173,15 @@ try
         ads_vals[2] = 0.0025;
         // ads_vals[1] = 0.0;
         // ads_vals[2] = 0.0;
+        std::vector<double> water_vel_vals(2, -1e100);
+        water_vel_vals[0] = 0.0;
+        water_vel_vals[1] = 10.0;
+        std::vector<double> shear_vrf_vals(2, -1e100);
+        shear_vrf_vals[0] = 1.0;
+        shear_vrf_vals[1] = 1.0;
         poly_props.set(c_max, mix_param, rock_density, dead_pore_vol, res_factor, c_max_ads,
                        static_cast<Opm::PolymerProperties::AdsorptionBehaviour>(ads_index),
-                       c_vals_visc,  visc_mult_vals, c_vals_ads, ads_vals);
+                       c_vals_visc,  visc_mult_vals, c_vals_ads, ads_vals, water_vel_vals, shear_vrf_vals);
     }
 
     // Warn if gravity but no density difference.

examples/test_singlecellsolves.cpp

@@ -129,9 +129,15 @@ try
     ads_vals[2] = 0.0025;
     // ads_vals[1] = 0.0;
     // ads_vals[2] = 0.0;
+    std::vector<double> water_vel_vals(2, -1e100);
+    water_vel_vals[0] = 0.0;
+    water_vel_vals[1] = 10.0;
+    std::vector<double> shear_vrf_vals(2, -1e100);
+    shear_vrf_vals[0] = 1.0;
+    shear_vrf_vals[1] = 1.0;
     poly_props.set(c_max, mix_param, rock_density, dead_pore_vol, res_factor, c_max_ads,
                    static_cast<Opm::PolymerProperties::AdsorptionBehaviour>(ads_index),
-                   c_vals_visc,  visc_mult_vals, c_vals_ads, ads_vals);
+                   c_vals_visc,  visc_mult_vals, c_vals_ads, ads_vals, water_vel_vals, shear_vrf_vals);
 
     // Initialising src
     int num_cells = grid->c_grid()->number_of_cells;

opm/polymer/PolymerProperties.cpp

@@ -59,6 +59,18 @@ namespace Opm
         return ads_index_;
     }
     
+    std::vector<double>
+    PolymerProperties::shearWaterVelocity() const
+    {
+        return water_vel_vals_;
+    } 
+
+    std::vector<double>
+    PolymerProperties::shearVrf() const
+    {
+        return shear_vrf_vals_;
+    }
+
     double PolymerProperties::viscMult(double c) const
     {
         return Opm::linearInterpolation(c_vals_visc_, visc_mult_vals_, c);

opm/polymer/PolymerProperties.hpp

@@ -49,6 +49,8 @@ namespace Opm
         /// \param[in] visc_mult_vals Array of effective vicosity multiplier
         /// \param[in] c_vals_ads     Array of concentration for adsorption values
         /// \param[in] ads_vals       Array of adsorption values
+        /// \param[in] water_vel_vals_ Array of water phase velocity for shear
+        /// \param[in] shear_vrf_vals_ Array of viscosity reduction factor
         PolymerProperties(double c_max,
                           double mix_param,
                           double rock_density,
@@ -59,7 +61,9 @@ namespace Opm
                           const std::vector<double>& c_vals_visc,
                           const std::vector<double>& visc_mult_vals,
                           const std::vector<double>& c_vals_ads,
-                          const std::vector<double>& ads_vals
+                          const std::vector<double>& ads_vals,
+                          const std::vector<double>& water_vel_vals,
+                          const std::vector<double>& shear_vrf_vals
                           )
             : c_max_(c_max),
               mix_param_(mix_param),
@@ -71,7 +75,9 @@ namespace Opm
               c_vals_visc_(c_vals_visc),
               visc_mult_vals_(visc_mult_vals),
               c_vals_ads_(c_vals_ads),
-              ads_vals_(ads_vals)
+              ads_vals_(ads_vals),
+              water_vel_vals_(water_vel_vals),
+              shear_vrf_vals_(shear_vrf_vals)
         {
         }
 
@@ -90,7 +96,9 @@ namespace Opm
                  const std::vector<double>& c_vals_visc,
                  const std::vector<double>& visc_mult_vals,
                  const std::vector<double>& c_vals_ads,
-                 const std::vector<double>& ads_vals
+                 const std::vector<double>& ads_vals,
+                 const std::vector<double>& water_vel_vals,
+                 const std::vector<double>& shear_vrf_vals
                  )
         {
             c_max_ = c_max;
@@ -104,6 +112,8 @@ namespace Opm
             c_vals_ads_ = c_vals_ads;
             ads_vals_ = ads_vals;
             ads_index_ = ads_index;
+            water_vel_vals_ = water_vel_vals;
+            shear_vrf_vals_ = shear_vrf_vals;
         }
 
         void readFromDeck(const EclipseGridParser& gridparser)
@@ -134,6 +144,10 @@ namespace Opm
             c_vals_ads_ = plyads.local_concentration_;
             ads_vals_ = plyads.adsorbed_concentration_;
             
+            // We assume NTPVT=1
+            const PLYSHEAR& plyshear = gridparser.getPLYSHEAR();
+            water_vel_vals_ = plyshear.water_velocity_;
+            shear_vrf_vals_ = plyshear.vrf_;
         }
 
         double cMax() const;
@@ -150,6 +164,10 @@ namespace Opm
 
         int adsIndex() const;
         
+        std::vector<double> shearWaterVelocity() const;
+
+        std::vector<double> shearVrf() const;
+
         double viscMult(double c) const;
 
         double viscMultWithDer(double c, double* der) const;
@@ -259,6 +277,8 @@ namespace Opm
         std::vector<double> visc_mult_vals_;
         std::vector<double> c_vals_ads_;
         std::vector<double> ads_vals_;
+        std::vector<double> water_vel_vals_;
+        std::vector<double> shear_vrf_vals_;
         void simpleAdsorptionBoth(double c, double& c_ads,
                                   double& dc_ads_dc, bool if_with_der) const;
         void adsorptionBoth(double c, double cmax,