Release notes of Version 1.1

This version inherits from all the features introduced in:

Highlights

  • NonLinearEvolutionProblem are now able to make a prediction of the solution, removing major convergence issues when imposed displacements are imposed.

  • The \(\bar{F}\) method, or FBar formulation, has been implemented for finite strain behaviours to handle nearly incompressible materials.

  • The regularization proposed by Faltus et al. in the context of the third medium contact has been implemented for plane strain, plane stress and tridmensional hypotheses.

  • MGIS’s contexts now handle gathering computation time information to create the performance table instead of the previously used CatchTimeSection.

  • Many methods have been deprecated to have a consistent error handling scheme based on MGIS’s one. As such, many methods and functions now take a MGIS’s Context as their first argument. The deprecated methods of AbstractNonLinearEvolutionProblem have been removed, see the known incompatibilities below.

  • The name of the materials and boundaries are automatically retrieved from mfem’s mesh.

Known incompatibilites

  • In previous versions, the failure of Hypre’s linear iterative solvers were discarded by MFEM/MGIS’s Newton solver due to the lack of methods to test their convergence in mfem’s version prior to 4.10. The parameter DiscardLinearSolverFailure can be passed to MFEM/MGIS’s Newton solver to recover the behavior of previous versions.

  • The deprecated methods of AbstractNonLinearEvolutionProblem and of its derived classes which did not take a Context have been removed, as well as the deprecated functions getGradient, getThermodynamicForce and getInternalStateVariable which did not take a Context. The overloads taking a Context as first argument shall be used instead.

  • Partial quadrature spaces are now associated with a location, i.e. a material or a boundary (see Issue 466):

    • the getId method of PartialQuadratureSpace has been removed, the getLocation method shall be used instead.

    • the constructors of PartialQuadratureSpace taking a material identifier have been removed. A LocationIdentifier, built from a MaterialIdentifier or a BoundaryIdentifier, shall be given instead, e.g. MaterialIdentifier{.id = m}.

    • a partial quadrature space defined on a boundary using a finite element discretization is now built on the boundary elements of the mesh rather than on the submesh of the boundary.

    • LocationIdentifier is no longer an aggregate and can not be assigned anymore.

    • QPDependency::getMaterialIdentifier, PartialQuadratureFunctionsSet::getMaterialIdentifiers and the overload of PartialQuadratureFunctionsSet::get taking a material identifier have been replaced by getLocation, getLocations and an overload taking a LocationIdentifier.

    • the constructor of QPDependency taking a material identifier now takes a LocationIdentifier.

    • the static buildQuadratureSpace methods of the behaviour integrators now take attributes::Throwing as first argument.

    • PartialQuadratureSpace::treatInvalidElementIndex now takes a LocationIdentifier instead of a material identifier.

    • PartialQuadratureSpace::getLocationName now returns a description of the location, e.g. material 'BEAM' (12). It previously returned the name alone for named materials and boundaries. The output of getInformation has changed accordingly.

    • the overloads of updateGridFunction taking a mesh or a submesh have been removed. updateGridFunction requires a grid function whose finite element space is managed by the finite element spaces manager of the partial quadrature functions, as the ones created by makeGridFunction.

    • updateL2Projection now reports an error if a function is defined on a material or a boundary which is not part of the mesh of the result. Such functions were previously ignored.

    • the setMaterialProperty and setExternalStateVariable methods of behaviour integrators, called with an evaluator, require an evaluator defined on an equivalent partial quadrature space, i.e. built on a finite element discretization sharing the same finite element spaces manager. Those calls are now collective.

  • updateGridFunction now takes a Context as first argument and returns a boolean (see Issue 475).

  • assign_values now reports an error if both views are not defined on the same mesh discretization, even if both discretizations are built from the same mesh (see Issue 476). Two problems exchanging values must be built on the same mesh discretization, for instance by sharing a finite element spaces manager.

  • the deprecated functions getMPIrank and getMPIsize have been removed. getMPIRankForDefaultCommunicator and getMPISizeForDefaultCommunicator shall be used instead (see Issue 474).

New features

Prediction of the solution

By default, a nonlinear evolution problem uses the solution at the beginning of the time step, modified by applying Dirichlet boundary conditions, as the initial guess of the solution at the end of the time step, see below for details.

This can be changed by using the setPredictionPolicy method, as follows:

// use the elastic operator by default
mechanics.setPredictionPolicy(
   {.strategy = PredictionStrategy::BEGINNING_OF_TIME_STEP_PREDICTION});

Available strategies

Default prediction (PredictionStrategy::DEFAULT_PREDICTION)

By default, a nonlinear evolution problem uses the solution at the beginning of the time step, modified by applying Dirichlet boundary conditions, as the initial guess of the solution at the end of the time step.

Warning

In mechanics, this may lead to very high increments of the deformation gradients or the strain in the neighboring elements of boundaries where evolving displacements are imposed.

Prediction for the state at the beginning of the time step (PredictionStrategy::BEGINNING_OF_TIME_STEP_PREDICTION)

The prediction for the state at the beginning of the time step strategy determines the increment of the displacement \(\Delta\,\mathbb{u}\) by solving the following linear system:

\[\mathbb{K}\,\cdot\,\Delta\,\mathbb{u} = \ets{\mathbb{F}_{e}}-\bts{\mathbb{F}_{i}}\]

where:

  • \(\mathbb{K}_{e}\) denotes one of the prediction operator (see below).

  • \(\bts{\mathbb{F}_{e}}\) denotes the external forces at the beginning of the time step.

  • \(\bts{\mathbb{F}_{i}}\) denotes the inner forces at the beginning of the time step.

  • \(\Delta\,\mathbb{u}\) is submitted to the increment of the imposed Dirichlet boundary conditions.

Note

Although the wording explicitly refers to mechanics, this equation applies to all physics.

The following prediction operators can be chosen:

  • PredictionOperator::ELASTIC: the elastic operator

  • PredictionOperator::SECANT: the secant operator is typically defined by the elastic operator

  • PredictionOperator::TANGENT_PREDICTION: the tangent operator, defined by the time-continuous derivative of the thermodynamic force with respect to the gradients.

  • PredictionOperator::LAST_ITERATE_OPERATOR: this operator reuses the one computed at the last iteration of the previous time step. At the first time step, the elastic operator is used.

Prediction based on a behaviour integration with constant gradients (PredictionStrategy::CONSTANT_GRADIENTS_INTEGRATION_PREDICTION)

The CONSTANT_GRADIENTS_INTEGRATION strategy determines the increment of the unknown \(\Delta\,{u}\) by solving the following linear system:

\[\tilde{\mathbb{K}}\,\cdot\,\Delta\,\mathbb{u} = \ets{\mathbb{F}_{e}}-\ets{\tilde{\mathbb{F}}_{i}}\]

where:

  • \(\tilde{\mathbb{K}}\) denotes the operator computed at the end of the behaviour integration.

  • \(\bts{\mathbb{F}_{e}}\) denotes the external forces at the beginning of the time step.

  • \(\ets{\tilde{\mathbb{F}}_{i}}\) denotes an approximation inner forces at the end of the time step computed by assuming that the gradients are constant over the time (and thus equal to their values at the beginning of the time step).

  • \(\Delta\,\mathbb{u}\) is submitted to the increment of the imposed Dirichlet boundary conditions.

The behaviour integration allows taking into account:

  • the evolution of stress-free strain over the time step (thermal expansion, swelling, etc..),

  • the viscoplastic relaxation of the stress. This relaxation can be discarded by integrating the behaviour with a null time step.

The following operators are available:

  • IntegrationOperator::ELASTIC: the elastic operator,

  • IntegrationOperator::SECANT: the secant operator is typically defined by the elastic operator affected by damage,

  • IntegrationOperator::TANGENT: the tangent operator, defined by the time-continuous derivative of the thermodynamic force with respect to the gradients,

  • IntegrationOperator::CONSISTENT_TANGENT: the consistent tangent operator, defined by the derivative of the thermodynamic force with respect to the gradients at the end of the time step. See [ST85] for details.

FBar formulation

The \(\bar{F}\) method, or FBar formulation, is implemented following [dSNPDO96]. This formulation is designed to handle nearly incompressible materials in large strain analysis by using a modified deformation gradient \(\bar{\underline{F}}\) that separates volumetric and deviatoric responses.

The method replaces the standard deformation gradient \(\underline{F}\) with an assumed modified counterpart \(\bar{\underline{F}}\) in the computation of stresses. This modification is based on a multiplicative split into volumetric and deviatoric parts:

\[\underline{F} = J^{1/3}\, \bar{\underline{F}}\]

where \(J = \det{\underline{F}}\) is the Jacobian of the deformation gradient. This formulation effectively avoids locking issues in nearly incompressible materials while maintaining accuracy.

The \(\bar{F}\) formulation is particularly suited for low-order finite elements and is applicable to arbitrary material models. It ensures quadratic rates of convergence in Newton-Raphson schemes.

This formulation is enabled by passing an additional parameter to the Mechanics behaviour integrator, as follows:

const auto fbar_parameters = mfem_mgis::dict{
    {"Regularization", mfem_mgis::dict{{{"FBar", mfem_mgis::list{}}}}};
mechanics.addBehaviourIntegrator(ctx, "Mechanics", library,
                                 behaviour, fbar_parameters) | or_die;

Faltus 2026 regularization

The regularization proposed by Faltus et al. in the context of contact mechanics using a third medium is implemented here [FAH]. This regularization only applies to finite strain behaviours. Currently, this regularization is only available for isotropic behaviours.

This regularization adds a contribution to the standard variational operator in finite strain and can be derived from an energy \(W\) which penalizes the difference between the deformation gradient \(\underline{F}\) at a given quadrature point and its value \(\bar{\underline{F}}\) at the centroid of the element:

\[W\left(\underline{F}, \bar{\underline{F}}\right) = \alpha\,\left(\underline{F}-\bar{\underline{F}}\right)\,\colon\, \left(\underline{F}-\bar{\underline{F}}\right)\]

where \(\alpha\) is a penalization coefficient.

This regularization is enabled by passing an additional parameter to the Mechanics behaviour integrator, as follows:

const auto faltus_parameters = mfem_mgis::Parameters{
    {"Regularization",
     mfem_mgis::Parameters{
         {"Faltus2026",
          mfem_mgis::Parameters{{{"PenalizationCoefficient", 1e11}}}}}};
mechanics.addBehaviourIntegrator(ctx, "Mechanics", "ThirdMedium", library,
                                 behaviour2, faltus_parameters) | or_die;

The info function

The info function allows displaying information about an object in an output stream.

Retrieving information on a finite element discretization

Example of usage
const auto& fed = problem.getFiniteElementDiscretization();
const auto success = mfem_mgis::info(ctx, fed);

Retrieving information on a partial quadrature space

The PartialQuadratureSpaceInformation structure contains some relevant information about a partial quadrature space:

  • the identifier and a description of the underlying material or boundary, including its name, if defined,

  • the total number of elements,

  • the total number of integration points,

  • the number of elements points per geometric type.

  • the number of quadrature points per geometric type.

This structure is created by:

  • getLocalInformation, which returns the information relative to the current process.

  • getInformation, which returns the information gathered from all processes.

The PartialQuadratureSpaceInformation structure can be printed to an output stream using the info function.

Example of usage
const auto& qspace =
   problem.getBehaviourIntegrator(1).getPartialQuadratureSpace();
const auto success = mfem_mgis::info(ctx, std::cout, qspace);

Evaluators of quantities at integration points

Overview

The setMaterialProperty and setExternalStateVariable methods of behaviour integrators

Resolution of dependencies of a nonlinear evolution problem using another

The resolveBehaviourIntegratorsDependencies function resolves the material properties and external state variables of the behaviour integrators of a nonlinear evolution problem using the gradients, thermodynamic forces and internal state variables of another one.

Both problems must exchange data defined on equivalent partial quadrature spaces. Hence, the finite element discretizations of both problems must share the same finite element spaces manager, i.e. the finite element discretization of one problem must be built from the finite element spaces manager of the other.

The Simulation class

Physical system, coupling schemes and models

Line-search-like handling of behaviour integration failures

Profiling Toolkit

Timers are now accessible by MGIS contexts, offering greater flexibility:

  • Eliminates the need for a single global object

  • Allows multiple contexts and the ability to distinguish different parts of the code managed through separate objects

  • Prevents potential negative interactions with other code using the same timer type.

Example of usage:

CatchTimeSection(ctx, "Class::FunctionName");

Refactoring of submesh and finite element space creation

The creation of submeshes and associated finite element spaces has been refactored to be handled through the MeshDiscretization and FiniteElementSpacesManager classes. This refactoring introduces the following new classes:

  • MeshDiscretization: A fundamental class that handles the lifetime of meshes and provides utilities for managing mesh-related operations, including submesh creation based on material or boundary identifiers.

  • FiniteElementDiscretization: Extends MeshDiscretization to handle the lifetime of finite element collections and spaces, providing a high-level interface for creating and managing finite element spaces.

  • FiniteElementSpacesManager: Manages similar finite element spaces (siblings) that share the same mesh and finite element collection but may have different vectorial dimensions, with automatic reuse of existing spaces.

This refactoring improves code organization, reduces duplication, and provides a more consistent and flexible API for working with submeshes and finite element spaces.

Partial quadrature functions defined on another mesh

A partial quadrature function can be updated from a grid function defined on another mesh of the same discretization, i.e. the whole mesh, a submesh of materials or a submesh of boundaries (see Issue 466). The correspondence between the elements of both meshes is given by MeshDiscretization::getElementsMapping.

A partial quadrature space defined on a boundary can be built on the boundary elements of the whole mesh or of a submesh of materials. Boundary elements oriented differently in a submesh of materials and in the whole mesh are refused.

The values of a partial quadrature function can be copied by assign_values into a function defined on the same location of another mesh of the same discretization. The values are then copied element by element, and the call is collective in parallel (see Issue 476).

Search of the behaviour integrator providing a quantity

The functions hasGradientProvider, hasThermodynamicForceProvider and hasInternalStateVariableProvider search, among the behaviour integrators defined on a given location, the one providing a given quantity. The ParaviewExportIntegrationPointResultsAtNodes post-processing relies on them, so that a result can be exported on a material carrying several behaviour integrators.

MPI utilities

The isMPIMainProcess, getMPIRank and getMPISize functions return the information for a given communicator or for the communicator of a mesh discretization. The isMPIMainProcessForDefaultCommunicator, getMPIRankForDefaultCommunicator and getMPISizeForDefaultCommunicator functions use the default communicator.

The collective synchronize function gathers the error messages of an execution context on all the processes of a communicator or of a mesh discretization. It does nothing for a sequential mesh discretization. A message shared by all the processes having one is kept unchanged. Otherwise, each message is preceded by the rank of its process:

- Error on rank 1: <error message>
- Error on rank 3: <error message>

Messages with several lines are preceded by a # Error on rank r title.

Issues fixed

  • Issue 491: Simulation: the message of a failed time step describes its first sub-step

  • Issue 482: MPI functions on a MeshDiscretization and updateGridFunction with no function

  • Issue 477: synchronize Context to merge error messages on MPI processes

  • Issue 476: Fix assign_values

  • Issue 475: Error handling in updateGridFunction

  • Issue 474: Add function isMPIRootProcess and remove tests getMPIrank() == 0

  • Issue 466: Refactor partial quadrature space to better support spaces on boundaries

  • Issue 429: ex1 and ex3 sources shall only live in mm

  • Issue 427: satoh, rve and mox examples are duplicates from mm-examples

  • Issue 404: Remove deprecated methods of AbstractNonLinearEvolutionProblem

  • Issue 292: Remove deprecated usage of getMaterial in ParaviewExportIntegrationPointResultsAtNodesBase::getPartialQuadratureFunctionViews and ParaviewExportIntegrationPointResultsAtNodesBase::getResultDescription

  • Issue 290: Add support for partial quadratures spaces on boundaries

  • Issue 288: master doesnt build without MPI

  • Issue 286: updateGridFunction dilutes interface values when the functions do not cover all the materials of the mesh

  • Issue 284: NonLinearModel does not export the initial state

  • Issue 281: Nodal export of integration point results uses the node index as an integration point index

  • Issue 279: ParaviewExportResults silently ignores unknown parameters

  • Issue 277: NonLinearModel calls the deprecated solve and loses the caller’s context

  • Issue 275: UniformHeatSourceBoundaryCondition: only the last integration point contributes to the residual 

  • Issue 273: Fix forgotten change while parallel execution bug in behaviour integration 

  • Issue 270: FirstIterationConvergenceCriterion never triggers coupling iterations + wrong post-processing time

  • Issue 268: Parallel execution bug in behaviour integration

  • Issue 264: Installation from the doc fails on a fresh setup

  • Issue 262: Refactor creation of submeshes and associated finite element spaces to handle it through MeshDiscretization and FiniteElementSpacesManager enhancement 

  • Issue 260: Refactor commented examples documentation

  • Issue 257: Allow NewtonSolver to discard linear solver failures

  • Issue 254: Missing parameter option for GMRESSolver

  • Issue 253: Incomplete linear solver convergence checks in NonLinearEvolutionProblemImplementation.cxx and NewtonSolver.cxx

  • Issue 248: Improve PartialQuadratureFunction interface

  • Issue 245: Check the size of the unknown when defining bricks

  • Issue 240: Small bug in LinearSolverFactory.cxx

  • Issue 237: [cmake] Add a build-tests target

  • Issue 218: [performance] synchronize success of the setup methods at a higher level to minimize collective communications enhancement

  • Issue 213:  Add a simple way to resolve dependencies (material properties, external state variables) of a NonLinearEvolutionProblem using the gradients, thermodynamic forces and internal state variables of another one.

  • Issue 211: Add coupling schemes and models

  • Issue 209: Introduce the Simulation class

  • Issue 206: Line-search-like handling of behaviour integration failures

  • Issue 200: automatically assign materials and boundaries’s names from mfem’s attributes 

  • Issue 198: Add the ability to define multiple behaviour integrators on the same material

  • Issue 193: Add support for other types of search operators for computing a prediction of the solution at the end of the time step enhancement

  • Issue 192: Take external forces into account when computing the prediction of the solution

  • Issue 188: retrieve information about a quadrature space

  • Issue 149: work on the prediction of the solution