Gas/water relperm test model not running

[EdmundStephens]

[v2025_10] Attached models designed to test the hysteresis scanning curve for gas-water systems. The model has 6 cells initialised with different initial SGAS, a forced water displacement through each cell and reporting from the block summary keywords (BWKR, BGKR, etc.). These are using EHYSTR mode 4.

Test models

• Store2r-F.zip ... CO2STORE model with saturation data in GSF, WSF tables ... this model runs successfully although still shows incorrect scanning curve interpolation, see discussion in Gas/water hysteresis test model run does not time step #6383.
• Swg2r-F.zip ... WATER/GAS model with saturation data in SGWFN tables ... this model does not time step.

So both these models show things are not yet quite working.

Output from the model Store2r-F.data

Output from the model Swg2r-F.data

================ Starting main simulation loop ===============

Report step  0/150 at day 0/100000, date = 01-Jan-2000
Using Newton nonlinear solver.

Starting time step 0, stepsize 1e-06 days, at day 0/0.001, date = 01-Jan-2000
Restart file written for report step   0/150, date = 01-Jan-2000 00:00:00

Problem: [./opm/simulators/wells/StandardWellPrimaryVariables.cpp:755] Infinite primary variable after update from wellState, well: P1

Warning: A zero well potential is returned for output purposes.

Problem: [./opm/simulators/wells/StandardWellPrimaryVariables.cpp:755] Infinite primary variable after update from wellState, well: P1

Warning: Inner well iterations failed for well P1 Treat the well as unconverged.

Problem: [./opm/simulators/wells/StandardWellPrimaryVariables.cpp:755] Infinite primary variable after update from wellState, well: P1

Warning: Inner well iterations failed for well P1 Treat the well as unconverged.

Problem: [./opm/simulators/wells/StandardWellPrimaryVariables.cpp:755] Infinite primary variable after update from wellState, well: P1

Warning: Inner well iterations failed for well P1 Treat the well as unconverged.
flow: ./opm/simulators/flow/FlowProblem.hpp:996: void Opm::FlowProblem<TypeTag>::source(RateVector&, unsigned int, unsigned int) const [with TypeTag = Opm::Properties::TTag::FlowGasWaterProblem; RateVector = Opm::BlackOilRateVector<Opm::Properties::TTag::FlowGasWaterProblem>]: Assertion `isfinite(rate[eqIdx])' failed.
[EdLaptop:00326] *** Process received signal ***
[EdLaptop:00326] Signal: Aborted (6)
[EdLaptop:00326] Signal code:  (-6)
flow: ./opm/simulators/flow/FlowProblem.hpp:996: void Opm::FlowProblem<TypeTag>::source(RateVector&, unsigned int, unsigned int) const [with TypeTag = Opm::Properties::TTag::FlowGasWaterProblem; RateVector = Opm::BlackOilRateVector<Opm::Properties::TTag::FlowGasWaterProblem>]: Assertion `isfinite(rate[eqIdx])' failed.

[totto82]

Thanks for reporting the issue. For me Swg2r-F.data runs but not Store2r-F.data. Is it the PC data that is wrong in the case that runs? You need to enable the PC fix by setting item 14 to 1 in EHYSTR. I hope to change the defaults before the release, but have not done so yet. I will investigate why the failing case fails and report back.

[totto82]

With item 14 = 0 (default)

With item 14 = 1

[totto82]

The Store2r case runs fine without DISGASW. For relperm testing, that may be sufficient. I don't know what's causing the problem with dissolved CO2 yet, but I hope to figure it out.

[EdmundStephens]

Interesting

For me Swg2r-F.data runs but not Store2r-F.data

I get the other way round, including when EHYSTR item 14 is set to 1. There should be no difference in the data as such, coming from the same spreadsheet, but these are loaded differently:

• Swg2r ... data in a SGWFN table
• Store2r ... data WSF and GSF tables

The Swg2r case fails and the Store2r runs (including with the DISGASW keyword).

I tried the same data with SWFN and SGFN tables, where these have the Pc data tabulated on the SWFN and SGFN tables respectively:

• SfnW2r ... data in SWFN and SGFN tables with Pc data in SWFN table
• SfnG2r ... data in SWFN and SGFN tables with Pc data in SGFN table

SfnW2r-F.zip
SfnG2r-F.zip

These fail to run on my laptop, e.g.

Warning: Inner well iterations failed for well P1 Treat the well as unconverged.

Problem: [./opm/simulators/wells/StandardWellPrimaryVariables.cpp:755] Infinite primary variable after update from wellState, well: P1

Warning: Inner well iterations failed for well P1 Treat the well as unconverged.
flow: ./opm/simulators/flow/FlowProblem.hpp:996: void Opm::FlowProblem<TypeTag>::source(RateVector&, unsigned int, unsigned int) const [with TypeTag = Opm::Properties::TTag::FlowGasWaterProblem; RateVector = Opm::BlackOilRateVector<Opm::Properties::TTag::FlowGasWaterProblem>]: Assertion `isfinite(rate[eqIdx])' failed.
flow: ./opm/simulators/flow/FlowProblem.hpp:996: void Opm::FlowProblem<TypeTag>::source(RateVector&, unsigned int, unsigned int) const [with TypeTag = Opm::Properties::TTag::FlowGasWaterProblem; RateVector = Opm::BlackOilRateVector<Opm::Properties::TTag::FlowGasWaterProblem>]: Assertion `isfinite(rate[eqIdx])' failed.

Overall quite confusing, hopefully these additional cases clarify in the debugger where differences might be emerging.

Since you get different run outcomes, possibly there is a compiler library issue ... I have a MS Surface Laptop with the ARM architecture build .. some weird difference there perhaps. On my laptop the Store2r case runs fine, and also when I remove DISGASW and set EHYSTR item 14 = 1:

[totto82]

Thanks for the update. I have dug down and found that the issue is related to a cell that is first filled with gas, then, gradually, the gas disappears, and then the dissolved CO2 disappears as pure water fills the cell from the injectors. In this cell, the Newton iterations stagnate and do not update, so the simulator does not converge. I have not found any fundamental issues that explain why the simulator fails when the amount of dissolved CO2 goes to zero yet, but I will continue to look, as this is not an uncommon situation in CO2 storage. Note that I have not seen this issue in other cases where dissolved CO2 goes to zero.

Also, are the results still wrong for the cases that run when items 13 and 14 are enabled?

[totto82]

As I see it, the key difference between Store2r-F.data and Swg2r-F.data is that the former uses CO2STORE with DISGASW, whereas the latter uses PVTW and PVDG to define the properties. The first will allow CO2 to dissolve in water, while the second will not.

The convergence problems I observe are related to cells where the gas saturation goes below the critical saturation due to dissolution. With master this could lead to negative mobilities. I added a fix (see OPM/opm-common#5073). The PR fixes the issue for the relperm, but I still have convergence issues related to PC hysteresis. I.e. setting item 5 to KR helps. I don't see any obvious issues in the code, but I can understand that the test case is challenging for the simulator when dPc goes towards -10 bars. Note also that I need to set item 13 and item 14 for the simulator to converge on these cases. I have also fixed the output of BGPC and BWPC for gas-water models.

[EdmundStephens]

Update for version 2026.40

The gas-water hysteresis updates are almost correct, except if the starting water saturation is below 0.275 in this case. (This happens to be exactly at reduced saturation Swn = 0.15, i.e. Swc + 0.15*(1 - Swc - Sgr)). So if there is a step change in the results of cells with Swi = 0.274 and 0.275.

These cases were run in version 2026.04, single processor and arguments: --parsing-strictness="low" --enable-esmry=false --enable-opm-rst-file=true.

GAS-WATER test case:
Swg2r-F.zip

The results are correct for capillary and relperm result except for the Pc in cell [1,1,3] in pink ... the BWPC is anomalous.

CO2STORE test case:
Store2r-F.zip

This run stalls at around t=4042 days where timesteps drop to 1E-12 days. The results are correct for capillary and relperm result except for the Pc in cell [1,1,3] in pink ... the BWPC is anomalous. This change in behavious is again triggered below Swi 0.275 (although the case still fails to run so there may be some other issue involved)

[EdmundStephens]

Hi, Good news: I was looking into what might be the significance of Sw = 0.275, and it looks related to the TOLCRIT default 1E-6, in this case triggers at that saturation.

This addition remedies the situation in the gas-water case Swg2r and the CO2 store case Store2r. (The slight departure of the Swi=0.20 scanning curve from the original imbibition table data is a bit weird but seems to happen in all simulators).

TOLCRIT
	1.0E-24 /

Short term I think you only need a note in the release notes documentation to this effect. I think TOLCRIT was used mainly in the past to immobilise fluid in cells with low saturations to help solver performance, so sometimes TOLCRIT will be set to a higher value and a patch needs to be robust to those cases.

TOLCRIT should only really be affecting the relperm scaling "Critical fluid saturations are determined from the relative permeability tables, that is the last saturation in a relative permeability table where the relative permeability of a phase is set equal to zero". Therefore I would suggest to ensure capillary end point scaling is always effected from the ISWLPC (=ISWL by default) to ISWU which should not be subject to the TOLCRIT value. The relperm scaling should still be effected from the critical saturations (ISWCR etc.) subject to the TOLCRIT value.

Thinking of robustness, there should still be a sensible result if SWLPC or ISWLPC are set to a higher values, when Sw is between SWL and SWLPC or ISWL and ISWLPC. I notice also that there is only one assignment array for SWLPC and ISWLPC, whereas directional assignment is possible for SWL (SWLX etc.) and ISWL (ISWLX etc.). That seems weird, I can understand how there can be a directional critical saturation (SWCRX, ISWCRX etc.) but actually do not for the tabulated minimum. Perhaps in the saturation end points checking code such situations are eliminated from ever occurring ?

[totto82]

Thanks for following up. I just started looking into this and can confirm that TOLCRIT seems to fix the issue.

TOLCRIT should only really be affecting the relperm scaling "Critical fluid saturations are determined from the relative permeability tables, that is the last saturation in a relative permeability table where the relative permeability of a phase is set equal to zero". Therefore I would suggest to ensure capillary end point scaling is always effected from the ISWLPC (=ISWL by default) to ISWU which should not be subject to the TOLCRIT value. The relperm scaling should still be effected from the critical saturations (ISWCR etc.) subject to the TOLCRIT value.

I assume you mean "ISWLPC (=ISWL by default) to ISWL which should not be subject to the TOLCRIT value."

Do I understand you correctly that ISWL used as default for ISWLPC should not honor TOLCRIT even though native ISWL does?

[EdmundStephens]

I meant that the scaling interval should be from ISWLPC to ISWU, i.e. the imb Pc table saturations scaled to saturation range Sw_min=ISWLPC and Sw_max=ISWU. (And of course the drainage curve equivalently scaled from Sw_min=SWLPC to Sw_max=SWU).

Regarding SWLPC vs. SWL issue, there seems to be a distinction here that is getting deeper than I understand. Where is Dave Ponting when we need him ?! My best guess is that SWLPC was introduced so that old-style SWATINIT models could have a single table having SWL=0, Sw(Pc) fitted to a single saturation-height, then introducing SWLPC and ISWLPC would somewhat artificially enable the dynamic capillary pressure to be modelled. But in practice I don't understand the difference, since I think by definition the end scaling ranges should be:

• Capillary: Sw_min = SWL, ISWL and Sw_max = SWU, ISWU
• Relperm: Sw_min = SWCR, ISWCR and Sw_max = SWU, ISWU

... which is to say I don't really understand what SWLPC, ISWLPC are supposed to do other than add a degree of freedom that these might be used in some simplification strategies. Then again I don't get how SWL and SWU can be directional. Given these keywords are in place, probably best to amend the above to:

• Capillary: Sw_min = SWLPC, ISWLPC and Sw_max = SWU, ISWU

where these default in almost all situations to (non-directional) SWL, ISWL.

---

By the way there is still a snafu with the Store2r model when the Pc curve reaches its minimum the gas saturation is allowed to go to zero even though the gas relperm is zero. I'll post this in a moment in this thread.

[EdmundStephens]

I notice that with the TOLCRIT fix in place the result is perfect for the gas water case Swg2r, but in the Store2r model when the Pc curve reaches its minimum the gas saturation is allowed to go to zero even though the gas relperm is zero.

Could this be due to gas dissolving into the water therefore allowing reduction of gas saturation to zero ? In that case this behaviour is probably correct but it would be nice to have some control on this:

• DRSDT to control dissolution rate, similarly DRVDT to control evaporation rate
• ... or maybe bespoke keywords DRSWDT and DRVWDT for the aqueous phase
• WCONINJE argument 10 in a gas-water model to set the injectant Rsw (or Rvw for gas injector).

The gas-water case Swg2r is correct: cells reduce their Pc monotonically, maxxing around - 3.3 bar.

Swg2r-F v2.zip

The CO2store case Store2r ... while you are into this issue, can you have a look to make sure you are happy with this outcome. These effects are occurring at extreme negative imbibition with gas near residual saturation, which might be the case around injector wells.

Store2r-F v2.zip

[EdmundStephens]

Confirm I think this is related to gas dissolution. The CO2store model without DISGASW gives the same result as the gas-water case.

I would guess that the experience of a gas bubble where its volume is gradually dissolving in water should have capillary force remaining high and negative, gas relperm zero and water relperm probably trending towards 1. So the output as it is is probably good enough, especially if the imbibition table is extended to zero gas saturation. See below for the zipped model and plot output for this variant.

--    Sg              Krg i           Krwg i          Pcw i
      0.0000          0.0000E+00      1.00000000      -10.0000
      0.3000          0.0000E+00      0.50000000      -10.0000
      0.3005          8.0000E-10      0.49883980      -6.8147
      ...

So: this looks Ok as a result. Also I think this situation is going to be quite rare, the gas needs to be at extreme negative imbibition with gas near residual saturation, so only occurs at the base of a very extended plume (column exceeding the table Pc end point), or where water is being injected into a gas zone ... that can happen during wash treatments but a rather recherche modelling case where solids precipitation and plugging would be involved.

Store2t-F.zip