Handle degenerate PTFlash phase splits - #5365
Conversation
GitPaean
added
the
manual:irrelevant
|
jekins build this opm-simulators=7403 please |
1 similar comment
|
jekins build this opm-simulators=7403 please |
|
jenkins build this opm-simulators=7403 please |
372b3d6 to
12c3560
Compare
|
jenkins build this opm-simulators=7403 please |
There was a problem hiding this comment.
🟡 Changes recommended
The newly added singular-Jacobian fallback remains untested.
Get a fresh assessment by requesting another Copilot review.
Pull request overview
Prevents PTFlash Newton iterations from accepting degenerate, non-physical two-phase solutions and enables hybrid SSI recovery.
Changes:
- Detects trivial Newton roots and validates converged phase bounds.
- Falls back to SSI after Newton or singular-Jacobian failures.
- Adds regression coverage for the recorded CO₂ state.
File summaries
| File | Description |
|---|---|
opm/material/constraintsolvers/PTFlash.hpp |
Adds guards, validation, and fallback handling. |
tests/material/test_ptflash_trivial_root.cpp |
Tests the recorded failure and recovery. |
CMakeLists_files.cmake |
Registers the new test. |
Review details
- Files reviewed: 3/3 changed files
- Comments generated: 1
- Review effort level: Balanced
💡 Add a code-review agent skill or configure MCP servers for context-aware, tailored reviews. Learn more in the docs.
There was a problem hiding this comment.
🔵 Needs a closer look
The trivial-root check incorrectly rejects otherwise valid mixtures containing an absent component.
Review details
Suppressed comments (1)
Previously missed (1) — in code that hasn't changed since the last review.
opm/material/constraintsolvers/PTFlash.hpp:827
- An exactly absent component has
x_i == y_i == 0, sonewton_state.K(i)was formed as0/0above. This condition then classifies every such otherwise nontrivial mixture as the trivial root and makes the explicit Newton method fail before assembling its valid equations. Ignore components absent from both phases when accumulating the K-based measure (while retaining the rejection for invalid ratios of active components).
- Files reviewed: 4/4 changed files
- Comments generated: 0 new
- Review effort level: Balanced
There was a problem hiding this comment.
🔵 Needs a closer look
The implementation is coherent and well tested, but changes sensitive numerical-solver convergence behavior requiring final human validation.
Review details
- Files reviewed: 4/4 changed files
- Comments generated: 0 new
- Review effort level: Balanced
|
jenkins build this opm-simulators=7403 please |
1 similar comment
|
jenkins build this opm-simulators=7403 please |
|
jenkins build this opm-simulators=7403 please |
1 similar comment
|
jenkins build this opm-simulators=7403 please |
|
jenkins build this opm-simulators=7403 please |
There was a problem hiding this comment.
🔵 Needs a closer look
The thermodynamic classification and derivative-path changes warrant final human validation despite strong regression coverage.
Review details
- Files reviewed: 4/4 changed files
- Comments generated: 0 new
- Review effort level: Balanced
ee73b66 to
2aa50df
Compare
|
jenkins build this opm-simulators=7403 please https://ci.opm-project.org/job/opm-common-PR-builder/10382/console |
|
jenkins build this please https://ci.opm-project.org/job/opm-common-PR-builder/10385/console The failure does not look like related to this PR. |
|
jenkins build this please |
8e258a8 to
44629b5
Compare
|
jenkins build this please |
44629b5 to
8c5ec60
Compare
Guard the two-phase Newton solve against coincident phases and non-physical composition roots. In the hybrid method, return to successive substitution when Newton fails or encounters a singular FieldMatrix Jacobian. Reassess rejected and coincident warm-started splits through stability analysis. Retry one unstable split from the resulting equilibrium ratios, and use the single-phase derivative path after reclassification. Preserve finite negative-flash fractions as phase-classification results. Exclude absent components from stability and two-phase SSI fugacity ratios. Use the fluid state's Wilson estimates for recovery, and report recoverable numerical failures without immediate error logging. Add regression coverage for the recorded near-pure CO2 false root, phase classification and derivatives, negative flashes, stability-based retry, absent components, composition bounds, and singular-Jacobian fallback.
Centralize the split operation and its two recoverable exception handlers in one local helper. Use the returned success flag to drive reclassification and stability recovery.
|
jenkins build this please |
There was a problem hiding this comment.
🔵 Needs a closer look
It makes broad numerical recovery and phase-classification changes in core thermodynamic logic that warrant final human validation.
Review details
- Files reviewed: 4/4 changed files
- Comments generated: 0 new
- Review effort level: Balanced
| * \exp\left[ 5.3727 (1 + \omega_i) \left( 1 - \frac{T_{c,i}}{T} \right) \right]. \f] | ||
| */ | ||
| template <class FlashFluidState> | ||
| static typename FlashFluidState::ValueType wilsonK_(const FlashFluidState& fluid_state, int compIdx) |
There was a problem hiding this comment.
The original function is an unused function before this PR, there is a function ValueType wilsonK_() in FluidStateCompositionModules.hpp, and in this PR, that wilsonK_() is used.
The function trivialSolutionMeasure_() is basically a new function with totally different purpose.
|
This PR was created to solve a corner situation with The PR OPM/opm-simulators#7403 was due to a failed solve in an equilibration case, which could not solve the flash below the convergence tolerance 1.e-8 . (stalling between 2-5e-8). The PR OPM/opm-simulators#7403 makes it converged without relaxing the tolerance. Furthermore, since SSI+Newton will fall back to SSI, I think it is good to use SSI+Newton as the default flash solver. I am marking the PR as ready for review and inviting reviewers. |
bska
left a comment
bska
left a comment
There was a problem hiding this comment.
There's nothing here that jumps out to me so I'll merge into master. After all, this feature is still in development so we can afford to redo certain aspects if experience suggests that revised approaches are needed.
The two-phase Newton equations admit a degenerate
x = yroot where thecomponent balance no longer constrains the liquid fraction. For the affected
near-pure CO2 state, Newton reported convergence with an enormous negative
liquid fraction and a negative mole fraction, which led to non-finite mobility.
Main changes
singular Jacobian and after a two-phase solve completes.
ssi+newton, fall back to SSI after recoverable Newton or singularFieldMatrixfailures, then reassess a rejected split through stabilityanalysis and retry it once from the resulting non-trivial equilibrium ratios.
single-phase, preserve the liquid/vapour direction, and use the single-phase
derivative path.
measures, validate final compositions, and add focused regression coverage
for each recovery path.
Flash recovery flow
flowchart TB A(["Flash input<br/>z, p, T, EOS, K, L"]) --> B{"Cold start?<br/>L <= 0 or L >= 1"} B -- yes --> C["Phase-stability test"] B -- no --> D["Solve Rachford-Rice equation"] C --> C1{"Stable?"} C1 -- yes --> SP["Single-phase result<br/>x = y = z<br/>phase labelled using Li's method"] C1 -- no --> D D -- failure --> R["Reset K using<br/>Wilson estimates"] D -- success --> E["Solve phase compositions<br/>SSI, Newton, or SSI + Newton"] E --> F{"Composition-solver outcome?"} F -- "Newton failure in SSI + Newton" --> G["Restart with plain SSI"] F -- "failure in Newton or SSI" --> R F -- success --> H{"Finite, physical compositions?"} G -- success --> H G -- failure --> R H -- no --> R H -- yes --> I{"Negative flash?<br/>L <= 0 or L >= 1"} I -- yes --> NF["Single phase<br/>clamp L to 0 or 1<br/>and preserve direction"] I -- no --> J{"Liquid and vapour<br/>phases coincide?"} J -- no --> TP["Valid two-phase result"] J -- "yes, warm start" --> R J -- "yes, cold start" --> ERR["Report numerical failure"] R --> K["Repeat phase-stability test"] K --> L{"Stable?"} L -- yes --> SP L -- no --> M["Retry the split once"] M --> N{"Valid retry?"} N -- "distinct phases" --> TP N -- "negative flash" --> NF N -- "failure or coincident phases" --> ERR SP --> SD(["Single-phase derivative path"]) NF --> SD TP --> TD(["Two-phase derivative path"])Recorded failure scenario
The regression was recorded in
co2_ptflash_ecfvfor the three-componentfluid system
[CO2, C1, C10]at:14981770.14236617 Pa;423.25 K;[0.9900000005514261, 0.009000003165640235, 0.0009999962829336465];0.2252754135811573; and[6.260223748149762, 12.126317087860864, 0.0071260520516577735].Because the incoming liquid fraction was strictly between zero and one, the
warm-start path skipped stability analysis. Before this fix, Newton reported a
small equation residual after reaching approximately
L = -2.5e16and a molefraction of approximately
-0.78. The result subsequently produced anon-finite mobility.
Plain successive substitution avoided the unbounded Newton result, but
converged to coincident normalized phases,
x = y = z, withL = 5.82077e-11while still reporting the state as two-phase. Reconstructingthe two-phase derivatives at that degenerate state produced liquid-fraction
sensitivities with magnitudes from about
7.1e8to1.9e11. In contrast, acold start ran stability analysis and classified the same state as
single-phase vapour with
L = 0.With this fix, Newton rejects a coincident starting state before assembling its
singular Jacobian, and
ssi+newtonfalls back to SSI when Newton or a singularFieldMatrixsolve fails. A rejected split is reassessed from the fluidstate's Wilson estimates and may be retried once with the non-trivial
equilibrium ratios returned by stability analysis. A warm-started split that
converges to coincident phases is reassessed in the same way. The recorded
state is consequently classified like the cold start, and uses the
single-phase derivative path with
dL = 0.Finite negative-flash roots outside
0 <= L <= 1are phase-classificationresults, not composition-solver failures. This matters for well flashes. For a
C1/n-C10 mixture at standard conditions (
1.01325 bar,288.71 K) withz(C1) = 0.001, SSI returns aboutL = 1.005; this is now classified assingle-phase liquid (
L = 1) instead of throwing. Likewise, forz = [0.5, 0.3, 0.2]at350 K, representative warm starts at180 bar(
L = 1.275, distinct phases) and220 bar(L = 1.9e7, coincident phases)both agree with the cold-start single-phase liquid result.
Additional implementation details:
L <= 0andL >= 1results as single-phase and reservesfailure for a non-finite liquid fraction or invalid phase compositions;
stability result;
NumericalProblemwithout immediate error logging for failures that anenclosing flash layer may recover;
ratios, and requires at least two active components before declaring
normalized phases coincident;
The coincident-phase measure is
sum(log(y_i / x_i)^2) < 1e-5over components with a defined positive ratio.A separate final-state check validates phase compositions using a bound slack
independent of the fugacity-residual tolerance. Computing the measure from the
normalized phase compositions also detects the uniformly scaled
Kvectorthat successive substitution can retain.
The regression tests cover the recorded false root, phase reclassification,
single-phase derivatives, cold and warm negative-flash results, well surface
flashes, stability-based retry, absent components, composition bounds, and the
exactly singular Jacobian fallback. The stored-reference
1D_COMPtests forboth
flow_compandflowexp_comp, andSIMPLE_COMP_SSHIFT, remain intact.The derivative change is intentionally limited to states reclassified by this
flash logic. General conditioning checks for a distinct, near-critical
two-phase derivative Jacobian are outside this PR.
This fixes the regression exposed by OPM/opm-simulators#7403.