Define user problems

Fedoo can be extended at two complementary levels:

• a weak form defines the equation assembled on the mesh;

• a constitutive law defines the material response used by a mechanical weak form.

For a linear equation, a DiffOp expression and fedoo.WeakForm are generally sufficient. For a nonlinear or history-dependent equation, derive from WeakFormBase. For a mechanical material, MechanicalUMAT is the recommended extension point; derive directly from Mechanical3D only when full control over the Fedoo constitutive lifecycle is required.

Define a weak form from differential operators

A modeling space owns the unknown fields and creates their differential operators. The virtual property selects the test-function side of a weak term. Ordinary Python arithmetic combines these objects into a DiffOp expression.

The following example defines the weak form of the Poisson equation

\[\int_\Omega \nabla \delta u \mathbin{\cdot} \nabla u\,d\Omega = 0.\]

import fedoo as fd

space = fd.ModelingSpace("2D")
space.new_variable("U")

dU_dX = space.derivative("U", "X")
dU_dY = space.derivative("U", "Y")

equation = dU_dX.virtual * dU_dX + dU_dY.virtual * dU_dY
weakform = fd.WeakForm(equation)

mesh = fd.mesh.rectangle_mesh()
assembly = fd.Assembly.create(weakform, mesh)
problem = fd.problem.Linear(assembly)

problem.bc.add("Dirichlet", "left", "U", 0)
problem.bc.add("Dirichlet", "right", "U", 1)
problem.solve()

The same construction applies to systems of equations. Declare every scalar unknown with fedoo.ModelingSpace.new_variable(), group components with fedoo.ModelingSpace.new_vector() when useful, and add the weak terms. Mixed formulations may assign different element interpolations to different variables.

Complete examples are available in poison_eq.py and fluid_mechanics.py.

State-dependent weak forms

When coefficients or residual terms depend on the current solution, derive a class from WeakFormBase. Its principal methods are:

initialize(assembly, problem)

Allocate values required by the formulation in assembly.sv.

update(assembly, problem)

Evaluate the current fields and update the residual state before assembly.

get_weak_equation(assembly, problem)

Return the linearized DiffOp, including the tangent and current residual terms.

set_start(assembly, problem)

Commit the state at the beginning of a new time increment.

to_start(assembly, problem)

Restore the beginning-of-increment state after a rejected increment.

reset()

Clear the complete problem history when necessary.

The problem may have problem.dtime == 0 during initialization. A custom weak form must not divide by dtime in that state. Time-dependent terms can be omitted until a positive increment is supplied. See navier_stokes.py for a complete nonlinear example.

Define a mechanical constitutive law

Use MechanicalUMAT for a conventional material-point update. It handles the Fedoo lifecycle, state allocation, temperature lookup, material local frames, finite-rotation increments, and the transformations of strains, stresses, and tangents. The callback only integrates the constitutive equation in the supplied corotational material coordinates.

Recommended: MechanicalUMAT

The callback has the following signature:

def umat(
    strain,
    dstrain,
    F0,
    F1,
    stress,
    DR,
    props,
    statev_start,
    time,
    dtime,
    wm_start,
    temperature,
    *,
    ndi,
    tangent_mode,
    start,
    corate,
):
    ...
    return stress_new, statev_new, wm_new, tangent

Its arguments follow these conventions:

Argument	Meaning
strain	Total strain at the beginning of the increment, shape (6, n_points). In a finite-strain analysis this is the logarithmic strain in the corotational material frame, and the infinitesimal strain otherwise.
dstrain	Current strain increment, shape (6, n_points).
F0, F1	Beginning and current deformation gradients, shape (3, 3, n_points) in a finite-strain analysis. They are empty in a small-strain analysis.
stress	Cauchy stress at the beginning of the increment, shape (6, n_points).
DR	Objective incremental rotation, shape (3, 3, n_points).
props	Material properties stored as a Fortran-contiguous column.
statev_start	Private state at the beginning of the increment, shape (n_statev, n_points).
time, dtime	Current problem time and time increment.
wm_start	Beginning-of-increment work variables, shape (4, n_points).
temperature	Gauss-point temperatures, or None when no temperature field is available.
ndi	Number of direct components: 3 normally and 2 for the plane-stress interface.
tangent_mode	User-defined tangent selector forwarded unchanged by Fedoo.
start	True during initialization; False for every material increment, including corrections at time zero.
corate	Objective-rate name as a string, as selected by weakform.corate in finite strain (e.g. "log_r", "log", "jaumann", or "green_naghdi"). In small strain it is "log_r".

Stress and strain use engineering Voigt notation. The component ordering, including the axisymmetric convention, is documented by Mechanical3D.

The callback returns the updated stress, state variables, work variables, and the material tangent. A tangent can have shape (6, 6) when it is uniform, or (6, 6, n_points) when it varies between material points. Do not modify statev_start or wm_start in place: they represent the committed state and must remain available if the global increment is rejected.

Finite-strain stress and tangent convention

Fedoo follows the convention used by the Simcoon UMAT interface. A solid constitutive law returns the true Cauchy stress

\[\boldsymbol{\sigma},\]

and an unnormalized corotational Kirchhoff (box) tangent

\[\overset{\Box}{\boldsymbol{\tau}} = \mathbb{C}^{\Box}:\boldsymbol{D}, \qquad \boldsymbol{\tau}=J\boldsymbol{\sigma}.\]

Thus, a law written in terms of Kirchhoff stress must divide the updated stress by J = det(F1) before returning it. At small strain J tends to one, so this convention reduces to the usual stress–strain tangent.

The constitutive law does not transform its tangent for the selected finite element formulation. fedoo.weakform.StressEquilibrium performs the conversion centrally:

• in TL, TangentMatrix is converted to the material tangent dS/dE;

• in UL, it is converted to the spatial Lie tangent consistent with the Cauchy-stress residual;

• the corresponding initial-stress geometric stiffness is enabled by default in finite strain.

This convention applies to continuum solid laws used by StressEquilibrium. Beam and shell laws use their own generalized forces, moments, strains, and curvatures.

The low-level convert_tangent=False option disables these transformations. It is intended only for a custom weak-form implementation that already provides dS/dE in TL or the spatial Lie tangent in UL. When conversion is disabled, geometric stiffness is not enabled automatically.

fedoo.constitutivelaw.ElastoPlasticity is a readable Python reference for this convention: it integrates its radial return in the Kirchhoff measure and divides by J before returning the stress.

See the Simcoon UMAT documentation for the material-update interface from which this convention is taken.

Initialization call

Fedoo calls the UMAT once during initialization with:

• start=True;

• time == 0 and dtime == 0;

• zero strain, strain increment, stress, and work arrays;

• the assembly’s initial state-variable values, or zeros when none were supplied;

• identity DR;

• identity F0 and F1 for finite strain, or empty deformation-gradient arrays for small strain.

The dtime == 0 call is not a physical material increment. It must return the initial elastic tangent without advancing irreversible state, accumulating work, or evaluating expressions that divide by dtime. This tangent is stored as ElasticMatrix when use_elastic_tangent=True and is restored at the beginning of each increment for the elastic predictor.

Use start to distinguish initialization from integration. All subsequent calls pass start=False, including Newton corrections of the first increment at time == 0. Return the material’s initialized state in statev_new; Fedoo stores it, with components set through set_initial_statev() taking precedence.

Initial state variables

State-variable initial conditions belong to an assembly rather than to the material, allowing one material instance to be reused by several problems with different histories. Set labeled values after creating the assembly and before initializing the problem:

import fedoo as fd

material = MechanicalUMAT(
    umat,
    props,
    n_statev=7,
    statev_label={"P": 0, "EP": slice(1, 7)},
)
weakform = fd.weakform.StressEquilibrium(material)
assembly = fd.Assembly.create(weakform, mesh)

material.set_initial_statev(assembly, "P", 0.02)
material.set_initial_statev(assembly, "EP", initial_plastic_strain)

A scalar is broadcast over the selected components and all Gauss points. A component vector is broadcast over the Gauss points, while a complete field has shape (n_components, assembly.n_gauss_points). For a scalar label, a one-dimensional array of length assembly.n_gauss_points supplies spatially varying values. Node and element fields are not converted implicitly.

Use set_initial_statev() rather than assigning assembly.sv["Statev"] directly: it records which components must be preserved during material initialization. The initial state is passed to the start=True UMAT call, so it may affect the initial tangent.

Calling assembly.reset() clears runtime state, including these initial values. Call set_initial_statev again before reinitializing a reset assembly.

Minimal elastic callback

This small 3D example illustrates the callback contract. Production laws should also validate their parameters and supported dimensional assumptions.

import numpy as np
from fedoo.core import MechanicalUMAT


def elastic_umat(
    strain, dstrain, F0, F1, stress, DR, props,
    statev_start, time, dtime, wm_start, temperature,
    *, ndi, tangent_mode, start, corate,
):
    E, nu = props[:, 0]
    shear = E / (2 * (1 + nu))
    lame = E * nu / ((1 + nu) * (1 - 2 * nu))

    tangent = np.zeros((6, 6))
    tangent[:3, :3] = lame
    tangent[0, 0] += 2 * shear
    tangent[1, 1] += 2 * shear
    tangent[2, 2] += 2 * shear
    tangent[3, 3] = shear
    tangent[4, 4] = shear
    tangent[5, 5] = shear

    if F1.size:
        J0 = np.linalg.det(F0.transpose(2, 0, 1))
        J1 = np.linalg.det(F1.transpose(2, 0, 1))
        tau = J0 * stress + tangent @ dstrain
        stress_new = tau / J1
    else:
        stress_new = stress + tangent @ dstrain
    return (
        stress_new,
        statev_start.copy(),
        wm_start.copy(),
        tangent,
    )


material = MechanicalUMAT(
    elastic_umat,
    props=[210_000.0, 0.3],
    props_label={"E": 0, "nu": 1},
    is_isotropic=True,
)

The constructor also accepts:

n_statev and statev_label

Allocate private history and expose selected components to post-processing.

density

Define the mass density for dynamic mechanical problems.

tangent_from_F

Mark a finite-strain law whose tangent is derived from deformation-gradient kinematics instead of the default corotational strain measure. Such a law requires a geometrically nonlinear formulation; F0 and F1 remain empty in small strain.

use_elastic_tangent

Control whether the initialization tangent is restored at the beginning of each increment.

tangent_mode

Select the tangent requested from the callback.

is_isotropic

Avoid unnecessary frame transformations for an invariant response.

Local frames and finite rotations

Assigning a frame with material.set_local_frame(...) defines the initial material coordinates. MechanicalUMAT transforms the strain, deformation gradient, starting stress, and incremental rotation into those coordinates, then transforms the returned stress and tangent back to the global frame. Uniform, nodal, elemental, and Gauss-point material frames are supported as described in Constitutive Law (fedoo.constitutivelaw).

Private state variables are opaque to Fedoo because their tensorial meaning is law-dependent. The callback must apply DR to any tensor-valued history that requires objective transport. Scalar history variables require no rotation.

The objective-rate formulation is selected on fedoo.weakform.StressEquilibrium through weakform.corate. A UMAT subclass with a strict requirement can set, for example:

class MyUMAT(MechanicalUMAT):
    required_corate = ("log_r", "log_r_inc")

Fedoo then validates the selected formulation when the material is initialized. The comparison is case-insensitive and applies only in finite strain. The default required_corate = None accepts any formulation. This attribute does not change weakform.corate; the selected name is passed as the callback’s corate argument regardless of this validation setting.

Advanced: derive directly from Mechanical3D

Deriving from Mechanical3D provides complete control but also transfers the constitutive lifecycle to the implementation. This route is appropriate when the response cannot be expressed by the UMAT contract, needs custom assembly state, or uses a specialized tangent representation.

A derived law is generally responsible for:

• calling super().__init__(name=name, density=density);

• declaring is_isotropic correctly;

• allocating its state and initial assembly.sv["TangentMatrix"] in initialize;

• reading assembly.sv["Strain"] and assembly.sv["DStrain"] and writing assembly.sv["Stress"] and assembly.sv["TangentMatrix"] in update;

• implementing commit, rollback, and reset behavior through set_start, to_start, and reset when it owns history variables;

• handling plane stress explicitly when that assumption is supported;

• using the frame-transformation helpers supplied by Mechanical3D when the law is anisotropic;

• setting _Lt_from_F only when the tangent follows deformation-gradient rather than strain-increment kinematics.

assembly.sv_start is a shallow beginning-of-increment snapshot. Do not modify arrays shared with it in place during a trial iteration. Replace trial arrays in assembly.sv or make explicit copies so that to_start can restore the committed state.

Direct Mechanical3D laws must also tolerate initialization before a positive time increment exists. They should always provide the initial tangent required to assemble the first global system. The implementations of fedoo.constitutivelaw.ElasticAnisotropic and the pedagogical fedoo.constitutivelaw.ElastoPlasticity provide useful references; the latter uses MechanicalUMAT for its lifecycle.