BlackoilModel: move to has-a instead of is-a modelling of penalty cards

CMakeLists_files.cmake

@@ -88,6 +88,7 @@ list (APPEND MAIN_SOURCE_FILES
   opm/simulators/flow/ActionHandler.cpp
   opm/simulators/flow/Banners.cpp
   opm/simulators/flow/BlackoilModelParameters.cpp
+  opm/simulators/flow/BlackoilModelConvergenceMonitor.cpp
   opm/simulators/flow/CollectDataOnIORank.cpp
   opm/simulators/flow/ConvergenceOutputConfiguration.cpp
   opm/simulators/flow/EclGenericWriter.cpp
@@ -792,6 +793,7 @@ list (APPEND PUBLIC_HEADER_FILES
   opm/simulators/flow/Banners.hpp
   opm/simulators/flow/BaseAquiferModel.hpp
   opm/simulators/flow/BlackoilModel.hpp
+  opm/simulators/flow/BlackoilModelConvergenceMonitor.hpp
   opm/simulators/flow/BlackoilModelNldd.hpp
   opm/simulators/flow/BlackoilModelParameters.hpp
   opm/simulators/flow/BlackoilModelProperties.hpp

opm/simulators/flow/BlackoilModel.hpp

@@ -36,6 +36,7 @@
 #include <opm/simulators/aquifers/AquiferGridUtils.hpp>
 #include <opm/simulators/aquifers/BlackoilAquiferModel.hpp>
 
+#include <opm/simulators/flow/BlackoilModelConvergenceMonitor.hpp>
 #include <opm/simulators/flow/BlackoilModelNldd.hpp>
 #include <opm/simulators/flow/BlackoilModelParameters.hpp>
 #include <opm/simulators/flow/BlackoilModelProperties.hpp>
@@ -160,6 +161,7 @@ namespace Opm {
         , terminal_output_ (terminal_output)
         , current_relaxation_(1.0)
         , dx_old_(simulator_.model().numGridDof())
+        , conv_monitor_(param_.monitor_params_)
         {
             // compute global sum of number of cells
             global_nc_ = detail::countGlobalCells(grid_);
@@ -298,7 +300,9 @@ namespace Opm {
                     OPM_THROW_NOLOG(NumericalProblem, "Too large residual found!");
                 } else if (severity == ConvergenceReport::Severity::ConvergenceMonitorFailure) {
                     failureReport_ += report;
-                    OPM_THROW_PROBLEM(ConvergenceMonitorFailure, "Total penalty count exceeded cut-off-limit of " + std::to_string(param_.convergence_monitoring_cutoff_));
+                    OPM_THROW_PROBLEM(ConvergenceMonitorFailure,
+                                      "Total penalty count exceeded cut-off-limit of " +
+                                          std::to_string(param_.monitor_params_.cutoff_));
                 }
             }
             report.update_time += perfTimer.stop();
@@ -322,9 +326,7 @@ namespace Opm {
                 // For each iteration we store in a vector the norms of the residual of
                 // the mass balance for each active phase, the well flux and the well equations.
                 residual_norms_history_.clear();
-                total_penaltyCard_.reset();
-                prev_above_tolerance_ = 0;
-                prev_distance_ = std::numeric_limits<double>::infinity();
+                conv_monitor_.reset();
                 current_relaxation_ = 1.0;
                 dx_old_ = 0.0;
                 convergence_reports_.push_back({timer.reportStepNum(), timer.currentStepNum(), {}});
@@ -1040,51 +1042,6 @@ namespace Opm {
             return report;
         }
 
-    void checkCardPenalty(ConvergenceReport& report, int iteration)
-    {
-
-        const auto& current_metrics = report.reservoirConvergence();
-        auto distances = std::vector<double>(current_metrics.size(), 0.0);
-        int current_above_tolerance = 0;
-
-        for (size_t i = 0; i < current_metrics.size(); ++i) {
-            distances[i] = std::max(std::log10(current_metrics[i].value()/current_metrics[i].tolerance()), 0.0);
-                // Count number of metrics above tolerance
-                if (current_metrics[i].value() > current_metrics[i].tolerance()) {
-                    current_above_tolerance++;
-                }
-            }
-
-        // use L1 norm of the distances vector
-        double current_distance = std::accumulate(distances.begin(), distances.end(), 0.0);
-
-        if (iteration > 0) {
-            // Add penalty if number of metrics above tolerance has increased
-            if (current_above_tolerance > prev_above_tolerance_) {
-                report.addNonConvergedPenalty();
-            }
-
-            if (current_distance > param_.convergence_monitoring_decay_factor_ * prev_distance_) {
-                report.addDistanceDecayPenalty();
-            }
-        }
-
-        prev_distance_ = current_distance;
-        prev_above_tolerance_ = current_above_tolerance;
-
-        if (report.wellFailures().size() > 0) {
-            report.addLargeWellResidualsPenalty();
-        }
-
-        total_penaltyCard_ += report.getPenaltyCard();
-
-        if (param_.convergence_monitoring_ && (total_penaltyCard_.total() > param_.convergence_monitoring_cutoff_)) {
-            report.setReservoirFailed({ConvergenceReport::ReservoirFailure::Type::ConvergenceMonitorFailure,
-                                       ConvergenceReport::Severity::ConvergenceMonitorFailure,
-                                       -1}); // -1 indicates it's not specific to any component
-        }
-    }
-
         /// Compute convergence based on total mass balance (tol_mb) and maximum
         /// residual mass balance (tol_cnv).
         /// \param[in]   timer       simulation timer
@@ -1107,7 +1064,7 @@ namespace Opm {
                 report += wellModel().getWellConvergence(B_avg, /*checkWellGroupControls*/report.converged());
             }
 
-            checkCardPenalty(report, iteration);
+            conv_monitor_.checkPenaltyCard(report, iteration);
 
             return report;
         }
@@ -1374,10 +1331,8 @@ namespace Opm {
         Scalar dsMax() const { return param_.ds_max_; }
         Scalar drMaxRel() const { return param_.dr_max_rel_; }
         Scalar maxResidualAllowed() const { return param_.max_residual_allowed_; }
-        double linear_solve_setup_time_;
-        ConvergenceReport::PenaltyCard total_penaltyCard_;
-        double prev_distance_ = std::numeric_limits<double>::infinity();
-        int prev_above_tolerance_ = 0;
+        double linear_solve_setup_time_{};
+        BlackoilModelConvergenceMonitor<Scalar> conv_monitor_;
         std::vector<bool> wasSwitched_;
     };
 

opm/simulators/flow/BlackoilModelConvergenceMonitor.cpp

@@ -0,0 +1,104 @@
+/*
+  Copyright 2013, 2015 SINTEF ICT, Applied Mathematics.
+  Copyright 2014, 2015 Dr. Blatt - HPC-Simulation-Software & Services
+  Copyright 2014, 2015 Statoil ASA.
+  Copyright 2015 NTNU
+  Copyright 2015, 2016, 2017 IRIS AS
+
+  This file is part of the Open Porous Media project (OPM).
+
+  OPM is free software: you can redistribute it and/or modify
+  it under the terms of the GNU General Public License as published by
+  the Free Software Foundation, either version 3 of the License, or
+  (at your option) any later version.
+
+  OPM is distributed in the hope that it will be useful,
+  but WITHOUT ANY WARRANTY; without even the implied warranty of
+  MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
+  GNU General Public License for more details.
+
+  You should have received a copy of the GNU General Public License
+  along with OPM.  If not, see <http://www.gnu.org/licenses/>.
+*/
+#include <config.h>
+#include <opm/simulators/flow/BlackoilModelConvergenceMonitor.hpp>
+
+#include <algorithm>
+#include <cmath>
+#include <limits>
+
+namespace Opm {
+
+template<class Scalar>
+BlackoilModelConvergenceMonitor<Scalar>::
+BlackoilModelConvergenceMonitor(const MonitorParams& param)
+    : param_(param)
+    , prev_distance_(std::numeric_limits<double>::infinity())
+    , prev_above_tolerance_(0)
+{
+}
+
+template<class Scalar>
+void BlackoilModelConvergenceMonitor<Scalar>::
+checkPenaltyCard(ConvergenceReport& report, int iteration)
+{
+    const auto& current_metrics = report.reservoirConvergence();
+    auto distances = std::vector<double>(current_metrics.size(), 0.0);
+    int current_above_tolerance = 0;
+
+    for (size_t i = 0; i < current_metrics.size(); ++i) {
+        distances[i] = std::max(std::log10(current_metrics[i].value() /
+                                           current_metrics[i].tolerance()), 0.0);
+            // Count number of metrics above tolerance
+            if (current_metrics[i].value() > current_metrics[i].tolerance()) {
+                current_above_tolerance++;
+            }
+        }
+
+    // use L1 norm of the distances vector
+    double current_distance = std::accumulate(distances.begin(), distances.end(), 0.0);
+
+    if (iteration > 0) {
+        // Add penalty if number of metrics above tolerance has increased
+        if (current_above_tolerance > prev_above_tolerance_) {
+            report.addNonConvergedPenalty();
+        }
+
+        if (current_distance > param_.decay_factor_ * prev_distance_) {
+            report.addDistanceDecayPenalty();
+        }
+    }
+
+    prev_distance_ = current_distance;
+    prev_above_tolerance_ = current_above_tolerance;
+
+    if (report.wellFailures().size() > 0) {
+        report.addLargeWellResidualsPenalty();
+    }
+
+    total_penaltyCard_ += report.getPenaltyCard();
+
+    if (param_.enabled_ && (total_penaltyCard_.total() > param_.cutoff_)) {
+        report.setReservoirFailed(
+            {ConvergenceReport::ReservoirFailure::Type::ConvergenceMonitorFailure,
+             ConvergenceReport::Severity::ConvergenceMonitorFailure,
+             -1}); // -1 indicates it's not specific to any component
+    }
+}
+
+template<class Scalar>
+void BlackoilModelConvergenceMonitor<Scalar>::
+reset()
+{
+    total_penaltyCard_.reset();
+    prev_above_tolerance_ = 0;
+    prev_distance_ = std::numeric_limits<double>::infinity();
+}
+
+template class BlackoilModelConvergenceMonitor<double>;
+
+#if FLOW_INSTANTIATE_FLOAT
+template class BlackoilModelConvergenceMonitor<float>;
+#endif
+
+} // namespace Opm

opm/simulators/flow/BlackoilModelConvergenceMonitor.hpp

@@ -0,0 +1,53 @@
+/*
+  Copyright 2013, 2015 SINTEF ICT, Applied Mathematics.
+  Copyright 2014, 2015 Dr. Blatt - HPC-Simulation-Software & Services
+  Copyright 2014, 2015 Statoil ASA.
+  Copyright 2015 NTNU
+  Copyright 2015, 2016, 2017 IRIS AS
+
+  This file is part of the Open Porous Media project (OPM).
+
+  OPM is free software: you can redistribute it and/or modify
+  it under the terms of the GNU General Public License as published by
+  the Free Software Foundation, either version 3 of the License, or
+  (at your option) any later version.
+
+  OPM is distributed in the hope that it will be useful,
+  but WITHOUT ANY WARRANTY; without even the implied warranty of
+  MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
+  GNU General Public License for more details.
+
+  You should have received a copy of the GNU General Public License
+  along with OPM.  If not, see <http://www.gnu.org/licenses/>.
+*/
+
+#ifndef OPM_BLACKOILMODEL_CONV_MONITOR_HEADER_INCLUDED
+#define OPM_BLACKOILMODEL_CONV_MONITOR_HEADER_INCLUDED
+
+#include <opm/simulators/flow/BlackoilModelParameters.hpp>
+#include <opm/simulators/timestepping/ConvergenceReport.hpp>
+
+namespace Opm {
+
+/// Implementation of penalty cards for three-phase black oil.
+template <class Scalar>
+class BlackoilModelConvergenceMonitor
+{
+public:
+    using MonitorParams = typename BlackoilModelParameters<Scalar>::ConvergenceMonitorParams;
+    explicit BlackoilModelConvergenceMonitor(const MonitorParams& param);
+
+    void checkPenaltyCard(ConvergenceReport& report, int iteration);
+
+    void reset();
+
+private:
+    const MonitorParams& param_;
+    ConvergenceReport::PenaltyCard total_penaltyCard_;
+    double prev_distance_;
+    int prev_above_tolerance_;
+};
+
+} // namespace Opm
+
+#endif

opm/simulators/flow/BlackoilModelParameters.cpp

@@ -100,9 +100,9 @@ BlackoilModelParameters<Scalar>::BlackoilModelParameters()
     local_domain_ordering_ = domainOrderingMeasureFromString(Parameters::Get<Parameters::LocalDomainsOrderingMeasure>());
     write_partitions_ = Parameters::Get<Parameters::DebugEmitCellPartition>();
 
-    convergence_monitoring_ = Parameters::Get<Parameters::ConvergenceMonitoring>();
-    convergence_monitoring_cutoff_ = Parameters::Get<Parameters::ConvergenceMonitoringCutOff>();
-    convergence_monitoring_decay_factor_ = Parameters::Get<Parameters::ConvergenceMonitoringDecayFactor<Scalar>>();
+    monitor_params_.enabled_ = Parameters::Get<Parameters::ConvergenceMonitoring>();
+    monitor_params_.cutoff_ = Parameters::Get<Parameters::ConvergenceMonitoringCutOff>();
+    monitor_params_.decay_factor_ = Parameters::Get<Parameters::ConvergenceMonitoringDecayFactor<Scalar>>();
 
     nupcol_group_rate_tolerance_ = Parameters::Get<Parameters::NupcolGroupRateTolerance<Scalar>>();
 }

opm/simulators/flow/BlackoilModelParameters.hpp

@@ -315,12 +315,18 @@ public:
 
     bool write_partitions_{false};
 
-    /// Whether to enable convergence monitoring
-    bool convergence_monitoring_;
-    /// Cut-off limit for convergence monitoring
-    int convergence_monitoring_cutoff_;
-    /// Decay factor used in convergence monitoring
-    Scalar convergence_monitoring_decay_factor_;
+    /// Struct holding convergence monitor params
+    struct ConvergenceMonitorParams
+    {
+        /// Whether to enable convergence monitoring
+        bool enabled_;
+        /// Cut-off limit for convergence monitoring
+        int cutoff_;
+        /// Decay factor used in convergence monitoring
+        Scalar decay_factor_;
+    };
+
+    ConvergenceMonitorParams monitor_params_; //!< Convergence monitoring parameters
 
     // Relative tolerance of group rates (VREP, REIN)
     // If violated the nupcol wellstate is updated