Two-point stress approximation (TPSA) geomechanics #6650

svenn-t · 2025-12-02T14:38:01Z

This PR introduces the implementation of TPSA geomechanics solver with coupling to Flow. Most important features are:

TPSA-Flow coupling done at the nonlinear solver step level:
- Lagged coupling -> TPSA lags behind one time step
- Fixed-stress -> Iterate until Flow and TPSA converge
Generic TPSA solver using Newton with associated linearizer
Linear elasticity equations implemented
Specialized version of ISTL linear solver for TPSA
Output to restart files and VTK

There are blackoil, gas-water dissolution (mostly for CO2-/H2STORE) and one-phase simulators implemented, but in principle any other simulator can be generated as well.

Reference to implementation paper can be found here.

Putting this PR in draft mode since it is rather large, but any feedback/comments/discussion on the implementation is welcome at this point!

Companion PR: OPM/opm-common#4846

flow/flow_blackoil_tpsa.cpp

bska · 2025-12-03T10:39:38Z

I'll just make one general comment here and then I'll be quiet. Adding 6,600 lines of new code in a single PR, especially since there's hardly any unit tests, is usually not the best way to go about making the PR easy to review.

Is there no way to split this up into smaller parts, ideally with documentation and unit test coverage for each part?

akva2 · 2025-12-04T13:01:04Z

Since it's entirely a one way dependency, how about we extract the block size 7 changes and do those first?

svenn-t · 2025-12-04T13:46:25Z

Since it's entirely a one way dependency, how about we extract the block size 7 changes and do those first?

I have made a new PR #6660 with only the block size increase. Thanks for the suggestion @akva2! The number of GPU-related files that needed to be modified was more than I thought initially! I will rebase this PR once the block size PR is merged.

I'll just make one general comment here and then I'll be quiet. Adding 6,600 lines of new code in a single PR, especially since there's hardly any unit tests, is usually not the best way to go about making the PR easy to review.

Is there no way to split this up into smaller parts, ideally with documentation and unit test coverage for each part?

Thanks for the comment @bska. I have added some additional unit tests, including the PR in opm-tests OPM/opm-tests#1440. I hope this covers and shows the workings of the TPSA code better. I can also add that TPSA, structurally speaking, is (almost) a straightforward copy of TPFA in Flow. I hope this will make it easier to review.

bska · 2025-12-05T08:34:24Z

I have added some additional unit tests, including the PR in opm-tests

Thanks for adding tests, that certainly aids in the review.

Tangentially, my own personal goal when I'm adding new code is to cover the public interface of each new class and public function with a comprehensive set of unit tests. I also try to make adding each class/public function into a separate PR and then having an integration PR at the end that ties everything together and activates the new logic. I don't always achieve this goal but it feels good when I do.

opm/models/blackoil/blackoilintensivequantities.hh

opm/models/discretization/common/tpsalinearizer.hpp

svenn-t · 2025-12-12T13:58:23Z

An issue I have found is that the flow executable does not show the TPSA-related parameters (i.e., command-line options starting with "--tpsa"). The TPSA-specific simulators, f.ex. flow_blackoil_tpsa, show the parameters just fine. I suspect it's related to the FlowEarlyBird type tag in Main.hpp not knowing about TPSA-related definitions, and thereby the TPSA parameters will not be registered for flow (?).

If anyone has any ideas how to get the TPSA parameters in flow or some hints/suggestions what I should try, it will be much appreciated!

akva2 · 2025-12-15T09:47:47Z

Sadly this is a limitation of the system. The same applies to other model specific options, such as solvent.
We never put in the effort to make this behave properly, nor made a decision on what is proper behavior.

To make it work, you have to register the parameters manually, since, as you allude to, they are not registered through the FlowEarlyBird typetag since it has the standard blackoil model. One way of achieving this conceptually is by adding a call to the registration function in

SimulatorFullyImplicitBlackoil::registerParameters()

hnil · 2025-12-15T11:06:43Z

I had a brief look. This seems very good new contributions. It is an alternative to part of opm-flowgeomechanics, which do vem discretization + some fracturing. You have done a very good job at getting things into the system with models +++, which is partly avoided in opm-flowgeomechanics. Also maybe the work in opm-flowgeomechanics should be changed a bit to be more like "opm" interface as the tpsa uses.

Some comments:

maybe even more could have been done with separate classes which inherits from main classes
it would be usefull to sync part of the classes with opm-flowgeomechanics so solvers can be interchanged. I see you had to adapt to some output which was already merged into master, probably the same is the case for input in particular boundary conditions.
Also I will recommand to look at using hypre based solver. It seems to be much more stable in particular non default partition is required which normally is better for mechanics.

svenn-t · 2025-12-15T13:07:59Z

Sadly this is a limitation of the system. The same applies to other model specific options, such as solvent. We never put in the effort to make this behave properly, nor made a decision on what is proper behavior.

To make it work, you have to register the parameters manually, since, as you allude to, they are not registered through the FlowEarlyBird typetag since it has the standard blackoil model. One way of achieving this conceptually is by adding a call to the registration function in
SimulatorFullyImplicitBlackoil::registerParameters()

Thanks for the hint! I have included linear solver and newton parameters for TPSA directly in SimulatorFullyImplicitBlackoil::registerParameters() in the latest update. It might not be the most elegant solution since these are very specific parameters for TPSA included in the general simulator class, but it seems to work at least!

svenn-t · 2025-12-15T13:29:37Z

I had a brief look. This seems very good new contributions. It is an alternative to part of opm-flowgeomechanics, which do vem discretization + some fracturing. You have done a very good job at getting things into the system with models +++, which is partly avoided in opm-flowgeomechanics. Also maybe the work in opm-flowgeomechanics should be changed a bit to be more like "opm" interface as the tpsa uses.

Thanks for the feedback, @hnil. I appreciate you taking the time to have a look at this.

Some comments:

maybe even more could have been done with separate classes which inherits from main classes

Yes, there may be some classes which I have name somethingTPSA (like FlowProblemTPSA, BlackoilModelTPSA, ISTLSolverTPSA) which could be generalized to fit any geomechanics model.

it would be usefull to sync part of the classes with opm-flowgeomechanics so solvers can be interchanged. I see you had to adapt to some output which was already merged into master, probably the same is the case for input in particular boundary conditions.

We are still figuring out some of the outputs (particularly stress and strain) and boundary conditions (particularly Robin conditions). It would be worth syncing up with opm-flowgeomechanics, in general, to avoid duplicating code!

Also I will recommand to look at using hypre based solver. It seems to be much more stable in particular non default partition is required which normally is better for mechanics.

Thanks, I will look into hypre-based solvers! ISTL solvers have worked great, at least on the examples we have run so far, which I think is due to the nature of TPSA being closely related to TPFA.

totto82

Thanks for the work. Most of my comments are minor. My main concern is related to code duplication (also raised by @hnil) . There are many files that are essentially copies of other files for TPFA. Can some of these copies be avoided? Let's discuss it file by file (in a meeting?)

totto82 · 2025-12-17T13:18:48Z

opm/models/tpsa/elasticitylocalresidualtpsa.hpp

+
+/*!
+* \brief Calculation of (linear) elasticity model terms for the residual
+*/


Add the equations that are solved. (+ ref to paper?)

totto82 · 2025-12-17T13:21:53Z

opm/models/discretization/common/tpsalinearizer.hpp

+  it under the terms of the GNU General Public License as published by
+  the Free Software Foundation, either version 2 of the License, or
+  (at your option) any later version.
+


This is essentially a copy of tpfalinearizer.hpp. Can the TPFA code be re-used?

Unfortunately there are a few fundamental differences with the linearization, so I could not make it work with tpfalinearizer.

totto82 · 2025-12-17T13:22:28Z

opm/models/tpsa/elasticityprimaryvariables.hpp

+
+/*!
+* \brief Primary variables in (linear) elasticity equations
+*/


Add what the primary variables are

totto82 · 2025-12-17T13:24:59Z

opm/models/tpsa/tpsamodel.hpp

+    {
+        SymTensor val;
+        return val;
+    }


Why is this not implemented? (and below)

As far as I understand, many of the variables are outputted as cell center variables, and computing these with TPSA is a bit of challenge as of now. But we should definitely add, at least stress and strain, when we figure out how to compute them.

It maybe a better way, but a least square approximation based on forces on faces should at least be one way.

totto82 · 2025-12-17T13:26:17Z

opm/models/tpsa/tpsanewtonconvergencewriter.hpp

+*/
+template <class TypeTag>
+class TpsaNewtonConvergenceWriter
+{


Why is this needed?

It was unnecessary, so it is now removed

totto82 · 2025-12-17T13:28:08Z

opm/models/tpsa/tpsanewtonmethod.hpp

+*/
+template <class TypeTag>
+class TpsaNewtonMethod
+{


This is just Newton. Why can't you use newtonmethod.hpp directly (or at least inherit from it?

Unfortunately, and as far as my understanding, I could not inherit the standard NewtonMethod class (from opm/models/nonlinear) because I needed to override the "using" statements with TPSA related properties. At least, when I tried, the Newton method called Flow specific linearizer, linear solver, etc, instead of the TPSA specific functions and classes.

I have, however, done a refactoring of the TPSA Newton code to remove all unecessary and unused code.

In principle you could make a separate typetag just for the newton related types. You can splice in a 'default model' of this to the default typetag, so nothing changes in "regular" code, but replace it in the instance of the newton method for the tpsa model. Whether it is worth it, I make no judgement on, just a how it could be done.

Thanks for the tip! I will at least give it a try. Are there any examples of splicing typetags in the code already that I can look at for reference?

You can for instance have a look at fvbaseproperties.hh and how the linearsolver settings is spliced in.
The gist of it is that you have a splice, you assign types to the splice, then this is spliced into the main typetag so you can request those types from the main typetag.

totto82 · 2025-12-17T13:31:32Z

opm/simulators/flow/FacePropertiesTPSA_impl.hpp

+    outside.faceCenter = inside.faceCenter = geometry.center();
+
+    // OBS: Have not checked if this points from cell with lower to higher index!
+    faceNormal = intersection.centerUnitOuterNormal();


Is this comment still valid?

Yes, I have not tested on grids other than CpGrid. I replaced the function content with a throw, since I do not know if TPSA works properly with for example polyhedral grid.

totto82 · 2025-12-17T13:51:04Z

opm/simulators/linalg/ISTLSolverTPSA.hpp

+    * \param simulator Simulator object
+    */
+    ISTLSolverTPSA(const Simulator& simulator)
+        : simulator_(simulator)


Can you inherit / use the standard ISTLSolver here?

I could not inherit ISTLSolver class (for the same reasons as with Newton), but I changed the code to inherit AbstractISTLSolver. That required a small change in templates for AbstractISTLSolver, which I hope does not break anything!

I have also removed the setupPropertyTreeTPSA function, and just use setupPropertyTree instead with a small change in input.

For a way to separate the code almost as much as possible see https://github.com/hnil/opm-flowgeomechanics

- TPSA-Flow coupling done at the nonlinear solver step level: - -- Lagged coupling -> TPSA lags behind one time step - -- Fixed-stress -> Iterate until Flow and TPSA converge - Generic TPSA solver using Newton with associated linearizer. - Linear elasticity equations implemented - Specialized version of ISTL linear solver for TPSA - Output to restart files and VTK

Addition of a cross product term for rotation equations.

Template change in AbstractISTLSolver necessary to generalize class.

Linear solver setup used setupPropertyTree, so new function unecesary. New boolean to avoid Flow-specific linear solver presets.

- Removed unnecessary (null) convergence writer - Removed constraints - Removed unused/empty functions - Removed auxiliary module calculations - Added verbosity level with runtime parameter

svenn-t mentioned this pull request Dec 2, 2025

Keywords and material state container for TPSA OPM/opm-common#4846

Draft

totto82 added the manual:new-feature This is a new feature and should be described in the manual label Dec 2, 2025

akva2 reviewed Dec 3, 2025

View reviewed changes

flow/flow_blackoil_tpsa.cpp Outdated Show resolved Hide resolved

svenn-t mentioned this pull request Dec 4, 2025

Test cases for TPSA geomechanics OPM/opm-tests#1440

Draft

svenn-t mentioned this pull request Dec 4, 2025

Increase block sizes to 7 #6660

Merged

svenn-t force-pushed the tpsa_pr branch from d73c3f9 to 33a8196 Compare December 5, 2025 08:29

akva2 reviewed Dec 8, 2025

View reviewed changes

opm/models/blackoil/blackoilintensivequantities.hh Outdated Show resolved Hide resolved

opm/models/blackoil/blackoilintensivequantities.hh Outdated Show resolved Hide resolved

opm/models/discretization/common/tpsalinearizer.hpp Outdated Show resolved Hide resolved

svenn-t force-pushed the tpsa_pr branch from 7defe30 to 999eaf4 Compare December 15, 2025 13:01

totto82 reviewed Dec 17, 2025

View reviewed changes

svenn-t added 12 commits December 19, 2025 11:58

Fix for simulators using properties from blackoilmodel.hh

c605586

Fix unused parameter warnings

7f635ce

Fix EOL at EOF

bf02ada

Fix sign compare warnings

1d1da04

Tests for TPSA primary variables

3f3b051

Tests for TPSA local residual computations

e7e5304

Rename and use check with tolerance for some quantities

ea55602

Added regression tests for TPSA

0d9bcfc

Doxygen updates

b0089ea

const fixes

573ad35

Fix for FREE boundary conditions.

e7fd3cb

Addition of a cross product term for rotation equations.

svenn-t added 6 commits December 19, 2025 11:58

Register TPSA runtime parameters for flow executable

605d769

Add docstring explanations

9db6d23

Inherit AbstractISTLSolver in ISTLSolverTPSA.

276f80b

Template change in AbstractISTLSolver necessary to generalize class.

Remove setupPropertyTreeTPSA.

f33c26f

Linear solver setup used setupPropertyTree, so new function unecesary. New boolean to avoid Flow-specific linear solver presets.

Major clean up of TPSA Newton method.

fe7b24a

- Removed unnecessary (null) convergence writer - Removed constraints - Removed unused/empty functions - Removed auxiliary module calculations - Added verbosity level with runtime parameter

Throw for DUNE grids other than CpGrid

d95ec01

svenn-t force-pushed the tpsa_pr branch from 999eaf4 to d95ec01 Compare December 19, 2025 11:17

Fix unused parameter warning

3556b43

Two-point stress approximation (TPSA) geomechanics #6650

Are you sure you want to change the base?

Two-point stress approximation (TPSA) geomechanics #6650

Conversation

svenn-t commented Dec 2, 2025 • edited Loading Uh oh! There was an error while loading. Please reload this page.

Uh oh!

Uh oh!

Uh oh!

bska commented Dec 3, 2025

Uh oh!

akva2 commented Dec 4, 2025

Uh oh!

svenn-t commented Dec 4, 2025 • edited Loading Uh oh! There was an error while loading. Please reload this page.

Uh oh!

Uh oh!

bska commented Dec 5, 2025

Uh oh!

Uh oh!

Uh oh!

Uh oh!

svenn-t commented Dec 12, 2025 • edited Loading Uh oh! There was an error while loading. Please reload this page.

Uh oh!

Uh oh!

akva2 commented Dec 15, 2025 • edited Loading Uh oh! There was an error while loading. Please reload this page.

Uh oh!

Uh oh!

hnil commented Dec 15, 2025

Uh oh!

svenn-t commented Dec 15, 2025

Uh oh!

svenn-t commented Dec 15, 2025

Uh oh!

totto82 left a comment

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

svenn-t Dec 19, 2025 • edited Loading Uh oh! There was an error while loading. Please reload this page.

Uh oh!

Choose a reason for hiding this comment

Uh oh!

akva2 Dec 19, 2025 • edited Loading Uh oh! There was an error while loading. Please reload this page.

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

svenn-t Dec 19, 2025 • edited Loading Uh oh! There was an error while loading. Please reload this page.

Uh oh!

Choose a reason for hiding this comment

svenn-t commented Dec 2, 2025 •

edited

Loading

svenn-t commented Dec 4, 2025 •

edited

Loading

svenn-t commented Dec 12, 2025 •

edited

Loading

akva2 commented Dec 15, 2025 •

edited

Loading

svenn-t Dec 19, 2025 •

edited

Loading

akva2 Dec 19, 2025 •

edited

Loading

svenn-t Dec 19, 2025 •

edited

Loading