Initialization by jwboth · Pull Request #1686 · pmgbergen/porepy

jwboth · 2026-06-01T04:43:03Z

Three-fold extension:

Distinction between reference, relative and absolute variables (explicitly defined for each physics - only considering HM cases so far). If reference states are not set/updated the simulations should run as initially.
Constitutive laws in terms of relative variables.
Standardized initialization strategy as part of the simulation preparation - seamless transition from initialization to simulation.

The latter is defined for initialization of quasi-static initial conditions as well as reference states.

Types of changes

What types of changes does this PR introduce to PorePy?
Put an x in the boxes that apply.

Minor change (e.g., dependency bumps, broken links).
Bugfix (non-breaking change which fixes an issue).
New feature (non-breaking change which adds functionality).
Breaking change (fix or feature that would cause existing functionality to not work as expected).
Testing (contribution related to testing of existing or new functionality).
Documentation (contribution related to adding, improving, or fixing documentation).
Maintenance (e.g., improve logic and performance, remove obsolete code).
Other:

Checklist

Put an x in the boxes that apply or explain briefly why the box is not relevant.

The documentation is up-to-date.
Static typing is included in the update.
This PR does not duplicate existing functionality.
The update is covered by the test suite (including tests added in the PR).
If new skipped tests have been introduced in this PR, pytest was run with the --run-skipped flag.

…nce (not introduced yet).

jwboth

@IvarStefansson @keileg I tried highlighting the pieces in the code that I suspect to require a critical look. Hope it helps for reviewing the suggested concepts.

jwboth · 2026-06-02T04:43:25Z

+    has_momentum_balance_reference_state: bool = False
+    """Flag whether reference states have been defined."""


The changes in this class are representative for changes in other core model classes. Be critical about:

explicit introduction of essentially three versions of a variable (the absolute one derives from the others by addition)

the concept behind updating and shifting values including the flag above - it feels like one should have maybe the possibility to do this in a more global way and not model by model. Is it possible to utilize some decorators or collect variables in respective containers? Bit like with previous_timestep(). Argument against a global approach would possibly be memory. One could also argue for (like @IvarStefansson did yesterday) that we could ignore reference states in the remaining physics (pressure, temperature etc) since displacements are the only singular object and require correct treatment in terms of modelling, while the remaining once may instead only impact the numerical properties like conditioning.

the imposed DataSavingMixin that always exports all three variants for the variables (maybe too much?). Reference variables should for instance also not change over time and could be handled as constants (see the not so much used functionality of the Exporter)

I agree that there seems to be unnecessary amount of boilerplate.

jwboth · 2026-06-02T04:44:30Z

+class MechanicalAperture:
+    """Mechanical aperture of fractures wrt reference configuration."""
+
+    normal_component: Callable[[list[pp.Grid]], pp.ad.Operator]
+    """Operator giving the normal component of a vector on fractures."""
+
+    displacement_jump: Callable[[list[pp.Grid]], pp.ad.Operator]
+    """Displacement jump on fractures wrt reference ."""


Be critical about the handling of the difference apertures and the now explicit use of relative displacements.

jwboth · 2026-06-02T04:45:35Z

+    NOTE: The implementation of the mechanical stress below is in absolute terms.
+    The stress could be formulated wrt a reference configuration, in which case
+    a reference (background) stress would be required, only depending on the reference
+    displacement. Under the assumption of linearity, the summation of the reference and
+    the linear mechanical stress below cancels the reference contributions. For
+    simplicity, background stress is not supported in the current implementation.
+


Should one think more about including explicit background stress measurements or possibly a homogeneous background stress and then formulate the mechanical stress in terms of relative quantities?

jwboth · 2026-06-02T04:50:42Z

@@ -4549,7 +4619,7 @@ def porosity(self, subdomains: list[pp.Grid]) -> pp.ad.Operator:


 class PoroMechanicsPorosity(pp.PorePyModel):


Note that I did not cover all constitutive laws, e.g. code related to BartonBandis , damage etc.

What should be the right strategy? Coverage or sanity checks throwing an error/warning if reference states are used e.g. for unsupported constitutive laws or models.

jwboth · 2026-06-02T04:53:53Z

+    def relative_contact_traction(self, subdomains: list[pp.Grid]) -> pp.ad.Operator:
+        """Fracture contact traction increment [-].


Be critical here: Do we really need reference and relative versions of all variables? There is no constitutive law that will benefit from contact traction changes, or? Same would apply for saturation-type variables in multiphase flow settings. If there should really be an impact on the linear/nonlinear conditioning of the problem by solving for increments, one should maybe consider solving for increments wrt the previous time step?

IvarStefansson · 2026-06-08T08:06:48Z

+        """
+        return Scalar(0.0, "reference_mechanical_aperture")
+
+    def mechanical_aperture(self, subdomains: list[pp.Grid]) -> pp.ad.Operator:


Consider if this makes reference fracture gap redundant/is already covered by that.

IvarStefansson · 2026-06-08T08:08:43Z



-class DisplacementJumpAperture(DimensionReduction):
+class DisplacementJumpAperture(DimensionReduction, MechanicalAperture):


We need to be very careful about what this does to the interpretation of aperture (hydraulic vs mechanical)

IvarStefansson · 2026-06-08T08:16:20Z

    """See :class:`VariablesSinglePhaseFlow`."""

-    interface_darcy_flux: Callable[
+    relative_interface_darcy_flux: Callable[


Do we need relative fluxes?

IvarStefansson · 2026-06-08T08:21:16Z

+    has_momentum_balance_reference_state: bool = False
+    """Flag whether reference states have been defined."""


I agree that there seems to be unnecessary amount of boilerplate.

keileg · 2026-06-30T06:11:23Z

@jwboth Have I understood correctly that this PR is being replaced by the series of newer PRs on the same topic?

jwboth · 2026-06-30T11:49:20Z

@keileg yes, it was kept for reference, but it was hopeless to base a solution for the identified weak points based on the herein presented implementation.

@jwboth Have I understood correctly that this PR is being replaced by the series of newer PRs on the same topic?

keileg · 2026-07-06T11:37:40Z

Superseded by #1696, #1707 #1708 etc.

jwboth added 11 commits May 28, 2026 14:05

MAINT: Make constitutive laws depending on variations wrt some refere…

4da8de4

…nce (not introduced yet).

FEAT: Enable reference states for single phase flow

5278876

FEAT: Enable reference states for momentum balance

08f92e2

FEAT: Enable reference states for contact mechanics

0a82cd2

FEAT: Enable reference states for poromechanics

c9beff1

FEAT: Provide iterative initializtion routine during prepare_simulation

ebfd10f

MAINT: Move initialization strategy to separate file

646097f

MAINT: Reduce code through reuse

42faa5b

MAINT: Provide user control to parameters in initialization

dc6b2bf

STY: ruff, isort

093ec45

MAINT: Refactor initialization.

01887d5

jwboth requested review from IvarStefansson and keileg as code owners June 1, 2026 04:43

jwboth marked this pull request as draft June 1, 2026 04:43

jwboth commented Jun 2, 2026

View reviewed changes

IvarStefansson reviewed Jun 8, 2026

View reviewed changes

jwboth added 2 commits June 11, 2026 08:38

MAINT: Remove perturbation from reference for stresses.

e6d3762

MAINT: Use default font color for logging

2ee6fa8

jwboth mentioned this pull request Jun 19, 2026

Reference states #1696

Merged

13 tasks

keileg closed this Jul 6, 2026

keileg deleted the initialization branch September 3, 2026 06:20

		has_momentum_balance_reference_state: bool = False
		"""Flag whether reference states have been defined."""

		@@ -4549,7 +4619,7 @@ def porosity(self, subdomains: list[pp.Grid]) -> pp.ad.Operator:


		class PoroMechanicsPorosity(pp.PorePyModel):

		def relative_contact_traction(self, subdomains: list[pp.Grid]) -> pp.ad.Operator:
		"""Fracture contact traction increment [-].



		class DisplacementJumpAperture(DimensionReduction):
		class DisplacementJumpAperture(DimensionReduction, MechanicalAperture):

Uh oh!

Conversation

jwboth commented Jun 1, 2026 • edited Loading Uh oh! There was an error while loading. Please reload this page.

Uh oh!

Types of changes

Checklist

Uh oh!

jwboth 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!

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!

Uh oh!

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

keileg commented Jun 30, 2026

Uh oh!

jwboth commented Jun 30, 2026

Uh oh!

keileg commented Jul 6, 2026

Uh oh!

Reviewers

Assignees

Labels

Projects

Milestone

Development

Uh oh!

3 participants

jwboth commented Jun 1, 2026 •

edited

Loading