Fixed implicit transport solver interface.

There were significant lifetime issues, now handled by moving
objects inside the class. Work in progress.

opm/core/simulator/SimulatorIncompTwophase.cpp

@@ -343,29 +343,19 @@ namespace Opm
                                                            gravity, wells_manager.c_wells(), src, bcs));
                                                            */
         } else {
-            //Opm::ImplicitTransportDetails::NRReport  rpt;
+            if (rock_comp_props->isActive()) {
-            Opm::ImplicitTransportDetails::NRControl ctrl;
+                THROW("The implicit pressure solver cannot handle rock compressibility.");
-            ctrl.max_it = param.getDefault("max_it", 20);
+            }
-            ctrl.verbosity = param.getDefault("verbosity", 0);
-            ctrl.max_it_ls = param.getDefault("max_it_ls", 5);
-            const bool guess_old_solution = param.getDefault("guess_old_solution", false);
-            Opm::SimpleFluid2pWrappingProps fluid(props);
             std::vector<double> porevol;
-            //if (rock_comp->isActive()) {
+            computePorevolume(grid, props.porosity(), porevol);
-            //    computePorevolume(grid, props->porosity(), rock_comp, state.pressure(), porevol);
-            //} else {
-                computePorevolume(grid, props.porosity(), porevol);
-            //}
-            SinglePointUpwindTwoPhase<Opm::SimpleFluid2pWrappingProps>
-                    model(fluid, grid, porevol, gravity, guess_old_solution);
-            model.initGravityTrans(grid_, psolver_.getHalfTrans());
             tsolver_.reset(new Opm::TransportSolverTwophaseImplicit(
                                wells_manager,
                                *rock_comp_props,
-                               ctrl,
-                               model,
                                grid,
                                props,
+                               porevol,
+                               gravity,
+                               psolver_.getHalfTrans(),
                                param));
         }
 

opm/core/transport/TransportSolverTwophaseImplicit.cpp

@@ -37,18 +37,24 @@ namespace Opm
     TransportSolverTwophaseImplicit::TransportSolverTwophaseImplicit(
             const Opm::WellsManager& wells,
             const Opm::RockCompressibility& rock_comp,
-            const ImplicitTransportDetails::NRControl& ctrl,
-            SinglePointUpwindTwoPhase<Opm::SimpleFluid2pWrappingProps>& model,
             const UnstructuredGrid& grid,
             const Opm::IncompPropertiesInterface& props,
+            const std::vector<double>& porevol,
+            const double* gravity,
+            const std::vector<double>& half_trans,
             const parameter::ParameterGroup& param)
-        : tsolver_(model),
+        : fluid_(props),
+          model_(fluid_, grid, porevol, gravity, param.getDefault("guess_old_solution", false)),
+          tsolver_(model_),
           grid_(grid),
-          ctrl_(ctrl),
           props_(props),
           rock_comp_(rock_comp),
           wells_(wells)
     {
+        ctrl_.max_it = param.getDefault("max_it", 20);
+        ctrl_.verbosity = param.getDefault("verbosity", 0);
+        ctrl_.max_it_ls = param.getDefault("max_it_ls", 5);
+        model_.initGravityTrans(grid_, half_trans);
         tsrc_ = create_transport_source(2, 2);
     }
 
@@ -75,11 +81,16 @@ namespace Opm
         double ssrc[] = { 1.0, 0.0 };
         double dummy[] = { 0.0, 0.0 };
         clear_transport_source(tsrc_);
+        const int num_phases = 2;
         for (int cell = 0; cell < grid_.number_of_cells; ++cell) {
+            int success = 1;
             if (source[cell] > 0.0) {
-                append_transport_source(cell, 2, 0, source[cell], ssrc, dummy, tsrc_);
+                success = append_transport_source(cell, num_phases, state.pressure()[cell], source[cell], ssrc, dummy, tsrc_);
             } else if (source[cell] < 0.0) {
-                append_transport_source(cell, 2, 0, source[cell], dummy, dummy, tsrc_);
+                success = append_transport_source(cell, num_phases, state.pressure()[cell], source[cell], dummy, dummy, tsrc_);
+            }
+            if (!success) {
+                THROW("Failed building TransportSource struct.");
             }
         }
         Opm::ImplicitTransportDetails::NRReport  rpt;

opm/core/transport/TransportSolverTwophaseImplicit.hpp

@@ -62,10 +62,11 @@ namespace Opm
         /// \param[in] maxit     Maximum number of non-linear iterations used.
         TransportSolverTwophaseImplicit(const Opm::WellsManager& wells,
                                         const Opm::RockCompressibility& rock_comp,
-                                        const ImplicitTransportDetails::NRControl& ctrl,
-                                        SinglePointUpwindTwoPhase<Opm::SimpleFluid2pWrappingProps>& model,
                                         const UnstructuredGrid& grid,
                                         const Opm::IncompPropertiesInterface& props,
+                                        const std::vector<double>& porevol,
+                                        const double* gravity,
+                                        const std::vector<double>& half_trans,
                                         const parameter::ParameterGroup& param);
 
         virtual ~TransportSolverTwophaseImplicit();
@@ -109,9 +110,11 @@ namespace Opm
 
         // Data members.
         Opm::ImplicitTransportLinAlgSupport::CSRMatrixUmfpackSolver linsolver_;
+        Opm::SimpleFluid2pWrappingProps fluid_;
+        SinglePointUpwindTwoPhase<Opm::SimpleFluid2pWrappingProps> model_;
         TransportSolver tsolver_;
         const UnstructuredGrid& grid_;
-        const Opm::ImplicitTransportDetails::NRControl& ctrl_;
+        Opm::ImplicitTransportDetails::NRControl ctrl_;
         const Opm::IncompPropertiesInterface& props_;
         const Opm::RockCompressibility& rock_comp_;
         const Opm::WellsManager& wells_;