From 745d8b5a97c1fcd40b34fd3808cf7efb5f6b51cf Mon Sep 17 00:00:00 2001
From: Svenn Tveit <svtv@norceresearch.no>
Date: Tue, 13 Jun 2023 12:25:47 +0200
Subject: [PATCH] Accommodate  for PPCWMAX in SWATINIT equilibration

---
 ebos/equil/initstateequil.hh                  |  4 +--
 ebos/equil/initstateequil_impl.hh             | 34 ++++++++++++++-----
 .../utils/UnsupportedFlowKeywords.cpp         |  1 -
 3 files changed, 27 insertions(+), 12 deletions(-)

diff --git a/ebos/equil/initstateequil.hh b/ebos/equil/initstateequil.hh
index d0df7e10d..b247e11d0 100644
--- a/ebos/equil/initstateequil.hh
+++ b/ebos/equil/initstateequil.hh
@@ -544,7 +544,7 @@ private:
     /// \param[in] pcow O/W capillary pressure value (Po - Pw).
     ///
     /// \return Water saturation value.
-    double applySwatInit(const double pcow);
+    std::tuple<double, bool> applySwatInit(const double pcow);
 
     /// Derive water saturation from SWATINIT data.
     ///
@@ -558,7 +558,7 @@ private:
     ///
     /// \return Water saturation value.  Input value, possibly mollified by
     ///    current set of material laws.
-    double applySwatInit(const double pc, const double sw);
+    std::tuple<double, bool> applySwatInit(const double pc, const double sw);
 
     /// Invoke material law container's capillary pressure calculator on
     /// current fluid state.
diff --git a/ebos/equil/initstateequil_impl.hh b/ebos/equil/initstateequil_impl.hh
index b2f3cf6bd..1ab8a5ef5 100644
--- a/ebos/equil/initstateequil_impl.hh
+++ b/ebos/equil/initstateequil_impl.hh
@@ -694,9 +694,16 @@ void PhaseSaturations<MaterialLawManager, FluidSystem, Region, CellID>::deriveWa
         // cell.
         const auto pcow = this->press_.oil - this->press_.water;
 
-        sw = this->swatInit_.empty()
-            ? this->invertCapPress(pcow, this->waterPos(), isIncr)
-            : this->applySwatInit(pcow);
+        if (this->swatInit_.empty()) {
+            sw = this->invertCapPress(pcow, this->waterPos(), isIncr);
+        }
+        else {
+            auto [swout, newSwatInit] = this->applySwatInit(pcow);
+            if (newSwatInit == true)
+                sw = this->invertCapPress(pcow, this->waterPos(), isIncr);
+            else
+                sw = swout;
+        }
     }
 }
 
@@ -718,7 +725,13 @@ fixUnphysicalTransition()
         // Re-scale Pc to reflect imposed sw for vanishing oil phase.  This
         // seems consistent with ECLIPSE, but fails to honour SWATINIT in
         // case of non-trivial gas/oil capillary pressure.
-        sw = this->applySwatInit(pcgw, sw);
+        auto [swout, newSwatInit] = this->applySwatInit(pcgw, sw);
+        if (newSwatInit){
+            const auto isIncr = false; // dPcow/dSw <= 0 for all Sw.
+            sw = this->invertCapPress(pcgw, this->waterPos(), isIncr);
+        }
+        else
+            sw = swout;
     }
 
     sw = satFromSumOfPcs<FluidSystem>
@@ -846,18 +859,21 @@ accountForScaledSaturations()
 }
 
 template <class MaterialLawManager, class FluidSystem, class Region, typename CellID>
-double PhaseSaturations<MaterialLawManager, FluidSystem, Region, CellID>::
+std::tuple<double, bool>
+PhaseSaturations<MaterialLawManager, FluidSystem, Region, CellID>::
 applySwatInit(const double pcow)
 {
-    return this->applySwatInit(pcow, this->swatInit_[this->evalPt_.position->cell]);
+    auto [swout, newSwatInit] = this->applySwatInit(pcow, this->swatInit_[this->evalPt_.position->cell]);
+    return {swout, newSwatInit};
 }
 
 template <class MaterialLawManager, class FluidSystem, class Region, typename CellID>
-double PhaseSaturations<MaterialLawManager, FluidSystem, Region, CellID>::
+std::tuple<double, bool>
+PhaseSaturations<MaterialLawManager, FluidSystem, Region, CellID>::
 applySwatInit(const double pcow, const double sw)
 {
-    return this->matLawMgr_
-        .applySwatinit(this->evalPt_.position->cell, pcow, sw);
+    auto [swout, newSwatInit] =this->matLawMgr_.applySwatinit(this->evalPt_.position->cell, pcow, sw);
+    return {swout, newSwatInit};
 }
 
 template <class MaterialLawManager, class FluidSystem, class Region, typename CellID>
diff --git a/opm/simulators/utils/UnsupportedFlowKeywords.cpp b/opm/simulators/utils/UnsupportedFlowKeywords.cpp
index 82c02edf6..cf6c619a4 100644
--- a/opm/simulators/utils/UnsupportedFlowKeywords.cpp
+++ b/opm/simulators/utils/UnsupportedFlowKeywords.cpp
@@ -490,7 +490,6 @@ const KeywordValidation::UnsupportedKeywords& unsupportedKeywords()
         {"PLYATEMP", {true, std::nullopt}},
         {"PLYCAMAX", {true, std::nullopt}},
         {"PLYDHFLF", {true, std::nullopt}},
-        {"PPCWMAX", {true, std::nullopt}},
         {"PRORDER", {true, std::nullopt}},
         {"PRVD", {true, std::nullopt}},
         {"PVTGWO", {true, std::nullopt}},