Solvers

Currently, EdelweissFE provides

• a nonlinear implicit static solver (NIST),

• a nonlinear explicit static solver (NEST),

• a nonlinear explicit dynamic solver (NED),

• a nonlinear implicit dynamic solver (NID, Newmark-beta time integration),

• a parallel nonlinear implicit static solver (NISTParallel),

• a parallel nonlinear explicit static solver (NESTParallel),

• a parallel nonlinear explicit dynamic solver (NEDParallel),

• a parallel nonlinear implicit dynamic solver (NIDParallel),

• and a parallel arc length solver (NISTPArcLength).

Choose the solver in the *solver definition:

*solver, name=mySolver, solver=NISTParallel

edelweissfe.config.solvers.getSolverByName(name)

Get the class type of the requested solver.

Resolved through the registry (solver category) rather than through a hand-maintained table private to this module. Such a table could only ever list solvers living inside this package, so an external package – EdelweissMeshfree, a plugin – had no way to contribute one; going through the registry means a built-in, an entry point and an in-process register() call are all equally reachable here. An unknown name now raises RegistryLookupError naming the available solvers, instead of a KeyError.

Solver names are now case-insensitive, deliberately. This resolver was the one case-sensitive registry in the codebase: it indexed a table with CamelCase keys and then read the class off the module under the same string, so the name doubled as the class attribute name and e.g. "nist" failed twice over, while 12 of the 13 legacy config/*.py registries already casefolded the name at the resolver. That audit amends rule (c) to sanction this: a name must not resolve differently depending on which front-end it arrived through, and the registry is reached by callers with no .inp parser in the loop. The change is strictly more permissive, so no existing input file changes meaning.

Parameters:

name (str) – The name of the solver to load (case insensitive).

Returns:

The solver class type.

Return type:

type

The residual and the tangent

Every entity that enters the equation system – an element, a load, a constraint – contributes two things: a vector \(\boldsymbol{P}\) to the residual, and its tangent \(\boldsymbol{K} = \partial \boldsymbol{P} / \partial \boldsymbol{U}\). The rule is the same for all of them:

• An element adds its internal forces, \(\boldsymbol{P} = \boldsymbol{f}_\mathrm{int} = \int_\Omega \boldsymbol{B}^T \boldsymbol{\sigma} \, \mathrm{d}V\), and its stiffness matrix.

• A load adds the negative of the external load it applies: a node force \(\boldsymbol{F}\) adds \(-\boldsymbol{F}\), a body force \(-\int_\Omega \boldsymbol{N}^T \boldsymbol{b} \, \mathrm{d}V\). A dead load has no tangent; a load that follows the deformation, such as a pressure on the current surface, has one.

• A constraint adds the negative of the force it exerts on its nodes – just what an element of the same mechanics would add as its internal force. A penalty spring of stiffness \(k\) that pulls a degree of freedom \(u\) towards \(\bar{u}\) adds \(P = k\,(u - \bar{u})\) and \(K = k\); a Lagrange multiplier \(\boldsymbol{\lambda}\) enforcing \(\boldsymbol{g}(\boldsymbol{u}) = \boldsymbol{0}\) adds \(\boldsymbol{\lambda}^T \partial \boldsymbol{g} / \partial \boldsymbol{u}\) to the rows of the nodes and \(\boldsymbol{g}\) to the rows of the multipliers. See Constraints for writing one.

A solver sums the contributions in three loops of every iteration: computeElements over the elements, assembleLoads over the loads – every load assembles its own contribution, wherever it acts (see LoadBase) – and assembleConstraints over the constraints that are active (see Constraints). Summed over all entities,

\[\boldsymbol{G}(\boldsymbol{U}) = \sum \boldsymbol{P} = \boldsymbol{f}_\mathrm{int} - \boldsymbol{f}_\mathrm{ext}\]

is the out-of-balance force, which equilibrium makes zero, and the summed \(\boldsymbol{K}\) is its derivative. Newton’s method replaces \(\boldsymbol{G}(\boldsymbol{U} + \Delta \boldsymbol{U}) = \boldsymbol{0}\) by its linearization \(\boldsymbol{G} + \boldsymbol{K}\, \Delta \boldsymbol{U} = \boldsymbol{0}\), so the implicit solvers solve

\[\boldsymbol{K}\, \Delta \boldsymbol{U} = \boldsymbol{R}, \qquad \boldsymbol{R} = -\boldsymbol{G} = \boldsymbol{f}_\mathrm{ext} - \boldsymbol{f}_\mathrm{int}\]

in every iteration. The right-hand side \(\boldsymbol{R}\) is what the convergence criteria check, and what the multi-point-constraint condensation and the Dirichlet boundary conditions act on. The other solvers use the same sums:

• NID subtracts the inertia and damping forces from \(\boldsymbol{R}\), and adds their tangents to \(\boldsymbol{K}\).

• NED accelerates the nodes with the net force \(\boldsymbol{R}\), and needs no tangent.

• NISTPArcLength sums the loads twice, at the start of the step (the dead load, \(\boldsymbol{P}_0\)) and at its end; their difference is the reference load, \(\boldsymbol{P}_\mathrm{ref}\). Scaled by the load parameter \(\lambda\), it enters \(\boldsymbol{R} = -(\boldsymbol{P}_0 + \lambda \boldsymbol{P}_\mathrm{ref} + \boldsymbol{P}_\mathrm{elements} + \boldsymbol{P}_\mathrm{constraints})\).

Why one rule? Because then \(\boldsymbol{K}\) is always the derivative of exactly what was added to \(\boldsymbol{P}\), for every entity alike. That can be checked by finite differences, entity by entity, and the same check catches a contribution whose sign alone is wrong: \(\boldsymbol{K}\) then comes out as \(-\partial \boldsymbol{P} / \partial \boldsymbol{U}\). A contribution negated together with its tangent passes it, so tests/test_residual_tangent_convention.py checks every constraint this way and also checks the direction of the forces it exerts (see Constraints); a dead load has no tangent, so there the sign of each element load is checked against the load integrated independently.

The elements’ sum alone, \(\boldsymbol{f}_\mathrm{int}\), is what the quasi-static solvers report as the nodal result P of a field output, and the dynamic solvers add the inertial and damping forces to it (Node field results): at a degree of freedom with a Dirichlet boundary condition and no load, it is the reaction force. NED’s net force \(\boldsymbol{R}\) is its own state, and is not reported.

What the solver asks the model

A solver solves; what the equation system is made of, the model says. NIST, NISTParallel and NISTPArcLength ask the model three questions, and know nothing of their own about what a model consists of:

1. Has the coupling of the degrees of freedom changed? At the start of every increment, and of every retry after a cutback, the topology update calls updateConnectivity(), before the increment is proposed. If a coupling changed, the solver rebuilds its equation system before it solves the increment. A change of the coupling is not a change of the mesh: no node moves or appears, so the solver does not re-equilibrate after it (equilibrateAfterModelChange), and it keeps its prediction if the rebuilt system numbers the degrees of freedom as before (extrapolateAcrossConnectivityChange).

2. What is the equation system made of? Whenever the solver (re)builds its equation system, it calls buildEquationSystem(), which returns the layout of the system, an EquationSystemLayout. The solver keeps it as long as the numbering of the degrees of freedom it holds – until the next rebuild – and takes from it:

   • the numbering of the degrees of freedom (the DofManager). It also locates the degrees of freedom of a field on a node set, those a Dirichlet boundary condition prescribes or a node force loads (locateFieldOnNodeSet());

   • the elements, a list of everything the solver computes like an element (computeElements);

   • the state the solution vector (and the velocity and acceleration of a dynamic solver) starts from after a rebuild, readStateFromModel(), and where the results of an accepted increment go, writeResultsToModel() (see Node field results);

   • the model the linear solver is set up for.

3. Is anything to be updated before the elements are computed? In every evaluation of the model – every Newton iteration, every trial state of a line search – the solver assembles the internal forces (assembleInternalForces()): it calls updateMaterialPoints(), right before it computes the elements.

The loops over the entities stay with the solver: NIST computes the elements one after the other, NISTParallel on all threads, each element into its own part of the residual and its own block of the system matrix (see System matrix assembly). A layout that updates material points loops over them through the solver’s forEachChunk(), so that they are evaluated on the same threads.

The other solvers ask a part of these questions. NID and NEST derive from NIST, and build, fill and write back their equation systems through the layout; NEST builds its system once per step and has no topology update, and NID assembles its mass matrix from the model’s finite elements. NED takes the numbering of the degrees of freedom from the layout and passes the state and the results through it, keeps the net force that drives its velocity update as its own state, and evaluates the model’s finite elements in its own loop. A step of it starts from the velocity and the acceleration the model holds, and leaves the velocity at its end time (advanceVelocity()), so that the steps of a dynamic analysis follow each other whichever solver solves them.

A finite element model answers as follows: its constraints update the nodes they couple (a contact constraint searches for the surface it touches); its equation system numbers all nodes and the active constraints with their scalar variables; its elements are the finite elements; every node carries its total solution from increment to increment; and nothing is to be updated before the elements, which evaluate the material at their quadrature points themselves.

A model of EdelweissMeshfree answers differently, with the same solvers:

• Its material points move through a fixed background grid of cells, and its particles carry their kernel functions along. When asked for the coupling, it first assigns the material points to the cells they lie in and the particles to the kernel functions they reach, then lets its constraints follow.

• Only the part of the grid that holds material points, and the nodes in reach of the particles, form the equation system (the active domain). So a node set of the grid may reach beyond the equation system: a Dirichlet boundary condition then prescribes, and a node force loads, the nodes of the set in the active domain.

• Its elements are its finite elements, the active cells, and the particles.

• Only finite element nodes carry a total solution from increment to increment. The material point method and RKPM are incremental on the current configuration – the cells are reset and the kernel functions move with their particles, while the state lives in the material points and particles – so the values at the nodes of the grid and of the kernel functions are coefficients of one increment, and start every rebuilt equation system at zero.

• Before the cells are computed, each cell interpolates the solution increment to its material points, and the material points evaluate their material; the cells then integrate the result.

NIST - Nonlinear Implicit Static

class edelweissfe.solvers.nonlinearimplicitstatic.NIST(jobInfo, journal, **kwargs)

This is the Nonlinear Implicit STatic – solver.

Parameters:

• jobInfo – A dictionary containing the job information.

• journal – The journal instance for logging.

Attributes:

convergenceTable

equationSystemLayout

incrementStatus

mpcTransformation

Methods

acceptIncrement(step, model, timeStep)	Commit the converged increment to the model, and keep the predictor's state; see acceptIncrement().
applyCorrection(dU, ddU, stepLength)	Move the iterate by the correction, scaled by the step length \(\alpha\) a line search chose (1 without a line search): \(\Delta \boldsymbol{U} \leftarrow \Delta \boldsymbol{U} + \alpha \, \Delta\Delta \boldsymbol{U}\).
applyDirichletToResidual(timeStep, R, dirichlets)	Impose the Dirichlet BCs on the residual.
applyDirichletToStiffness(K, dirichlets[, rhs])	Impose the Dirichlet BCs on the system matrix -- and, given the right-hand side(s) rhs (modified in place), eliminate the constrained DOFs from the columns as well; see edelweissfe.solvers.base.dirichlet.applyDirichletToStiffness().
applyOptionsOverride(fieldValues)	Apply an override of this solver's options, see applyOptionsOverride(); an override that sets an option where it means nothing is refused before anything of it applies (see refuseOptionsThatMeanNothing()).
applyStepActionsAtIncrementStart(model, ...)	Called when all step actions should be applied at the start of a step.
applyStepActionsAtStepEnd(model, stepActions)	Called when all step actions should finish a step.
applyStepActionsAtStepStart(model, step)	Apply every step action's step-start part (initial conditions, material initialization, prescribed fields, ...).
assembleAdditionalTerms(dU, R, F, K, timeStep)	Assemble terms of the equation system beyond the internal and external forces.
assembleConstraintForces(constraints, U_np, ...)	Evaluate every active constraint and add its residual contribution to the given vector; a constraint that is switched off (see active) is skipped.
assembleConstraints(constraints, U_np, dU, ...)	Loop over all constraints, and evaluate the active ones: each adds its residual contribution and the tangent of it, see computeKernels().
assembleInternalForces(U_np, dU, P, K, F, ...)	Assemble the internal forces of the model at a trial state, and their tangent: those of the elements of its equation system, see computeElements().
assembleLoads(loads, U_np, dU, P, K, timeStep)	Assemble the residual contributions of all loads: each adds the negative of the external load it applies, and the tangent of that contribution, which is nonzero only for a load that follows the deformation; see The residual and the tangent.
assembleLoadsAtIncrementEnd(stepActions, ...)	Assemble the residual contributions of the loads acting at the end of the increment, and their tangents; called by assembleResidualAndTangent().
assembleResidualAndTangent(U_n, dU, U_np, R, ...)	Evaluate the model at the trial solution \(U_{n+1} = U_n + \Delta U\): the right-hand side \(R = f_\mathrm{ext} - f_\mathrm{int}\) (see The residual and the tangent) with the terms of assembleAdditionalTerms(), condensed by the multi-point constraints, and the tangent.
assembleStiffnessCSR(K)	The assembled system matrix in CSR format, see assembledMatrix().
assembleSystemMatrix(K, quasiNewton)	The system matrix of this Newton iteration, before the Dirichlet boundary conditions are imposed: the tangent in CSR layout, condensed by the multi-point constraints -- or, once the increment has fallen back to the quasi-Newton method, the tangent of its first iteration.
attemptIncrement(step, model, timeStep)	Solve the increment by Newton's method, rebuilding the equation system first if the topology update changed the model, and re-equilibrating a changed model first if so configured; see attemptIncrement().
beginStep(step, model, ...)	Set up the linear solver, and start the step; see beginStep().
buildMPCTransformation(model[, stepActions])	Collect the linear dependency records from all multi-point constraints of the model and assemble the master-slave condensation operator for the current equation system.
checkConvergence(R, ddU, dU, F, ...)	Check the convergence, individually for each field and for the scalar variables, based on the current flux residual and the correction, using the convergence criterion of the solver (see edelweissfe.solvers.base.convergencecriteria).
checkDivergingSolution(...)	Check if the iterative solution scheme is diverging: if the flux residual of a field grew, from one checked iteration to the next, more than maxGrowingIter times in this increment.
checkDivergingSolutionOrFallBack(...)	Cut an increment back whose residual grew more often than maxGrowingIter allows (see checkDivergingSolution()) -- or, with the option fallBackToQuasiNewton, let it continue with the tangent of its first iteration (a quasi-Newton method), see assembleSystemMatrix().
checkMPCDirichletConflicts(transformation, ...)	Raise if any Dirichlet boundary condition of the step prescribes a DOF that is a slave DOF of a multi-point constraint.
collectLoads(stepActions)	The loads among a step's actions: every action that is a LoadBase, whatever step action module it is defined in, so that a new kind of load needs no change to any solver.
computeElements(elements, U_np, dU, P, K, F, ...)	Loop over all elements, and evalute them: each adds its internal forces and their tangent, see computeKernels().
createConvergenceCriterion()	Create the convergence criterion selected by the option convergenceCriterion.
createSystemMatrixAssembly(assemblyType)	The assembly of the system matrix of the equation system just built (option systemMatrixAssembly), on the complete CSR pattern of csrGenerator, with numberOfAssemblyWorkers() workers.
endStep(step, model)	Report the step's performance timing.
extrapolateLastIncrement(extrapolation, ...)	Predict the solution increment, the point Newton's method starts from, by extrapolating the last accepted increments of the step.
extrapolates(extrapolation, prevTimeStep)	Whether an increment that follows the accepted increment prevTimeStep starts from an extrapolated prediction.
extrapolationWeights(extrapolation, ...)	The weights of the last accepted increment and of the one before it in the prediction of extrapolateLastIncrement(): \(\Delta U = w_n \, \Delta U_n + w_{n-1} \, \Delta U_{n-1}\), with
finalizeIncrement(model)	Store what the converged increment leaves behind.
firstIterateNeedsCorrection(dirichlets)	Whether the first iterate of an increment without a prediction -- the converged state of the last increment -- needs a correction before it can be checked for convergence (see solveIncrement()): if Dirichlet boundary conditions act, whose prescribed increment the first correction makes.
forEachChunk(entities, work)	Call work(chunk) for chunks that together hold all entities: the solver's loop over independent entities, other than the elements, that the model's equation system layout evaluates, e.g. the material points of a meshfree model (see updateMaterialPoints()).
forgetEquationSystem()	Drop the equation system of the last build -- the model's layout of it with the DOF manager, the system matrix assembly and its CSR pattern, the multi-point-constraint condensation, the DOF vectors -- so that it is gone before the next one is built, never both at once.
forgetPredictionAfterFailure()	A failed increment drops the prediction, which may have caused the failure, unless retryWithoutExtrapolation has answered that already: then the accepted history is kept, and the smaller increment that follows the cutback extrapolates it.
getRestartData()	The state this solver carries from one increment to the next, as it is; see checkpointedState.
initializeIncrement(U_n, stepActions, model, ...)	Prepare what assembleAdditionalTerms() needs during this increment's Newton iterations.
isOutputIncrement(timeStep)	Whether field outputs, output managers and restart checkpoints are written after this increment.
lineSearchIsDue(iterationCounter)	Whether the correction of this iteration is scaled by a line search (lineSearch): after lineSearchAfterIterations iterations, in every iteration whose number is a multiple of lineSearchEveryIterations.
lineSearchStepLengths()	The trial step lengths of the line search, from the option lineSearchStepLengths.
linearSolve(A, b)	Solve the linear equation system.
locateConstrainedDofs(dirichlets)	Determine, up front, which global DOFs each Dirichlet BC constrains.
numberOfAssemblyWorkers()	How many threads write blocks of the system matrix at the same time: one, the elements are computed in the calling thread.
predictionSurvivesRebuild(modelHasChanged, ...)	Whether this increment may extrapolate the last accepted increments after the equation system was rebuilt.
prepareIncrement(step, model, isRetry)	The topology update, at the start of every increment; see prepareIncrement().
printResidualOutlierNodes(residualOutliers)	Print which nodes have the largest residuals.
quadraticLineSearch(largestResidualAt)	The step length \(\alpha\) by which a Newton correction is scaled, for a Newton iteration that overshoots or oscillates (e.g., materials switching between loading and unloading).
readStateFromModel(U)	Read the state the next increment starts from out of the model, into the vectors of the equation system just (re)built: here, the solution, through the model's layout of the equation system (see readStateFromModel()).
refuseOptionsThatMeanNothing(options, ...)	Refuse an option that a declaration of options sets where it means nothing, with the options in force after the declaration.
schema	Option schema for this solver, per OptionSchemaProvider.
selectLinearSolverForStep()	Make linSolver the linear solver the step about to begin names (linsolver): the one of the previous step if it names the same, else a new one.
setRestartData(data)	Restore the state getRestartData() returned, so that the resumed step continues exactly where the uninterrupted run would have.
setTrialState(dU, ddU, stepLength)	Set the trial state of a line search at a step length \(\alpha\) along the correction, and return its solution increment, \(\Delta \boldsymbol{U} + \alpha \, \Delta\Delta \boldsymbol{U}\), at which the model is evaluated next; the iterate stays where it is.
solveForCorrection(K_, R, dirichlets)	Solve the linearized system of a Newton iteration, \(\boldsymbol{K} \, \Delta\Delta \boldsymbol{U} = \boldsymbol{R}\), for the correction \(\Delta\Delta \boldsymbol{U}\).
solveIncrement(U_n, dU, P, K, stepActions, ...)	Newton's method for the solution increment of an increment.
updateRigidBodies(model, timeStep)	Refresh the kinematics of all rigid bodies in the model after a converged increment.
updateTopologyAndConnectivity(model, step[, ...])	Run the topology update, then let every mesh-dependent consumer catch up on it.
validateModelCapabilities(model)	Validate whether the solver supports the active features/constraints of the model.
writeResultsToModel(U, dU, P[, V, A, ...])	Write the results of a converged increment to the model, through the model's layout of the equation system: the incremental degrees of freedom are set to zero first, from where the next increment starts them (see zeroIncrementalCoefficients()), then the solution, its increment and the internal force are written (see writeResultsToModel()).

schema

Option schema for this solver, per OptionSchemaProvider.

Attributes:

pruneCondensedMatrixZeros

zeroFluxThreshold

zeroIncrementThreshold

alias of NISTSchema

acceptIncrement(step, model, timeStep)

Commit the converged increment to the model, and keep the predictor’s state; see acceptIncrement().

Parameters:

• step – The step being solved.

• model (FEModel) – The model tree.

• timeStep (TimeStep) – The converged increment.

applyCorrection(dU, ddU, stepLength)

Move the iterate by the correction, scaled by the step length \(\alpha\) a line search chose (1 without a line search): \(\Delta \boldsymbol{U} \leftarrow \Delta \boldsymbol{U} + \alpha \, \Delta\Delta \boldsymbol{U}\). The arc length method moves the load parameter, too, see applyCorrection().

Parameters:

• dU (DofVector) – The solution increment of the iterate, moved in place.

• ddU (DofVector) – The correction of the iterate, see solveForCorrection().

• stepLength (float) – The step length.

Returns:

The correction the iterate was moved by, which the convergence check of the next iteration measures.

Return type:

DofVector

applyDirichletToStiffness(K, dirichlets, rhs=None)

Impose the Dirichlet BCs on the system matrix – and, given the right-hand side(s) rhs (modified in place), eliminate the constrained DOFs from the columns as well; see edelweissfe.solvers.base.dirichlet.applyDirichletToStiffness().

Parameters:

• K (csr_matrix)

• dirichlets (list[StepActionBase])

Return type:

csr_matrix

applyOptionsOverride(fieldValues)

Apply an override of this solver’s options, see applyOptionsOverride(); an override that sets an option where it means nothing is refused before anything of it applies (see refuseOptionsThatMeanNothing()).

Parameters:

fieldValues (dict) – Maps schema field name to its new, already-coerced value.

Return type:

None

assembleAdditionalTerms(dU, R, F, K, timeStep)

Assemble terms of the equation system beyond the internal and external forces. Called by assembleResidualAndTangent() in every evaluation of the model, next to computeElements(), assembleLoadsAtIncrementEnd() and assembleConstraints(): right after the right-hand side \(R = f_\mathrm{ext} - f_\mathrm{int}\) is formed (see The residual and the tangent), and before the multi-point-constraint condensation, the Dirichlet handling, the convergence test and the linear solve – which therefore all act on what is added here.

A static analysis has no such terms. A dynamic one, for example, subtracts the inertia force \(M \ddot{u}\) from the residual and adds the corresponding mass term to the stiffness – see NonlinearImplicitDynamic.

Parameters:

• dU (DofVector) – The current trial solution increment.

• R (DofVector) – The residual, to be augmented in place.

• F (DofVector) – The reference flux scale of the convergence test, to be augmented in place.

• K (SystemMatrixAssembly) – The tangent being assembled, to be augmented.

• timeStep (TimeStep) – The time step.

assembleConstraints(constraints, U_np, dU, P, K, timeStep, FConstraints=None)

Loop over all constraints, and evaluate the active ones: each adds its residual contribution and the tangent of it, see computeKernels(). A constraint that is switched off (see active) is not part of the equation system and is skipped. Is called by solveIncrement() in each iteration.

Parameters:

• constraints (list[ConstraintBase]) – The list of constraints.

• U_np (DofVector) – The current solution vector.

• dU (DofVector) – The current solution increment vector.

• P (DofVector) – The residual contribution vector to assemble into.

• K (SystemMatrixAssembly) – The system matrix.

• timeStep (TimeStep) – The time step.

• FConstraints (DofVector) – The vector of accumulated absolute constraint fluxes for convergence checks, if given.

Returns:

• The modified residual contribution vector.

• The modified system matrix.

Return type:

tuple[DofVector,SystemMatrixAssembly]

assembleInternalForces(U_np, dU, P, K, F, timeStep)

Assemble the internal forces of the model at a trial state, and their tangent: those of the elements of its equation system, see computeElements(). First, the model updates the material points that belong to no element, whose state the elements integrate (see updateMaterialPoints(); nothing to do in a finite element model).

The one evaluation of the elements of the solvers that derive from NIST: in every Newton iteration and for every trial state of a line search (assembleResidualAndTangent()), in every Runge-Kutta stage of NEST, and for the initial acceleration of NonlinearImplicitDynamic.

Parameters:

• U_np (DofVector) – The trial solution.

• dU (DofVector) – The trial solution increment.

• P (DofVector) – The internal forces, augmented.

• K (SystemMatrixAssembly) – The tangent being assembled.

• F (DofVector) – The accumulated absolute element fluxes, the reference of the convergence test, augmented.

• timeStep (TimeStep) – The time step.

Returns:

The internal forces, the tangent and the accumulated element fluxes.

Return type:

tuple[DofVector,SystemMatrixAssembly,DofVector]

assembleLoadsAtIncrementEnd(stepActions, U_np, dU, P, K, timeStep)

Assemble the residual contributions of the loads acting at the end of the increment, and their tangents; called by assembleResidualAndTangent(). Here, these are the loads of the step at the time of the end of the increment (see assembleLoads()). The arc length method replaces them by the dead load plus the reference load scaled by its load parameter, see NISTPArcLength.

Parameters:

• stepActions (dict) – The active step actions.

• U_np (DofVector) – The trial solution.

• dU (DofVector) – The trial solution increment.

• P (DofVector) – The residual contribution vector to assemble into.

• K (SystemMatrixAssembly) – The tangent being assembled.

• timeStep (TimeStep) – The time step.

Returns:

The updated residual contribution vector and tangent.

Return type:

tuple[DofVector,SystemMatrixAssembly]

assembleResidualAndTangent(U_n, dU, U_np, R, P, PLoadsAndConstraints, F, FConstraints, K, stepActions, model, timeStep)

Evaluate the model at the trial solution \(U_{n+1} = U_n + \Delta U\): the right-hand side \(R = f_\mathrm{ext} - f_\mathrm{int}\) (see The residual and the tangent) with the terms of assembleAdditionalTerms(), condensed by the multi-point constraints, and the tangent. The one place where the model is evaluated during the Newton iterations: called by solveIncrement() once per iteration, and by its line search once per trial step length. The loads are those at the end of the increment, see assembleLoadsAtIncrementEnd(). The Dirichlet boundary conditions are imposed by the caller.

The internal forces are those of assembleInternalForces(), which lets the model update the material points that belong to no element first.

Parameters:

• U_n (DofVector) – The solution at the start of the increment.

• dU (DofVector) – The trial solution increment.

• U_np (DofVector) – The trial solution, overwritten.

• R (DofVector) – The residual, overwritten.

• P (DofVector) – The internal forces (the reaction), overwritten.

• PLoadsAndConstraints (DofVector) – The residual contributions of the loads and the constraints, overwritten.

• F (DofVector) – The accumulated absolute element fluxes, the reference of the convergence test, overwritten.

• FConstraints (DofVector) – The accumulated absolute constraint fluxes, overwritten.

• K (SystemMatrixAssembly) – The tangent being assembled, assembled anew.

• stepActions (dict) – The active step actions.

• model (FEModel) – The model tree.

• timeStep (TimeStep) – The time step.

Returns:

The internal forces, the tangent and the accumulated element fluxes.

Return type:

tuple[DofVector,SystemMatrixAssembly,DofVector]

assembleSystemMatrix(K, quasiNewton)

The system matrix of this Newton iteration, before the Dirichlet boundary conditions are imposed: the tangent in CSR layout, condensed by the multi-point constraints – or, once the increment has fallen back to the quasi-Newton method, the tangent of its first iteration.

Parameters:

• K (SystemMatrixAssembly) – The tangent of this iteration, assembled.

• quasiNewton (QuasiNewtonFallback) – The quasi-Newton fallback of this attempt of the increment, which keeps the tangent of the first iteration.

Returns:

The system matrix.

Return type:

csr_matrix

attemptIncrement(step, model, timeStep)

Solve the increment by Newton’s method, rebuilding the equation system first if the topology update changed the model, and re-equilibrating a changed model first if so configured; see attemptIncrement().

Parameters:

• step – The step being solved.

• model (FEModel) – The model tree.

• timeStep (TimeStep) – The increment.

Raises:

IncrementFailed – If Newton’s method does not converge, or a material requests a cutback.

beginStep(step, model, fieldOutputController, outputmanagers)

Set up the linear solver, and start the step; see beginStep().

Parameters:

• step – The step to be solved.

• model (FEModel) – The model tree.

• fieldOutputController (FieldOutputController) – The field output controller.

• outputmanagers (dict[str, OutputManagerBase]) – The output managers.

checkDivergingSolution(incrementResidualHistory, maxGrowingIter)

Check if the iterative solution scheme is diverging: if the flux residual of a field grew, from one checked iteration to the next, more than maxGrowingIter times in this increment.

Two solver options decide what counts as a growth:

• countFirstResidualAsGrowth: the residual history of a field starts from zero, so the first checked residual of an increment is a growth already (EdelweissFE’s original counting); otherwise, counting starts with the second checked iteration.

• countGrowthOfConvergedFields: a field that is converged, while other fields keep the Newton iterations going, may still grow at the level of its round-off; such growths count too (EdelweissFE’s original counting), or are ignored.

Parameters:

• incrementResidualHistory (dict) – The convergence state of every field in each checked iteration of the increment, see checkConvergence().

• maxGrowingIter (int) – The maximum allowed number of growths of a residual during the iterative solution scheme.

Returns:

True if solution is diverging.

Return type:

bool

checkDivergingSolutionOrFallBack(quasiNewton, incrementResidualHistory, maxGrowingIter, nodesWithLargestResidual)

Cut an increment back whose residual grew more often than maxGrowingIter allows (see checkDivergingSolution()) – or, with the option fallBackToQuasiNewton, let it continue with the tangent of its first iteration (a quasi-Newton method), see assembleSystemMatrix(). Called by the Newton loop after every check that does not converge.

Parameters:

• quasiNewton (QuasiNewtonFallback) – The quasi-Newton fallback of this attempt of the increment.

• incrementResidualHistory (dict) – The convergence state of every field in each checked iteration of the increment.

• maxGrowingIter (int) – The maximum allowed number of growths of a residual.

• nodesWithLargestResidual (dict) – The nodes with the largest residual, per field, reported if the increment diverges.

Raises:

DivergingSolution – If the residual grew too often, and the increment does not fall back to the quasi-Newton method.

computeElements(elements, U_np, dU, P, K, F, timeStep)

Loop over all elements, and evalute them: each adds its internal forces and their tangent, see computeKernels(). Is called by assembleResidualAndTangent() in each iteration.

The elements are those of the model’s equation system (elements): the finite elements, and whatever else contributes like an element – in a meshfree model, its cells and particles.

Parameters:

• elements (list) – The list of elements.

• U_np (DofVector) – The current solution vector.

• dU (DofVector) – The current solution increment vector.

• P (DofVector) – The reaction vector.

• K (SystemMatrixAssembly) – The system matrix.

• F (DofVector) – The vector of accumulated fluxes for convergence checks.

• timeStep (TimeStep) – The time step.

Returns:

• The modified reaction vector.

• The modified system matrix.

• The modified accumulated flux vector.

Return type:

tuple[DofVector,SystemMatrixAssembly,DofVector]

createConvergenceCriterion()

Create the convergence criterion selected by the option convergenceCriterion.

Returns:

The convergence criterion.

Return type:

ConvergenceCriterionBase

createSystemMatrixAssembly(assemblyType)

The assembly of the system matrix of the equation system just built (option systemMatrixAssembly), on the complete CSR pattern of csrGenerator, with numberOfAssemblyWorkers() workers.

Parameters:

assemblyType (type) – The type of the assembly, one of systemMatrixAssemblies.

Returns:

The assembly.

Return type:

SystemMatrixAssembly

endStep(step, model)

Report the step’s performance timing.

Parameters:

• step – The step that was solved.

• model (FEModel) – The model tree.

extrapolateLastIncrement(extrapolation, timeStep, dU, dirichlets, prevTimeStep, model)

Predict the solution increment, the point Newton’s method starts from, by extrapolating the last accepted increments of the step. With \(h\) the size (step progress) of the new increment and \(h_n, \Delta U_n\) the size and the solution increment of the last accepted one:

• off: no prediction, \(\Delta U = 0\).

• linear: the last increment, scaled to the new size, \(\Delta U = \Delta U_n \, h / h_n\).

• quadratic: the parabola through the last three accepted states, evaluated at the end of the new increment. With the increment before the last one, \(h_{n-1}, \Delta U_{n-1}\),

  \[\Delta U = \Delta U_n \frac{h}{h_n} + \left( \frac{\Delta U_n}{h_n} - \frac{\Delta U_{n-1}}{h_{n-1}} \right) \frac{h \, (h + h_n)}{h_{n-1} + h_n},\]

  the linear prediction plus the change of the slope of the path. It reproduces a quadratic path exactly. While only one accepted increment of the step of nonzero size is available, the prediction is linear.

The prescribed degrees of freedom take their prescribed increment.

Parameters:

• extrapolation (str) – The type of extrapolation.

• timeStep (TimeStep) – The new increment.

• dU (DofVector) – The solution increment to predict, overwritten.

• dirichlets (list) – The active dirichlet boundary conditions.

• prevTimeStep (TimeStep) – The last accepted increment, or None for no prediction.

• model – The model tree.

Returns:

• The predicted solution increment.

• True if an extrapolation was performed.

Return type:

tuple[DofVector,bool]

extrapolates(extrapolation, prevTimeStep)

Whether an increment that follows the accepted increment prevTimeStep starts from an extrapolated prediction.

Parameters:

• extrapolation (str) – The type of extrapolation, off, linear or quadratic.

• prevTimeStep (TimeStep | None) – The last accepted increment, or None if there is none to extrapolate.

Returns:

True if the increment is extrapolated.

Return type:

bool

extrapolationWeights(extrapolation, timeStep, prevTimeStep)

The weights of the last accepted increment and of the one before it in the prediction of extrapolateLastIncrement(): \(\Delta U = w_n \, \Delta U_n + w_{n-1} \, \Delta U_{n-1}\), with

\[w_n = \frac{h}{h_n} + \frac{c}{h_n}, \qquad w_{n-1} = -\frac{c}{h_{n-1}}, \qquad c = \frac{h \, (h + h_n)}{h_{n-1} + h_n}\]

for a quadratic extrapolation, and \(w_n = h / h_n\) alone for a linear one. Every quantity that follows the path of the solution is predicted with these weights, e.g. the load parameter of the arc length method, too.

Parameters:

• extrapolation (str) – The type of extrapolation, linear or quadratic.

• timeStep (TimeStep) – The new increment.

• prevTimeStep (TimeStep) – The last accepted increment, which the prediction extrapolates.

Returns:

• The weight of the last accepted increment.

• The weight of the accepted increment before it, None for a linear prediction.

Return type:

tuple[float,float|None]

finalizeIncrement(model)

Store what the converged increment leaves behind. Called by solveIncrement() once its Newton loop has converged – only for an accepted increment, since a failing one leaves by exception before it.

A static analysis stores nothing here. A dynamic one, for example, stores the velocity and acceleration of the converged increment, as the starting point of the next one.

Parameters:

model (FEModel) – The model tree.

firstIterateNeedsCorrection(dirichlets)

Whether the first iterate of an increment without a prediction – the converged state of the last increment – needs a correction before it can be checked for convergence (see solveIncrement()): if Dirichlet boundary conditions act, whose prescribed increment the first correction makes. Without them, the first iterate is checked like every later one, and is converged already if the loads have not changed. The arc length method always corrects it first, see firstIterateNeedsCorrection().

Parameters:

dirichlets (list) – The Dirichlet boundary conditions of the step.

Returns:

True if the first iterate is corrected before the first convergence check.

Return type:

bool

forEachChunk(entities, work)

Call work(chunk) for chunks that together hold all entities: the solver’s loop over independent entities, other than the elements, that the model’s equation system layout evaluates, e.g. the material points of a meshfree model (see updateMaterialPoints()). Here, one chunk of all entities, in the calling thread, as the elements are computed.

Parameters:

• entities (list) – The entities, a list.

• work (Callable[[list], None]) – Called for each chunk; it must write only to what belongs to the entities of its chunk.

forgetEquationSystem()

Drop the equation system of the last build – the model’s layout of it with the DOF manager, the system matrix assembly and its CSR pattern, the multi-point-constraint condensation, the DOF vectors – so that it is gone before the next one is built, never both at once.

Return type:

None

forgetPredictionAfterFailure()

A failed increment drops the prediction, which may have caused the failure, unless retryWithoutExtrapolation has answered that already: then the accepted history is kept, and the smaller increment that follows the cutback extrapolates it.

initializeIncrement(U_n, stepActions, model, timeStep)

Prepare what assembleAdditionalTerms() needs during this increment’s Newton iterations. Called by solveIncrement() once per increment, after the step actions have been applied at the increment’s start and before the increment is extrapolated.

A static analysis needs nothing here. A dynamic one, for example, assembles its mass matrix here and computes the terms the time integration adds to the stiffness for this time increment – see NonlinearImplicitDynamic.

Parameters:

• U_n (DofVector) – The solution at the start of the increment.

• stepActions (dict) – The active step actions.

• model (FEModel) – The model tree.

• timeStep (TimeStep) – The time step.

lineSearchIsDue(iterationCounter)

Whether the correction of this iteration is scaled by a line search (lineSearch): after lineSearchAfterIterations iterations, in every iteration whose number is a multiple of lineSearchEveryIterations.

Parameters:

iterationCounter (int) – The current iteration number.

Returns:

True if a line search is due.

Return type:

bool

lineSearchStepLengths()

The trial step lengths of the line search, from the option lineSearchStepLengths.

Returns:

The trial step lengths.

Return type:

list[float]

numberOfAssemblyWorkers()

How many threads write blocks of the system matrix at the same time: one, the elements are computed in the calling thread.

Returns:

The number of workers.

Return type:

int

predictionSurvivesRebuild(modelHasChanged, connectivityHasChanged)

Whether this increment may extrapolate the last accepted increments after the equation system was rebuilt. They are vectors of the equation system they were solved in, so they fit the rebuilt one only if it numbers the degrees of freedom alike, which the digests of the two numberings tell (dofNumberingDigest; the one of the last accepted increments is dofNumberingOfLastIncrement):

• after a change of the topology (e.g., a mesh refinement), never;

• after a change of the connectivity alone (e.g., a contact search), only with extrapolateAcrossConnectivityChange, and only if the numbering is unchanged;

• at the first increment of a step that continues from a restart checkpoint – resumed, or fallen back to – only if the rebuilt system numbers the degrees of freedom as the checkpoint’s last accepted increments are numbered: e.g., the active cells of a material point model may have changed since. A prediction without the digest of its numbering has nothing to be compared with, and is dropped.

Parameters:

• modelHasChanged (bool) – Whether the topology update changed the model.

• connectivityHasChanged (bool) – Whether it changed the connectivity alone.

Returns:

True if the prediction may extrapolate the last accepted increments.

Return type:

bool

prepareIncrement(step, model, isRetry)

The topology update, at the start of every increment; see prepareIncrement().

Parameters:

• step – The step being solved.

• model (FEModel) – The model tree.

• isRetry (bool) – True if the increment is retried after a cutback: the model modifiers already decided on the state the retry starts from.

quadraticLineSearch(largestResidualAt)

The step length \(\alpha\) by which a Newton correction is scaled, for a Newton iteration that overshoots or oscillates (e.g., materials switching between loading and unloading).

The largest residual along the correction, \(\| R(\alpha) \|_\infty\) with the Dirichlet degrees of freedom left out, is evaluated at the trial step lengths lineSearchStepLengths, and a parabola is fitted to these values (least squares for more than three). The step length minimizes the parabola on lineSearchStepLengthRange: it is the vertex of a parabola that opens upwards, limited to that range, and otherwise the end of the range where the parabola is smaller. Each trial costs one evaluation of the model; the next Newton iteration evaluates the chosen one.

Two kinds of trials fit no parabola:

• A trial state at which the residual is not finite, which the model could not evaluate: the trials decide themselves, and the step length is the trial step length of the smallest finite residual (of equal ones, the shortest), limited to the range. If no trial residual is finite, the increment is cut back.

• Trial residuals that are equal to round-off, which do not tell the step lengths apart: the correction is taken in full, the step length is 1.

The trial states are those setTrialState() sets: the solution increment moved along the correction \(\Delta\Delta U\), and, under arc length control, the load parameter moved with it, such that every trial state satisfies the controller’s equation (see setTrialState()).

Parameters:

largestResidualAt (Callable[[float], float]) – Evaluates the model at the trial state of a step length, and returns its largest residual, without the Dirichlet degrees of freedom.

Returns:

The step length.

Return type:

float

Raises:

DivergingSolution – If the residual is not finite at any trial step length.

readStateFromModel(U)

Read the state the next increment starts from out of the model, into the vectors of the equation system just (re)built: here, the solution, through the model’s layout of the equation system (see readStateFromModel()). A solver that integrates in time reads its velocity and acceleration as well.

Parameters:

U (DofVector) – The solution vector of the equation system, zero, overwritten.

refuseOptionsThatMeanNothing(options, declaredOptions)

Refuse an option that a declaration of options sets where it means nothing, with the options in force after the declaration. A declaration is the definition of the solver (its *solver datalines, or the arguments of its constructor), or an >>options block of a step. Refused are:

• retryWithoutExtrapolation and extrapolateAcrossConnectivityChange switched on with extrapolation=off: there is no prediction to retry without, nor one to keep across a rebuild;

• an option of the abaqus criterion (abaqus...) with convergenceCriterion=legacy.

An option the declaration does not set is not refused: an option set for some steps keeps its value in a step that switches off what it acts on, and rests there.

Parameters:

• options (dict) – The options in force after the declaration.

• declaredOptions – The names of the options the declaration sets.

Raises:

ValueError – If the declaration sets an option where it means nothing.

Return type:

None

selectLinearSolverForStep()

Make linSolver the linear solver the step about to begin names (linsolver): the one of the previous step if it names the same, else a new one.

A kept solver builds anew for the step’s equation system (setModel) and frees what it built before. A new solver per step would leave the previous one to Python’s cyclic garbage collector whenever its tree refers back to itself (a function in a solver script, whose globals hold the tree), and the steps’ matrices would pile up.

Return type:

None

setTrialState(dU, ddU, stepLength)

Set the trial state of a line search at a step length \(\alpha\) along the correction, and return its solution increment, \(\Delta \boldsymbol{U} + \alpha \, \Delta\Delta \boldsymbol{U}\), at which the model is evaluated next; the iterate stays where it is. The arc length method sets the load parameter of the trial state, too, see setTrialState().

Parameters:

• dU (DofVector) – The solution increment of the iterate.

• ddU (DofVector) – The correction of the iterate, see solveForCorrection().

• stepLength (float) – The step length.

Returns:

The solution increment of the trial state.

Return type:

DofVector

solveForCorrection(K_, R, dirichlets)

Solve the linearized system of a Newton iteration, \(\boldsymbol{K} \, \Delta\Delta \boldsymbol{U} = \boldsymbol{R}\), for the correction \(\Delta\Delta \boldsymbol{U}\). The Dirichlet boundary conditions replace the rows of their degrees of freedom by identity rows, and eliminate their columns into \(\boldsymbol{R}\) (see applyDirichletToStiffness()); with multi-point constraints, the slave degrees of freedom take the part their masters interpolate. The arc length method solves for two right-hand sides here, see solveForCorrection().

Parameters:

• K – The system matrix of the iteration, see assembleSystemMatrix().

• R (DofVector) – The right-hand side, the out-of-balance force with the Dirichlet boundary conditions imposed; the columns of the Dirichlet degrees of freedom are eliminated into it.

• dirichlets (list) – The Dirichlet boundary conditions of the step.

• K_ (csr_matrix)

Returns:

The correction.

Return type:

DofVector

solveIncrement(U_n, dU, P, K, stepActions, model, timeStep, prevTimeStep, extrapolation, maxIter, maxGrowingIter)

Newton’s method for the solution increment of an increment. It starts from the prediction (extrapolateLastIncrement()), and every iteration

1. evaluates the model at the iterate (assembleResidualAndTangent()) and imposes the Dirichlet boundary conditions on the residual,

2. checks the iterate for convergence – except the first iterate of an increment without a prediction, if it needs a correction first (firstIterateNeedsCorrection()),

3. solves the linearized system for the correction (solveForCorrection()),

4. and moves the iterate by the correction (applyCorrection()), scaled by the step length of a line search, if one is due, which evaluates the trial states setTrialState() sets.

This is the Newton loop of every implicit solver. The dynamics enter it through initializeIncrement(), assembleAdditionalTerms() and finalizeIncrement() (NonlinearImplicitDynamic). The arc length method makes the load parameter an unknown of the iterations, through the predictor and the methods of steps 2 to 4 (solveIncrement()).

Parameters:

• U_n (DofVector) – The solution at the start of the increment.

• dU (DofVector) – The solution increment, overwritten.

• P (DofVector) – The internal forces (the reaction), overwritten.

• K (SystemMatrixAssembly) – The tangent being assembled, assembled anew.

• stepActions (dict) – The active step actions.

• model (FEModel) – The model tree.

• timeStep (TimeStep) – The time step.

• prevTimeStep (TimeStep) – The last accepted increment, from which the prediction extrapolates, or None.

• extrapolation (str) – The type of extrapolation.

• maxIter (int) – The maximum number of iterations.

• maxGrowingIter (int) – The maximum number of residual growths before the increment is cut back.

Returns:

A tuple containing

• the new solution vector

• the solution increment

• the new reaction vector

• the number of required iterations

• the history of residuals per field

Return type:

tuple[DofVector,DofVector,DofVector,int,dict]

writeResultsToModel(U, dU, P, V=None, A=None, dynamicFields=None)

Write the results of a converged increment to the model, through the model’s layout of the equation system: the incremental degrees of freedom are set to zero first, from where the next increment starts them (see zeroIncrementalCoefficients()), then the solution, its increment and the internal force are written (see writeResultsToModel()). A solver that integrates in time passes its velocity and acceleration as well, and the internal force with its inertial and damping forces.

Parameters:

• U (DofVector) – The solution.

• dU (DofVector) – The solution increment.

• P (DofVector) – The internal force (with the inertial and damping forces of a solver that integrates in time); None to leave the model’s as it is.

• V (DofVector) – The velocity of a solver that integrates in time; None for a static solver.

• A (DofVector) – The acceleration of a solver that integrates in time; None for a static solver.

• dynamicFields (list) – The names of the fields the solver integrates in time; None for a static solver.

checkpointedState: dict | None = {'convergenceCriterion': <class 'edelweissfe.solvers.base.convergencecriteria.ConvergenceCriterionBase'>, 'dUOfIncrementBeforeLast': <class 'numpy.ndarray'>, 'dUOfLastIncrement': <class 'numpy.ndarray'>, 'dofNumberingOfLastIncrement': <class 'str'>, 'prevTimeStep': <class 'edelweissfe.timesteppers.timestep.TimeStep'>, 'timeStepOfIncrementBeforeLast': <class 'edelweissfe.timesteppers.timestep.TimeStep'>}

the last accepted increment and its dU, for a quadratic extrapolation the accepted increment before it and its dU, and the digest of the numbering of the degrees of freedom these dU are vectors of; and the state of the convergence criterion.

Type:

The predictor’s state between increments

convergenceCriterion

The convergence criterion of the current step, created when the step begins (None before).

dUOfIncrementBeforeLast

The solution increment of timeStepOfIncrementBeforeLast.

dUOfLastIncrement

a copy, since the working increment of the equation system changes in every iteration, and the model’s entry dU is a result, which no solver reads back.

Type:

The solution increment of the last accepted increment, kept for the predictor

dofNumberingOfLastIncrement

The digest of the numbering of the degrees of freedom of the equation system the last increment was accepted in (dofNumberingDigest): dUOfLastIncrement is a vector of it, and so is dUOfIncrementBeforeLast, which is kept only as long as the prediction survives the rebuilds of the equation system (see predictionSurvivesRebuild()).

equationSystemLayout = None

What the model’s current equation system is made of (see buildEquationSystem()), rebuilt with it; None until it is built.

linSolver

The linear solver, and the linsolver it was made for (see selectLinearSolverForStep()).

lineSearchStepLengthRange = (0.1, 1.5)

a shorter step would hardly move the iterate, a longer one would trust the fitted parabola far beyond the trial step lengths.

Type:

The step lengths a line search may choose from

linearSolverDefinitions

The *linearSolver blocks of the input file by name, which linsolver= may name (none without an input file).

prevTimeStep

no extrapolation).

Type:

The last accepted increment, from which the predictor extrapolates (None

supportsMPC = True

Whether this solver supports master-slave condensation / multi-point constraints (e.g. surface ties). Subclasses supporting MPCs must set this to True.

supportsModelModifiers = True

Whether this solver runs the topology update (e.g. h-adaptivity) at all; without it a modifier never runs. A solver that runs it only before its increment loop sets this too, and refuses the modifiers that would need it later (see ModelModifierBase.actsOnlyAtSimulationStart and NED.validateModelCapabilities).

systemMatrixAssemblies = {'direct': <class 'edelweissfe.numerics.systemmatrixassembly.DirectCSRSystemMatrixAssembly'>, 'vij': <class 'edelweissfe.numerics.systemmatrixassembly.VIJSystemMatrixAssembly'>}

The ways of summing the system matrix (option systemMatrixAssembly) this solver supports, see createSystemMatrixAssembly().

timeStepOfIncrementBeforeLast

The accepted increment before prevTimeStep, kept for a quadratic extrapolation only (else None).

class edelweissfe.solvers.nonlinearimplicitstatic.NISTSchema(defaultMaxIter=10, defaultCriticalIter=5, defaultMaxGrowingIter=10, extrapolation='linear', retryWithoutExtrapolation=False, extrapolateAcrossConnectivityChange=False, extrapolateAfterModelChange=True, equilibrateAfterModelChange=False, linsolver='pardiso', pruneCondensedMatrixZeros=None, systemMatrixAssembly='vij', useThreadedMPCCondensation=False, reportPerformanceFrequency=0, lineSearch=False, lineSearchAfterIterations=5, lineSearchEveryIterations=2, lineSearchStepLengths='0.8, 1.0, 1.2', fallBackToQuasiNewton=False, convergenceCriterion='legacy', iterationsForAlternativeTolerance=15, zeroFluxThreshold=None, zeroIncrementThreshold=None, minimumIterations=0, countFirstResidualAsGrowth=True, countGrowthOfConvergedFields=True, abaqusResidualTolerance=0.005, abaqusCorrectionTolerance=0.01, abaqusResidualToleranceAlternative=0.02, abaqusIterationsForAlternativeTolerance=9, abaqusZeroFluxThreshold=1e-05, abaqusZeroFluxCorrectionTolerance=0.001, abaqusLinearResidualTolerance=1e-08, abaqusInitialTimeAveragedFlux=0.01)

The options of the *solver datalines and of an >>options block routed to this solver, owned by this module and never mutated from outside it.

Mirrors NIST.SolverSpecificOptions one-for-one; that dict remains the actual source of truth consulted at runtime (self.options, a plain mutable dict) – this schema exists so the registry and the name-based >>options override mechanism have a typed description of what this solver accepts, without requiring every internal self.options[...] access to become a dataclass attribute access.

Attributes:

pruneCondensedMatrixZeros

zeroFluxThreshold

zeroIncrementThreshold

Parameters:

• defaultMaxIter (int | None)

• defaultCriticalIter (int | None)

• defaultMaxGrowingIter (int | None)

• extrapolation (str | None)

• retryWithoutExtrapolation (bool | None)

• extrapolateAcrossConnectivityChange (bool | None)

• extrapolateAfterModelChange (bool | None)

• equilibrateAfterModelChange (bool | None)

• linsolver (str | None)

• pruneCondensedMatrixZeros (bool | None)

• systemMatrixAssembly (str | None)

• useThreadedMPCCondensation (bool | None)

• reportPerformanceFrequency (int | None)

• lineSearch (bool | None)

• lineSearchAfterIterations (int | None)

• lineSearchEveryIterations (int | None)

• lineSearchStepLengths (str | None)

• fallBackToQuasiNewton (bool | None)

• convergenceCriterion (str | None)

• iterationsForAlternativeTolerance (int | None)

• zeroFluxThreshold (float | None)

• zeroIncrementThreshold (float | None)

• minimumIterations (int | None)

• countFirstResidualAsGrowth (bool | None)

• countGrowthOfConvergedFields (bool | None)

• abaqusResidualTolerance (float | None)

• abaqusCorrectionTolerance (float | None)

• abaqusResidualToleranceAlternative (float | None)

• abaqusIterationsForAlternativeTolerance (int | None)

• abaqusZeroFluxThreshold (float | None)

• abaqusZeroFluxCorrectionTolerance (float | None)

• abaqusLinearResidualTolerance (float | None)

• abaqusInitialTimeAveragedFlux (float | None)

class edelweissfe.solvers.nonlinearimplicitstatic.QuasiNewtonFallback(isEnabled)

The quasi-Newton fallback of one attempt of an increment (solver option fallBackToQuasiNewton): once the residual has grown more often than the step option maxGrowIter allows, the increment continues with the tangent of its first iteration instead of being cut back. A tangent that is not updated converges only linearly, but it cannot switch abruptly either, as a tangent of materials switching between loading and unloading does.

Parameters:

isEnabled (bool) – Whether the option is on; only then is the tangent of the first iteration kept.

isEngaged

True once the increment continues with the tangent of its first iteration.

tangentOfFirstIteration

The tangent of the first iteration, in CSR layout and condensed, without the Dirichlet rows and columns.

NIST

Option	Type / Default	Description
defaultMaxIter	int, default=10	The default maximum number of iterations.
defaultCriticalIter	int, default=5	The default number of critical iterations.
defaultMaxGrowingIter	int, default=10	The default number of allowed residual growths.
extrapolation	str, default=’linear’	The prediction a new increment starts its Newton iterations from (off|linear|quadratic): none, the last accepted increment scaled to the new size, or the parabola through the last three accepted states (linear while only one accepted increment of the step of nonzero size is available).
retryWithoutExtrapolation	bool, default=False	If an extrapolated increment fails, attempt it once more from a zero prediction before cutting it back. A failed increment then also keeps the accepted history the next attempt extrapolates from.
extrapolateAcrossConnectivityChange	bool, default=False	Keep extrapolating across a rebuild of the equation system caused by a change of the connectivity alone (e.g., a contact search), as long as the rebuilt system numbers the degrees of freedom as before. Off: every rebuild starts the increment from a zero prediction.
extrapolateAfterModelChange	bool, default=True	Whether to extrapolate the predictor on the increment FOLLOWING a model change (adaptive mesh refinement). Set False to start that increment from a zero predictor, avoiding extrapolation of the one-off warm-start/remesh settling transient.
equilibrateAfterModelChange	bool, default=False	Whether to insert one constant-load, zero-time re-equilibration increment immediately after an adaptive mesh refinement, before advancing the load. When True, the warm-started refined mesh is first settled to equilibrium at the last converged load level (no load advance, no Dirichlet increment, zero time increment) so the subsequent load-advancing increment starts from an equilibrated state. Intended for softening problems where remeshing near the process zone otherwise couples the load advance with the warm-start settling transient in one solve. Under arc length control (NISTPArcLength with a controller), the re-equilibration is an arc length increment of size zero: the controlled quantity is held where the last accepted increment left it, and the load parameter is solved for, and with it the loads and the prescribed displacements it scales. Note: the equilibration solve integrates materials with dT=0, which suits rate-independent models; rate-dependent materials see no time advance during it (by design).
linsolver	str, default=’pardiso’	The linear solver: the name of a *linearSolver block, or a type of linear solver with its defaults (e.g. pardiso, superlu, blockamg).
pruneCondensedMatrixZeros	bool, default=None	Compact explicitly stored zeros out of the multi-point-constraint-condensed system matrix before solving. Pruning removes whichever entries happen to be exactly zero on this particular iteration, so the pattern changes on almost every Newton iteration – which defeats any linear solver that relies on a stable pattern across iterations (e.g. blockamg reusing a per-field AMG hierarchy built for the old one). Left unset (the default), the choice is automatic: pruned if the linear solver in force reports that it needs it (:attr:`~edelweissfe.linsolve.base.LinearSolver.needsPrunedSystemMatrix`), which a direct factorization does – PARDISO’s reordering is sensitive to the extra structural entries on these path-dependent condensed systems, and keeping them has been observed to drift from the converged reference path – and a Krylov solver with an AMG preconditioner does not: it gives the same iterations with and without the stored zeros, and pruning costs it a copy of the matrix and a pass over it in every Newton iteration. An explicit True/False here always wins over the automatic choice; turn pruning off for a direct solver only after verifying the load path against a reference run.
systemMatrixAssembly	str, default=’vij’	How the blocks of the elements, constraints and loads are summed into the system matrix (vij|direct): staged in a value array of the VIJ pattern’s length and gathered into CSR, or scattered into CSR right away, without that array. See the documentation page on the system matrix assembly.
useThreadedMPCCondensation	bool, default=False	Condense the multi-point-constraint system matrix, T^T K T + C, in place, with an OpenMP-threaded kernel, instead of with SciPy’s single-threaded CSR routines, which build a second matrix. The CSR pattern of K gains, once per equation system, the entries where a master collects the couplings of its slaves; the result does not depend on the number of threads. Explicit zeros (the slave columns among them) are pruned from a copy according to pruneCondensedMatrixZeros. Off by default.
report-performance-frequency	int, default=0	Print the cumulative performance table (totals since the start of the step) every this many increments, in addition to the one printed at the end of the step. 0 (default) disables the periodic report.
lineSearch	bool, default=False	Scale the Newton correction by the step length at which a parabola, fitted to the largest residual evaluated at the trial step lengths lineSearchStepLengths, is smallest (a quadratic line search).
lineSearchAfterIterations	int, default=5	The line search starts after this many iterations of an increment, at least 0: the first correction of an increment, which makes the prescribed increments, is never scaled.
lineSearchEveryIterations	int, default=2	From then on, the line search is done in every iteration whose number is a multiple of this.
lineSearchStepLengths	str, default=’0.8, 1.0, 1.2’	The trial step lengths of the line search, comma-separated, at least three.
fallBackToQuasiNewton	bool, default=False	Instead of cutting an increment back once the residual grew more often than maxGrowIter allows, continue it with the tangent of its first iteration (a quasi-Newton method); maxIter still applies.
convergenceCriterion	str, default=’legacy’	The convergence criterion of the Newton iterations, legacy (the original EdelweissFE one) or abaqus (the one of Abaqus/Standard); see edelweissfe.solvers.base.convergencecriteria.
iterationsForAlternativeTolerance	int, default=15	legacy criterion: I, the iteration from which on the flux residual is checked against fluxResidualToleranceAlternative instead of fluxResidualTolerance.
zeroFluxThreshold	float, default=None	legacy criterion: a field whose spatial average flux does not exceed this value carries no flux its residual could be compared with, so its flux residual counts as converged; needs minimumIterations >= 1, and refuses a smaller one. Off if not set.
zeroIncrementThreshold	float, default=None	legacy criterion: a field whose largest increment does not exceed this value has not moved, so there is no increment to relate its correction to, and its correction counts as converged. Off if not set.
minimumIterations	int, default=0	legacy criterion: no flux residual counts as converged before this iteration, so an increment takes at least this many Newton iterations; 0 accepts a prediction that already satisfies the criterion.
countFirstResidualAsGrowth	bool, default=True	Residual growths (the step option maxGrowIter): the residual history of a field starts from zero, so the first checked residual of an increment counts as a growth; if False, counting starts with the second.
countGrowthOfConvergedFields	bool, default=True	Residual growths (the step option maxGrowIter): count a growth of the flux residual of a field also in an iteration in which this field is converged and only other fields keep the iterations going.
abaqusResidualTolerance	float, default=0.005	abaqus criterion: R_n, the flux residual tolerance relative to the time average flux.
abaqusCorrectionTolerance	float, default=0.01	abaqus criterion: C_n, the correction tolerance relative to the largest increment of the field.
abaqusResidualToleranceAlternative	float, default=0.02	abaqus criterion: R_p, the flux residual tolerance used from iteration I_p on.
abaqusIterationsForAlternativeTolerance	int, default=9	abaqus criterion: I_p, the iteration from which on R_p is used.
abaqusZeroFluxThreshold	float, default=1e-05	abaqus criterion: epsilon, below which (relative to the time average flux) a field carries no flux.
abaqusZeroFluxCorrectionTolerance	float, default=0.001	abaqus criterion: C_epsilon, the relative correction tolerance of a field without flux.
abaqusLinearResidualTolerance	float, default=1e-08	abaqus criterion: R_l, below which a field is converged irrespective of its correction.
abaqusInitialTimeAveragedFlux	float, default=0.01	abaqus criterion: q~0, the time average flux of a field that has not carried flux yet (first step).

NISTParallel - Nonlinear Implicit Static (parallel)

Parallel implementation of the NIST solver.

class edelweissfe.solvers.nonlinearimplicitstaticparallel.NISTParallel(jobInfo, journal, **kwargs)

NIST, with the elements computed on the available threads (computeElements()).

Parameters:

• jobInfo – A dictionary containing the job information.

• journal – The journal instance for logging.

Attributes:

convergenceTable

equationSystemLayout

incrementStatus

mpcTransformation

Methods

acceptIncrement(step, model, timeStep)	Commit the converged increment to the model, and keep the predictor's state; see acceptIncrement().
applyCorrection(dU, ddU, stepLength)	Move the iterate by the correction, scaled by the step length \(\alpha\) a line search chose (1 without a line search): \(\Delta \boldsymbol{U} \leftarrow \Delta \boldsymbol{U} + \alpha \, \Delta\Delta \boldsymbol{U}\).
applyDirichletToResidual(timeStep, R, dirichlets)	Impose the Dirichlet BCs on the residual.
applyDirichletToStiffness(K, dirichlets[, rhs])	Impose the Dirichlet BCs on the system matrix -- and, given the right-hand side(s) rhs (modified in place), eliminate the constrained DOFs from the columns as well; see edelweissfe.solvers.base.dirichlet.applyDirichletToStiffness().
applyOptionsOverride(fieldValues)	Apply an override of this solver's options, see applyOptionsOverride(); an override that sets an option where it means nothing is refused before anything of it applies (see refuseOptionsThatMeanNothing()).
applyStepActionsAtIncrementStart(model, ...)	Called when all step actions should be applied at the start of a step.
applyStepActionsAtStepEnd(model, stepActions)	Called when all step actions should finish a step.
applyStepActionsAtStepStart(model, step)	Apply every step action's step-start part (initial conditions, material initialization, prescribed fields, ...).
assembleAdditionalTerms(dU, R, F, K, timeStep)	Assemble terms of the equation system beyond the internal and external forces.
assembleConstraintForces(constraints, U_np, ...)	Evaluate every active constraint and add its residual contribution to the given vector; a constraint that is switched off (see active) is skipped.
assembleConstraints(constraints, U_np, dU, ...)	Loop over all constraints, and evaluate the active ones: each adds its residual contribution and the tangent of it, see computeKernels().
assembleInternalForces(U_np, dU, P, K, F, ...)	Assemble the internal forces of the model at a trial state, and their tangent: those of the elements of its equation system, see computeElements().
assembleLoads(loads, U_np, dU, P, K, timeStep)	Assemble the residual contributions of all loads: each adds the negative of the external load it applies, and the tangent of that contribution, which is nonzero only for a load that follows the deformation; see The residual and the tangent.
assembleLoadsAtIncrementEnd(stepActions, ...)	Assemble the residual contributions of the loads acting at the end of the increment, and their tangents; called by assembleResidualAndTangent().
assembleResidualAndTangent(U_n, dU, U_np, R, ...)	Evaluate the model at the trial solution \(U_{n+1} = U_n + \Delta U\): the right-hand side \(R = f_\mathrm{ext} - f_\mathrm{int}\) (see The residual and the tangent) with the terms of assembleAdditionalTerms(), condensed by the multi-point constraints, and the tangent.
assembleStiffnessCSR(K)	The assembled system matrix in CSR format, see assembledMatrix().
assembleSystemMatrix(K, quasiNewton)	The system matrix of this Newton iteration, before the Dirichlet boundary conditions are imposed: the tangent in CSR layout, condensed by the multi-point constraints -- or, once the increment has fallen back to the quasi-Newton method, the tangent of its first iteration.
attemptIncrement(step, model, timeStep)	Solve the increment by Newton's method, rebuilding the equation system first if the topology update changed the model, and re-equilibrating a changed model first if so configured; see attemptIncrement().
beginStep(step, model, ...)	Set up the linear solver, and start the step; see beginStep().
buildMPCTransformation(model[, stepActions])	Collect the linear dependency records from all multi-point constraints of the model and assemble the master-slave condensation operator for the current equation system.
checkConvergence(R, ddU, dU, F, ...)	Check the convergence, individually for each field and for the scalar variables, based on the current flux residual and the correction, using the convergence criterion of the solver (see edelweissfe.solvers.base.convergencecriteria).
checkDivergingSolution(...)	Check if the iterative solution scheme is diverging: if the flux residual of a field grew, from one checked iteration to the next, more than maxGrowingIter times in this increment.
checkDivergingSolutionOrFallBack(...)	Cut an increment back whose residual grew more often than maxGrowingIter allows (see checkDivergingSolution()) -- or, with the option fallBackToQuasiNewton, let it continue with the tangent of its first iteration (a quasi-Newton method), see assembleSystemMatrix().
checkMPCDirichletConflicts(transformation, ...)	Raise if any Dirichlet boundary condition of the step prescribes a DOF that is a slave DOF of a multi-point constraint.
collectLoads(stepActions)	The loads among a step's actions: every action that is a LoadBase, whatever step action module it is defined in, so that a new kind of load needs no change to any solver.
computeElements(elements, Un1, dU, P, K, F, ...)	Compute the elements on the available threads, see computeElementsInParallel(); otherwise as computeElements().
createConvergenceCriterion()	Create the convergence criterion selected by the option convergenceCriterion.
createSystemMatrixAssembly(assemblyType)	The assembly of the system matrix of the equation system just built (option systemMatrixAssembly), on the complete CSR pattern of csrGenerator, with numberOfAssemblyWorkers() workers.
endStep(step, model)	Report the step's performance timing.
extrapolateLastIncrement(extrapolation, ...)	Predict the solution increment, the point Newton's method starts from, by extrapolating the last accepted increments of the step.
extrapolates(extrapolation, prevTimeStep)	Whether an increment that follows the accepted increment prevTimeStep starts from an extrapolated prediction.
extrapolationWeights(extrapolation, ...)	The weights of the last accepted increment and of the one before it in the prediction of extrapolateLastIncrement(): \(\Delta U = w_n \, \Delta U_n + w_{n-1} \, \Delta U_{n-1}\), with
finalizeIncrement(model)	Store what the converged increment leaves behind.
firstIterateNeedsCorrection(dirichlets)	Whether the first iterate of an increment without a prediction -- the converged state of the last increment -- needs a correction before it can be checked for convergence (see solveIncrement()): if Dirichlet boundary conditions act, whose prescribed increment the first correction makes.
forEachChunk(entities, work)	Call work(chunk) for chunks of entities on the available threads, as the elements are computed (see mapChunksInParallel()); otherwise as forEachChunk().
forgetEquationSystem()	Drop the equation system of the last build -- the model's layout of it with the DOF manager, the system matrix assembly and its CSR pattern, the multi-point-constraint condensation, the DOF vectors -- so that it is gone before the next one is built, never both at once.
forgetPredictionAfterFailure()	A failed increment drops the prediction, which may have caused the failure, unless retryWithoutExtrapolation has answered that already: then the accepted history is kept, and the smaller increment that follows the cutback extrapolates it.
getRestartData()	The state this solver carries from one increment to the next, as it is; see checkpointedState.
initializeIncrement(U_n, stepActions, model, ...)	Prepare what assembleAdditionalTerms() needs during this increment's Newton iterations.
isOutputIncrement(timeStep)	Whether field outputs, output managers and restart checkpoints are written after this increment.
lineSearchIsDue(iterationCounter)	Whether the correction of this iteration is scaled by a line search (lineSearch): after lineSearchAfterIterations iterations, in every iteration whose number is a multiple of lineSearchEveryIterations.
lineSearchStepLengths()	The trial step lengths of the line search, from the option lineSearchStepLengths.
linearSolve(A, b)	Solve the linear equation system.
locateConstrainedDofs(dirichlets)	Determine, up front, which global DOFs each Dirichlet BC constrains.
numberOfAssemblyWorkers()	How many threads write blocks of the system matrix at the same time: those that compute the elements.
predictionSurvivesRebuild(modelHasChanged, ...)	Whether this increment may extrapolate the last accepted increments after the equation system was rebuilt.
prepareIncrement(step, model, isRetry)	The topology update, at the start of every increment; see prepareIncrement().
printResidualOutlierNodes(residualOutliers)	Print which nodes have the largest residuals.
quadraticLineSearch(largestResidualAt)	The step length \(\alpha\) by which a Newton correction is scaled, for a Newton iteration that overshoots or oscillates (e.g., materials switching between loading and unloading).
readStateFromModel(U)	Read the state the next increment starts from out of the model, into the vectors of the equation system just (re)built: here, the solution, through the model's layout of the equation system (see readStateFromModel()).
refuseOptionsThatMeanNothing(options, ...)	Refuse an option that a declaration of options sets where it means nothing, with the options in force after the declaration.
schema
selectLinearSolverForStep()	Make linSolver the linear solver the step about to begin names (linsolver): the one of the previous step if it names the same, else a new one.
setRestartData(data)	Restore the state getRestartData() returned, so that the resumed step continues exactly where the uninterrupted run would have.
setTrialState(dU, ddU, stepLength)	Set the trial state of a line search at a step length \(\alpha\) along the correction, and return its solution increment, \(\Delta \boldsymbol{U} + \alpha \, \Delta\Delta \boldsymbol{U}\), at which the model is evaluated next; the iterate stays where it is.
solveForCorrection(K_, R, dirichlets)	Solve the linearized system of a Newton iteration, \(\boldsymbol{K} \, \Delta\Delta \boldsymbol{U} = \boldsymbol{R}\), for the correction \(\Delta\Delta \boldsymbol{U}\).
solveIncrement(U_n, dU, P, K, stepActions, ...)	Newton's method for the solution increment of an increment.
updateRigidBodies(model, timeStep)	Refresh the kinematics of all rigid bodies in the model after a converged increment.
updateTopologyAndConnectivity(model, step[, ...])	Run the topology update, then let every mesh-dependent consumer catch up on it.
validateModelCapabilities(model)	Validate whether the solver supports the active features/constraints of the model.
writeResultsToModel(U, dU, P[, V, A, ...])	Write the results of a converged increment to the model, through the model's layout of the equation system: the incremental degrees of freedom are set to zero first, from where the next increment starts them (see zeroIncrementalCoefficients()), then the solution, its increment and the internal force are written (see writeResultsToModel()).

beginStep(step, model, fieldOutputController, outputmanagers)

Set up the linear solver, and start the step; see beginStep().

Parameters:

• step – The step to be solved.

• model – The model tree.

• fieldOutputController – The field output controller.

• outputmanagers – The output managers.

computeElements(elements, Un1, dU, P, K, F, timeStep)

Compute the elements on the available threads, see computeElementsInParallel(); otherwise as computeElements().

forEachChunk(entities, work)

Call work(chunk) for chunks of entities on the available threads, as the elements are computed (see mapChunksInParallel()); otherwise as forEachChunk().

Parameters:

entities (list)

numberOfAssemblyWorkers()

How many threads write blocks of the system matrix at the same time: those that compute the elements.

Returns:

The number of workers.

Return type:

int

NISTParallel

Option	Type / Default	Description
defaultMaxIter	int, default=10	The default maximum number of iterations.
defaultCriticalIter	int, default=5	The default number of critical iterations.
defaultMaxGrowingIter	int, default=10	The default number of allowed residual growths.
extrapolation	str, default=’linear’	The prediction a new increment starts its Newton iterations from (off|linear|quadratic): none, the last accepted increment scaled to the new size, or the parabola through the last three accepted states (linear while only one accepted increment of the step of nonzero size is available).
retryWithoutExtrapolation	bool, default=False	If an extrapolated increment fails, attempt it once more from a zero prediction before cutting it back. A failed increment then also keeps the accepted history the next attempt extrapolates from.
extrapolateAcrossConnectivityChange	bool, default=False	Keep extrapolating across a rebuild of the equation system caused by a change of the connectivity alone (e.g., a contact search), as long as the rebuilt system numbers the degrees of freedom as before. Off: every rebuild starts the increment from a zero prediction.
extrapolateAfterModelChange	bool, default=True	Whether to extrapolate the predictor on the increment FOLLOWING a model change (adaptive mesh refinement). Set False to start that increment from a zero predictor, avoiding extrapolation of the one-off warm-start/remesh settling transient.
equilibrateAfterModelChange	bool, default=False	Whether to insert one constant-load, zero-time re-equilibration increment immediately after an adaptive mesh refinement, before advancing the load. When True, the warm-started refined mesh is first settled to equilibrium at the last converged load level (no load advance, no Dirichlet increment, zero time increment) so the subsequent load-advancing increment starts from an equilibrated state. Intended for softening problems where remeshing near the process zone otherwise couples the load advance with the warm-start settling transient in one solve. Under arc length control (NISTPArcLength with a controller), the re-equilibration is an arc length increment of size zero: the controlled quantity is held where the last accepted increment left it, and the load parameter is solved for, and with it the loads and the prescribed displacements it scales. Note: the equilibration solve integrates materials with dT=0, which suits rate-independent models; rate-dependent materials see no time advance during it (by design).
linsolver	str, default=’pardiso’	The linear solver: the name of a *linearSolver block, or a type of linear solver with its defaults (e.g. pardiso, superlu, blockamg).
pruneCondensedMatrixZeros	bool, default=None	Compact explicitly stored zeros out of the multi-point-constraint-condensed system matrix before solving. Pruning removes whichever entries happen to be exactly zero on this particular iteration, so the pattern changes on almost every Newton iteration – which defeats any linear solver that relies on a stable pattern across iterations (e.g. blockamg reusing a per-field AMG hierarchy built for the old one). Left unset (the default), the choice is automatic: pruned if the linear solver in force reports that it needs it (:attr:`~edelweissfe.linsolve.base.LinearSolver.needsPrunedSystemMatrix`), which a direct factorization does – PARDISO’s reordering is sensitive to the extra structural entries on these path-dependent condensed systems, and keeping them has been observed to drift from the converged reference path – and a Krylov solver with an AMG preconditioner does not: it gives the same iterations with and without the stored zeros, and pruning costs it a copy of the matrix and a pass over it in every Newton iteration. An explicit True/False here always wins over the automatic choice; turn pruning off for a direct solver only after verifying the load path against a reference run.
systemMatrixAssembly	str, default=’vij’	How the blocks of the elements, constraints and loads are summed into the system matrix (vij|direct): staged in a value array of the VIJ pattern’s length and gathered into CSR, or scattered into CSR right away, without that array. See the documentation page on the system matrix assembly.
useThreadedMPCCondensation	bool, default=False	Condense the multi-point-constraint system matrix, T^T K T + C, in place, with an OpenMP-threaded kernel, instead of with SciPy’s single-threaded CSR routines, which build a second matrix. The CSR pattern of K gains, once per equation system, the entries where a master collects the couplings of its slaves; the result does not depend on the number of threads. Explicit zeros (the slave columns among them) are pruned from a copy according to pruneCondensedMatrixZeros. Off by default.
report-performance-frequency	int, default=0	Print the cumulative performance table (totals since the start of the step) every this many increments, in addition to the one printed at the end of the step. 0 (default) disables the periodic report.
lineSearch	bool, default=False	Scale the Newton correction by the step length at which a parabola, fitted to the largest residual evaluated at the trial step lengths lineSearchStepLengths, is smallest (a quadratic line search).
lineSearchAfterIterations	int, default=5	The line search starts after this many iterations of an increment, at least 0: the first correction of an increment, which makes the prescribed increments, is never scaled.
lineSearchEveryIterations	int, default=2	From then on, the line search is done in every iteration whose number is a multiple of this.
lineSearchStepLengths	str, default=’0.8, 1.0, 1.2’	The trial step lengths of the line search, comma-separated, at least three.
fallBackToQuasiNewton	bool, default=False	Instead of cutting an increment back once the residual grew more often than maxGrowIter allows, continue it with the tangent of its first iteration (a quasi-Newton method); maxIter still applies.
convergenceCriterion	str, default=’legacy’	The convergence criterion of the Newton iterations, legacy (the original EdelweissFE one) or abaqus (the one of Abaqus/Standard); see edelweissfe.solvers.base.convergencecriteria.
iterationsForAlternativeTolerance	int, default=15	legacy criterion: I, the iteration from which on the flux residual is checked against fluxResidualToleranceAlternative instead of fluxResidualTolerance.
zeroFluxThreshold	float, default=None	legacy criterion: a field whose spatial average flux does not exceed this value carries no flux its residual could be compared with, so its flux residual counts as converged; needs minimumIterations >= 1, and refuses a smaller one. Off if not set.
zeroIncrementThreshold	float, default=None	legacy criterion: a field whose largest increment does not exceed this value has not moved, so there is no increment to relate its correction to, and its correction counts as converged. Off if not set.
minimumIterations	int, default=0	legacy criterion: no flux residual counts as converged before this iteration, so an increment takes at least this many Newton iterations; 0 accepts a prediction that already satisfies the criterion.
countFirstResidualAsGrowth	bool, default=True	Residual growths (the step option maxGrowIter): the residual history of a field starts from zero, so the first checked residual of an increment counts as a growth; if False, counting starts with the second.
countGrowthOfConvergedFields	bool, default=True	Residual growths (the step option maxGrowIter): count a growth of the flux residual of a field also in an iteration in which this field is converged and only other fields keep the iterations going.
abaqusResidualTolerance	float, default=0.005	abaqus criterion: R_n, the flux residual tolerance relative to the time average flux.
abaqusCorrectionTolerance	float, default=0.01	abaqus criterion: C_n, the correction tolerance relative to the largest increment of the field.
abaqusResidualToleranceAlternative	float, default=0.02	abaqus criterion: R_p, the flux residual tolerance used from iteration I_p on.
abaqusIterationsForAlternativeTolerance	int, default=9	abaqus criterion: I_p, the iteration from which on R_p is used.
abaqusZeroFluxThreshold	float, default=1e-05	abaqus criterion: epsilon, below which (relative to the time average flux) a field carries no flux.
abaqusZeroFluxCorrectionTolerance	float, default=0.001	abaqus criterion: C_epsilon, the relative correction tolerance of a field without flux.
abaqusLinearResidualTolerance	float, default=1e-08	abaqus criterion: R_l, below which a field is converged irrespective of its correction.
abaqusInitialTimeAveragedFlux	float, default=0.01	abaqus criterion: q~0, the time average flux of a field that has not carried flux yet (first step).

NISTPArcLength - Nonlinear Implicit Static - Arc length

(Parallel) Arc length solver, based on the approach of Jirásek/Bažant 2001: the Newton iterations of the NISTParallel solver, with the load parameter as an unknown of every increment besides the solution.

class edelweissfe.solvers.nonlinearimplicitstaticparallelarclength.ArcLengthIncrement(loadParameterAtStart, timeStep, referenceTimeStep, deadLoad, referenceLoad, rightHandSides, referenceCorrection=None, loadParameterIncrement=0.0, loadParameter=0.0)

What the arc length method adds to the Newton iterations of an increment: the load parameter \(\lambda\), an unknown besides the solution, and the external load it scales, \(\boldsymbol{f}_\mathrm{ext} = \boldsymbol{f}_0 + \lambda \, \boldsymbol{f}_\mathrm{ref}\), with the residual contributions of the dead load (the loads at the start of the step, \(-\boldsymbol{f}_0\)) and of the reference load (their change over the step, \(-\boldsymbol{f}_\mathrm{ref}\)). Created for every attempt of an increment by NISTPArcLength.solveIncrement(); the loads are assembled by NISTPArcLength.assembleLoadsAtIncrementEnd(), which adds their tangents to the tangent of the increment.

Attributes:

referenceCorrection

Parameters:

• loadParameterAtStart (float)

• timeStep (TimeStep)

• referenceTimeStep (TimeStep)

• deadLoad (DofVector)

• referenceLoad (DofVector)

• rightHandSides (ndarray)

• referenceCorrection (ndarray | None)

• loadParameterIncrement (float)

• loadParameter (float)

deadLoad: DofVector

The residual contribution of the dead load.

loadParameter: float = 0.0

of the iterate, or of a trial state of a line search.

Type:

The load parameter of the state being evaluated

loadParameterAtStart: float

The load parameter at the start of the increment, the one the last accepted increment reached.

loadParameterIncrement: float = 0.0

The increment of the load parameter of the iterate, \(\\Delta\\lambda\).

referenceCorrection: ndarray | None = None

The correction of the solution \(\\boldsymbol{R}_f\) causes, \(\\Delta\\Delta \\boldsymbol{U}_f\), of the last linear solve.

referenceLoad: DofVector

The residual contribution of the reference load.

referenceTimeStep: TimeStep

at it, the Dirichlet boundary conditions give the change of their prescribed values over the step, which the load parameter scales like the reference load.

Type:

The whole step as one time step

rightHandSides: ndarray

the out-of-balance force \(\\boldsymbol{R}_0\), of the dead load and of the reference load scaled by the load parameter of the iterate, and the reference load alone, \(\\boldsymbol{R}_f\).

Type:

The two right-hand sides of the linearized system, solved at once as its columns

timeStep: TimeStep

The increment being solved.

class edelweissfe.solvers.nonlinearimplicitstaticparallelarclength.NISTPArcLength(jobInfo, journal, **kwargs)

The arc length method (Jirásek/Bažant 2001): a step’s reference load is scaled by a load parameter \(\lambda\), which is an unknown of every increment besides the solution, and which a controller step action determines, e.g. such that a crack opening grows steadily. Paths with snap-back or snap-through, which neither a prescribed load nor a prescribed displacement can follow, are traced this way.

The controller is the step action of the step that the option arcLengthController names (see ArcLengthControllerBase); without one (arcLengthController=off), the step is solved as by NISTParallel.

The Newton iterations are those of solveIncrement(); the load parameter enters them through its hooks: the prediction (extrapolateLastIncrement()), the loads of the evaluated state (assembleLoadsAtIncrementEnd()), the first iterate (firstIterateNeedsCorrection()), the linear solve for two right-hand sides (solveForCorrection()), and the corrections of the pair of the solution and the load parameter (setTrialState(), applyCorrection()).

Every arc length step scales its own reference load – the change of its loads and prescribed displacements over the step – from \(\lambda = 0\), so it starts from the loads the steps before it ended with; at its end, every load keeps the value the step applied (see applyStepActionsAtStepEnd()). The load parameter of the step being solved is written to the model at the acceptance of every increment, as model.additionalParameters["arc length parameter"], e.g. for a load-displacement curve.

The load parameter is the solver’s state, as the solution is the model’s: it advances only once an increment has converged, together with the solution it belongs to (at the end of solveIncrement(), which a failing attempt leaves by an exception before). A cut back increment therefore starts again from the load parameter of the last accepted one.

The re-equilibration after a model change (option equilibrateAfterModelChange) is an arc length increment of size zero, solved by solveIncrement() like any other: the controller’s equation asks for a change of zero, so the controlled quantity stays where the last accepted increment left it, and the load parameter is solved for, and with it the loads and the prescribed degrees of freedom it scales. This is what a zero-time re-equilibration holds under displacement control, and, unlike a constant load, it has a solution past a limit point. As for any increment, the stop condition is checked before it, and the controller’s finishIncrement is called once it has converged; the load parameter it reaches, committed with its state, is the one the next increment starts from.

Parameters:

• jobInfo – A dictionary containing the job information.

• journal – The journal instance for logging.

Attributes:

convergenceTable

equationSystemLayout

incrementStatus

mpcTransformation

Methods

acceptIncrement(step, model, timeStep)	Accept the increment as acceptIncrement() does, write the load parameter to the model with the results, and keep the increment of the load parameter for the predictor alongside the solution increment.
applyCorrection(dU, ddU_0, stepLength)	Move the iterate, the pair of the solution increment and the increment of the load parameter, by their corrections at the step length a line search chose (1 without one), see correctionAtStepLength().
applyDirichletToResidual(timeStep, R, dirichlets)	Impose the Dirichlet BCs on the residual.
applyDirichletToStiffness(K, dirichlets[, rhs])	Impose the Dirichlet BCs on the system matrix -- and, given the right-hand side(s) rhs (modified in place), eliminate the constrained DOFs from the columns as well; see edelweissfe.solvers.base.dirichlet.applyDirichletToStiffness().
applyOptionsOverride(fieldValues)	Apply an override of this solver's options, see applyOptionsOverride().
applyStepActionsAtIncrementStart(model, ...)	Called when all step actions should be applied at the start of a step.
applyStepActionsAtStepEnd(model, stepActions)	Finish the step actions at the end of a step.
applyStepActionsAtStepStart(model, step)	Apply every step action's step-start part (initial conditions, material initialization, prescribed fields, ...).
arcLengthControllerOfStep(step)	The arc length controller the option arcLengthController names among the step actions of the step.
assembleAdditionalTerms(dU, R, F, K, timeStep)	Assemble terms of the equation system beyond the internal and external forces.
assembleConstraintForces(constraints, U_np, ...)	Evaluate every active constraint and add its residual contribution to the given vector; a constraint that is switched off (see active) is skipped.
assembleConstraints(constraints, U_np, dU, ...)	Loop over all constraints, and evaluate the active ones: each adds its residual contribution and the tangent of it, see computeKernels().
assembleInternalForces(U_np, dU, P, K, F, ...)	Assemble the internal forces of the model at a trial state, and their tangent: those of the elements of its equation system, see computeElements().
assembleLoads(loads, U_np, dU, P, K, timeStep)	Assemble the residual contributions of all loads: each adds the negative of the external load it applies, and the tangent of that contribution, which is nonzero only for a load that follows the deformation; see The residual and the tangent.
assembleLoadsAtIncrementEnd(stepActions, ...)	Assemble the loads at the end of the increment: under arc length control, the dead load plus the reference load scaled by the load parameter of the evaluated state, \(-(\boldsymbol{f}_0 + \lambda \, \boldsymbol{f}_\mathrm{ref})\), with their tangents.
assembleResidualAndTangent(U_n, dU, U_np, R, ...)	Evaluate the model at the trial solution \(U_{n+1} = U_n + \Delta U\): the right-hand side \(R = f_\mathrm{ext} - f_\mathrm{int}\) (see The residual and the tangent) with the terms of assembleAdditionalTerms(), condensed by the multi-point constraints, and the tangent.
assembleStiffnessCSR(K)	The assembled system matrix in CSR format, see assembledMatrix().
assembleSystemMatrix(K, quasiNewton)	The system matrix of this Newton iteration, before the Dirichlet boundary conditions are imposed: the tangent in CSR layout, condensed by the multi-point constraints -- or, once the increment has fallen back to the quasi-Newton method, the tangent of its first iteration.
attemptIncrement(step, model, timeStep)	Solve the increment by Newton's method, rebuilding the equation system first if the topology update changed the model, and re-equilibrating a changed model first if so configured; see attemptIncrement().
beginStep(step, model, ...)	Take the controller of the step, and start the step at the load parameter zero; see beginStep().
buildMPCTransformation(model[, stepActions])	Collect the linear dependency records from all multi-point constraints of the model and assemble the master-slave condensation operator for the current equation system.
checkConvergence(R, ddU, dU, F, ...)	Check the convergence, individually for each field and for the scalar variables, based on the current flux residual and the correction, using the convergence criterion of the solver (see edelweissfe.solvers.base.convergencecriteria).
checkDivergingSolution(...)	Check if the iterative solution scheme is diverging: if the flux residual of a field grew, from one checked iteration to the next, more than maxGrowingIter times in this increment.
checkDivergingSolutionOrFallBack(...)	Cut an increment back whose residual grew more often than maxGrowingIter allows (see checkDivergingSolution()) -- or, with the option fallBackToQuasiNewton, let it continue with the tangent of its first iteration (a quasi-Newton method), see assembleSystemMatrix().
checkMPCDirichletConflicts(transformation, ...)	Raise if any Dirichlet boundary condition of the step prescribes a DOF that is a slave DOF of a multi-point constraint.
collectLoads(stepActions)	The loads among a step's actions: every action that is a LoadBase, whatever step action module it is defined in, so that a new kind of load needs no change to any solver.
computeElements(elements, Un1, dU, P, K, F, ...)	Compute the elements on the available threads, see computeElementsInParallel(); otherwise as computeElements().
correctionAtStepLength(dU, ddU_0, stepLength)	The correction of the pair of the solution and the load parameter at a step length \(\alpha\): the correction the out-of-balance force causes, scaled by \(\alpha\), which the controller completes with the correction the reference load causes, \(\Delta\Delta U = \alpha \, \Delta\Delta U_0 + \Delta\Delta\lambda \, \Delta\Delta U_f\), by the correction of the load parameter \(\Delta\Delta\lambda\) its equation asks for.
createConvergenceCriterion()	Create the convergence criterion selected by the option convergenceCriterion.
createSystemMatrixAssembly(assemblyType)	The assembly of the system matrix of the equation system just built (option systemMatrixAssembly), on the complete CSR pattern of csrGenerator, with numberOfAssemblyWorkers() workers.
endStep(step, model)	Report the step's performance timing.
extrapolateLastIncrement(extrapolation, ...)	Predict the pair of the solution increment and the increment of the load parameter: the solution increment as extrapolateLastIncrement() predicts it, the increment of the load parameter with the same weights (extrapolateLoadParameterIncrement()), and the prescribed degrees of freedom where the predicted load parameter puts them: their prescribed change over the step, scaled by it.
extrapolateLoadParameterIncrement(...)	Predict the increment of the load parameter, as the solution increment is predicted, see extrapolateLastIncrement(): from the same accepted increments, with the same weights (extrapolationWeights()), so that the prediction of the pair stays on one extrapolated path.
extrapolates(extrapolation, prevTimeStep)	Whether an increment that follows the accepted increment prevTimeStep starts from an extrapolated prediction.
extrapolationWeights(extrapolation, ...)	The weights of the last accepted increment and of the one before it in the prediction of extrapolateLastIncrement(): \(\Delta U = w_n \, \Delta U_n + w_{n-1} \, \Delta U_{n-1}\), with
finalizeIncrement(model)	Store what the converged increment leaves behind.
firstIterateNeedsCorrection(dirichlets)	Under arc length control, the first iterate of an increment without a prediction always needs a correction: it is the converged state of the last increment, at the load parameter that increment reached, and the first correction moves the load parameter as the controller prescribes.
forEachChunk(entities, work)	Call work(chunk) for chunks of entities on the available threads, as the elements are computed (see mapChunksInParallel()); otherwise as forEachChunk().
forgetEquationSystem()	Drop the equation system of the last build -- the model's layout of it with the DOF manager, the system matrix assembly and its CSR pattern, the multi-point-constraint condensation, the DOF vectors -- so that it is gone before the next one is built, never both at once.
forgetPredictionAfterFailure()	A failed increment drops the prediction, which may have caused the failure, unless retryWithoutExtrapolation has answered that already: then the accepted history is kept, and the smaller increment that follows the cutback extrapolates it.
getRestartData()	The state this solver carries from one increment to the next, as it is; see checkpointedState.
initializeIncrement(U_n, stepActions, model, ...)	Prepare what assembleAdditionalTerms() needs during this increment's Newton iterations.
isOutputIncrement(timeStep)	Whether field outputs, output managers and restart checkpoints are written after this increment.
lineSearchIsDue(iterationCounter)	Whether the correction of this iteration is scaled by a line search (lineSearch): after lineSearchAfterIterations iterations, in every iteration whose number is a multiple of lineSearchEveryIterations.
lineSearchStepLengths()	The trial step lengths of the line search, from the option lineSearchStepLengths.
linearSolve(A, b)	Solve the linear equation system.
locateConstrainedDofs(dirichlets)	Determine, up front, which global DOFs each Dirichlet BC constrains.
numberOfAssemblyWorkers()	How many threads write blocks of the system matrix at the same time: those that compute the elements.
predictionSurvivesRebuild(modelHasChanged, ...)	Whether this increment may extrapolate the last accepted increments after the equation system was rebuilt.
prepareIncrement(step, model, isRetry)	The topology update, at the start of every increment; see prepareIncrement().
printResidualOutlierNodes(residualOutliers)	Print which nodes have the largest residuals.
quadraticLineSearch(largestResidualAt)	The step length \(\alpha\) by which a Newton correction is scaled, for a Newton iteration that overshoots or oscillates (e.g., materials switching between loading and unloading).
readStateFromModel(U)	Read the state the next increment starts from out of the model, into the vectors of the equation system just (re)built: here, the solution, through the model's layout of the equation system (see readStateFromModel()).
refuseOptionsThatMeanNothing(options, ...)	Refuse an option that a declaration of options sets where it means nothing, with the options in force after the declaration.
schema	Option schema for this solver, per OptionSchemaProvider.
selectLinearSolverForStep()	Make linSolver the linear solver the step about to begin names (linsolver): the one of the previous step if it names the same, else a new one.
setRestartData(data)	Restore the state getRestartData() returned, so that the resumed step continues exactly where the uninterrupted run would have.
setTrialState(dU, ddU_0, stepLength)	Set the trial state of a line search at a step length along the correction of the pair: the solution increment and the load parameter, both moved by their corrections at that step length, which satisfy the controller's equation (see correctionAtStepLength()); the iterate stays where it is.
solveForCorrection(K_, R, dirichlets)	Solve the linearized system for two right-hand sides at once: the out-of-balance force \(R_0 = R\), and the reference load \(R_f\), for the corrections \(\Delta\Delta U_0\) and \(\Delta\Delta U_f\).
solveIncrement(U_n, dU, P, K, stepActions, ...)	Newton's method for the solution increment and the increment of the load parameter together: the iterations of solveIncrement(), into which the load parameter enters through the hooks of this class.
updateRigidBodies(model, timeStep)	Refresh the kinematics of all rigid bodies in the model after a converged increment.
updateTopologyAndConnectivity(model, step[, ...])	Run the topology update, then let every mesh-dependent consumer catch up on it.
validateModelCapabilities(model)	Validate whether the solver supports the active features/constraints of the model.
writeResultsToModel(U, dU, P[, V, A, ...])	Write the results of a converged increment to the model, through the model's layout of the equation system: the incremental degrees of freedom are set to zero first, from where the next increment starts them (see zeroIncrementalCoefficients()), then the solution, its increment and the internal force are written (see writeResultsToModel()).

schema

Option schema for this solver, per OptionSchemaProvider.

Attributes:

pruneCondensedMatrixZeros

stopCondition

zeroFluxThreshold

zeroIncrementThreshold

alias of NISTPArcLengthSchema

acceptIncrement(step, model, timeStep)

Accept the increment as acceptIncrement() does, write the load parameter to the model with the results, and keep the increment of the load parameter for the predictor alongside the solution increment.

Parameters:

• step – The step being solved.

• model (FEModel) – The model tree.

• timeStep (TimeStep) – The converged increment.

applyCorrection(dU, ddU_0, stepLength)

Move the iterate, the pair of the solution increment and the increment of the load parameter, by their corrections at the step length a line search chose (1 without one), see correctionAtStepLength(). Without a controller, as applyCorrection().

Parameters:

• dU (DofVector) – The solution increment of the iterate, moved in place.

• ddU_0 (DofVector) – The correction the out-of-balance force causes, see solveForCorrection().

• stepLength (float) – The step length.

Returns:

The correction the solution increment was moved by.

Return type:

DofVector

applyOptionsOverride(fieldValues)

Apply an override of this solver’s options, see applyOptionsOverride().

A stop condition belongs to the arc length control it is declared with: an override that names the controller (arcLengthController, off too) without a stopCondition returns to the stop condition of the definition of the solver (stopConditionOfSolverDefinition), and drops the one of an earlier step. Otherwise a step would stop by the condition that ended the step before it – which typically holds already – at its first check. An override that sets a stopCondition for a step without arc length control is refused before anything of it applies: it would never be checked. (The stop condition of the definition of the solver may stand without a controller: it applies to the steps that name one.)

Parameters:

fieldValues (dict) – Maps schema field name to its new, already-coerced value.

Raises:

ValueError – If the override sets a stop condition for a step without arc length control.

Return type:

None

applyStepActionsAtStepEnd(model, stepActions)

Finish the step actions at the end of a step. Under arc length control, every load keeps the value the step applied at the load parameter it ended at, \(f(0) + \lambda \, (f(1) - f(0))\) with its own amplitude (see applyAtStepEnd()), before the other step actions finish as in any step: the next step starts from these loads.

Parameters:

• model – The model tree.

• stepActions (dict[str, StepActionBase]) – The list of active step actions.

arcLengthControllerOfStep(step)

The arc length controller the option arcLengthController names among the step actions of the step.

Parameters:

step – The step to be solved.

Returns:

The controller, or None for arcLengthController=off.

Return type:

ArcLengthControllerBase | None

Raises:

ValueError – If the step has no arc length controller of this name, or more than one.

assembleLoadsAtIncrementEnd(stepActions, U_np, dU, P, K, timeStep)

Assemble the loads at the end of the increment: under arc length control, the dead load plus the reference load scaled by the load parameter of the evaluated state, \(-(\boldsymbol{f}_0 + \lambda \, \boldsymbol{f}_\mathrm{ref})\), with their tangents. The loads are evaluated at the start of the step (the dead load) and at its end, and their difference is the reference load (see ArcLengthIncrement). Without a controller, the loads of the step at the time of the increment, as in assembleLoadsAtIncrementEnd().

Parameters:

• stepActions (dict) – The active step actions.

• U_np (DofVector) – The trial solution.

• dU (DofVector) – The trial solution increment.

• P (DofVector) – The residual contribution vector to assemble into.

• K (SystemMatrixAssembly) – The tangent being assembled.

• timeStep (TimeStep) – The time step.

Returns:

The updated residual contribution vector and tangent.

Return type:

tuple[DofVector,SystemMatrixAssembly]

beginStep(step, model, fieldOutputController, outputmanagers)

Take the controller of the step, and start the step at the load parameter zero; see beginStep().

The load parameter scales the change of the loads over this step, so it starts from zero: the external load is then the dead load, the loads the steps before ended with. A resumed step takes over the checkpointed load parameter instead.

Parameters:

• step – The step to be solved.

• model (FEModel) – The model tree.

• fieldOutputController (FieldOutputController) – The field output controller.

• outputmanagers (dict[str, OutputManagerBase]) – The output managers.

correctionAtStepLength(dU, ddU_0, stepLength)

The correction of the pair of the solution and the load parameter at a step length \(\alpha\): the correction the out-of-balance force causes, scaled by \(\alpha\), which the controller completes with the correction the reference load causes, \(\Delta\Delta U = \alpha \, \Delta\Delta U_0 + \Delta\Delta\lambda \, \Delta\Delta U_f\), by the correction of the load parameter \(\Delta\Delta\lambda\) its equation asks for. So the controller’s equation holds at every step length, as for the full correction (\(\alpha = 1\)): a line search cannot leave it unsatisfied, also not from an iterate that does not satisfy it yet.

Parameters:

• dU (DofVector) – The solution increment of the iterate.

• ddU_0 (DofVector) – The correction the out-of-balance force causes, see solveForCorrection().

• stepLength (float) – The step length.

Returns:

• The correction of the solution.

• The correction of the load parameter.

Return type:

tuple[DofVector,float]

extrapolateLastIncrement(extrapolation, timeStep, dU, dirichlets, prevTimeStep, model)

Predict the pair of the solution increment and the increment of the load parameter: the solution increment as extrapolateLastIncrement() predicts it, the increment of the load parameter with the same weights (extrapolateLoadParameterIncrement()), and the prescribed degrees of freedom where the predicted load parameter puts them: their prescribed change over the step, scaled by it. Without a controller, as NIST.

Parameters:

• extrapolation (str) – The type of extrapolation.

• timeStep (TimeStep) – The new increment.

• dU (DofVector) – The solution increment to predict, overwritten.

• dirichlets (list) – The active dirichlet boundary conditions.

• prevTimeStep (TimeStep) – The last accepted increment, or None for no prediction.

• model (FEModel) – The model tree.

Returns:

• The predicted solution increment.

• True if an extrapolation was performed.

Return type:

tuple[DofVector,bool]

extrapolateLoadParameterIncrement(extrapolation, timeStep, prevTimeStep)

Predict the increment of the load parameter, as the solution increment is predicted, see extrapolateLastIncrement(): from the same accepted increments, with the same weights (extrapolationWeights()), so that the prediction of the pair stays on one extrapolated path.

Parameters:

• extrapolation (str) – The type of extrapolation.

• timeStep (TimeStep) – The new increment.

• prevTimeStep (TimeStep | None) – The last accepted increment, or None for no prediction.

Returns:

The predicted increment of the load parameter.

Return type:

float

firstIterateNeedsCorrection(dirichlets)

Under arc length control, the first iterate of an increment without a prediction always needs a correction: it is the converged state of the last increment, at the load parameter that increment reached, and the first correction moves the load parameter as the controller prescribes. Without a controller, as firstIterateNeedsCorrection().

Parameters:

dirichlets (list) – The Dirichlet boundary conditions of the step.

Returns:

True if the first iterate is corrected before the first convergence check.

Return type:

bool

setTrialState(dU, ddU_0, stepLength)

Set the trial state of a line search at a step length along the correction of the pair: the solution increment and the load parameter, both moved by their corrections at that step length, which satisfy the controller’s equation (see correctionAtStepLength()); the iterate stays where it is. Without a controller, as setTrialState().

Parameters:

• dU (DofVector) – The solution increment of the iterate.

• ddU_0 (DofVector) – The correction the out-of-balance force causes, see solveForCorrection().

• stepLength (float) – The step length.

Returns:

The solution increment of the trial state.

Return type:

DofVector

solveForCorrection(K_, R, dirichlets)

Solve the linearized system for two right-hand sides at once: the out-of-balance force \(R_0 = R\), and the reference load \(R_f\), for the corrections \(\Delta\Delta U_0\) and \(\Delta\Delta U_f\). Without a controller, as solveForCorrection().

The prescribed degrees of freedom move with the load parameter, and only with it: the reference solve carries their prescribed change over the step, so a correction \(\Delta\Delta\lambda\) moves them by \(\Delta\Delta\lambda\) times that change, as it scales the reference load. The out-of-balance solve leaves them where they are, at the place the load parameter of the iterate puts them, in every iteration – also in the first iteration of an increment without a prediction, where NIST’s would move them by their prescribed increment.

The controller completes \(\Delta\Delta U_0\) with \(\Delta\Delta U_f\) when the iterate is moved (setTrialState(), applyCorrection()).

Parameters:

• K – The system matrix of the iteration.

• R (DofVector) – The out-of-balance force.

• dirichlets (list) – The Dirichlet boundary conditions of the step.

• K_ (csr_matrix)

Returns:

The correction the out-of-balance force causes, \(\Delta\Delta U_0\).

Return type:

DofVector

solveIncrement(U_n, dU, P, K, stepActions, model, timeStep, prevTimeStep, extrapolation, maxIter, maxGrowingIter)

Newton’s method for the solution increment and the increment of the load parameter together: the iterations of solveIncrement(), into which the load parameter enters through the hooks of this class.

Each iteration evaluates the model at the iterate (with the loads of assembleLoadsAtIncrementEnd()), and solves the linearized system for two right-hand sides at once: the out-of-balance force \(R_0\), and the reference load \(R_f\) (solveForCorrection()). The controller combines the two corrections, \(\Delta\Delta U = \Delta\Delta U_0 + \Delta\Delta\lambda \, \Delta\Delta U_f\), and returns the correction of the load parameter \(\Delta\Delta\lambda\) (applyCorrection()). The prescribed degrees of freedom follow the load parameter, like the loads: the reference solve carries their prescribed change over the step, so every correction of the load parameter moves them along, while the out-of-balance solve leaves them where they are – they are no unknowns, and their rows are left out of the convergence check. The options of Newton’s method apply to the pair \((\Delta U, \Delta\lambda)\): the prediction extrapolates both with the same weights (extrapolateLastIncrement()), a line search scales the correction the out-of-balance force causes, which the controller completes at every step length, so that every trial state satisfies its equation (setTrialState()), and the quasi-Newton fallback keeps the tangent of the first iteration for both right-hand sides.

Before the iterations, the step’s stop condition is checked, if an increment of nonzero size has been accepted before (see the option stopCondition); once they have converged, the load parameter advances with the solution it belongs to.

Parameters:

• U_n (DofVector) – The solution at the start of the increment.

• dU (DofVector) – The solution increment, overwritten.

• P (DofVector) – The internal forces (the reaction), overwritten.

• K (SystemMatrixAssembly) – The tangent being assembled, assembled anew.

• stepActions (dict) – The active step actions.

• model (FEModel) – The model tree.

• timeStep (TimeStep) – The time step.

• prevTimeStep (TimeStep) – The last accepted increment, from which the prediction extrapolates, or None.

• extrapolation (str) – The type of extrapolation.

• maxIter (int) – The maximum number of iterations.

• maxGrowingIter (int) – The maximum number of residual growths before the increment is cut back.

Returns:

A tuple containing

• the new solution vector

• the solution increment

• the new reaction vector

• the number of required iterations

• the history of residuals per field

Return type:

tuple[DofVector,DofVector,DofVector,int,dict]

Lambda

The load parameter of the current step, reached by its last converged increment.

arcLengthController

The controller of the current step, or None.

arcLengthIncrement

What the arc length method adds to the iterations of the increment being solved; None between increments.

checkpointedState: dict | None = {'Lambda': <class 'float'>, 'convergenceCriterion': <class 'edelweissfe.solvers.base.convergencecriteria.ConvergenceCriterionBase'>, 'dLambdaOfIncrementBeforeLast': <class 'float'>, 'dLambdaOfLastIncrement': <class 'float'>, 'dUOfIncrementBeforeLast': <class 'numpy.ndarray'>, 'dUOfLastIncrement': <class 'numpy.ndarray'>, 'dofNumberingOfLastIncrement': <class 'str'>, 'prevTimeStep': <class 'edelweissfe.timesteppers.timestep.TimeStep'>, 'timeStepOfIncrementBeforeLast': <class 'edelweissfe.timesteppers.timestep.TimeStep'>}

the load parameter the step has reached, and the predictor’s state, NIST’s plus the increments of the load parameter that belong to its accepted increments.

Type:

The state between increments

dLambdaOfIncrementBeforeLast

The increment of the load parameter in the accepted increment before, for a quadratic extrapolation.

dLambdaOfLastIncrement

The increment of the load parameter in the last accepted increment, from which the predictor extrapolates.

stopConditionOfSolverDefinition

The stop condition of the definition of the solver (its *solver datalines or constructor arguments), or None: the one of every step whose >>options name the controller without a stop condition of their own.

supportsMPC = False

Whether this solver supports master-slave condensation / multi-point constraints (e.g. surface ties). Subclasses supporting MPCs must set this to True.

class edelweissfe.solvers.nonlinearimplicitstaticparallelarclength.NISTPArcLengthSchema(defaultMaxIter=10, defaultCriticalIter=5, defaultMaxGrowingIter=10, extrapolation='linear', retryWithoutExtrapolation=False, extrapolateAcrossConnectivityChange=False, extrapolateAfterModelChange=True, equilibrateAfterModelChange=False, linsolver='pardiso', pruneCondensedMatrixZeros=None, systemMatrixAssembly='vij', useThreadedMPCCondensation=False, reportPerformanceFrequency=0, lineSearch=False, lineSearchAfterIterations=5, lineSearchEveryIterations=2, lineSearchStepLengths='0.8, 1.0, 1.2', fallBackToQuasiNewton=False, convergenceCriterion='legacy', iterationsForAlternativeTolerance=15, zeroFluxThreshold=None, zeroIncrementThreshold=None, minimumIterations=0, countFirstResidualAsGrowth=True, countGrowthOfConvergedFields=True, abaqusResidualTolerance=0.005, abaqusCorrectionTolerance=0.01, abaqusResidualToleranceAlternative=0.02, abaqusIterationsForAlternativeTolerance=9, abaqusZeroFluxThreshold=1e-05, abaqusZeroFluxCorrectionTolerance=0.001, abaqusLinearResidualTolerance=1e-08, abaqusInitialTimeAveragedFlux=0.01, arcLengthController='off', stopCondition=None)

NISTSchema’s options, plus the two of the arc length method; like every option of this solver, they are *solver datalines, >>options and constructor keyword arguments alike.

Attributes:

pruneCondensedMatrixZeros

stopCondition

zeroFluxThreshold

zeroIncrementThreshold

Parameters:

• defaultMaxIter (int | None)

• defaultCriticalIter (int | None)

• defaultMaxGrowingIter (int | None)

• extrapolation (str | None)

• retryWithoutExtrapolation (bool | None)

• extrapolateAcrossConnectivityChange (bool | None)

• extrapolateAfterModelChange (bool | None)

• equilibrateAfterModelChange (bool | None)

• linsolver (str | None)

• pruneCondensedMatrixZeros (bool | None)

• systemMatrixAssembly (str | None)

• useThreadedMPCCondensation (bool | None)

• reportPerformanceFrequency (int | None)

• lineSearch (bool | None)

• lineSearchAfterIterations (int | None)

• lineSearchEveryIterations (int | None)

• lineSearchStepLengths (str | None)

• fallBackToQuasiNewton (bool | None)

• convergenceCriterion (str | None)

• iterationsForAlternativeTolerance (int | None)

• zeroFluxThreshold (float | None)

• zeroIncrementThreshold (float | None)

• minimumIterations (int | None)

• countFirstResidualAsGrowth (bool | None)

• countGrowthOfConvergedFields (bool | None)

• abaqusResidualTolerance (float | None)

• abaqusCorrectionTolerance (float | None)

• abaqusResidualToleranceAlternative (float | None)

• abaqusIterationsForAlternativeTolerance (int | None)

• abaqusZeroFluxThreshold (float | None)

• abaqusZeroFluxCorrectionTolerance (float | None)

• abaqusLinearResidualTolerance (float | None)

• abaqusInitialTimeAveragedFlux (float | None)

• arcLengthController (str | None)

• stopCondition (str | None)

NISTPArcLength

Option	Type / Default	Description
defaultMaxIter	int, default=10	The default maximum number of iterations.
defaultCriticalIter	int, default=5	The default number of critical iterations.
defaultMaxGrowingIter	int, default=10	The default number of allowed residual growths.
extrapolation	str, default=’linear’	The prediction a new increment starts its Newton iterations from (off|linear|quadratic): none, the last accepted increment scaled to the new size, or the parabola through the last three accepted states (linear while only one accepted increment of the step of nonzero size is available).
retryWithoutExtrapolation	bool, default=False	If an extrapolated increment fails, attempt it once more from a zero prediction before cutting it back. A failed increment then also keeps the accepted history the next attempt extrapolates from.
extrapolateAcrossConnectivityChange	bool, default=False	Keep extrapolating across a rebuild of the equation system caused by a change of the connectivity alone (e.g., a contact search), as long as the rebuilt system numbers the degrees of freedom as before. Off: every rebuild starts the increment from a zero prediction.
extrapolateAfterModelChange	bool, default=True	Whether to extrapolate the predictor on the increment FOLLOWING a model change (adaptive mesh refinement). Set False to start that increment from a zero predictor, avoiding extrapolation of the one-off warm-start/remesh settling transient.
equilibrateAfterModelChange	bool, default=False	Whether to insert one constant-load, zero-time re-equilibration increment immediately after an adaptive mesh refinement, before advancing the load. When True, the warm-started refined mesh is first settled to equilibrium at the last converged load level (no load advance, no Dirichlet increment, zero time increment) so the subsequent load-advancing increment starts from an equilibrated state. Intended for softening problems where remeshing near the process zone otherwise couples the load advance with the warm-start settling transient in one solve. Under arc length control (NISTPArcLength with a controller), the re-equilibration is an arc length increment of size zero: the controlled quantity is held where the last accepted increment left it, and the load parameter is solved for, and with it the loads and the prescribed displacements it scales. Note: the equilibration solve integrates materials with dT=0, which suits rate-independent models; rate-dependent materials see no time advance during it (by design).
linsolver	str, default=’pardiso’	The linear solver: the name of a *linearSolver block, or a type of linear solver with its defaults (e.g. pardiso, superlu, blockamg).
pruneCondensedMatrixZeros	bool, default=None	Compact explicitly stored zeros out of the multi-point-constraint-condensed system matrix before solving. Pruning removes whichever entries happen to be exactly zero on this particular iteration, so the pattern changes on almost every Newton iteration – which defeats any linear solver that relies on a stable pattern across iterations (e.g. blockamg reusing a per-field AMG hierarchy built for the old one). Left unset (the default), the choice is automatic: pruned if the linear solver in force reports that it needs it (:attr:`~edelweissfe.linsolve.base.LinearSolver.needsPrunedSystemMatrix`), which a direct factorization does – PARDISO’s reordering is sensitive to the extra structural entries on these path-dependent condensed systems, and keeping them has been observed to drift from the converged reference path – and a Krylov solver with an AMG preconditioner does not: it gives the same iterations with and without the stored zeros, and pruning costs it a copy of the matrix and a pass over it in every Newton iteration. An explicit True/False here always wins over the automatic choice; turn pruning off for a direct solver only after verifying the load path against a reference run.
systemMatrixAssembly	str, default=’vij’	How the blocks of the elements, constraints and loads are summed into the system matrix (vij|direct): staged in a value array of the VIJ pattern’s length and gathered into CSR, or scattered into CSR right away, without that array. See the documentation page on the system matrix assembly.
useThreadedMPCCondensation	bool, default=False	Condense the multi-point-constraint system matrix, T^T K T + C, in place, with an OpenMP-threaded kernel, instead of with SciPy’s single-threaded CSR routines, which build a second matrix. The CSR pattern of K gains, once per equation system, the entries where a master collects the couplings of its slaves; the result does not depend on the number of threads. Explicit zeros (the slave columns among them) are pruned from a copy according to pruneCondensedMatrixZeros. Off by default.
report-performance-frequency	int, default=0	Print the cumulative performance table (totals since the start of the step) every this many increments, in addition to the one printed at the end of the step. 0 (default) disables the periodic report.
lineSearch	bool, default=False	Scale the Newton correction by the step length at which a parabola, fitted to the largest residual evaluated at the trial step lengths lineSearchStepLengths, is smallest (a quadratic line search).
lineSearchAfterIterations	int, default=5	The line search starts after this many iterations of an increment, at least 0: the first correction of an increment, which makes the prescribed increments, is never scaled.
lineSearchEveryIterations	int, default=2	From then on, the line search is done in every iteration whose number is a multiple of this.
lineSearchStepLengths	str, default=’0.8, 1.0, 1.2’	The trial step lengths of the line search, comma-separated, at least three.
fallBackToQuasiNewton	bool, default=False	Instead of cutting an increment back once the residual grew more often than maxGrowIter allows, continue it with the tangent of its first iteration (a quasi-Newton method); maxIter still applies.
convergenceCriterion	str, default=’legacy’	The convergence criterion of the Newton iterations, legacy (the original EdelweissFE one) or abaqus (the one of Abaqus/Standard); see edelweissfe.solvers.base.convergencecriteria.
iterationsForAlternativeTolerance	int, default=15	legacy criterion: I, the iteration from which on the flux residual is checked against fluxResidualToleranceAlternative instead of fluxResidualTolerance.
zeroFluxThreshold	float, default=None	legacy criterion: a field whose spatial average flux does not exceed this value carries no flux its residual could be compared with, so its flux residual counts as converged; needs minimumIterations >= 1, and refuses a smaller one. Off if not set.
zeroIncrementThreshold	float, default=None	legacy criterion: a field whose largest increment does not exceed this value has not moved, so there is no increment to relate its correction to, and its correction counts as converged. Off if not set.
minimumIterations	int, default=0	legacy criterion: no flux residual counts as converged before this iteration, so an increment takes at least this many Newton iterations; 0 accepts a prediction that already satisfies the criterion.
countFirstResidualAsGrowth	bool, default=True	Residual growths (the step option maxGrowIter): the residual history of a field starts from zero, so the first checked residual of an increment counts as a growth; if False, counting starts with the second.
countGrowthOfConvergedFields	bool, default=True	Residual growths (the step option maxGrowIter): count a growth of the flux residual of a field also in an iteration in which this field is converged and only other fields keep the iterations going.
abaqusResidualTolerance	float, default=0.005	abaqus criterion: R_n, the flux residual tolerance relative to the time average flux.
abaqusCorrectionTolerance	float, default=0.01	abaqus criterion: C_n, the correction tolerance relative to the largest increment of the field.
abaqusResidualToleranceAlternative	float, default=0.02	abaqus criterion: R_p, the flux residual tolerance used from iteration I_p on.
abaqusIterationsForAlternativeTolerance	int, default=9	abaqus criterion: I_p, the iteration from which on R_p is used.
abaqusZeroFluxThreshold	float, default=1e-05	abaqus criterion: epsilon, below which (relative to the time average flux) a field carries no flux.
abaqusZeroFluxCorrectionTolerance	float, default=0.001	abaqus criterion: C_epsilon, the relative correction tolerance of a field without flux.
abaqusLinearResidualTolerance	float, default=1e-08	abaqus criterion: R_l, below which a field is converged irrespective of its correction.
abaqusInitialTimeAveragedFlux	float, default=0.01	abaqus criterion: q~0, the time average flux of a field that has not carried flux yet (first step).
arcLengthController	str, default=’off’	The name of the step action of the step that controls the load parameter (an arc length controller, e.g. an indirectcontrol), or off: the step is solved without arc length control.
stopCondition	str, default=None	A model accessible expression; the step ends once it is true, checked after every accepted increment of nonzero size (under arc length control). It belongs to the controller it is declared with: an >>options block that names arcLengthController without a stopCondition returns to the stop condition of the *solver datalines (none by default), and one that sets a stopCondition without arc length control is refused.

Arc length control

A load-controlled analysis cannot pass a limit point: beyond the peak of the load, no larger load is in equilibrium. A displacement-controlled one cannot follow a path that snaps back, where the controlled displacement decreases again. The arc length method (here the indirect displacement control of Jirásek and Bažant) controls neither the load nor a displacement directly. The external load of a step becomes

\[\boldsymbol{f}_\mathrm{ext} = \boldsymbol{f}_0 + \lambda \, \boldsymbol{f}_\mathrm{ref},\]

the loads at the start of the step (the dead load \(\boldsymbol{f}_0\)) plus the change of the loads over the step (the reference load \(\boldsymbol{f}_\mathrm{ref}\)) scaled by a load parameter \(\lambda\), which is an unknown of every increment besides the solution. A controller adds the equation that determines it, typically that a combination of displacements – a crack opening, the mean contraction of a tunnel ring – grows steadily over the step. Every Newton iteration solves the linearized system for two right-hand sides at once, the out-of-balance force and the reference load, and the controller combines the two corrections into the correction of the solution and of the load parameter. The Newton iterations are those of NIST, into which the load parameter enters through the hooks of the solver (see its class documentation above). Every arc length step scales its own reference load: its load parameter starts from zero, where the external load is the dead load – the loads the steps before it ended with. An analysis that continues under arc length control in a later step declares there the change of the loads it applies on top of them.

Prescribed displacements belong to the reference loading of the step, like the loads: a Dirichlet boundary condition that prescribes a change \(\Delta\bar{\boldsymbol{u}}\) over the step moves its degrees of freedom by \(\lambda \, \Delta\bar{\boldsymbol{u}}\). They move with the load parameter, and only with it. The reference solve carries \(\Delta\bar{\boldsymbol{u}}\) in their rows, so a correction \(\Delta\Delta\lambda\) moves them by \(\Delta\Delta\lambda \, \Delta\bar{\boldsymbol{u}}\); eliminating their columns turns \(\Delta\bar{\boldsymbol{u}}\) into the forces it takes to impose it, which join the reference load. The out-of-balance solve leaves them where they are. So the prescribed degrees of freedom are no unknowns of Newton’s method, and their rows, which hold no equation, are left out of the convergence check, as in NIST. A predicted increment places them where the predicted load parameter puts them.

At the end of an arc length step, every load keeps the value the step applied, \(f(0) + \lambda \, (f(1) - f(0))\), from which the next step starts. With its amplitude \(a(t)\), a load keeps the change it declared for the step scaled by \(a(0) + \lambda \, (a(1) - a(0))\) (see LoadBase): a load ramped linearly over the step keeps \(\lambda\) times its change, and a load the step holds at a constant amplitude, such as a self-weight declared with f(t)='1', is part of the dead load and keeps its full value, whatever the load parameter.

The controller is a step action of the step, e.g. an indirectcontrol, and the solver option arcLengthController names it (off, the default, solves the step without arc length control):

*step, solver=theSolver
maxInc=1e-2, startInc=1e-2, minInc=1e-6, maxNumInc=1000, maxIter=25
>>options, name=theSolver, arcLengthController=opening
>>indirectcontrol, name=opening, dof1='model.nodeSets["IDC1"][0].fields["displacement"]', cVector1='-1, 0', dof2='model.nodeSets["IDC2"][0].fields["displacement"]', cVector2='1, 0', L=0.01

From Python, without an input file, the controller is a step action built through its typed constructor, and the option is a constructor argument of the solver:

stepActions["indirectcontrol"]["opening"] = IndirectControl("opening", dof1, dof2, [-1.0, 0.0], [1.0, 0.0], 0.01, model, journal)
solver = NISTPArcLength(jobInfo, journal, arcLengthController="opening")

A step under arc length control may end before its time is up: the solver option stopCondition, an expression of the model and its field outputs, e.g. stopCondition='fieldOutputs["UTop"] < -1.0', is checked after every accepted increment of nonzero size, and the step ends once it is true. A stop condition belongs to the controller it is declared with: an >>options block that names arcLengthController (off too) without a stopCondition drops the one of an earlier step, so that a later step does not stop by the condition that ended the step before it:

>>options, name=theSolver, arcLengthController=opening, stopCondition='fieldOutputs["UTop"] < -1.0'

A stop condition in the *solver datalines belongs to the solver: every step whose >>options name a controller without a stop condition of their own returns to it. An >>options block that sets a stopCondition for a step without arc length control is refused, since the condition would never be checked.

Every controller derives from ArcLengthControllerBase, which states what the solver asks of it:

The step actions that control the load parameter of the arc length method, NISTPArcLength.

class edelweissfe.stepactions.base.arclengthcontrollerbase.ArcLengthControllerBase(name, definition, jobInfo, model, fieldOutputController, journal)

The base of the step actions that control the load parameter of the arc length method, NISTPArcLength. The solver takes the controller of a step by its name, from its option arcLengthController.

The arc length method scales the reference load of a step by a load parameter \(\lambda\), which becomes an unknown of the increment besides the solution. The controller adds the equation that determines it, typically that a chosen combination of displacements – a crack opening, say – grows by a prescribed amount over the step. In every Newton iteration, the solver solves for the correction of the solution caused by the out-of-balance force, \(\Delta\Delta U_0\), and for the one caused by the reference load alone, \(\Delta\Delta U_f\). The controller returns the correction of the load parameter \(\Delta\Delta\lambda\) such that the corrected increment \(\Delta U + \Delta\Delta U_0 + \Delta\Delta\lambda \, \Delta\Delta U_f\) satisfies its equation. A line search asks it the same for every trial state, with \(\Delta\Delta U_0\) scaled by the trial step length, so that no trial state leaves the equation unsatisfied; the controller is asked several times per iteration, and the answers must not change it.

Attributes:

checkpointedState

schema

Parameters:

• name (str)

• definition (dict)

• jobInfo (dict)

• model (FEModel)

• fieldOutputController (FieldOutputController)

• journal (Journal)

Methods

applyAtIncrementStart(model, timeStep)	Is called when a step increment starts.
applyAtStepEnd(model)	Is called when a step successfully finished.
applyAtStepStart(model)	Is called when a step starts.
computeDDLambda(dU, ddU_0, ddU_f, timeStep, ...)	Compute the correction of the load parameter in the current Newton iteration.
finishIncrement(U, dU, dLambda, timeStep, ...)	Called once an increment has converged, e.g. to report or to store the controlled quantity it reached.
fromStepActionDefinition(name, definition, ...)	Create this step action from a parsed .inp step action definition.
getRestartData()	The state this step action carries; by default the attributes declared in checkpointedState.
setRestartData(data)	Restore the state getRestartData() returned.
updateStepAction(definition, jobInfo, model, ...)	Is called when an updated definition is present for a new step.
updateStepActionFromDefinition(definition, ...)	Update this step action from a parsed .inp step action definition.

abstractmethod computeDDLambda(dU, ddU_0, ddU_f, timeStep, dofManager)

Compute the correction of the load parameter in the current Newton iteration.

Parameters:

• dU (DofVector) – The current increment of the solution.

• ddU_0 (DofVector) – The correction of the solution caused by the out-of-balance force, of the dead load and the reference load scaled by the current load parameter; for a trial state of a line search, scaled by its step length.

• ddU_f (DofVector) – The correction of the solution caused by the reference load alone.

• timeStep (TimeStep) – The increment being solved: timeStep.stepProgress is the progress of the step at its end, and timeStep.stepProgressIncrement its size, so it starts at their difference.

• dofManager (DofManager) – The dof manager of the current equation system.

Returns:

The correction of the load parameter.

Return type:

float

abstractmethod finishIncrement(U, dU, dLambda, timeStep, dofManager)

Called once an increment has converged, e.g. to report or to store the controlled quantity it reached.

Parameters:

• U (DofVector) – The solution at the start of the increment.

• dU (DofVector) – The converged increment of the solution.

• dLambda (float) – The converged increment of the load parameter.

• timeStep (TimeStep) – The converged time step.

• dofManager (DofManager) – The dof manager of the current equation system.

name: str

The name of this step action, by which the solver option arcLengthController selects it.

The controllers of EdelweissFE, indirectcontrol and indirectcontractioncontrol, are indirect displacement controls: their controlled quantity is a linear combination of degrees of freedom, and their target L is, with absolute (the default), the value it reaches at the end of the step, counted from its value in the model when the controller is declared – whether for the first time, in any step, or anew in a later one. Without absolute, L is the change over the step:

The base of the arc length controllers of the indirect displacement control, whose controlled quantity is a linear combination of degrees of freedom.

class edelweissfe.stepactions.base.indirectdisplacementcontrolbase.IndirectDisplacementControlBase(name, absolute)

The base of the arc length controllers of the indirect displacement control (Jirásek and Bažant): the controlled quantity is a linear combination of degrees of freedom,

\[q = \boldsymbol{c} \cdot \boldsymbol{U},\]

with the weights \(\boldsymbol{c}\) of the controlled degrees of freedom (cVector), e.g. the opening of a crack between two nodes, or the mean contraction of a ring of nodes. Over its step, \(q\) changes by \(\Delta L\), and every increment by its share, \(\Delta p \, \Delta L\) for an increment of the step progress \(\Delta p\). The controller’s equation for the corrected increment of the solution,

\[\boldsymbol{c} \cdot (\Delta \boldsymbol{U} + \Delta\Delta \boldsymbol{U}_0 + \Delta\Delta\lambda \, \Delta\Delta \boldsymbol{U}_f) = \Delta p \, \Delta L,\]

is linear in the correction of the load parameter \(\Delta\Delta\lambda\), which computeDDLambda() solves it for.

The target L of a declaration – the construction of the controller, or its re-declaration in a later step, which updates it – is, with absolute, the value of \(q\) at the end of the step: the change of the step is \(\Delta L = L - q\), with \(q\) of the model when the controller is declared, i.e. of the state the step starts from (declareTarget()). A controller first declared in a later step aims at L as well, and so does one whose quantity changed over steps it did not control. Without absolute, L is the change \(\Delta L\) itself. A step that does not re-declare the controller changes \(q\) by the \(\Delta L\) of its last declaration once more.

A controller of this kind names its controlled degrees of freedom and their weights: cVector, controlledDofIndices() and controlledValuesOfModel().

Parameters:

• name (str) – The name of this step action.

• absolute (bool) – Whether the target L of a declaration is the value of the controlled quantity at the end of the step, or its change over the step. Fixed at construction: a re-declaration cannot change it.

Attributes:

schema

Methods

applyAtIncrementStart(model, timeStep)	Is called when a step increment starts.
applyAtStepEnd(model)	Is called when a step successfully finished.
applyAtStepStart(model)	Is called when a step starts.
computeDDLambda(dU, ddU_0, ddU_f, timeStep, ...)	The correction of the load parameter that makes the corrected increment change the controlled quantity by the share of the increment, \(\Delta p \, \Delta L\); see computeDDLambda().
controlledDofIndices(dofManager)	The indices of the controlled degrees of freedom in the dof vector of an equation system.
controlledQuantity()	The controlled quantity of the model, \(q = \boldsymbol{c} \cdot \boldsymbol{U}\).
controlledValuesOfModel()	The values of the controlled degrees of freedom in the model, as its node fields hold them.
declareTarget(L)	Declare the target of the step, from which the change of the controlled quantity over the step, changeOverStep, follows: \(\Delta L = L - q\), with the controlled quantity \(q\) of the model now, if absolute is set, else \(\Delta L = L\).
finishIncrement(U, dU, dLambda, timeStep, ...)	Called once an increment has converged, e.g. to report or to store the controlled quantity it reached.
fromStepActionDefinition(name, definition, ...)	Create this step action from a parsed .inp step action definition.
getRestartData()	The state this step action carries; by default the attributes declared in checkpointedState.
setRestartData(data)	Restore the state getRestartData() returned.
updateStepAction(definition, jobInfo, model, ...)	Is called when an updated definition is present for a new step.
updateStepActionFromDefinition(definition, ...)	Update this step action from a parsed .inp step action definition.

computeDDLambda(dU, ddU_0, ddU_f, timeStep, dofManager)

The correction of the load parameter that makes the corrected increment change the controlled quantity by the share of the increment, \(\Delta p \, \Delta L\); see computeDDLambda().

Parameters:

• dU (DofVector) – The current increment of the solution.

• ddU_0 (DofVector) – The correction of the solution caused by the out-of-balance force.

• ddU_f (DofVector) – The correction of the solution caused by the reference load alone.

• timeStep (TimeStep) – The increment being solved, whose size \(\Delta p\) is timeStep.stepProgressIncrement.

• dofManager (DofManager) – The dof manager of the current equation system.

Returns:

The correction of the load parameter.

Return type:

float

abstractmethod controlledDofIndices(dofManager)

The indices of the controlled degrees of freedom in the dof vector of an equation system.

Parameters:

dofManager (DofManager) – The dof manager of the equation system.

Returns:

The indices, in the order of cVector.

Return type:

np.ndarray

controlledQuantity()

The controlled quantity of the model, \(q = \boldsymbol{c} \cdot \boldsymbol{U}\).

Returns:

The controlled quantity.

Return type:

float

abstractmethod controlledValuesOfModel()

The values of the controlled degrees of freedom in the model, as its node fields hold them.

Returns:

The values, in the order of cVector.

Return type:

np.ndarray

declareTarget(L)

Declare the target of the step, from which the change of the controlled quantity over the step, changeOverStep, follows: \(\Delta L = L - q\), with the controlled quantity \(q\) of the model now, if absolute is set, else \(\Delta L = L\). Called with the controlled degrees of freedom and their weights in place.

Parameters:

L (float) – The target: the value of the controlled quantity at the end of the step with absolute, else its change over the step.

absolute

Whether the target L of a declaration is the value of the controlled quantity at the end of the step.

cVector: ndarray

The weights \(\\boldsymbol{c}\) of the controlled degrees of freedom, in the order of controlledDofIndices() and controlledValuesOfModel().

changeOverStep

The change of the controlled quantity over the step, \(\\Delta L\); see declareTarget().

checkpointedState: dict | None = {'changeOverStep': <class 'float'>}

The change of the controlled quantity over the step, \(\\Delta L\), which a declaration derives from its target and the model; a resumed run declares the controller anew on the restored model, so a checkpoint holds it.

The options of Newton’s method act on the pair of the solution increment and the increment of the load parameter \((\Delta U, \Delta\lambda)\):

• extrapolation: both are predicted from the same accepted increments with the same weights, linearly or quadratically, so that the prediction stays on one extrapolated path; retryWithoutExtrapolation retries a failed increment from a zero prediction of both.

• lineSearch: a line search scales the correction that the out-of-balance force causes, and the controller completes every trial state with the correction of the load parameter its equation asks for at that step length, as it completes the full correction. So every trial state satisfies the controller’s equation, and so does the iterate the line search chooses: the line search cannot leave the equation unsatisfied, also not from an iterate that does not satisfy it yet.

• fallBackToQuasiNewton: the tangent of the first iteration serves both right-hand sides.

Every accepted increment writes the load parameter of its step to the model, with its results, as model.additionalParameters["arc length parameter"], e.g. for a load-displacement curve. The solver checkpoints the load parameter of the step and the increments of it its predictor extrapolates, the loads what they keep from step to step, and the controller the change of its controlled quantity over the step, so an arc length controlled run resumes exactly where it stopped. The load parameter advances only when an increment has converged, with its solution, so a cut back increment starts again from the load parameter of the last accepted one.

The re-equilibration after a model change (equilibrateAfterModelChange, see Model modifiers) is, under arc length control, an arc length increment of size zero: the controller asks for a change of zero, so the controlled quantity stays where the last accepted increment left it, and the load parameter is solved for, and with it the loads and the prescribed displacements it scales – the counterpart of a zero-time increment under displacement control, which, unlike a constant load, has a solution past a limit point. As for any increment, the stop condition is checked before it, and the controller’s finishIncrement is called once it has converged. The load parameter it reaches is committed with its state, and the next increment starts from it.

NID - Nonlinear Implicit Dynamic (Newmark-beta)

The nonlinear implicit dynamic solver: Newmark-beta time integration on top of the Newton loop of NIST.

An increment solves the semi-discrete equation of motion

\[\boldsymbol{M} \ddot{\boldsymbol{u}}_{n+1} + \boldsymbol{C} \dot{\boldsymbol{u}}_{n+1} + \boldsymbol{P}_\mathrm{int}(\boldsymbol{u}_{n+1}) = \boldsymbol{P}_\mathrm{ext}(t_{n+1})\]

for the displacement at the end of the increment, with the velocity and the acceleration expressed through the displacement increment \(\Delta \boldsymbol{u} = \boldsymbol{u}_{n+1} - \boldsymbol{u}_n\) by the two Newmark relations (Newmark 1959; Hughes, The Finite Element Method, ch. 9)

\[\begin{split}\boldsymbol{u}_{n+1} &= \boldsymbol{u}_n + \Delta t \, \dot{\boldsymbol{u}}_n + \tfrac{\Delta t^2}{2} \left[ (1 - 2\beta) \, \ddot{\boldsymbol{u}}_n + 2\beta \, \ddot{\boldsymbol{u}}_{n+1} \right], \\ \dot{\boldsymbol{u}}_{n+1} &= \dot{\boldsymbol{u}}_n + \Delta t \left[ (1 - \gamma) \, \ddot{\boldsymbol{u}}_n + \gamma \, \ddot{\boldsymbol{u}}_{n+1} \right].\end{split}\]

Solved for the end-of-increment kinematics these read

\[\begin{split}\ddot{\boldsymbol{u}}_{n+1} &= \frac{\Delta \boldsymbol{u} - \Delta t \, \dot{\boldsymbol{u}}_n - \Delta t^2 \left( \tfrac{1}{2} - \beta \right) \ddot{\boldsymbol{u}}_n}{\beta \, \Delta t^2}, \\ \dot{\boldsymbol{u}}_{n+1} &= \dot{\boldsymbol{u}}_n + \Delta t \left[ (1 - \gamma) \, \ddot{\boldsymbol{u}}_n + \gamma \, \ddot{\boldsymbol{u}}_{n+1} \right],\end{split}\]

so that the effective residual and tangent handed to the Newton loop are

\[\begin{split}\boldsymbol{R} &= \boldsymbol{f}_\mathrm{ext} - \boldsymbol{f}_\mathrm{int} - \boldsymbol{M} \ddot{\boldsymbol{u}}_{n+1} - \boldsymbol{C} \dot{\boldsymbol{u}}_{n+1}, \\ \boldsymbol{K}_\mathrm{eff} &= \boldsymbol{K} + \frac{1}{\beta \, \Delta t^2} \boldsymbol{M} + \frac{\gamma}{\beta \, \Delta t} \boldsymbol{C}.\end{split}\]

Everything else – element evaluation, loads, constraints, multi-point-constraint condensation, Dirichlet handling, the convergence test, cutbacks and the linear solver – is the parent’s: this solver runs solveIncrement() itself and enters it only through its three hooks – NonlinearImplicitDynamic.initializeIncrement() (mass, damping, dynamic stiffness, initial acceleration), NonlinearImplicitDynamic.assembleAdditionalTerms() (two terms in the residual and two in the tangent, where the parent forms R = f_ext - f_int and before the tangent is converted to CSR) and NonlinearImplicitDynamic.finalizeIncrement() (the converged velocity and acceleration become the state).

Parameters. The default \(\beta = 1/4\), \(\gamma = 1/2\) is the average-acceleration (trapezoidal) rule: second-order accurate, unconditionally stable for linear problems and free of numerical damping, so the discrete energy of an undamped linear system is conserved exactly and any drift in it points at a bug rather than at the scheme. \(\gamma > 1/2\) adds algorithmic damping at the price of first-order accuracy; \(\gamma \geq 1/2\) together with \(\beta \geq (\gamma + 1/2)^2 / 4\) keeps the scheme unconditionally stable. A choice outside that region is accepted with a warning, \(\beta = 0\) (the explicit central-difference member of the family) is refused, because the formulation above divides by it.

Mass. The mass is the CONSISTENT one, \(\int_\Omega \rho N^T N \, dV\), assembled from computeConsistentInertia() into the very sparsity pattern the stiffness uses – both are scattered through the same VIJ layout of the DofManager, so adding \(\boldsymbol{M}/(\beta \Delta t^2)\) to the tangent is an entry-wise scaled addition of two value vectors, before the parent’s in-place CSR update, its multi-point-constraint condensation and its Dirichlet row replacement, all of which therefore act on the effective matrix without knowing it is one. For the same reason the mass is summed into CSR form by the stiffness’ own CSR generator. The damping \(\boldsymbol{C}\) is the diagonal each element reports through computeLumpedDamping; being diagonal, it is kept only as its nonzero entries and their positions in the same layout, not as a full value vector. Both are assembled in one loop over the elements (in chunks of elements on the thread pool, see Parallelization) when the equation system is (re)built – at a step’s start and after a topology change – not per Newton iteration. A rebuild caused only by a constraint changing its connectivity (a contact candidate list, which can change on every increment) leaves every element, and with it every element’s slot in the layout, where it was: the operators are then reused in the new layout rather than reassembled – see NonlinearImplicitDynamic._reuseMassAndDamping(). The one exception is a step that changes a material property mid-step (>>changematerialproperty), which invalidates the density and makes them reassemble every increment.

Checks on the mass. The assembled mass and damping are checked at runtime for non-finite entries, the damping for negative entries, and every time-integrated degree of freedom for a nonzero mass; the total mass of each field is reported. The symmetry of the mass is not checked at runtime – on the assembled matrix that costs a transposed copy of it, gigabytes on a large model – but pinned once per element type in tests/test_nid_element_mass.py, together with the total mass \(\rho V\). It does matter: a non-symmetric mass still converges, to a wrong response. That test covers only the element types it lists – the Marmot hexahedra C3D8, C3D20, C3D20R and the gradient-enhanced GC3D8, GC3D20R; any other element type (the plane elements, for one) is not guarded anywhere, so add it to that list before using it with this solver.

Which fields. Only the fields whose inertia is a mass (carriesLinearMomentum(), i.e. the displacement) are integrated in time. Rows and columns of the mass and the damping belonging to any other field are zeroed, and those fields keep the parent’s quasi-static treatment. A nonzero inertia an element reports there – a rotational inertia, or the micro-inertia of a gradient-enhanced element – is discarded with a warning, once per field. A model with no such field is refused.

Zero-length and negligible increments. The time stepper proposes a zero-length increment before the first real one of every step. A quasi-static solver equilibrates it; this solver skips it with the state kept, because in zero time nothing can move: displacement and velocity are continuous, and a load appearing at that instant is answered by the acceleration below, not by a displacement. Solving it statically would put a suddenly loaded model at its static deflection before the first real increment – the whole dynamic response skipped. An increment shorter than \(10^{-10}\) of the time elapsed in the step is skipped for the same reason: it is the round-off remainder the stepper’s progress accumulation leaves at the end of a step (400 increments of 0.0025 leave 2e-14), and the Newmark update would divide by its square. A skipped increment evaluates no force, so it leaves the node field result P as the last accepted increment wrote it.

Initial acceleration. At the start of a step the acceleration is computed from equilibrium, \(\boldsymbol{M} \ddot{\boldsymbol{u}}_0 = \boldsymbol{P}_\mathrm{ext}(t_0) - \boldsymbol{P}_\mathrm{int}(\boldsymbol{u}_0) - \boldsymbol{C} \dot{\boldsymbol{u}}_0\), with the loads evaluated at the step’s own start (an f(t) amplitude at \(t = 0\)). Without it, a load applied suddenly at the start of a step would enter the trapezoidal rule as though it were zero at \(t_0\), which costs one first-order error in the momentum, and the scheme’s second order is lost. It is computed in the first real increment of the step. A step that continues from a restart checkpoint – a resumed run, or a step that falls back to one – therefore computes it only if the checkpoint was written before that, after the zero increment: every later checkpoint carries the acceleration, and recomputing it would replace a consistent state with one that agrees only to solver tolerance. Switch it off with computeInitialAcceleration=False to continue a step from whatever acceleration was carried over.

State and restart. Velocity and acceleration are ordinary node-field entries, V and A, read from the model whenever the equation system is (re)built and written at the acceptance of every increment alongside the parent’s U/P/dU, both through the model’s layout of the equation system – the same mechanism the explicit solver uses for its velocity. That makes them part of every checkpoint the *output, type=restart manager writes and of every *fieldOutput (result=V, result=A), with no solver-specific restart code beyond the marker of whether the initial acceleration is still to be computed. The force written as P is the internal force with the inertial and damping forces of the converged state, \(\boldsymbol{f}_\mathrm{int} + \boldsymbol{M} \ddot{\boldsymbol{u}}_{n+1} + \boldsymbol{C} \dot{\boldsymbol{u}}_{n+1}\), as the documentation page Node field results defines it: at a support without load, the reaction.

Live h-adaptivity. A *modelModifier, type=hAdaptivity refining the mesh mid-step is supported. Three things make that work, and only the third is this solver’s own:

• Being node-field entries, V and A are carried onto the nodes a refinement creates by the same isoparametric interpolation from the parent element that already carries U (WARM_STARTED_NODE_FIELD_ENTRIES). Nothing here interpolates anything itself.

• The mass and the damping are not transferred at all: they are reassembled from the elements that now exist, because the parent rebuilds its DofManager on a topology change and everything bound to it is rebuilt with it. There is therefore no re-lumping error of the kind an explicit solver has to account for – the total mass is an assembly of the new mesh, not a redistribution of the old.

• The interpolated acceleration is only as good as an interpolation of the previous mesh’s acceleration, and it is not in equilibrium with the operators that were just reassembled. So the increment after a topology change recomputes the acceleration from equilibrium, with the cold-start machinery a step start already uses (computeInitialAcceleration, which also switches this off). The displacement and the velocity are kept: they are the kinematic state, and interpolation is the best statement of them the new mesh admits. Only a change that the topology pipeline actually recorded re-arms it – a mere change of a constraint’s connectivity (a contact candidate list) rebuilds the equation system without moving a node, and must not restart the acceleration of every degree of freedom in the model.

What is exactly conserved across a refinement, and what is not, is reported at every event and follows from the same partition-of-unity argument the explicit solver’s own live-adaptivity notes make: the assembled total mass is a geometric identity (children tile their parent at the same density) and is asserted; the linear momentum and the kinetic energy are exact for a spatially uniform velocity field and differ at second order in the velocity gradient across a refined parent otherwise, which is discretisation error rather than a defect, and are reported only. Unlike an explicit run, the increment that follows re-establishes equilibrium by Newton iteration, so an imbalance the transfer leaves behind is absorbed rather than radiated.

Prescribed motion. On a degree of freedom with a Dirichlet condition the acceleration is the one the Newmark relations imply for the prescribed displacement history, starting from zero at the step’s start. For a prescribed displacement that is not smooth in time – a velocity ramp switched on suddenly – the average-acceleration rule answers with an acceleration that oscillates from increment to increment on those degrees of freedom, and the consistent mass, having off-diagonal entries, passes the corresponding inertia forces on to the free neighbours. A smooth amplitude, or some algorithmic damping (newmarkGamma > 1/2), keeps that out of the response.

Not supported in this version. Rayleigh or any other damping model beyond what the elements report as their lumped damping; elements without a consistent mass (the pure-Python element library), which are refused at the first increment.

class edelweissfe.solvers.nonlinearimplicitdynamic.NIDSchema(defaultMaxIter=10, defaultCriticalIter=5, defaultMaxGrowingIter=10, extrapolation='linear', retryWithoutExtrapolation=False, extrapolateAcrossConnectivityChange=False, extrapolateAfterModelChange=True, equilibrateAfterModelChange=False, linsolver='pardiso', pruneCondensedMatrixZeros=None, systemMatrixAssembly='vij', useThreadedMPCCondensation=False, reportPerformanceFrequency=0, lineSearch=False, lineSearchAfterIterations=5, lineSearchEveryIterations=2, lineSearchStepLengths='0.8, 1.0, 1.2', fallBackToQuasiNewton=False, convergenceCriterion='legacy', iterationsForAlternativeTolerance=15, zeroFluxThreshold=None, zeroIncrementThreshold=None, minimumIterations=0, countFirstResidualAsGrowth=True, countGrowthOfConvergedFields=True, abaqusResidualTolerance=0.005, abaqusCorrectionTolerance=0.01, abaqusResidualToleranceAlternative=0.02, abaqusIterationsForAlternativeTolerance=9, abaqusZeroFluxThreshold=1e-05, abaqusZeroFluxCorrectionTolerance=0.001, abaqusLinearResidualTolerance=1e-08, abaqusInitialTimeAveragedFlux=0.01, newmarkBeta=0.25, newmarkGamma=0.5, computeInitialAcceleration=True, massConservationTolerance=1e-06)

NISTSchema’s options, plus the Newmark parameters and the initial-acceleration switch. Mirrors NonlinearImplicitDynamic.SolverSpecificOptions one-for-one; that dict remains the runtime source of truth, for the reasons given on the parent schema.

Attributes:

pruneCondensedMatrixZeros

zeroFluxThreshold

zeroIncrementThreshold

Parameters:

• defaultMaxIter (int | None)

• defaultCriticalIter (int | None)

• defaultMaxGrowingIter (int | None)

• extrapolation (str | None)

• retryWithoutExtrapolation (bool | None)

• extrapolateAcrossConnectivityChange (bool | None)

• extrapolateAfterModelChange (bool | None)

• equilibrateAfterModelChange (bool | None)

• linsolver (str | None)

• pruneCondensedMatrixZeros (bool | None)

• systemMatrixAssembly (str | None)

• useThreadedMPCCondensation (bool | None)

• reportPerformanceFrequency (int | None)

• lineSearch (bool | None)

• lineSearchAfterIterations (int | None)

• lineSearchEveryIterations (int | None)

• lineSearchStepLengths (str | None)

• fallBackToQuasiNewton (bool | None)

• convergenceCriterion (str | None)

• iterationsForAlternativeTolerance (int | None)

• zeroFluxThreshold (float | None)

• zeroIncrementThreshold (float | None)

• minimumIterations (int | None)

• countFirstResidualAsGrowth (bool | None)

• countGrowthOfConvergedFields (bool | None)

• abaqusResidualTolerance (float | None)

• abaqusCorrectionTolerance (float | None)

• abaqusResidualToleranceAlternative (float | None)

• abaqusIterationsForAlternativeTolerance (int | None)

• abaqusZeroFluxThreshold (float | None)

• abaqusZeroFluxCorrectionTolerance (float | None)

• abaqusLinearResidualTolerance (float | None)

• abaqusInitialTimeAveragedFlux (float | None)

• newmarkBeta (float | None)

• newmarkGamma (float | None)

• computeInitialAcceleration (bool | None)

• massConservationTolerance (float | None)

class edelweissfe.solvers.nonlinearimplicitdynamic.NonlinearImplicitDynamic(jobInfo, journal, **kwargs)

This is the Nonlinear Implicit Dynamic – solver (NID), Newmark-beta time integration on top of the Newton loop of NIST.

Parameters:

• jobInfo – A dictionary containing the job information.

• journal – The journal instance for logging.

Attributes:

convergenceTable

equationSystemLayout

incrementStatus

mpcTransformation

Methods

acceptIncrement(step, model, timeStep)	Commit the converged increment to the model, and keep the predictor's state; see acceptIncrement().
applyCorrection(dU, ddU, stepLength)	Move the iterate by the correction, scaled by the step length \(\alpha\) a line search chose (1 without a line search): \(\Delta \boldsymbol{U} \leftarrow \Delta \boldsymbol{U} + \alpha \, \Delta\Delta \boldsymbol{U}\).
applyDirichletToResidual(timeStep, R, dirichlets)	Impose the Dirichlet BCs on the residual.
applyDirichletToStiffness(K, dirichlets[, rhs])	Impose the Dirichlet BCs on the system matrix -- and, given the right-hand side(s) rhs (modified in place), eliminate the constrained DOFs from the columns as well; see edelweissfe.solvers.base.dirichlet.applyDirichletToStiffness().
applyOptionsOverride(fieldValues)	Apply an override of this solver's options, see applyOptionsOverride(); an override that sets an option where it means nothing is refused before anything of it applies (see refuseOptionsThatMeanNothing()).
applyStepActionsAtIncrementStart(model, ...)	Called when all step actions should be applied at the start of a step.
applyStepActionsAtStepEnd(model, stepActions)	Called when all step actions should finish a step.
applyStepActionsAtStepStart(model, step)	Apply every step action's step-start part (initial conditions, material initialization, prescribed fields, ...).
assembleAdditionalTerms(dU, R, F, K, timeStep)	Add the inertia and damping forces to the residual and their derivatives to the tangent; see the parent's hook.
assembleConstraintForces(constraints, U_np, ...)	Evaluate every active constraint and add its residual contribution to the given vector; a constraint that is switched off (see active) is skipped.
assembleConstraints(constraints, U_np, dU, ...)	Loop over all constraints, and evaluate the active ones: each adds its residual contribution and the tangent of it, see computeKernels().
assembleInternalForces(U_np, dU, P, K, F, ...)	Assemble the internal forces of the model at a trial state, and their tangent: those of the elements of its equation system, see computeElements().
assembleLoads(loads, U_np, dU, P, K, timeStep)	Assemble the residual contributions of all loads: each adds the negative of the external load it applies, and the tangent of that contribution, which is nonzero only for a load that follows the deformation; see The residual and the tangent.
assembleLoadsAtIncrementEnd(stepActions, ...)	Assemble the residual contributions of the loads acting at the end of the increment, and their tangents; called by assembleResidualAndTangent().
assembleResidualAndTangent(U_n, dU, U_np, R, ...)	Evaluate the model at the trial solution \(U_{n+1} = U_n + \Delta U\): the right-hand side \(R = f_\mathrm{ext} - f_\mathrm{int}\) (see The residual and the tangent) with the terms of assembleAdditionalTerms(), condensed by the multi-point constraints, and the tangent.
assembleStiffnessCSR(K)	The assembled system matrix in CSR format, see assembledMatrix().
assembleSystemMatrix(K, quasiNewton)	The system matrix of this Newton iteration, before the Dirichlet boundary conditions are imposed: the tangent in CSR layout, condensed by the multi-point constraints -- or, once the increment has fallen back to the quasi-Newton method, the tangent of its first iteration.
attemptIncrement(step, model, timeStep)	Solve the increment by Newton's method, rebuilding the equation system first if the topology update changed the model, and re-equilibrating a changed model first if so configured; see attemptIncrement().
beginStep(step, model, ...)	Start a step; see the parent.
buildMPCTransformation(model[, stepActions])	Collect the linear dependency records from all multi-point constraints of the model and assemble the master-slave condensation operator for the current equation system.
checkConvergence(R, ddU, dU, F, ...)	Check the convergence, individually for each field and for the scalar variables, based on the current flux residual and the correction, using the convergence criterion of the solver (see edelweissfe.solvers.base.convergencecriteria).
checkDivergingSolution(...)	Check if the iterative solution scheme is diverging: if the flux residual of a field grew, from one checked iteration to the next, more than maxGrowingIter times in this increment.
checkDivergingSolutionOrFallBack(...)	Cut an increment back whose residual grew more often than maxGrowingIter allows (see checkDivergingSolution()) -- or, with the option fallBackToQuasiNewton, let it continue with the tangent of its first iteration (a quasi-Newton method), see assembleSystemMatrix().
checkMPCDirichletConflicts(transformation, ...)	Raise if any Dirichlet boundary condition of the step prescribes a DOF that is a slave DOF of a multi-point constraint.
collectLoads(stepActions)	The loads among a step's actions: every action that is a LoadBase, whatever step action module it is defined in, so that a new kind of load needs no change to any solver.
computeElements(elements, U_np, dU, P, K, F, ...)	Loop over all elements, and evalute them: each adds its internal forces and their tangent, see computeKernels().
createConvergenceCriterion()	Create the convergence criterion selected by the option convergenceCriterion.
createSystemMatrixAssembly(assemblyType)	The assembly of the system matrix of the equation system just built (option systemMatrixAssembly), on the complete CSR pattern of csrGenerator, with numberOfAssemblyWorkers() workers.
endStep(step, model)	Report the step's performance timing.
extrapolateLastIncrement(extrapolation, ...)	Predict the solution increment, the point Newton's method starts from, by extrapolating the last accepted increments of the step.
extrapolates(extrapolation, prevTimeStep)	Whether an increment that follows the accepted increment prevTimeStep starts from an extrapolated prediction.
extrapolationWeights(extrapolation, ...)	The weights of the last accepted increment and of the one before it in the prediction of extrapolateLastIncrement(): \(\Delta U = w_n \, \Delta U_n + w_{n-1} \, \Delta U_{n-1}\), with
finalizeIncrement(model)	The converged trial kinematics become the state of the solver, which the model receives at the acceptance of the increment (see writeResultsToModel()); see the parent's hook.
firstIterateNeedsCorrection(dirichlets)	Whether the first iterate of an increment without a prediction -- the converged state of the last increment -- needs a correction before it can be checked for convergence (see solveIncrement()): if Dirichlet boundary conditions act, whose prescribed increment the first correction makes.
forEachChunk(entities, work)	Call work(chunk) for chunks that together hold all entities: the solver's loop over independent entities, other than the elements, that the model's equation system layout evaluates, e.g. the material points of a meshfree model (see updateMaterialPoints()).
forgetEquationSystem()	Drop the equation system of the last build -- the model's layout of it with the DOF manager, the system matrix assembly and its CSR pattern, the multi-point-constraint condensation, the DOF vectors -- so that it is gone before the next one is built, never both at once.
forgetPredictionAfterFailure()	A failed increment drops the prediction, which may have caused the failure, unless retryWithoutExtrapolation has answered that already: then the accepted history is kept, and the smaller increment that follows the cutback extrapolates it.
getRestartData()	The state this solver carries from one increment to the next, as it is; see checkpointedState.
initializeIncrement(U_n, stepActions, model, ...)	Bring the Newmark operators up to date with the equation system, form the two terms the dynamics add to the tangent, and -- at a cold step start or after a topology change -- compute the acceleration from equilibrium; see the parent's hook.
isOutputIncrement(timeStep)	Whether field outputs, output managers and restart checkpoints are written after this increment.
lineSearchIsDue(iterationCounter)	Whether the correction of this iteration is scaled by a line search (lineSearch): after lineSearchAfterIterations iterations, in every iteration whose number is a multiple of lineSearchEveryIterations.
lineSearchStepLengths()	The trial step lengths of the line search, from the option lineSearchStepLengths.
linearSolve(A, b)	Solve the linear equation system.
locateConstrainedDofs(dirichlets)	Determine, up front, which global DOFs each Dirichlet BC constrains.
numberOfAssemblyWorkers()	How many threads write blocks of the system matrix at the same time: one, the elements are computed in the calling thread.
predictionSurvivesRebuild(modelHasChanged, ...)	Whether this increment may extrapolate the last accepted increments after the equation system was rebuilt.
prepareIncrement(step, model, isRetry)	The topology update, at the start of every increment; see prepareIncrement().
printResidualOutlierNodes(residualOutliers)	Print which nodes have the largest residuals.
quadraticLineSearch(largestResidualAt)	The step length \(\alpha\) by which a Newton correction is scaled, for a Newton iteration that overshoots or oscillates (e.g., materials switching between loading and unloading).
readStateFromModel(U)	Read the state the next increment starts from out of the model, into the vectors of the equation system just (re)built: the solution, as in the parent, and the velocity and the acceleration of the fields carrying a mass, where the last accepted increment (or the checkpoint, or nothing -- a start from rest) left them, and where a refinement has meanwhile interpolated them onto the nodes it created.
refuseOptionsThatMeanNothing(options, ...)	Refuse an option that a declaration of options sets where it means nothing, with the options in force after the declaration.
reportTopologyChangeConservation(before, after)	Report what a topology change did to the quantities that ought to survive it: every dynamic field's total mass is checked (and raises when violated), momentum and kinetic energy are reported -- see conservationchecks for which of them is exact when.
schema	Option schema for this solver, per OptionSchemaProvider.
selectLinearSolverForStep()	Make linSolver the linear solver the step about to begin names (linsolver): the one of the previous step if it names the same, else a new one.
setRestartData(data)	Restore the state getRestartData() returned, so that the resumed step continues exactly where the uninterrupted run would have.
setTrialState(dU, ddU, stepLength)	Set the trial state of a line search at a step length \(\alpha\) along the correction, and return its solution increment, \(\Delta \boldsymbol{U} + \alpha \, \Delta\Delta \boldsymbol{U}\), at which the model is evaluated next; the iterate stays where it is.
solveForCorrection(K_, R, dirichlets)	Solve the linearized system of a Newton iteration, \(\boldsymbol{K} \, \Delta\Delta \boldsymbol{U} = \boldsymbol{R}\), for the correction \(\Delta\Delta \boldsymbol{U}\).
solveIncrement(U_n, dU, P, K, stepActions, ...)	Newton-Raphson scheme on the Newmark-effective equation of motion of an increment.
updateRigidBodies(model, timeStep)	Refresh the kinematics of all rigid bodies in the model after a converged increment.
updateTopologyAndConnectivity(model, step[, ...])	Run the topology update, then let every mesh-dependent consumer catch up on it.
validateModelCapabilities(model)	Validate whether the solver supports the active features/constraints of the model.
writeResultsToModel(U, dU, P)	Write the results of a converged increment to the model, as the parent does, and the velocity and the acceleration of the fields carrying a mass with them: what field outputs read, what a checkpoint stores, and what the next equation system reads back.

schema

Option schema for this solver, per OptionSchemaProvider.

Attributes:

pruneCondensedMatrixZeros

zeroFluxThreshold

zeroIncrementThreshold

alias of NIDSchema

assembleAdditionalTerms(dU, R, F, K, timeStep)

Add the inertia and damping forces to the residual and their derivatives to the tangent; see the parent’s hook. The trial displacement increment fixes the trial kinematics.

Parameters:

• dU (DofVector) – The current trial solution increment.

• R (DofVector) – The residual, augmented in place.

• F (DofVector) – The reference flux scale of the convergence test, augmented in place.

• K (SystemMatrixAssembly) – The tangent being assembled; its VIJ value array is augmented in place.

• timeStep (TimeStep) – The time step.

beginStep(step, model, fieldOutputController, outputmanagers)

Start a step; see the parent.

Arms the initial-acceleration computation and drops the operators of the previous step – the parent rebuilds its equation system at the start of every step, and the operators follow that system. A step that continues from a restart checkpoint takes over, with the checkpointed state restored after this, whether the initial acceleration is still to be computed: the first real increment of the step computes it, and a checkpoint written after the zero increment comes before that. The velocity and acceleration are not dropped: they are re-read from the V/A node-field entries the previous step (or the checkpoint) left, which is what carries the kinematic state across steps.

Parameters:

• step – The step to be solved.

• model (FEModel) – The model tree.

• fieldOutputController (FieldOutputController) – The field output controller.

• outputmanagers (dict[str, OutputManagerBase]) – The output managers.

finalizeIncrement(model)

The converged trial kinematics become the state of the solver, which the model receives at the acceptance of the increment (see writeResultsToModel()); see the parent’s hook. Only a converged increment gets here – a failing one leaves by exception and keeps the previous velocity and acceleration for the cutback.

Parameters:

model (FEModel) – The model tree.

initializeIncrement(U_n, stepActions, model, timeStep)

Bring the Newmark operators up to date with the equation system, form the two terms the dynamics add to the tangent, and – at a cold step start or after a topology change – compute the acceleration from equilibrium; see the parent’s hook.

After the step actions, which the parent has applied just before: a >>changematerialproperty acting on this increment changes the density the mass is assembled from.

Parameters:

• U_n (DofVector) – The solution at the start of the increment.

• stepActions (dict) – The active step actions.

• model (FEModel) – The model tree.

• timeStep (TimeStep) – The time step.

readStateFromModel(U)

Read the state the next increment starts from out of the model, into the vectors of the equation system just (re)built: the solution, as in the parent, and the velocity and the acceleration of the fields carrying a mass, where the last accepted increment (or the checkpoint, or nothing – a start from rest) left them, and where a refinement has meanwhile interpolated them onto the nodes it created.

Parameters:

U (DofVector) – The solution vector of the equation system, zero, overwritten.

reportTopologyChangeConservation(before, after)

Report what a topology change did to the quantities that ought to survive it: every dynamic field’s total mass is checked (and raises when violated), momentum and kinetic energy are reported – see conservationchecks for which of them is exact when. The mass is reassembled here rather than transferred, so its conservation is an identity of the assembly and the cheapest correctness check on the whole refinement.

Parameters:

• before (_ConservedQuantities) – The quantities measured on the system before and on the one after the change.

• after (_ConservedQuantities) – The quantities measured on the system before and on the one after the change.

Raises:

RuntimeError – If any dynamic field’s total mass changed by more than massConservationTolerance.

solveIncrement(U_n, dU, P, K, stepActions, model, timeStep, prevTimeStep, extrapolation, maxIter, maxGrowingIter)

Newton-Raphson scheme on the Newmark-effective equation of motion of an increment.

The parent’s loop, into which the dynamics enter through its hooks (initializeIncrement(), assembleAdditionalTerms(), finalizeIncrement()). A zero-length increment – the one the time stepper yields before the first real increment of a step – has no equation of motion and is skipped with the state kept: displacement and velocity are continuous in time, and a load appearing at that instant is answered by the initial acceleration, not by a displacement.

Parameters:

• U_n (DofVector) – The old solution vector.

• dU (DofVector) – The old solution increment.

• P (DofVector) – The old reaction vector.

• K (SystemMatrixAssembly) – The system matrix to be used.

• stepActions (list) – The list of active step actions.

• model (FEModel) – The model tree.

• timeStep (TimeStep) – The time step.

• prevTimeStep (TimeStep) – The previous time step.

• extrapolation (str) – The type of extrapolation to be used.

• maxIter (int) – The maximum number of iterations to be used.

• maxGrowingIter (int) – The maximum number of growing residuals until the Newton-Raphson is terminated.

Returns:

A tuple containing

• the new solution vector

• the solution increment

• the new reaction vector

• the number of required iterations

• the history of residuals per field

Return type:

tuple[DofVector,DofVector,DofVector,int,dict]

writeResultsToModel(U, dU, P)

Write the results of a converged increment to the model, as the parent does, and the velocity and the acceleration of the fields carrying a mass with them: what field outputs read, what a checkpoint stores, and what the next equation system reads back.

The force written is the internal force with the inertial and damping forces of the converged state, \(\boldsymbol{f}_\mathrm{int} + \boldsymbol{M} \ddot{\boldsymbol{u}}_{n+1} + \boldsymbol{C} \dot{\boldsymbol{u}}_{n+1}\) (see the documentation page Node field results): at a support that no load acts on, the reaction, the inertia of the support’s own mass included. A skipped (zero-length or negligible) increment evaluates no force; the model keeps the force the last accepted increment wrote, of this step or of the step before.

Parameters:

• U (DofVector) – The solution.

• dU (DofVector) – The solution increment.

• P (DofVector) – The internal force, the elements’ residual contributions.

checkpointedState: dict | None = {'_conservationCheck': <class 'edelweissfe.solvers.base.conservationchecks.ConservationCheck'>, '_conservedAtLastAccept': <class 'edelweissfe.solvers.nonlinearimplicitdynamic._ConservedQuantities'>, '_initialAccelerationPending': <class 'bool'>, 'convergenceCriterion': <class 'edelweissfe.solvers.base.convergencecriteria.ConvergenceCriterionBase'>, 'dUOfIncrementBeforeLast': <class 'numpy.ndarray'>, 'dUOfLastIncrement': <class 'numpy.ndarray'>, 'dofNumberingOfLastIncrement': <class 'str'>, 'prevTimeStep': <class 'edelweissfe.timesteppers.timestep.TimeStep'>, 'timeStepOfIncrementBeforeLast': <class 'edelweissfe.timesteppers.timestep.TimeStep'>}

The predictor’s state, as in the parent, the drift the conservation checks accumulate over the step, the conserved totals after the last accepted increment, which a topology change is compared against, and whether the next increment computes the acceleration from equilibrium first, which a checkpoint written after the zero increment of a step still has to do. The velocity and acceleration are node-field entries, restored with the model.

supportsModelModifiers = True

Live h-adaptivity is supported; see the module docstring.

systemMatrixAssemblies = {'vij': <class 'edelweissfe.numerics.systemmatrixassembly.VIJSystemMatrixAssembly'>}

the effective tangent adds the mass and the damping entry by entry to the VIJ value array (see assembleAdditionalTerms()), and the mass is gathered into CSR by the generator’s gather map.

Type:

Only the VIJ assembly

NID

Option	Type / Default	Description
defaultMaxIter	int, default=10	The default maximum number of iterations.
defaultCriticalIter	int, default=5	The default number of critical iterations.
defaultMaxGrowingIter	int, default=10	The default number of allowed residual growths.
extrapolation	str, default=’linear’	The prediction a new increment starts its Newton iterations from (off|linear|quadratic): none, the last accepted increment scaled to the new size, or the parabola through the last three accepted states (linear while only one accepted increment of the step of nonzero size is available).
retryWithoutExtrapolation	bool, default=False	If an extrapolated increment fails, attempt it once more from a zero prediction before cutting it back. A failed increment then also keeps the accepted history the next attempt extrapolates from.
extrapolateAcrossConnectivityChange	bool, default=False	Keep extrapolating across a rebuild of the equation system caused by a change of the connectivity alone (e.g., a contact search), as long as the rebuilt system numbers the degrees of freedom as before. Off: every rebuild starts the increment from a zero prediction.
extrapolateAfterModelChange	bool, default=True	Whether to extrapolate the predictor on the increment FOLLOWING a model change (adaptive mesh refinement). Set False to start that increment from a zero predictor, avoiding extrapolation of the one-off warm-start/remesh settling transient.
equilibrateAfterModelChange	bool, default=False	Whether to insert one constant-load, zero-time re-equilibration increment immediately after an adaptive mesh refinement, before advancing the load. When True, the warm-started refined mesh is first settled to equilibrium at the last converged load level (no load advance, no Dirichlet increment, zero time increment) so the subsequent load-advancing increment starts from an equilibrated state. Intended for softening problems where remeshing near the process zone otherwise couples the load advance with the warm-start settling transient in one solve. Under arc length control (NISTPArcLength with a controller), the re-equilibration is an arc length increment of size zero: the controlled quantity is held where the last accepted increment left it, and the load parameter is solved for, and with it the loads and the prescribed displacements it scales. Note: the equilibration solve integrates materials with dT=0, which suits rate-independent models; rate-dependent materials see no time advance during it (by design).
linsolver	str, default=’pardiso’	The linear solver: the name of a *linearSolver block, or a type of linear solver with its defaults (e.g. pardiso, superlu, blockamg).
pruneCondensedMatrixZeros	bool, default=None	Compact explicitly stored zeros out of the multi-point-constraint-condensed system matrix before solving. Pruning removes whichever entries happen to be exactly zero on this particular iteration, so the pattern changes on almost every Newton iteration – which defeats any linear solver that relies on a stable pattern across iterations (e.g. blockamg reusing a per-field AMG hierarchy built for the old one). Left unset (the default), the choice is automatic: pruned if the linear solver in force reports that it needs it (:attr:`~edelweissfe.linsolve.base.LinearSolver.needsPrunedSystemMatrix`), which a direct factorization does – PARDISO’s reordering is sensitive to the extra structural entries on these path-dependent condensed systems, and keeping them has been observed to drift from the converged reference path – and a Krylov solver with an AMG preconditioner does not: it gives the same iterations with and without the stored zeros, and pruning costs it a copy of the matrix and a pass over it in every Newton iteration. An explicit True/False here always wins over the automatic choice; turn pruning off for a direct solver only after verifying the load path against a reference run.
systemMatrixAssembly	str, default=’vij’	How the blocks of the elements, constraints and loads are summed into the system matrix (vij|direct): staged in a value array of the VIJ pattern’s length and gathered into CSR, or scattered into CSR right away, without that array. See the documentation page on the system matrix assembly.
useThreadedMPCCondensation	bool, default=False	Condense the multi-point-constraint system matrix, T^T K T + C, in place, with an OpenMP-threaded kernel, instead of with SciPy’s single-threaded CSR routines, which build a second matrix. The CSR pattern of K gains, once per equation system, the entries where a master collects the couplings of its slaves; the result does not depend on the number of threads. Explicit zeros (the slave columns among them) are pruned from a copy according to pruneCondensedMatrixZeros. Off by default.
report-performance-frequency	int, default=0	Print the cumulative performance table (totals since the start of the step) every this many increments, in addition to the one printed at the end of the step. 0 (default) disables the periodic report.
lineSearch	bool, default=False	Scale the Newton correction by the step length at which a parabola, fitted to the largest residual evaluated at the trial step lengths lineSearchStepLengths, is smallest (a quadratic line search).
lineSearchAfterIterations	int, default=5	The line search starts after this many iterations of an increment, at least 0: the first correction of an increment, which makes the prescribed increments, is never scaled.
lineSearchEveryIterations	int, default=2	From then on, the line search is done in every iteration whose number is a multiple of this.
lineSearchStepLengths	str, default=’0.8, 1.0, 1.2’	The trial step lengths of the line search, comma-separated, at least three.
fallBackToQuasiNewton	bool, default=False	Instead of cutting an increment back once the residual grew more often than maxGrowIter allows, continue it with the tangent of its first iteration (a quasi-Newton method); maxIter still applies.
convergenceCriterion	str, default=’legacy’	The convergence criterion of the Newton iterations, legacy (the original EdelweissFE one) or abaqus (the one of Abaqus/Standard); see edelweissfe.solvers.base.convergencecriteria.
iterationsForAlternativeTolerance	int, default=15	legacy criterion: I, the iteration from which on the flux residual is checked against fluxResidualToleranceAlternative instead of fluxResidualTolerance.
zeroFluxThreshold	float, default=None	legacy criterion: a field whose spatial average flux does not exceed this value carries no flux its residual could be compared with, so its flux residual counts as converged; needs minimumIterations >= 1, and refuses a smaller one. Off if not set.
zeroIncrementThreshold	float, default=None	legacy criterion: a field whose largest increment does not exceed this value has not moved, so there is no increment to relate its correction to, and its correction counts as converged. Off if not set.
minimumIterations	int, default=0	legacy criterion: no flux residual counts as converged before this iteration, so an increment takes at least this many Newton iterations; 0 accepts a prediction that already satisfies the criterion.
countFirstResidualAsGrowth	bool, default=True	Residual growths (the step option maxGrowIter): the residual history of a field starts from zero, so the first checked residual of an increment counts as a growth; if False, counting starts with the second.
countGrowthOfConvergedFields	bool, default=True	Residual growths (the step option maxGrowIter): count a growth of the flux residual of a field also in an iteration in which this field is converged and only other fields keep the iterations going.
abaqusResidualTolerance	float, default=0.005	abaqus criterion: R_n, the flux residual tolerance relative to the time average flux.
abaqusCorrectionTolerance	float, default=0.01	abaqus criterion: C_n, the correction tolerance relative to the largest increment of the field.
abaqusResidualToleranceAlternative	float, default=0.02	abaqus criterion: R_p, the flux residual tolerance used from iteration I_p on.
abaqusIterationsForAlternativeTolerance	int, default=9	abaqus criterion: I_p, the iteration from which on R_p is used.
abaqusZeroFluxThreshold	float, default=1e-05	abaqus criterion: epsilon, below which (relative to the time average flux) a field carries no flux.
abaqusZeroFluxCorrectionTolerance	float, default=0.001	abaqus criterion: C_epsilon, the relative correction tolerance of a field without flux.
abaqusLinearResidualTolerance	float, default=1e-08	abaqus criterion: R_l, below which a field is converged irrespective of its correction.
abaqusInitialTimeAveragedFlux	float, default=0.01	abaqus criterion: q~0, the time average flux of a field that has not carried flux yet (first step).
newmarkBeta	float, default=0.25	Newmark parameter beta, weighting the end-of-increment acceleration in the displacement update. Must be > 0. 1/4 with gamma = 1/2 is the average-acceleration (trapezoidal) rule.
newmarkGamma	float, default=0.5	Newmark parameter gamma, weighting the end-of-increment acceleration in the velocity update. 1/2 gives second-order accuracy and no algorithmic damping; larger values damp.
computeInitialAcceleration	bool, default=True	Compute the acceleration at the start of every step from equilibrium with the loads at the step’s start, instead of continuing from the carried-over acceleration. Skipped automatically when the step resumes from a restart checkpoint. Also re-armed for the increment following a topology change (h-adaptivity), whose interpolated acceleration is not in equilibrium with the operators reassembled on the refined mesh.
massConservationTolerance	float, default=1e-06	Relative-change tolerance for the per-topology-change check on a dynamic field’s total assembled mass. Children of a refined element tile it and carry the same density, so this is conserved geometrically; the default is already generous relative to floating-point precision and exists to catch a genuinely wrong refinement, not to absorb ordinary quadrature noise. The counterpart of the explicit solver’s lumped-quantity-conservation-tolerance.

The relevant blocks of a deck, abridged from the regression case testfiles/marmot/NID – a bar under a suddenly applied end load, oscillating about its static deflection:

*solver, solver=NID, name=theSolver
newmarkBeta=0.25
newmarkGamma=0.5

*fieldOutput
>>perNode, nSet=gen_right, field=displacement, result=U, name=tipU, f(x)='np.mean(x[:,0])', saveHistory=True
>>perNode, nSet=gen_right, field=displacement, result=V, name=tipV, f(x)='np.mean(x[:,0])', saveHistory=True
>>perNode, nSet=gen_right, field=displacement, result=A, name=tipA, f(x)='np.mean(x[:,0])', saveHistory=True

*step, solver=theSolver
stepLength=2.0, startInc=0.01, maxInc=0.01, minInc=1e-4, maxNumInc=1000, maxIter=25
>>dirichlet, name=fixedEnd, nSet=gen_left, field=displacement, 1=0.0
>>nodeforces, name=endLoad, nSet=gen_right, field=displacement, 1=0.0987, f(t)='1'

The time increment is the step’s increment (stepLength times the step-progress increment), so a constant dt is obtained with startInc = maxInc; a cutback shrinks it like any other implicit increment. The velocity and acceleration are available to every *fieldOutput as result=V and result=A on the displacement field, and travel with the *output, type=restart checkpoints, which is what makes a resumed run continue the same trajectory rather than restart it from rest.

NIDParallel is the same solver with the element loop evaluated in parallel (see below). The mass and the damping are assembled at a step’s start and after every topology change; when a contact constraint merely changes its connectivity, they are reused in the rebuilt equation system instead, which matters for contact problems whose candidate lists change on nearly every increment.

Reduced-integration elements (C3D20R, GC3D20R, …) need a Marmot whose consistent mass is integrated with the full rule of the element shape (Marmot PR #101): with the element’s own reduced rule the mass is singular, and the equilibrium solve for the initial acceleration with it.

The material must report a density (for LinearElastic, the third material parameter): the mass is assembled from it, and a time-integrated degree of freedom that receives none is refused rather than integrated as though it had inertia.

A *modelModifier, type=hAdaptivity refining the mesh mid-step is supported: V and A are carried onto new nodes by the node-field warm start, the mass and damping are reassembled on the refined mesh, and the acceleration is re-solved from equilibrium on the increment that follows – see Model modifiers, “Refinement under a dynamic solver”, for what is conserved exactly across such an event and what is not.

Verified in tests/test_nid_newmark.py, tests/test_nid_restart.py and tests/test_nid_amr.py against the exact discrete Newmark recurrence of the same problem, against the continuous closed form at second order, for exact conservation of the discrete energy, for restart equivalence of U, V and A, and – under live refinement – for exact transfer of a velocity field of degree at most one onto new nodes together with exact conservation of mass, momentum and kinetic energy in that case.

NIDParallel - Nonlinear Implicit Dynamic (Newmark-beta, parallel)

Parallel implementation of the NID solver.

The dynamics add nothing to the element loop, so the parallel variant is the serial one with NISTParallel’s thread-parallel computeElements mixed in. The method resolution order NIDParallel -> NID -> NISTParallel -> NIST routes every super() call of NID through NISTParallel, so both the Newton loop and the initial-acceleration solve evaluate the elements in parallel.

class edelweissfe.solvers.nonlinearimplicitdynamicparallel.NIDParallel(jobInfo, journal, **kwargs)

This is the parallel Nonlinear Implicit Dynamic – solver (NIDParallel).

Parameters:

• jobInfo – A dictionary containing the job information.

• journal – The journal instance for logging.

Attributes:

convergenceTable

equationSystemLayout

incrementStatus

mpcTransformation

Methods

acceptIncrement(step, model, timeStep)	Commit the converged increment to the model, and keep the predictor's state; see acceptIncrement().
applyCorrection(dU, ddU, stepLength)	Move the iterate by the correction, scaled by the step length \(\alpha\) a line search chose (1 without a line search): \(\Delta \boldsymbol{U} \leftarrow \Delta \boldsymbol{U} + \alpha \, \Delta\Delta \boldsymbol{U}\).
applyDirichletToResidual(timeStep, R, dirichlets)	Impose the Dirichlet BCs on the residual.
applyDirichletToStiffness(K, dirichlets[, rhs])	Impose the Dirichlet BCs on the system matrix -- and, given the right-hand side(s) rhs (modified in place), eliminate the constrained DOFs from the columns as well; see edelweissfe.solvers.base.dirichlet.applyDirichletToStiffness().
applyOptionsOverride(fieldValues)	Apply an override of this solver's options, see applyOptionsOverride(); an override that sets an option where it means nothing is refused before anything of it applies (see refuseOptionsThatMeanNothing()).
applyStepActionsAtIncrementStart(model, ...)	Called when all step actions should be applied at the start of a step.
applyStepActionsAtStepEnd(model, stepActions)	Called when all step actions should finish a step.
applyStepActionsAtStepStart(model, step)	Apply every step action's step-start part (initial conditions, material initialization, prescribed fields, ...).
assembleAdditionalTerms(dU, R, F, K, timeStep)	Add the inertia and damping forces to the residual and their derivatives to the tangent; see the parent's hook.
assembleConstraintForces(constraints, U_np, ...)	Evaluate every active constraint and add its residual contribution to the given vector; a constraint that is switched off (see active) is skipped.
assembleConstraints(constraints, U_np, dU, ...)	Loop over all constraints, and evaluate the active ones: each adds its residual contribution and the tangent of it, see computeKernels().
assembleInternalForces(U_np, dU, P, K, F, ...)	Assemble the internal forces of the model at a trial state, and their tangent: those of the elements of its equation system, see computeElements().
assembleLoads(loads, U_np, dU, P, K, timeStep)	Assemble the residual contributions of all loads: each adds the negative of the external load it applies, and the tangent of that contribution, which is nonzero only for a load that follows the deformation; see The residual and the tangent.
assembleLoadsAtIncrementEnd(stepActions, ...)	Assemble the residual contributions of the loads acting at the end of the increment, and their tangents; called by assembleResidualAndTangent().
assembleResidualAndTangent(U_n, dU, U_np, R, ...)	Evaluate the model at the trial solution \(U_{n+1} = U_n + \Delta U\): the right-hand side \(R = f_\mathrm{ext} - f_\mathrm{int}\) (see The residual and the tangent) with the terms of assembleAdditionalTerms(), condensed by the multi-point constraints, and the tangent.
assembleStiffnessCSR(K)	The assembled system matrix in CSR format, see assembledMatrix().
assembleSystemMatrix(K, quasiNewton)	The system matrix of this Newton iteration, before the Dirichlet boundary conditions are imposed: the tangent in CSR layout, condensed by the multi-point constraints -- or, once the increment has fallen back to the quasi-Newton method, the tangent of its first iteration.
attemptIncrement(step, model, timeStep)	Solve the increment by Newton's method, rebuilding the equation system first if the topology update changed the model, and re-equilibrating a changed model first if so configured; see attemptIncrement().
beginStep(step, model, ...)	Start a step; see the parent.
buildMPCTransformation(model[, stepActions])	Collect the linear dependency records from all multi-point constraints of the model and assemble the master-slave condensation operator for the current equation system.
checkConvergence(R, ddU, dU, F, ...)	Check the convergence, individually for each field and for the scalar variables, based on the current flux residual and the correction, using the convergence criterion of the solver (see edelweissfe.solvers.base.convergencecriteria).
checkDivergingSolution(...)	Check if the iterative solution scheme is diverging: if the flux residual of a field grew, from one checked iteration to the next, more than maxGrowingIter times in this increment.
checkDivergingSolutionOrFallBack(...)	Cut an increment back whose residual grew more often than maxGrowingIter allows (see checkDivergingSolution()) -- or, with the option fallBackToQuasiNewton, let it continue with the tangent of its first iteration (a quasi-Newton method), see assembleSystemMatrix().
checkMPCDirichletConflicts(transformation, ...)	Raise if any Dirichlet boundary condition of the step prescribes a DOF that is a slave DOF of a multi-point constraint.
collectLoads(stepActions)	The loads among a step's actions: every action that is a LoadBase, whatever step action module it is defined in, so that a new kind of load needs no change to any solver.
computeElements(elements, Un1, dU, P, K, F, ...)	Compute the elements on the available threads, see computeElementsInParallel(); otherwise as computeElements().
createConvergenceCriterion()	Create the convergence criterion selected by the option convergenceCriterion.
createSystemMatrixAssembly(assemblyType)	The assembly of the system matrix of the equation system just built (option systemMatrixAssembly), on the complete CSR pattern of csrGenerator, with numberOfAssemblyWorkers() workers.
endStep(step, model)	Report the step's performance timing.
extrapolateLastIncrement(extrapolation, ...)	Predict the solution increment, the point Newton's method starts from, by extrapolating the last accepted increments of the step.
extrapolates(extrapolation, prevTimeStep)	Whether an increment that follows the accepted increment prevTimeStep starts from an extrapolated prediction.
extrapolationWeights(extrapolation, ...)	The weights of the last accepted increment and of the one before it in the prediction of extrapolateLastIncrement(): \(\Delta U = w_n \, \Delta U_n + w_{n-1} \, \Delta U_{n-1}\), with
finalizeIncrement(model)	The converged trial kinematics become the state of the solver, which the model receives at the acceptance of the increment (see writeResultsToModel()); see the parent's hook.
firstIterateNeedsCorrection(dirichlets)	Whether the first iterate of an increment without a prediction -- the converged state of the last increment -- needs a correction before it can be checked for convergence (see solveIncrement()): if Dirichlet boundary conditions act, whose prescribed increment the first correction makes.
forEachChunk(entities, work)	Call work(chunk) for chunks of entities on the available threads, as the elements are computed (see mapChunksInParallel()); otherwise as forEachChunk().
forgetEquationSystem()	Drop the equation system of the last build -- the model's layout of it with the DOF manager, the system matrix assembly and its CSR pattern, the multi-point-constraint condensation, the DOF vectors -- so that it is gone before the next one is built, never both at once.
forgetPredictionAfterFailure()	A failed increment drops the prediction, which may have caused the failure, unless retryWithoutExtrapolation has answered that already: then the accepted history is kept, and the smaller increment that follows the cutback extrapolates it.
getRestartData()	The state this solver carries from one increment to the next, as it is; see checkpointedState.
initializeIncrement(U_n, stepActions, model, ...)	Bring the Newmark operators up to date with the equation system, form the two terms the dynamics add to the tangent, and -- at a cold step start or after a topology change -- compute the acceleration from equilibrium; see the parent's hook.
isOutputIncrement(timeStep)	Whether field outputs, output managers and restart checkpoints are written after this increment.
lineSearchIsDue(iterationCounter)	Whether the correction of this iteration is scaled by a line search (lineSearch): after lineSearchAfterIterations iterations, in every iteration whose number is a multiple of lineSearchEveryIterations.
lineSearchStepLengths()	The trial step lengths of the line search, from the option lineSearchStepLengths.
linearSolve(A, b)	Solve the linear equation system.
locateConstrainedDofs(dirichlets)	Determine, up front, which global DOFs each Dirichlet BC constrains.
numberOfAssemblyWorkers()	How many threads write blocks of the system matrix at the same time: those that compute the elements.
predictionSurvivesRebuild(modelHasChanged, ...)	Whether this increment may extrapolate the last accepted increments after the equation system was rebuilt.
prepareIncrement(step, model, isRetry)	The topology update, at the start of every increment; see prepareIncrement().
printResidualOutlierNodes(residualOutliers)	Print which nodes have the largest residuals.
quadraticLineSearch(largestResidualAt)	The step length \(\alpha\) by which a Newton correction is scaled, for a Newton iteration that overshoots or oscillates (e.g., materials switching between loading and unloading).
readStateFromModel(U)	Read the state the next increment starts from out of the model, into the vectors of the equation system just (re)built: the solution, as in the parent, and the velocity and the acceleration of the fields carrying a mass, where the last accepted increment (or the checkpoint, or nothing -- a start from rest) left them, and where a refinement has meanwhile interpolated them onto the nodes it created.
refuseOptionsThatMeanNothing(options, ...)	Refuse an option that a declaration of options sets where it means nothing, with the options in force after the declaration.
reportTopologyChangeConservation(before, after)	Report what a topology change did to the quantities that ought to survive it: every dynamic field's total mass is checked (and raises when violated), momentum and kinetic energy are reported -- see conservationchecks for which of them is exact when.
schema
selectLinearSolverForStep()	Make linSolver the linear solver the step about to begin names (linsolver): the one of the previous step if it names the same, else a new one.
setRestartData(data)	Restore the state getRestartData() returned, so that the resumed step continues exactly where the uninterrupted run would have.
setTrialState(dU, ddU, stepLength)	Set the trial state of a line search at a step length \(\alpha\) along the correction, and return its solution increment, \(\Delta \boldsymbol{U} + \alpha \, \Delta\Delta \boldsymbol{U}\), at which the model is evaluated next; the iterate stays where it is.
solveForCorrection(K_, R, dirichlets)	Solve the linearized system of a Newton iteration, \(\boldsymbol{K} \, \Delta\Delta \boldsymbol{U} = \boldsymbol{R}\), for the correction \(\Delta\Delta \boldsymbol{U}\).
solveIncrement(U_n, dU, P, K, stepActions, ...)	Newton-Raphson scheme on the Newmark-effective equation of motion of an increment.
updateRigidBodies(model, timeStep)	Refresh the kinematics of all rigid bodies in the model after a converged increment.
updateTopologyAndConnectivity(model, step[, ...])	Run the topology update, then let every mesh-dependent consumer catch up on it.
validateModelCapabilities(model)	Validate whether the solver supports the active features/constraints of the model.
writeResultsToModel(U, dU, P)	Write the results of a converged increment to the model, as the parent does, and the velocity and the acceleration of the fields carrying a mass with them: what field outputs read, what a checkpoint stores, and what the next equation system reads back.

NIDParallel

Option	Type / Default	Description
defaultMaxIter	int, default=10	The default maximum number of iterations.
defaultCriticalIter	int, default=5	The default number of critical iterations.
defaultMaxGrowingIter	int, default=10	The default number of allowed residual growths.
extrapolation	str, default=’linear’	The prediction a new increment starts its Newton iterations from (off|linear|quadratic): none, the last accepted increment scaled to the new size, or the parabola through the last three accepted states (linear while only one accepted increment of the step of nonzero size is available).
retryWithoutExtrapolation	bool, default=False	If an extrapolated increment fails, attempt it once more from a zero prediction before cutting it back. A failed increment then also keeps the accepted history the next attempt extrapolates from.
extrapolateAcrossConnectivityChange	bool, default=False	Keep extrapolating across a rebuild of the equation system caused by a change of the connectivity alone (e.g., a contact search), as long as the rebuilt system numbers the degrees of freedom as before. Off: every rebuild starts the increment from a zero prediction.
extrapolateAfterModelChange	bool, default=True	Whether to extrapolate the predictor on the increment FOLLOWING a model change (adaptive mesh refinement). Set False to start that increment from a zero predictor, avoiding extrapolation of the one-off warm-start/remesh settling transient.
equilibrateAfterModelChange	bool, default=False	Whether to insert one constant-load, zero-time re-equilibration increment immediately after an adaptive mesh refinement, before advancing the load. When True, the warm-started refined mesh is first settled to equilibrium at the last converged load level (no load advance, no Dirichlet increment, zero time increment) so the subsequent load-advancing increment starts from an equilibrated state. Intended for softening problems where remeshing near the process zone otherwise couples the load advance with the warm-start settling transient in one solve. Under arc length control (NISTPArcLength with a controller), the re-equilibration is an arc length increment of size zero: the controlled quantity is held where the last accepted increment left it, and the load parameter is solved for, and with it the loads and the prescribed displacements it scales. Note: the equilibration solve integrates materials with dT=0, which suits rate-independent models; rate-dependent materials see no time advance during it (by design).
linsolver	str, default=’pardiso’	The linear solver: the name of a *linearSolver block, or a type of linear solver with its defaults (e.g. pardiso, superlu, blockamg).
pruneCondensedMatrixZeros	bool, default=None	Compact explicitly stored zeros out of the multi-point-constraint-condensed system matrix before solving. Pruning removes whichever entries happen to be exactly zero on this particular iteration, so the pattern changes on almost every Newton iteration – which defeats any linear solver that relies on a stable pattern across iterations (e.g. blockamg reusing a per-field AMG hierarchy built for the old one). Left unset (the default), the choice is automatic: pruned if the linear solver in force reports that it needs it (:attr:`~edelweissfe.linsolve.base.LinearSolver.needsPrunedSystemMatrix`), which a direct factorization does – PARDISO’s reordering is sensitive to the extra structural entries on these path-dependent condensed systems, and keeping them has been observed to drift from the converged reference path – and a Krylov solver with an AMG preconditioner does not: it gives the same iterations with and without the stored zeros, and pruning costs it a copy of the matrix and a pass over it in every Newton iteration. An explicit True/False here always wins over the automatic choice; turn pruning off for a direct solver only after verifying the load path against a reference run.
systemMatrixAssembly	str, default=’vij’	How the blocks of the elements, constraints and loads are summed into the system matrix (vij|direct): staged in a value array of the VIJ pattern’s length and gathered into CSR, or scattered into CSR right away, without that array. See the documentation page on the system matrix assembly.
useThreadedMPCCondensation	bool, default=False	Condense the multi-point-constraint system matrix, T^T K T + C, in place, with an OpenMP-threaded kernel, instead of with SciPy’s single-threaded CSR routines, which build a second matrix. The CSR pattern of K gains, once per equation system, the entries where a master collects the couplings of its slaves; the result does not depend on the number of threads. Explicit zeros (the slave columns among them) are pruned from a copy according to pruneCondensedMatrixZeros. Off by default.
report-performance-frequency	int, default=0	Print the cumulative performance table (totals since the start of the step) every this many increments, in addition to the one printed at the end of the step. 0 (default) disables the periodic report.
lineSearch	bool, default=False	Scale the Newton correction by the step length at which a parabola, fitted to the largest residual evaluated at the trial step lengths lineSearchStepLengths, is smallest (a quadratic line search).
lineSearchAfterIterations	int, default=5	The line search starts after this many iterations of an increment, at least 0: the first correction of an increment, which makes the prescribed increments, is never scaled.
lineSearchEveryIterations	int, default=2	From then on, the line search is done in every iteration whose number is a multiple of this.
lineSearchStepLengths	str, default=’0.8, 1.0, 1.2’	The trial step lengths of the line search, comma-separated, at least three.
fallBackToQuasiNewton	bool, default=False	Instead of cutting an increment back once the residual grew more often than maxGrowIter allows, continue it with the tangent of its first iteration (a quasi-Newton method); maxIter still applies.
convergenceCriterion	str, default=’legacy’	The convergence criterion of the Newton iterations, legacy (the original EdelweissFE one) or abaqus (the one of Abaqus/Standard); see edelweissfe.solvers.base.convergencecriteria.
iterationsForAlternativeTolerance	int, default=15	legacy criterion: I, the iteration from which on the flux residual is checked against fluxResidualToleranceAlternative instead of fluxResidualTolerance.
zeroFluxThreshold	float, default=None	legacy criterion: a field whose spatial average flux does not exceed this value carries no flux its residual could be compared with, so its flux residual counts as converged; needs minimumIterations >= 1, and refuses a smaller one. Off if not set.
zeroIncrementThreshold	float, default=None	legacy criterion: a field whose largest increment does not exceed this value has not moved, so there is no increment to relate its correction to, and its correction counts as converged. Off if not set.
minimumIterations	int, default=0	legacy criterion: no flux residual counts as converged before this iteration, so an increment takes at least this many Newton iterations; 0 accepts a prediction that already satisfies the criterion.
countFirstResidualAsGrowth	bool, default=True	Residual growths (the step option maxGrowIter): the residual history of a field starts from zero, so the first checked residual of an increment counts as a growth; if False, counting starts with the second.
countGrowthOfConvergedFields	bool, default=True	Residual growths (the step option maxGrowIter): count a growth of the flux residual of a field also in an iteration in which this field is converged and only other fields keep the iterations going.
abaqusResidualTolerance	float, default=0.005	abaqus criterion: R_n, the flux residual tolerance relative to the time average flux.
abaqusCorrectionTolerance	float, default=0.01	abaqus criterion: C_n, the correction tolerance relative to the largest increment of the field.
abaqusResidualToleranceAlternative	float, default=0.02	abaqus criterion: R_p, the flux residual tolerance used from iteration I_p on.
abaqusIterationsForAlternativeTolerance	int, default=9	abaqus criterion: I_p, the iteration from which on R_p is used.
abaqusZeroFluxThreshold	float, default=1e-05	abaqus criterion: epsilon, below which (relative to the time average flux) a field carries no flux.
abaqusZeroFluxCorrectionTolerance	float, default=0.001	abaqus criterion: C_epsilon, the relative correction tolerance of a field without flux.
abaqusLinearResidualTolerance	float, default=1e-08	abaqus criterion: R_l, below which a field is converged irrespective of its correction.
abaqusInitialTimeAveragedFlux	float, default=0.01	abaqus criterion: q~0, the time average flux of a field that has not carried flux yet (first step).
newmarkBeta	float, default=0.25	Newmark parameter beta, weighting the end-of-increment acceleration in the displacement update. Must be > 0. 1/4 with gamma = 1/2 is the average-acceleration (trapezoidal) rule.
newmarkGamma	float, default=0.5	Newmark parameter gamma, weighting the end-of-increment acceleration in the velocity update. 1/2 gives second-order accuracy and no algorithmic damping; larger values damp.
computeInitialAcceleration	bool, default=True	Compute the acceleration at the start of every step from equilibrium with the loads at the step’s start, instead of continuing from the carried-over acceleration. Skipped automatically when the step resumes from a restart checkpoint. Also re-armed for the increment following a topology change (h-adaptivity), whose interpolated acceleration is not in equilibrium with the operators reassembled on the refined mesh.
massConservationTolerance	float, default=1e-06	Relative-change tolerance for the per-topology-change check on a dynamic field’s total assembled mass. Children of a refined element tile it and carry the same density, so this is conserved geometrically; the default is already generous relative to floating-point precision and exists to catch a genuinely wrong refinement, not to absorb ordinary quadrature noise. The counterpart of the explicit solver’s lumped-quantity-conservation-tolerance.

Prediction of an increment in the implicit solvers

Newton’s method starts each increment from a prediction of its solution increment. The solver option extrapolation chooses it:

• off: from zero, i.e., from the last converged state.

• linear (default): the last accepted increment, scaled to the size of the new one. The path of the solution is continued along its last secant.

• quadratic: the parabola through the last three accepted states, evaluated at the end of the new increment, which continues the path along its curvature, too. For a solution that grows quadratically with the load, the prediction is exact, and the increment converges without a correction. While only one increment of the step of nonzero size has been accepted, the prediction is linear.

Sizes are measured in step progress. The prescribed degrees of freedom take their prescribed increment, so the first iteration of a predicted increment already checks convergence. Under arc length control, the increment of the load parameter is predicted from the same accepted increments with the same weights, see Arc length control.

A prediction can lead Newton’s method astray, e.g., when the path turns sharply. Two options decide what happens then and around rebuilds of the equation system, both off by default:

• retryWithoutExtrapolation: if an extrapolated increment fails, it is attempted once more from a zero prediction before it is cut back. A failed increment then keeps the accepted history, and the smaller increment after the cutback is extrapolated from it as well. Without this option, a failed increment drops the prediction, and the increments up to the next accepted one start from zero.

• extrapolateAcrossConnectivityChange: a change of the connectivity, e.g., a contact search that couples other nodes, rebuilds the equation system, and the increment starts from zero. With this option, it is extrapolated nevertheless, as long as the rebuilt system numbers the degrees of freedom exactly as the one the last increments were accepted in (compared by the digests of the two numberings, dofNumberingDigest): the last accepted increments are vectors of the same degrees of freedom then. A change of the topology, e.g., a mesh refinement, always starts from zero (see extrapolateAfterModelChange for the increment that follows it).

Both act on a prediction: a declaration of options – the *solver datalines, or an >>options block – that switches either on with extrapolation=off is refused.

The last accepted increments are checkpointed with the digest of their numbering. A step that continues from a restart checkpoint, resumed or fallen back to (see Restart), rebuilds its equation system at its first increment, and extrapolates them only if the rebuilt system numbers the degrees of freedom alike: the numbering of a meshfree model, e.g., follows its active cells, which may have changed since the checkpoint was written.

Line search and quasi-Newton fallback of the implicit solvers

Newton’s method converges quadratically close to the solution, but far from it, or where the tangent changes abruptly (e.g., material points switching between loading and unloading), a full correction can overshoot, and the iterations oscillate or diverge. Two remedies, solver options of the implicit solvers, are off by default. Under arc length control, both act on the pair of the solution and the load parameter, see Arc length control.

Line search (lineSearch). The correction \(\Delta\Delta U\) is scaled by a step length \(\alpha\): the largest residual \(\| R(\Delta U + \alpha \, \Delta\Delta U) \|_\infty\) (without the Dirichlet degrees of freedom) is evaluated at the trial step lengths lineSearchStepLengths (default 0.8, 1.0, 1.2), a parabola is fitted to these values, and \(\alpha\) minimizes the parabola on \([0.1, 1.5]\): its vertex, limited to that range, or, for a parabola that does not open upwards, the end of the range where it is smaller. A trial at which the residual is not finite fits no parabola: \(\alpha\) is then the trial step length of the smallest finite residual, and without one, the increment is cut back. Trial residuals that are equal to round-off do not tell the step lengths apart, and the correction is taken in full. Each trial costs one evaluation of the model, so the line search is done only after lineSearchAfterIterations iterations (default 5), and then in every iteration whose number is a multiple of lineSearchEveryIterations (default 2). Small trial step lengths, e.g. lineSearchStepLengths='0.15, 0.4, 0.7, 1.0', damp an iteration that cycles between two states.

Quasi-Newton fallback (fallBackToQuasiNewton). An increment whose residual grew more often than the step option maxGrowIter allows is cut back. With this option, it continues instead with the tangent of its first iteration, i.e., of the converged state at the start of the increment, or of the prediction if the increment is extrapolated. The convergence becomes linear, but a tangent that is not updated cannot switch abruptly either; maxIter still ends an increment that does not converge.

>>options, name=theSolver, lineSearch=True, lineSearchAfterIterations=2, lineSearchEveryIterations=1, fallBackToQuasiNewton=True

Convergence criteria of the implicit solvers

The implicit solvers (NIST, NISTParallel, NISTPArcLength, NID, NIDParallel) select the criterion of their Newton iterations with the solver option convergenceCriterion, e.g.:

>>options, name=theSolver, convergenceCriterion=abaqus

The options of the abaqus criterion are named abaqus...; a declaration of options – the *solver datalines, or an >>options block – that sets one of them with convergenceCriterion=legacy is refused, since it would mean nothing.

The legacy criterion is configured by tolerances per field and by solver options. EdelweissMeshfree’s parameter set, for instance, reads for a model with a displacement field

*updateConfiguration, configuration=fluxResidualTolerance
displacement=1e-4
*updateConfiguration, configuration=fluxResidualToleranceAlternative
displacement=1e-3
*updateConfiguration, configuration=absoluteFluxResidualTolerance
displacement=1e-14
*updateConfiguration, configuration=fieldCorrectionTolerance
displacement=1e-14
*updateConfiguration, configuration=relativeFieldCorrectionTolerance
displacement=1e-3

together with

>>options, name=theSolver, iterationsForAlternativeTolerance=5, zeroFluxThreshold=1e-9, zeroIncrementThreshold=1e-13, minimumIterations=1, countFirstResidualAsGrowth=False, countGrowthOfConvergedFields=False

The last two options belong to the test for a diverging Newton iteration, which holds for either criterion: an increment is abandoned, and retried smaller, once the flux residual of a field has grown more than maxGrowIter (a step option) times in it. A growth is a flux residual larger than that of the previous checked iteration, and the two options decide what else counts:

• countFirstResidualAsGrowth (default True): the residual history of a field starts at zero, so the first checked residual of an increment is a growth already. This is EdelweissFE’s original counting; with False, counting starts with the second checked iteration.

• countGrowthOfConvergedFields (default True): a growth counts even in an iteration in which the field is converged, and only another field keeps the iterations going. With False, such growths are ignored: the residual of a converged field fluctuates at the level of its round-off, which says nothing about divergence.

Independent of the criterion, a residual that contains a NaN ends the increment at once.

Field-wise convergence criteria for the Newton-Raphson iterations of the implicit solvers.

Two criteria are available, selected by the solver option convergenceCriterion:

legacy (default)

EdelweissFE’s own criterion. Field by field, it compares the largest flux residual \(r_{\max}\) with the spatial average flux \(\bar{q}\) of the current iteration, computed from the element fluxes only, and the largest correction \(c_{\max}\) with an absolute tolerance and with a tolerance relative to \(\Delta u_{\max}\), the largest increment of the field in the current increment. A field is converged if

\[\left( r_{\max} \leq R_a \;\text{or}\; r_{\max} \leq R \, \bar{q} \right) \quad \text{and} \quad \left( c_{\max} < C_a \;\text{or}\; c_{\max} < C \, \Delta u_{\max} \right),\]

with the tolerances of the field:

• \(R\), the flux residual tolerance relative to the spatial average flux (fluxResidualTolerance), replaced by fluxResidualToleranceAlternative from iteration \(I\) (iterationsForAlternativeTolerance) on, where a Newton iteration that converges only slowly may stop at a looser tolerance;

• \(R_a\), the absolute flux residual tolerance (absoluteFluxResidualTolerance): a floor for a field whose flux is so small that a relative tolerance would ask for the round-off of the residual;

• \(C_a\), the absolute correction tolerance (fieldCorrectionTolerance);

• \(C\), the correction tolerance relative to the largest increment (relativeFieldCorrectionTolerance); the default 0 accepts a correction by its absolute size only.

The tolerances of each field are set with *updateConfiguration, for instance

*updateConfiguration, configuration=absoluteFluxResidualTolerance
displacement=1e-10

and default to edelweissfe.config.phenomena. Three rules, set by solver options, are off by default:

• Zero flux (zeroFluxThreshold, \(\bar{q}_0\)): a field whose spatial average flux does not exceed \(\bar{q}_0\) carries no flux its residual could be compared with, so its flux residual counts as converged. An example is a nonlocal damage field before the onset of damage.

• Zero increment (zeroIncrementThreshold, \(\Delta u_0\)): a field whose largest increment does not exceed \(\Delta u_0\) has not moved, so there is no increment to relate its correction to; its correction counts as converged.

• Minimum iterations (minimumIterations, \(n\)): no flux residual counts as converged before iteration \(n\), so an increment takes at least \(n\) Newton iterations. With \(n = 0\), a prediction that already satisfies the criterion is accepted as it is. The zero flux rule needs \(n \geq 1\), and refuses a smaller \(n\): before the first correction, the elements of a model at rest report no flux at all.

The scalar variables (e.g., Lagrange multipliers) are checked like a field, with the tolerances of the entry scalar variables. No element flux acts on them, so their flux residual is compared with the tolerances directly, as if their spatial average flux were one, and the zero flux rule does not apply to them: it would accept any residual of theirs for a threshold of one or more.

Two parameter sets of this criterion are in use: EdelweissFE’s defaults, which are those of the original criterion, and the settings of EdelweissMeshfree’s solvers (for every field):

parameter	EdelweissFE	EdelweissMeshfree
fluxResidualTolerance	per field, see phenomena	1e-4
fluxResidualToleranceAlternative	per field, see phenomena	1e-3
iterationsForAlternativeTolerance	15	5
absoluteFluxResidualTolerance	1e-7, scalar variables 0	1e-14
fieldCorrectionTolerance	per field, see phenomena	1e-14
relativeFieldCorrectionTolerance	0	1e-3
zeroFluxThreshold	off	1e-9
zeroIncrementThreshold	off	1e-13
minimumIterations	0	1

EdelweissMeshfree also counts the residual growths of a diverging iteration differently: its solvers set the solver options countFirstResidualAsGrowth and countGrowthOfConvergedFields to False.

abaqus

The criterion of Abaqus/Standard (“Convergence criteria for nonlinear problems”, Abaqus Analysis User’s Guide), with its default parameters as printed in the .msg file of an Abaqus run without *CONTROLS. A field is converged if

\[r_{\max} \leq R_n \, \tilde{q} \quad \text{and} \quad c_{\max} \leq C_n \, \Delta u_{\max},\]

with \(R_n\) replaced by \(R_p\) from iteration \(I_p\) on. \(\Delta u_{\max}\) is the largest increment of the field in the current increment, so the correction criterion is relative. \(\tilde{q}\) is the time average flux: the spatial average flux \(\bar{q}\) averaged over all converged increments of the step in which the field carried flux, and the current iteration. The spatial average \(\bar{q}\) includes the fluxes of constraints (e.g., contact forces). Until a field has carried flux, \(\tilde{q}\) is the time average flux of the previous step or, in the first step, the initial value \(\tilde{q}^0\) (Abaqus’ “initial value of time average force”, 1e-2).

A field whose current spatial average flux is below \(\varepsilon \, \tilde{q}\) carries no flux; it is converged if \(c_{\max} \leq C_\varepsilon \, \Delta u_{\max}\). A field with \(r_{\max} \leq R_l \, \tilde{q}\) is converged irrespective of the correction (linear case).

The criterion is formulated for fields; the scalar variables are checked by the legacy criterion.

Not implemented: the estimate of the next correction from the convergence rate (Abaqus accepts an increment if it predicts a small enough correction), and the detection of inactive regions for the spatial average. Both make Abaqus’ criterion more permissive than this implementation.

class edelweissfe.solvers.base.convergencecriteria.AbaqusConvergenceCriterion(scalarVariableCriterion, residualTolerance=0.005, correctionTolerance=0.01, residualToleranceAlternative=0.02, iterationsForAlternativeTolerance=9, zeroFluxThreshold=1e-05, zeroFluxCorrectionTolerance=0.001, linearResidualTolerance=1e-08, initialTimeAveragedFlux=0.01)

The criterion of Abaqus/Standard; see the module documentation.

Parameters:

• scalarVariableCriterion (ConvergenceCriterionBase) – The criterion that checks the scalar variables, which Abaqus’ criterion does not cover: the legacy one.

• residualTolerance (float) – \(R_n\), the flux residual tolerance relative to the time average flux.

• correctionTolerance (float) – \(C_n\), the correction tolerance relative to the largest increment of the field.

• residualToleranceAlternative (float) – \(R_p\), the flux residual tolerance used from iteration \(I_p\) on.

• iterationsForAlternativeTolerance (int) – \(I_p\).

• zeroFluxThreshold (float) – \(\varepsilon\), below which a field is considered to carry no flux, relative to the time average flux.

• zeroFluxCorrectionTolerance (float) – \(C_\varepsilon\), the correction tolerance of a field without flux.

• linearResidualTolerance (float) – \(R_l\), below which a field is converged irrespective of its correction.

• initialTimeAveragedFlux (float) – \(\tilde{q}^0\), the time average flux of a field that has not carried flux yet, in the first step.

Methods

acceptIncrement()	Called once an increment is converged and accepted.
checkField(field, R, ddU, dU, ...)	Check the convergence of a single field.
checkScalarVariables(R, ddU, dU, ...)	Check the convergence of the scalar variables (e.g., Lagrange multipliers), together as one column.
computeSpatialAveragedFluxes(F, ...)	Compute the spatial average flux of every field.
startStep(previousCriterion)	Called at the beginning of each step.
timeAveragedFlux(field, spatialAveragedFlux)	The time average flux \(\tilde{q}\), including the current iteration if the field carries flux.
timeAveragedFluxesAtStepEnd()	The time average fluxes at the end of the step, used as reference by the criterion of the next step.

acceptIncrement()

Called once an increment is converged and accepted.

checkField(field, R, ddU, dU, spatialAveragedFlux, iterationCounter)

Check the convergence of a single field.

Parameters:

• field – The name of the field.

• R – The residual of the field.

• ddU – The current correction of the field, None if no correction was computed yet in this increment.

• dU – The current increment of the field.

• spatialAveragedFlux – The spatial average flux \(\bar{q}\) of the field in the current iteration.

• iterationCounter – The current iteration number.

Returns:

The convergence state of the field.

Return type:

FieldConvergence

checkScalarVariables(R, ddU, dU, iterationCounter)

Check the convergence of the scalar variables (e.g., Lagrange multipliers), together as one column.

Parameters:

• R – The residual of the scalar variables.

• ddU – The current correction of the scalar variables, None if no correction was computed yet in this increment.

• dU – The current increment of the scalar variables.

• iterationCounter – The current iteration number.

Returns:

The convergence state of the scalar variables.

Return type:

FieldConvergence

computeSpatialAveragedFluxes(F, FConstraints, theDofManager)

Compute the spatial average flux of every field.

Parameters:

• F – The accumulated absolute element fluxes.

• FConstraints – The accumulated absolute constraint fluxes.

• theDofManager – The DofManager.

Returns:

The spatial average flux of every field.

Return type:

dict[str,float]

startStep(previousCriterion)

Called at the beginning of each step.

Parameters:

previousCriterion – The criterion of the previous step, if any.

timeAveragedFlux(field, spatialAveragedFlux)

The time average flux \(\tilde{q}\), including the current iteration if the field carries flux.

Parameters:

• field (str) – The name of the field.

• spatialAveragedFlux (float) – The spatial average flux \(\bar{q}\) of the field in the current iteration.

Returns:

The time average flux.

Return type:

float

timeAveragedFluxesAtStepEnd()

The time average fluxes at the end of the step, used as reference by the criterion of the next step.

Returns:

The time average flux per field; empty for criteria without time averaging.

Return type:

dict[str,float]

checkpointedState: dict | None = {'_nIncrementsWithFlux': <class 'dict'>, '_sumOfConvergedSpatialAveragedFluxes': <class 'dict'>, '_timeAveragedFluxOfPreviousStep': <class 'dict'>}

The time average flux of every field, as the converged increments of the step and the previous step left it.

class edelweissfe.solvers.base.convergencecriteria.ConvergenceCriterionBase

Base class for the field-wise convergence criteria of the implicit solvers.

Attributes:

checkpointedState

Methods

acceptIncrement()	Called once an increment is converged and accepted.
checkField(field, R, ddU, dU, ...)	Check the convergence of a single field.
checkScalarVariables(R, ddU, dU, ...)	Check the convergence of the scalar variables (e.g., Lagrange multipliers), together as one column.
computeSpatialAveragedFluxes(F, ...)	Compute the spatial average flux of every field.
startStep(previousCriterion)	Called at the beginning of each step.
timeAveragedFluxesAtStepEnd()	The time average fluxes at the end of the step, used as reference by the criterion of the next step.

acceptIncrement()

Called once an increment is converged and accepted.

abstractmethod checkField(field, R, ddU, dU, spatialAveragedFlux, iterationCounter)

Check the convergence of a single field.

Parameters:

• field (str) – The name of the field.

• R (ndarray) – The residual of the field.

• ddU (ndarray | None) – The current correction of the field, None if no correction was computed yet in this increment.

• dU (ndarray) – The current increment of the field.

• spatialAveragedFlux (float) – The spatial average flux \(\bar{q}\) of the field in the current iteration.

• iterationCounter (int) – The current iteration number.

Returns:

The convergence state of the field.

Return type:

FieldConvergence

abstractmethod checkScalarVariables(R, ddU, dU, iterationCounter)

Check the convergence of the scalar variables (e.g., Lagrange multipliers), together as one column.

Parameters:

• R (ndarray) – The residual of the scalar variables.

• ddU (ndarray | None) – The current correction of the scalar variables, None if no correction was computed yet in this increment.

• dU (ndarray) – The current increment of the scalar variables.

• iterationCounter (int) – The current iteration number.

Returns:

The convergence state of the scalar variables.

Return type:

FieldConvergence

abstractmethod computeSpatialAveragedFluxes(F, FConstraints, theDofManager)

Compute the spatial average flux of every field.

Parameters:

• F (DofVector) – The accumulated absolute element fluxes.

• FConstraints (DofVector) – The accumulated absolute constraint fluxes.

• theDofManager (DofManager) – The DofManager.

Returns:

The spatial average flux of every field.

Return type:

dict[str,float]

startStep(previousCriterion)

Called at the beginning of each step.

Parameters:

previousCriterion (ConvergenceCriterionBase | None) – The criterion of the previous step, if any.

timeAveragedFluxesAtStepEnd()

The time average fluxes at the end of the step, used as reference by the criterion of the next step.

Returns:

The time average flux per field; empty for criteria without time averaging.

Return type:

dict[str,float]

checkpointedState: dict | None = None

The state this criterion carries from one increment to the next, checkpointed with the solver; see checkpointedstate. Every criterion declares it, an empty mapping if it carries nothing: one whose declaration is still None cannot be checkpointed.

class edelweissfe.solvers.base.convergencecriteria.FieldConvergence(fluxResidual, fluxResidualConverged, correction, correctionConverged, indexOfLargestResidual)

The convergence state of one field in one iteration.

Parameters:

• fluxResidual (float) – The largest absolute flux residual \(r_{\max}\).

• fluxResidualConverged (bool) – True if the flux residual criterion is satisfied.

• correction (float) – The largest absolute correction \(c_{\max}\).

• correctionConverged (bool) – True if the correction criterion is satisfied.

• indexOfLargestResidual (int) – The index of \(r_{\max}\) in the dof vector of the field.

class edelweissfe.solvers.base.convergencecriteria.LegacyConvergenceCriterion(fluxResidualTolerances, fluxResidualTolerancesAlternative, absoluteFluxResidualTolerances, fieldCorrectionTolerances, relativeFieldCorrectionTolerances, iterationsForAlternativeTolerance, zeroFluxThreshold, zeroIncrementThreshold, minimumIterations)

EdelweissFE’s own criterion; see the module documentation.

Parameters:

• fluxResidualTolerances (dict[str, float]) – \(R\) per field: the flux residual tolerances relative to the spatial average flux.

• fluxResidualTolerancesAlternative (dict[str, float]) – The relative flux residual tolerances per field from iteration \(I\) on.

• absoluteFluxResidualTolerances (dict[str, float]) – \(R_a\) per field: the absolute flux residual tolerances.

• fieldCorrectionTolerances (dict[str, float]) – \(C_a\) per field: the absolute correction tolerances.

• relativeFieldCorrectionTolerances (dict[str, float]) – \(C\) per field: the correction tolerances relative to the largest increment of the field.

• iterationsForAlternativeTolerance (int) – \(I\), the iteration from which on the alternative flux residual tolerances apply.

• zeroFluxThreshold (float | None) – \(\bar{q}_0\), the spatial average flux up to which a field carries no flux; None for no zero flux rule.

• zeroIncrementThreshold (float | None) – \(\Delta u_0\), the largest increment up to which a field has not moved; None for no zero increment rule.

• minimumIterations (int) – \(n\), the iteration before which no flux residual counts as converged.

Raises:

ValueError – If the zero flux rule is on with \(n < 1\).

Methods

acceptIncrement()	Called once an increment is converged and accepted.
checkField(field, R, ddU, dU, ...)	Check the convergence of a single field.
checkScalarVariables(R, ddU, dU, ...)	Check the convergence of the scalar variables (e.g., Lagrange multipliers), together as one column.
computeSpatialAveragedFluxes(F, ...)	Compute the spatial average flux of every field.
correctionConverged(field, correction, dU)	The correction part of the criterion.
fluxResidualConverged(field, fluxResidual, ...)	The flux residual part of the criterion.
startStep(previousCriterion)	Called at the beginning of each step.
timeAveragedFluxesAtStepEnd()	The time average fluxes at the end of the step, used as reference by the criterion of the next step.

checkField(field, R, ddU, dU, spatialAveragedFlux, iterationCounter)

Check the convergence of a single field.

Parameters:

• field – The name of the field.

• R – The residual of the field.

• ddU – The current correction of the field, None if no correction was computed yet in this increment.

• dU – The current increment of the field.

• spatialAveragedFlux – The spatial average flux \(\bar{q}\) of the field in the current iteration.

• iterationCounter – The current iteration number.

Returns:

The convergence state of the field.

Return type:

FieldConvergence

checkScalarVariables(R, ddU, dU, iterationCounter)

Check the convergence of the scalar variables (e.g., Lagrange multipliers), together as one column.

Parameters:

• R – The residual of the scalar variables.

• ddU – The current correction of the scalar variables, None if no correction was computed yet in this increment.

• dU – The current increment of the scalar variables.

• iterationCounter – The current iteration number.

Returns:

The convergence state of the scalar variables.

Return type:

FieldConvergence

computeSpatialAveragedFluxes(F, FConstraints, theDofManager)

Compute the spatial average flux of every field.

Parameters:

• F – The accumulated absolute element fluxes.

• FConstraints – The accumulated absolute constraint fluxes.

• theDofManager – The DofManager.

Returns:

The spatial average flux of every field.

Return type:

dict[str,float]

correctionConverged(field, correction, dU)

The correction part of the criterion.

Parameters:

• field (str) – The name of the field.

• correction (float) – The largest absolute correction \(c_{\max}\) of the field.

• dU (ndarray) – The current increment of the field.

Returns:

True if the correction is converged.

Return type:

bool

fluxResidualConverged(field, fluxResidual, spatialAveragedFlux, iterationCounter)

The flux residual part of the criterion.

Parameters:

• field (str) – The name of the field.

• fluxResidual (float) – The largest absolute flux residual \(r_{\max}\) of the field.

• spatialAveragedFlux (float) – The spatial average flux \(\bar{q}\) of the field; None for the scalar variables, on which no element flux acts: their flux residual is compared with the tolerances directly, as if \(\bar{q}\) were one, and the zero flux rule, which tells a field that carries no element flux, does not apply to them.

• iterationCounter (int) – The current iteration number.

Returns:

True if the flux residual is converged.

Return type:

bool

checkpointedState: dict | None = {}

every decision is taken from the current iteration alone.

Type:

Nothing

Iteration table of the implicit solvers

Every Newton iteration writes one row to the console: per field, the largest flux residual and the largest correction, each marked ✓ once the criterion is met. With a linear solver that reports per-solve diagnostics (e.g. blockamg at verbosity=info), a linear solve column follows with the solver’s iteration count and achieved residual. The rows are never wrapped, so a model with many fields gives long lines rather than broken ones.

   displacement            nonlocal damage              linear solve
||R||∞     ||ddU||∞      ||R||∞     ||ddU||∞         iters        ‖r‖
2.07e+01    8.35e-03     2.12e-02    8.15e-07          14        3.01e-04✓

The table of the Newton-Raphson iterations that the implicit solvers write to the journal.

One column per field (and one for the scalar variables, if any), each showing the largest flux residual \(\|R\|_\infty\) and the largest correction \(\|\Delta\Delta U\|_\infty\) of the iteration, marked ✓ once the convergence criterion is met. Linear solvers that report per-solve diagnostics (see reportsSolveSummary) get one more column with the solver’s iteration count and achieved residual.

The rows are written with preformatted(), so they are never wrapped, however many columns the model has.

class edelweissfe.solvers.base.convergencetable.ConvergenceTable(columnNames, showLinearSolve)

The column layout of the iteration table: the header lines, and one row per iteration.

Parameters:

• columnNames (list[str]) – The names of the field columns, in the order of the rows’ entries.

• showLinearSolve (bool) – Add the “linear solve” column.

Methods

headerLines()	The two header lines: the column names, and the labels of the values below them.
linearSolveDetailLine(text)	A line of linear solver detail, right-aligned inside the "linear solve" column, ending where the rows end.
row(fieldColumns, linearSolve)	The row of one iteration.

headerLines()

The two header lines: the column names, and the labels of the values below them.

Returns:

The column-name line and the label line.

Return type:

list[str]

linearSolveDetailLine(text)

A line of linear solver detail, right-aligned inside the “linear solve” column, ending where the rows end.

Parameters:

text (str) – The detail, at most one column less one character wide to stay inside the column.

Returns:

The line, blank in front of the “linear solve” column.

Return type:

str

row(fieldColumns, linearSolve)

The row of one iteration.

Parameters:

• fieldColumns (list[FieldConvergence]) – The convergence state of each column, in the order of the column names.

• linearSolve (LinearSolveSummary | None) – The diagnostics of the linear solve that produced this iteration’s correction, or None if there was none (the first iteration of an increment); the “linear solve” column is then left blank.

Returns:

The row.

Return type:

str

columnWidth = 25

two 12 characters wide cells (a value and its ✓ marker), and one space.

Type:

The width of one column

fieldCell = '{:11.2e}{:1}{:11.2e}{:1} '

The format of one field column in a row.

linearSolveGap = '   '

The space between the field columns and the “linear solve” column.

retriesSuperscript = {0: '', 1: '¹', 2: '²', 3: '³'}

The marker of a linear solve that needed extra internal attempts, by the number of attempts.

subHeaderCell = '{:^12}{:^12} '

The format of the two labels below a column name, each centered over its 12 characters wide cell.

Performance table of the implicit solvers

At the end of a step, and every report-performance-frequency increments if that option is set, the implicit solvers print where the step’s time went (performancetiming): one row per timed part, nested parts indented below their parent, each with its share of the step’s wall clock, its number of calls and its time per call. The table ends with the sum of the top-level rows, the wall clock since the start of the step, and the share left unaccounted, the time spent in code no timer covers. Every part of an increment is timed, so this share stays small; a large one points at work that has no row yet.

Besides the element loop and the linear solve, the rows include the parts that cost a lot on large models although no element is involved:

• topology update: the model modifiers (e.g. AMR), then the ties and the contact catching up;

• equation system rebuild, after a change of the mesh or of a constraint’s DOF footprint: dof manager, linear solver setup, system matrix pattern, mpc transformation build and mpc dirichlet conflict check;

• mpc transform system matrix: condensing the slave DOFs of the multi-point constraints out of the system;

• inertia and damping terms (NID, NIDParallel), per Newton iteration: newmark kinematics, inertia and damping forces and effective tangent;

• reset residual and tangent, per Newton iteration;

• accept increment, finalize increment and conservation report, per accepted increment;

• output, per output increment: field output and one row per output manager, under the name of its *output line.

The timers start from zero together with the step’s wall clock, so the rows and the wall clock of a step’s table cover the same time, whatever ran before the step. A resumed run (*restart, readFrom=...) whose checkpoint replays a topology history prints that replay in a table of its own (top-level AMR rows), before the step’s table.

NEST - Nonlinear Explicit Static

class edelweissfe.solvers.nonlinearexplicitstatic.NEST(jobInfo, journal, **kwargs)

This is the Nonlinear Explicit STatic – solver.

Parameters:

• jobInfo – A dictionary containing the job information.

• journal – The journal instance for logging.

Attributes:

checkpointedState

convergenceTable

equationSystemLayout

incrementStatus

mpcTransformation

Methods

acceptIncrement(step, model, timeStep)	Commit the increment to the model, and let the error estimate size the next one; see acceptIncrement().
applyCorrection(dU, ddU, stepLength)	Move the iterate by the correction, scaled by the step length \(\alpha\) a line search chose (1 without a line search): \(\Delta \boldsymbol{U} \leftarrow \Delta \boldsymbol{U} + \alpha \, \Delta\Delta \boldsymbol{U}\).
applyDirichletToResidual(timeStep, R, dirichlets)	Impose the Dirichlet BCs on the residual.
applyDirichletToStiffness(K, dirichlets[, rhs])	Impose the Dirichlet BCs on the system matrix -- and, given the right-hand side(s) rhs (modified in place), eliminate the constrained DOFs from the columns as well; see edelweissfe.solvers.base.dirichlet.applyDirichletToStiffness().
applyOptionsOverride(fieldValues)	Apply an override of this solver's options, see applyOptionsOverride(); an override that sets an option where it means nothing is refused before anything of it applies (see refuseOptionsThatMeanNothing()).
applyStepActionsAtIncrementStart(model, ...)	Called when all step actions should be applied at the start of a step.
applyStepActionsAtStepEnd(model, stepActions)	Called when all step actions should finish a step.
applyStepActionsAtStepStart(model, step)	Apply every step action's step-start part (initial conditions, material initialization, prescribed fields, ...).
assembleAdditionalTerms(dU, R, F, K, timeStep)	Assemble terms of the equation system beyond the internal and external forces.
assembleConstraintForces(constraints, U_np, ...)	Evaluate every active constraint and add its residual contribution to the given vector; a constraint that is switched off (see active) is skipped.
assembleConstraints(constraints, U_np, dU, ...)	Loop over all constraints, and evaluate the active ones: each adds its residual contribution and the tangent of it, see computeKernels().
assembleInternalForces(U_np, dU, P, K, F, ...)	Assemble the internal forces of the model at a trial state, and their tangent: those of the elements of its equation system, see computeElements().
assembleLoads(loads, U_np, dU, P, K, timeStep)	Assemble the residual contributions of all loads: each adds the negative of the external load it applies, and the tangent of that contribution, which is nonzero only for a load that follows the deformation; see The residual and the tangent.
assembleLoadsAtIncrementEnd(stepActions, ...)	Assemble the residual contributions of the loads acting at the end of the increment, and their tangents; called by assembleResidualAndTangent().
assembleResidualAndTangent(U_n, dU, U_np, R, ...)	Evaluate the model at the trial solution \(U_{n+1} = U_n + \Delta U\): the right-hand side \(R = f_\mathrm{ext} - f_\mathrm{int}\) (see The residual and the tangent) with the terms of assembleAdditionalTerms(), condensed by the multi-point constraints, and the tangent.
assembleStiffnessCSR(K)	The assembled system matrix in CSR format, see assembledMatrix().
assembleSystemMatrix(K, quasiNewton)	The system matrix of this Newton iteration, before the Dirichlet boundary conditions are imposed: the tangent in CSR layout, condensed by the multi-point constraints -- or, once the increment has fallen back to the quasi-Newton method, the tangent of its first iteration.
attemptIncrement(step, model, timeStep)	Integrate the increment with the embedded Runge-Kutta scheme; see attemptIncrement().
beginStep(step, model, ...)	Build the equation system, once per step, and start the step; see beginStep().
buildMPCTransformation(model[, stepActions])	Collect the linear dependency records from all multi-point constraints of the model and assemble the master-slave condensation operator for the current equation system.
checkConvergence(R, ddU, dU, F, ...)	Check the convergence, individually for each field and for the scalar variables, based on the current flux residual and the correction, using the convergence criterion of the solver (see edelweissfe.solvers.base.convergencecriteria).
checkDivergingSolution(...)	Check if the iterative solution scheme is diverging: if the flux residual of a field grew, from one checked iteration to the next, more than maxGrowingIter times in this increment.
checkDivergingSolutionOrFallBack(...)	Cut an increment back whose residual grew more often than maxGrowingIter allows (see checkDivergingSolution()) -- or, with the option fallBackToQuasiNewton, let it continue with the tangent of its first iteration (a quasi-Newton method), see assembleSystemMatrix().
checkMPCDirichletConflicts(transformation, ...)	Raise if any Dirichlet boundary condition of the step prescribes a DOF that is a slave DOF of a multi-point constraint.
collectLoads(stepActions)	The loads among a step's actions: every action that is a LoadBase, whatever step action module it is defined in, so that a new kind of load needs no change to any solver.
computeElements(elements, U_np, dU, P, K, F, ...)	Loop over all elements, and evalute them: each adds its internal forces and their tangent, see computeKernels().
createConvergenceCriterion()	Create the convergence criterion selected by the option convergenceCriterion.
createSystemMatrixAssembly(assemblyType)	The assembly of the system matrix of the equation system just built (option systemMatrixAssembly), on the complete CSR pattern of csrGenerator, with numberOfAssemblyWorkers() workers.
endStep(step, model)	Report the step's performance timing.
extrapolateLastIncrement(extrapolation, ...)	Predict the solution increment, the point Newton's method starts from, by extrapolating the last accepted increments of the step.
extrapolates(extrapolation, prevTimeStep)	Whether an increment that follows the accepted increment prevTimeStep starts from an extrapolated prediction.
extrapolationWeights(extrapolation, ...)	The weights of the last accepted increment and of the one before it in the prediction of extrapolateLastIncrement(): \(\Delta U = w_n \, \Delta U_n + w_{n-1} \, \Delta U_{n-1}\), with
finalizeIncrement(model)	Store what the converged increment leaves behind.
firstIterateNeedsCorrection(dirichlets)	Whether the first iterate of an increment without a prediction -- the converged state of the last increment -- needs a correction before it can be checked for convergence (see solveIncrement()): if Dirichlet boundary conditions act, whose prescribed increment the first correction makes.
forEachChunk(entities, work)	Call work(chunk) for chunks that together hold all entities: the solver's loop over independent entities, other than the elements, that the model's equation system layout evaluates, e.g. the material points of a meshfree model (see updateMaterialPoints()).
forgetEquationSystem()	Drop the equation system of the last build -- the model's layout of it with the DOF manager, the system matrix assembly and its CSR pattern, the multi-point-constraint condensation, the DOF vectors -- so that it is gone before the next one is built, never both at once.
forgetPredictionAfterFailure()	A failed increment drops the prediction, which may have caused the failure, unless retryWithoutExtrapolation has answered that already: then the accepted history is kept, and the smaller increment that follows the cutback extrapolates it.
getRestartData()	The state this solver carries from one increment to the next, as it is; see checkpointedState.
initializeIncrement(U_n, stepActions, model, ...)	Prepare what assembleAdditionalTerms() needs during this increment's Newton iterations.
isOutputIncrement(timeStep)	Whether field outputs, output managers and restart checkpoints are written after this increment.
lineSearchIsDue(iterationCounter)	Whether the correction of this iteration is scaled by a line search (lineSearch): after lineSearchAfterIterations iterations, in every iteration whose number is a multiple of lineSearchEveryIterations.
lineSearchStepLengths()	The trial step lengths of the line search, from the option lineSearchStepLengths.
linearSolve(A, b)	Solve the linear equation system.
locateConstrainedDofs(dirichlets)	Determine, up front, which global DOFs each Dirichlet BC constrains.
numberOfAssemblyWorkers()	How many threads write blocks of the system matrix at the same time: one, the elements are computed in the calling thread.
predictionSurvivesRebuild(modelHasChanged, ...)	Whether this increment may extrapolate the last accepted increments after the equation system was rebuilt.
prepareIncrement(step, model, isRetry)	Nothing to prepare: this solver has neither topology updates nor contact.
printResidualOutlierNodes(residualOutliers)	Print which nodes have the largest residuals.
quadraticLineSearch(largestResidualAt)	The step length \(\alpha\) by which a Newton correction is scaled, for a Newton iteration that overshoots or oscillates (e.g., materials switching between loading and unloading).
readStateFromModel(U)	Read the state the next increment starts from out of the model, into the vectors of the equation system just (re)built: here, the solution, through the model's layout of the equation system (see readStateFromModel()).
refuseOptionsThatMeanNothing(options, ...)	Refuse an option that a declaration of options sets where it means nothing, with the options in force after the declaration.
schema	Option schema for this solver, per OptionSchemaProvider.
selectLinearSolverForStep()	Make linSolver the linear solver the step about to begin names (linsolver): the one of the previous step if it names the same, else a new one.
setRestartData(data)	Restore the state getRestartData() returned, so that the resumed step continues exactly where the uninterrupted run would have.
setTrialState(dU, ddU, stepLength)	Set the trial state of a line search at a step length \(\alpha\) along the correction, and return its solution increment, \(\Delta \boldsymbol{U} + \alpha \, \Delta\Delta \boldsymbol{U}\), at which the model is evaluated next; the iterate stays where it is.
solveForCorrection(K_, R, dirichlets)	Solve the linearized system of a Newton iteration, \(\boldsymbol{K} \, \Delta\Delta \boldsymbol{U} = \boldsymbol{R}\), for the correction \(\Delta\Delta \boldsymbol{U}\).
solveIncrement(U_n, dU, P, K, stepActions, ...)	Explicit Runge-Kutta scheme to solve for an increment.
updateRigidBodies(model, timeStep)	Refresh the kinematics of all rigid bodies in the model after a converged increment.
updateTopologyAndConnectivity(model, step[, ...])	Run the topology update, then let every mesh-dependent consumer catch up on it.
validateModelCapabilities(model)	Refuse, on top of the base checks, a constraint that is switched on and needs its connectivity updated in every increment, i.e. a contact, which searches for what it touches: this solver never performs the update.
writeResultsToModel(U, dU, P[, V, A, ...])	Write the results of a converged increment to the model, through the model's layout of the equation system: the incremental degrees of freedom are set to zero first, from where the next increment starts them (see zeroIncrementalCoefficients()), then the solution, its increment and the internal force are written (see writeResultsToModel()).

schema

Option schema for this solver, per OptionSchemaProvider.

Attributes:

pruneCondensedMatrixZeros

alias of NESTSchema

acceptIncrement(step, model, timeStep)

Commit the increment to the model, and let the error estimate size the next one; see acceptIncrement().

Parameters:

• step – The step being solved.

• model (FEModel) – The model tree.

• timeStep (TimeStep) – The solved increment.

attemptIncrement(step, model, timeStep)

Integrate the increment with the embedded Runge-Kutta scheme; see attemptIncrement().

Parameters:

• step – The step being solved.

• model (FEModel) – The model tree.

• timeStep (TimeStep) – The increment.

Raises:

IncrementFailed – If the scheme requests a smaller increment.

beginStep(step, model, fieldOutputController, outputmanagers)

Build the equation system, once per step, and start the step; see beginStep().

Parameters:

• step – The step to be solved.

• model (FEModel) – The model tree.

• fieldOutputController (FieldOutputController) – The field output controller.

• outputmanagers (dict[str, OutputManagerBase]) – The output managers.

endStep(step, model)

Report the step’s performance timing.

Parameters:

• step – The step that was solved.

• model (FEModel) – The model tree.

prepareIncrement(step, model, isRetry)

Nothing to prepare: this solver has neither topology updates nor contact.

Parameters:

• model (FEModel)

• isRetry (bool)

solveIncrement(U_n, dU, P, K, stepActions, model, timeStep, prevTimeStep)

Explicit Runge-Kutta scheme to solve for an increment.

Parameters:

• Un – The old solution vector.

• dU (DofVector) – The old solution increment.

• P (DofVector) – The old reaction vector.

• K (SystemMatrixAssembly) – The system matrix to be used.

• elements – The dictionary containing all elements.

• stepActions (list) – The list of active step actions.

• model (FEModel) – The model tree.

• timeStep (TimeStep) – The current time step.

• prevTimeStep (TimeStep) – The previous time step.

• U_n (DofVector)

Returns:

A tuple containing

• the new solution vector

• the solution increment

• the new reaction vector

• the scale factor for the next increment size

Return type:

tuple[DofVector,DofVector,DofVector,int,dict]

validateModelCapabilities(model)

Refuse, on top of the base checks, a constraint that is switched on and needs its connectivity updated in every increment, i.e. a contact, which searches for what it touches: this solver never performs the update. A constraint switched off from outside (switchedOn) is not refused: it rests, and takes no part in the step.

Parameters:

model (FEModel) – The model tree.

checkpointedState: dict | None = None

no restart scenario covers this solver, so it does not claim to be.

Type:

Not restartable

linSolver

The linear solver, and the linsolver it was made for (see selectLinearSolverForStep()).

linearSolverDefinitions

The *linearSolver blocks of the input file by name, which linsolver= may name (none without an input file).

supportsMPC = False

Whether this solver supports master-slave condensation / multi-point constraints (e.g. surface ties). Subclasses supporting MPCs must set this to True.

supportsModelModifiers = False

Whether this solver runs the topology update (e.g. h-adaptivity) at all; without it a modifier never runs. A solver that runs it only before its increment loop sets this too, and refuses the modifiers that would need it later (see ModelModifierBase.actsOnlyAtSimulationStart and NED.validateModelCapabilities).

class edelweissfe.solvers.nonlinearexplicitstatic.NESTSchema(rungeKuttaStages=2, rungeKuttaErrorTolerance=0.001, rungeKuttaErrorControl='on', linsolver='pardiso', pruneCondensedMatrixZeros=None)

The options of the *solver datalines and of an >>options block routed to this solver, owned by this module and never mutated from outside it.

Mirrors NEST.SolverSpecificOptions one-for-one; the plain self.options dict remains the actual source of truth consulted at runtime (see NISTSchema for why). The runge-kutta-* option names are not valid Python identifiers, hence the optionName indirection.

Attributes:

pruneCondensedMatrixZeros

Parameters:

• rungeKuttaStages (int | None)

• rungeKuttaErrorTolerance (float | None)

• rungeKuttaErrorControl (str | None)

• linsolver (str | None)

• pruneCondensedMatrixZeros (bool | None)

NEST

Option	Type / Default	Description
runge-kutta-stages	int, default=2	The number of Runge-Kutta stages.
runge-kutta-error-tolerance	float, default=0.001	The error tolerance for the Runge-Kutta error control.
runge-kutta-error-control	str, default=’on’	Activate the Runge-Kutta error control (on|off).
linsolver	str, default=’pardiso’	The linear solver: the name of a *linearSolver block, or a type of linear solver with its defaults (e.g. pardiso, superlu, blockamg).
pruneCondensedMatrixZeros	bool, default=None	Compact explicitly stored zeros out of the multi-point-constraint-condensed system matrix before solving. See NISTSchema for the trade-off; unset (the default) decides automatically from the linear solver in force.

NESTParallel - Nonlinear Explicit Static (parallel)

Parallel implementation of the NEST solver.

class edelweissfe.solvers.nonlinearexplicitstaticparallel.NESTParallel(jobInfo, journal, **kwargs)

Attributes:

checkpointedState

convergenceTable

equationSystemLayout

incrementStatus

mpcTransformation

Methods

acceptIncrement(step, model, timeStep)	Commit the increment to the model, and let the error estimate size the next one; see acceptIncrement().
applyCorrection(dU, ddU, stepLength)	Move the iterate by the correction, scaled by the step length \(\alpha\) a line search chose (1 without a line search): \(\Delta \boldsymbol{U} \leftarrow \Delta \boldsymbol{U} + \alpha \, \Delta\Delta \boldsymbol{U}\).
applyDirichletToResidual(timeStep, R, dirichlets)	Impose the Dirichlet BCs on the residual.
applyDirichletToStiffness(K, dirichlets[, rhs])	Impose the Dirichlet BCs on the system matrix -- and, given the right-hand side(s) rhs (modified in place), eliminate the constrained DOFs from the columns as well; see edelweissfe.solvers.base.dirichlet.applyDirichletToStiffness().
applyOptionsOverride(fieldValues)	Apply an override of this solver's options, see applyOptionsOverride(); an override that sets an option where it means nothing is refused before anything of it applies (see refuseOptionsThatMeanNothing()).
applyStepActionsAtIncrementStart(model, ...)	Called when all step actions should be applied at the start of a step.
applyStepActionsAtStepEnd(model, stepActions)	Called when all step actions should finish a step.
applyStepActionsAtStepStart(model, step)	Apply every step action's step-start part (initial conditions, material initialization, prescribed fields, ...).
assembleAdditionalTerms(dU, R, F, K, timeStep)	Assemble terms of the equation system beyond the internal and external forces.
assembleConstraintForces(constraints, U_np, ...)	Evaluate every active constraint and add its residual contribution to the given vector; a constraint that is switched off (see active) is skipped.
assembleConstraints(constraints, U_np, dU, ...)	Loop over all constraints, and evaluate the active ones: each adds its residual contribution and the tangent of it, see computeKernels().
assembleInternalForces(U_np, dU, P, K, F, ...)	Assemble the internal forces of the model at a trial state, and their tangent: those of the elements of its equation system, see computeElements().
assembleLoads(loads, U_np, dU, P, K, timeStep)	Assemble the residual contributions of all loads: each adds the negative of the external load it applies, and the tangent of that contribution, which is nonzero only for a load that follows the deformation; see The residual and the tangent.
assembleLoadsAtIncrementEnd(stepActions, ...)	Assemble the residual contributions of the loads acting at the end of the increment, and their tangents; called by assembleResidualAndTangent().
assembleResidualAndTangent(U_n, dU, U_np, R, ...)	Evaluate the model at the trial solution \(U_{n+1} = U_n + \Delta U\): the right-hand side \(R = f_\mathrm{ext} - f_\mathrm{int}\) (see The residual and the tangent) with the terms of assembleAdditionalTerms(), condensed by the multi-point constraints, and the tangent.
assembleStiffnessCSR(K)	The assembled system matrix in CSR format, see assembledMatrix().
assembleSystemMatrix(K, quasiNewton)	The system matrix of this Newton iteration, before the Dirichlet boundary conditions are imposed: the tangent in CSR layout, condensed by the multi-point constraints -- or, once the increment has fallen back to the quasi-Newton method, the tangent of its first iteration.
attemptIncrement(step, model, timeStep)	Integrate the increment with the embedded Runge-Kutta scheme; see attemptIncrement().
beginStep(step, model, ...)	Build the equation system, once per step, and start the step; see beginStep().
buildMPCTransformation(model[, stepActions])	Collect the linear dependency records from all multi-point constraints of the model and assemble the master-slave condensation operator for the current equation system.
checkConvergence(R, ddU, dU, F, ...)	Check the convergence, individually for each field and for the scalar variables, based on the current flux residual and the correction, using the convergence criterion of the solver (see edelweissfe.solvers.base.convergencecriteria).
checkDivergingSolution(...)	Check if the iterative solution scheme is diverging: if the flux residual of a field grew, from one checked iteration to the next, more than maxGrowingIter times in this increment.
checkDivergingSolutionOrFallBack(...)	Cut an increment back whose residual grew more often than maxGrowingIter allows (see checkDivergingSolution()) -- or, with the option fallBackToQuasiNewton, let it continue with the tangent of its first iteration (a quasi-Newton method), see assembleSystemMatrix().
checkMPCDirichletConflicts(transformation, ...)	Raise if any Dirichlet boundary condition of the step prescribes a DOF that is a slave DOF of a multi-point constraint.
collectLoads(stepActions)	The loads among a step's actions: every action that is a LoadBase, whatever step action module it is defined in, so that a new kind of load needs no change to any solver.
createConvergenceCriterion()	Create the convergence criterion selected by the option convergenceCriterion.
createSystemMatrixAssembly(assemblyType)	The assembly of the system matrix of the equation system just built (option systemMatrixAssembly), on the complete CSR pattern of csrGenerator, with numberOfAssemblyWorkers() workers.
endStep(step, model)	Report the step's performance timing.
extrapolateLastIncrement(extrapolation, ...)	Predict the solution increment, the point Newton's method starts from, by extrapolating the last accepted increments of the step.
extrapolates(extrapolation, prevTimeStep)	Whether an increment that follows the accepted increment prevTimeStep starts from an extrapolated prediction.
extrapolationWeights(extrapolation, ...)	The weights of the last accepted increment and of the one before it in the prediction of extrapolateLastIncrement(): \(\Delta U = w_n \, \Delta U_n + w_{n-1} \, \Delta U_{n-1}\), with
finalizeIncrement(model)	Store what the converged increment leaves behind.
firstIterateNeedsCorrection(dirichlets)	Whether the first iterate of an increment without a prediction -- the converged state of the last increment -- needs a correction before it can be checked for convergence (see solveIncrement()): if Dirichlet boundary conditions act, whose prescribed increment the first correction makes.
forEachChunk(entities, work)	Call work(chunk) for chunks of entities on the available threads, as the elements are computed (see mapChunksInParallel()); otherwise as forEachChunk().
forgetEquationSystem()	Drop the equation system of the last build -- the model's layout of it with the DOF manager, the system matrix assembly and its CSR pattern, the multi-point-constraint condensation, the DOF vectors -- so that it is gone before the next one is built, never both at once.
forgetPredictionAfterFailure()	A failed increment drops the prediction, which may have caused the failure, unless retryWithoutExtrapolation has answered that already: then the accepted history is kept, and the smaller increment that follows the cutback extrapolates it.
getRestartData()	The state this solver carries from one increment to the next, as it is; see checkpointedState.
initializeIncrement(U_n, stepActions, model, ...)	Prepare what assembleAdditionalTerms() needs during this increment's Newton iterations.
isOutputIncrement(timeStep)	Whether field outputs, output managers and restart checkpoints are written after this increment.
lineSearchIsDue(iterationCounter)	Whether the correction of this iteration is scaled by a line search (lineSearch): after lineSearchAfterIterations iterations, in every iteration whose number is a multiple of lineSearchEveryIterations.
lineSearchStepLengths()	The trial step lengths of the line search, from the option lineSearchStepLengths.
linearSolve(A, b)	Solve the linear equation system.
locateConstrainedDofs(dirichlets)	Determine, up front, which global DOFs each Dirichlet BC constrains.
numberOfAssemblyWorkers()	How many threads write blocks of the system matrix at the same time: those that compute the elements.
predictionSurvivesRebuild(modelHasChanged, ...)	Whether this increment may extrapolate the last accepted increments after the equation system was rebuilt.
prepareIncrement(step, model, isRetry)	Nothing to prepare: this solver has neither topology updates nor contact.
printResidualOutlierNodes(residualOutliers)	Print which nodes have the largest residuals.
quadraticLineSearch(largestResidualAt)	The step length \(\alpha\) by which a Newton correction is scaled, for a Newton iteration that overshoots or oscillates (e.g., materials switching between loading and unloading).
readStateFromModel(U)	Read the state the next increment starts from out of the model, into the vectors of the equation system just (re)built: here, the solution, through the model's layout of the equation system (see readStateFromModel()).
refuseOptionsThatMeanNothing(options, ...)	Refuse an option that a declaration of options sets where it means nothing, with the options in force after the declaration.
schema
selectLinearSolverForStep()	Make linSolver the linear solver the step about to begin names (linsolver): the one of the previous step if it names the same, else a new one.
setRestartData(data)	Restore the state getRestartData() returned, so that the resumed step continues exactly where the uninterrupted run would have.
setTrialState(dU, ddU, stepLength)	Set the trial state of a line search at a step length \(\alpha\) along the correction, and return its solution increment, \(\Delta \boldsymbol{U} + \alpha \, \Delta\Delta \boldsymbol{U}\), at which the model is evaluated next; the iterate stays where it is.
solveForCorrection(K_, R, dirichlets)	Solve the linearized system of a Newton iteration, \(\boldsymbol{K} \, \Delta\Delta \boldsymbol{U} = \boldsymbol{R}\), for the correction \(\Delta\Delta \boldsymbol{U}\).
solveIncrement(U_n, dU, P, K, stepActions, ...)	Explicit Runge-Kutta scheme to solve for an increment.
updateRigidBodies(model, timeStep)	Refresh the kinematics of all rigid bodies in the model after a converged increment.
updateTopologyAndConnectivity(model, step[, ...])	Run the topology update, then let every mesh-dependent consumer catch up on it.
validateModelCapabilities(model)	Refuse, on top of the base checks, a constraint that is switched on and needs its connectivity updated in every increment, i.e. a contact, which searches for what it touches: this solver never performs the update.
writeResultsToModel(U, dU, P[, V, A, ...])	Write the results of a converged increment to the model, through the model's layout of the equation system: the incremental degrees of freedom are set to zero first, from where the next increment starts them (see zeroIncrementalCoefficients()), then the solution, its increment and the internal force are written (see writeResultsToModel()).

computeElements

beginStep(step, model, fieldOutputController, outputmanagers)

Build the equation system, once per step, and start the step; see beginStep().

Parameters:

• step – The step to be solved.

• model – The model tree.

• fieldOutputController – The field output controller.

• outputmanagers – The output managers.

computeElements(elements, Un1, dU, P, K, F, timeStep)

Loop over all elements, and evalute them: each adds its internal forces and their tangent, see computeKernels(). Is called by assembleResidualAndTangent() in each iteration.

The elements are those of the model’s equation system (elements): the finite elements, and whatever else contributes like an element – in a meshfree model, its cells and particles.

Parameters:

• elements – The list of elements.

• U_np – The current solution vector.

• dU – The current solution increment vector.

• P – The reaction vector.

• K – The system matrix.

• F – The vector of accumulated fluxes for convergence checks.

• timeStep – The time step.

Returns:

• The modified reaction vector.

• The modified system matrix.

• The modified accumulated flux vector.

Return type:

tuple[DofVector,SystemMatrixAssembly,DofVector]

forEachChunk(entities, work)

Call work(chunk) for chunks of entities on the available threads, as the elements are computed (see mapChunksInParallel()); otherwise as forEachChunk().

Parameters:

entities (list)

numberOfAssemblyWorkers()

How many threads write blocks of the system matrix at the same time: those that compute the elements.

Returns:

The number of workers.

Return type:

int

NESTParallel

Option	Type / Default	Description
runge-kutta-stages	int, default=2	The number of Runge-Kutta stages.
runge-kutta-error-tolerance	float, default=0.001	The error tolerance for the Runge-Kutta error control.
runge-kutta-error-control	str, default=’on’	Activate the Runge-Kutta error control (on|off).
linsolver	str, default=’pardiso’	The linear solver: the name of a *linearSolver block, or a type of linear solver with its defaults (e.g. pardiso, superlu, blockamg).
pruneCondensedMatrixZeros	bool, default=None	Compact explicitly stored zeros out of the multi-point-constraint-condensed system matrix before solving. See NISTSchema for the trade-off; unset (the default) decides automatically from the linear solver in force.

NED - Nonlinear Explicit Dynamic

The nonlinear explicit dynamic solver.

An increment solves no equation system. Second-order fields are advanced by central differences (a leapfrog, with the velocity staggered half an increment behind the displacement), first-order fields – a nonlocal damage field, for instance – by forward Euler, and both divide by a lumped mass, so the cost of an increment is one element pass plus vector arithmetic. There is no linear solver, no Newton loop, and no cutback: the time increment is set by stability, not by convergence, so a material that fails to integrate at the stable step is reported as an error rather than answered with a smaller step.

Inertia and damping. An element reports two diagonals: the coefficient of each field’s second time derivative (its inertia) and of its first (its damping). The integrator divides by the one and relaxes with the other, and which scheme a field gets is read off them – an inertia means central differences, a damping alone means forward Euler. The deck declares no field ordering, because an inertia is exactly what a central-difference update divides by. At zero damping the update is bit-identical to the undamped one.

The diagnostics need one thing more, and it is about units: which entries of the inertia diagonal may be summed. That belongs to the FIELD rather than to a deck, and is recorded in edelweissfe.config.phenomena as a mass, a rotational inertia, or neither. Linear momentum sums the first kind, the energy balance the first two; conservation across a topology change is checked per field regardless.

The case that forces it: a first-order field’s forward-Euler limit falls off with the SQUARE of the element size, so on a refined mesh it, not the mechanical field, bounds the increment. Giving it an inertia makes it second order – a damped wave equation, limit linear in h – with its old coefficient keeping its value and becoming the damping. Marmot supplies one through the nonlocal micro inertia element property; that inertia is a time squared, so 0.5 m_k rate^2 is a volume rather than an energy and the field is registered as carrying a non-mechanical inertia, which keeps it out of both balances.

first-order-fields, second-order-fields, first-order-scheme and second-order-scheme are removed rather than deprecated: a deck carrying any of them now fails with Invalid option, and deleting those lines is the whole of the migration.

The stable increment is computed once per mesh from the element wave speeds and scaled by courant-number. It is recomputed whenever the mesh changes and deliberately not recomputed when it has not: as a material softens the true limit only grows, so reusing it is the conservative direction, and re-deriving it costs a full element pass.

Three cadences, all in increments. output-frequency is how often progress is logged and field outputs are finalized (a divisor, so it must be at least 1). topology-check-frequency is how often model modifiers may change the mesh mid-run, and must be a multiple of output-frequency because a marker reads the last finalized output; 0 runs the topology update only once, before the increment loop. contact-update-frequency is how often constraints whose connectivity is the outcome of a search re-run it; 0 disables that. The reasoning behind the latter two, and what they cost in accuracy, is documented under Adaptive refinement: error estimation and marking and Contact mechanics: theory respectively.

Constraints. Multi-point constraints are enforced by folding a slave’s mass and force onto its masters (row-sum, mass-conserving) and slaving it kinematically, which leaves the critical time step untouched. Other constraints – contact above all – are evaluated through computeKernelsExplicit(), which asks for forces and no tangent, since a tangent assembled here would enter nothing. A constraint that can act only through its tangent contributes nothing to an explicit increment and is refused at validation rather than silently ignored.

Energy balance. Every reporting increment prints the energy the model held at the start of the step, kinetic and internal (E_0), the external work accumulated at prescribed degrees of freedom since, the kinetic energy, the internal (strain) energy, and what remains unaccounted. Its purpose is one check that nothing else performs: the kinetic energy cannot exceed the energy at the start of the step plus the external work, because the terms omitted from the balance are all non-negative, so a violation means energy is appearing from nowhere and the time step is above the true stability limit – which the critical time step does not see in full, as it never sees contact penalty stiffness, and sees the nonlocal field only where that field carries a non-mechanical inertia. Such a field’s own 0.5 m_k rate^2 is reported on a separate line and kept out of the balance: with a coefficient in seconds squared it has the units of a volume, not of an energy. E_0 is read off the model’s state, nothing evaluated: zero for a step from rest in a model without strain energy, and, for a step that continues another one, the energy that one left, its kinetic energy included. The continued step’s balance is therefore the one the two steps would report as one, less what the first of them dissipated.

Two honest limits, both of which the solver states in the log rather than leaving to be discovered. The external work accumulates only at prescribed degrees of freedom, so a step driven by body forces, node forces, a distributed load or an initial velocity has none, and the check cannot fire, whatever energy the step starts with; and the internal energy is whatever the materials publish as a strain energy, which not every material does. When either is identically zero – or the external work is negative, a structure returning work at its boundary – the balance is not usable as a quasi-static criterion, and the honest substitute is to integrate a reaction force against its prescribed displacement, both of which are available as saveHistory field outputs.

Diagnostics across a topology change. Because the lumped mass is the operator an increment divides by, a refinement’s effect on it is reported: total mass, per-component linear momentum and kinetic energy before and after, plus the relative mass drift accumulated over every such change, so that many individually-tolerable changes cannot silently add up to a meaningful one. The smallest coefficient the increment divides by is reported alongside – an inertia at a second-order degree of freedom, a damping at a first-order one – since that is what bounds the step.

The force it publishes. The node field result P is, as in every solver, the internal force with the inertial and damping forces, \(\boldsymbol{f}_\mathrm{int} + \boldsymbol{M} \boldsymbol{a} + \boldsymbol{C} \boldsymbol{v}\), without the loads and the constraints (see Node field results). The central difference makes the inertial and damping forces at the end of an increment exactly the net force \(\boldsymbol{f}_\mathrm{net}\) the next velocity update divides by – that update solves \(\boldsymbol{M} (\boldsymbol{v}_{n+3/2} - \boldsymbol{v}_{n+1/2}) / \Delta \bar{t} + \boldsymbol{C} (\boldsymbol{v}_{n+3/2} + \boldsymbol{v}_{n+1/2}) / 2 = \boldsymbol{f}_\mathrm{net}\) – so P is \(\boldsymbol{f}_\mathrm{int} + \boldsymbol{f}_\mathrm{net}\) at every degree of freedom with an equation of motion. Two kinds have none. At a degree of freedom with a prescribed motion, the acceleration at the end of the increment is set by the prescribed velocity of the next increment, which is not known when the increment is accepted, so P there is \(\boldsymbol{f}_\mathrm{int}\) alone: the reaction, without the inertia of the support’s own lumped mass. A multi-point-constraint slave carries no inertia of its own – it is folded onto its masters, whose P holds it – so P there is its \(\boldsymbol{f}_\mathrm{int}\) as well.

The velocity and the acceleration it publishes. Within a step, the velocity V the model receives at every accepted increment is the leapfrog’s own, \(\boldsymbol{v}_{n+1/2}\), half an increment behind the solution; on a field integrated by forward Euler, it is the rate of the increment just taken. The acceleration A is the one at the end of the increment: the one with which the velocity at that time, \(\boldsymbol{v}_{n+1} = \boldsymbol{v}_{n+1/2} + \tfrac{1}{2} \Delta t \, \boldsymbol{a}_{n+1}\), satisfies the equation of motion, \(\boldsymbol{M} \boldsymbol{a}_{n+1} + \boldsymbol{C} \boldsymbol{v}_{n+1} = \boldsymbol{f}_\mathrm{net}\), so

\[\left( \boldsymbol{M} + \tfrac{1}{2} \Delta t \, \boldsymbol{C} \right) \boldsymbol{a}_{n+1} = \boldsymbol{f}_\mathrm{net} - \boldsymbol{C} \boldsymbol{v}_{n+1/2},\]

which is \(\boldsymbol{M}^{-1} \boldsymbol{f}_\mathrm{net}\) without damping. On a field integrated by forward Euler, which has no inertia, it is the change of the rate from the increment just taken to the rate at the end of it. It is zero where the motion is prescribed, as the next increment’s prescribed velocity sets it, and a multi-point-constraint slave takes it from its masters, as it takes its velocity. At the end of the step, the velocity is advanced by the half increment it lags behind (advanceVelocity()), so that the step leaves the velocity at its end time, as every solver does (see Node field results).

The net force is solver state, and a step starts from the model’s state. The net force lives within one equation system, and is never read back from the model. A rebuild after a change of a constraint’s connectivity carries it over by numbering, and a resumed step takes it from its checkpoint. After a mesh change it starts at zero: re-evaluating it would run the constitutive law off-cycle and write into the material state, and zero is an O(dT) error for one increment. At the start of a step, it is the net force of the state the model holds, \(\boldsymbol{f}_\mathrm{net} = \boldsymbol{M} \boldsymbol{a}_0 + \boldsymbol{C} \boldsymbol{v}_0\), formed from the model’s velocity and acceleration with the lumped operators, without evaluating anything (and zero where the motion is prescribed), and the first increment is the leapfrog’s half step from \(\boldsymbol{v}_0\). After a quasi-static step, which leaves the model at rest, it is zero: the start state is taken as in equilibrium, exact after a quasi-static step that converged there (a gravity preload, for instance), and a lag of one increment for a load applied suddenly at the start. After an explicit step, it is the net force that step ended with, and the leapfrog continues as though the two steps were one (with damping, to second order in the time increment). After an implicit dynamic step, the step starts from the velocity and the acceleration that step ended with.

Restart. A resumed run must not repeat the half-step that starts a leapfrog: the velocity in the checkpoint already carries the half-increment offset, so applying the startup again would apply one half-impulse too few. The solver therefore takes the previous time increment and the net force from the restart state rather than synthesising them, and the velocity field is checkpointed in its own right, with the model – an implicit solver reconstructs everything it needs from the displacement, a central-difference scheme does not.

class edelweissfe.solvers.nonlinearexplicitdynamic.ExplicitSystem(Minv, U, dU, V, A, netForce, internalForce, criticalTimeStep)

Everything an explicit increment operates on that is sized by the current equation system.

None of these are independent state: each is indexed by the DofManager in force when it was created, so the moment that system changes – an h-adaptivity event creating nodes and elements, a contact constraint re-assigning its slave nodes to different master facets – every one of them has to be rebuilt together, and any that is not becomes either a length mismatch or, worse, a silently mis-indexed vector. Bundling them makes “rebuild the system” one assignment at the call site instead of nine, which is what keeps the build before the increment loop and a rebuild inside it the same code path rather than two that drift apart.

Parameters:

• Minv (edelweissfe.numerics.dofvector.DofVector) – The inverse lumped mass. Zero on multi-point-constraint slave DOFs, which integrate no equation of motion of their own. The lumped mass it was built from is not carried here: the diagnostics need the unfolded mass, which is held as _rawLumpedMass, and a second mass on this dataclass that no caller reads would only invite reading the wrong one.

• U (edelweissfe.numerics.dofvector.DofVector) – The solution vector.

• dU (edelweissfe.numerics.dofvector.DofVector) – The solution increment of the increment being computed.

• V (edelweissfe.numerics.dofvector.DofVector) – The velocity vector, staggered half an increment behind U.

• A (edelweissfe.numerics.dofvector.DofVector) – The acceleration at the end of the last increment (see NED.accelerationAtTheEndOfTheIncrement()), which the model receives with the velocity. At the start of a step, the one the model holds, from which the step’s first net force is formed (see NED.formNetForceOfTheModelsState()).

• netForce (edelweissfe.numerics.dofvector.DofVector) – The net nodal force – external minus internal, with the loads and the constraints – which drives the velocity update of the next increment.

• internalForce (edelweissfe.numerics.dofvector.DofVector) – The elements’ sum of the increment, \(\boldsymbol{f}_\mathrm{int}\), from which the published force is formed at acceptance (see NED.internalInertialAndDampingForce()).

• criticalTimeStep (float) – The stable time increment for the current mesh, already scaled by the courant number.

class edelweissfe.solvers.nonlinearexplicitdynamic.NED(jobInfo, journal, **kwargs)

This is the Nonlinear Explicit Dynamic – solver.

Parameters:

• jobInfo – A dictionary containing the job information.

• journal – The journal instance for logging.

Attributes:

convergenceTable

incrementStatus

mpcTransformation

Methods

accelerationAtTheEndOfTheIncrement(...)	The acceleration this solver publishes as the node field entry A of the increment just completed, the one at its end (see the module documentation, "The velocity and the acceleration it publishes"): with the velocity at that time, half an increment after the staggered one, it satisfies the equation of motion, \(\left( \boldsymbol{M} + \tfrac{1}{2} \Delta t \, \boldsymbol{C} \right) \boldsymbol{a} = \boldsymbol{f}_\mathrm{net} - \boldsymbol{C} \boldsymbol{v}_{n+1/2}\).
acceptIncrement(step, model, timeStep)	Commit the increment to the model; see acceptIncrement().
applyDirichletToResidual(timeStep, R, dirichlets)	Impose the Dirichlet BCs on the residual.
applyOptionsOverride(fieldValues)	Apply a partial override of this solver's own schema fields onto self.options.
applyStepActionsAtIncrementStart(model, ...)	Called when all step actions should be applied at the start of a step.
applyStepActionsAtStepEnd(model, stepActions)	Called when all step actions should finish a step.
applyStepActionsAtStepStart(model, step)	Apply every step action's step-start part (initial conditions, material initialization, prescribed fields, ...).
assembleConstraintForces(constraints, U_np, ...)	Evaluate every active constraint and add its residual contribution to the given vector; a constraint that is switched off (see active) is skipped.
assembleLoads(loads, U_np, dU, P, K, timeStep)	Assemble the residual contributions of all loads: each adds the negative of the external load it applies, and the tangent of that contribution, which is nonzero only for a load that follows the deformation; see The residual and the tangent.
assembleStiffnessCSR(K)	The assembled system matrix in CSR format, see assembledMatrix().
attemptIncrement(step, model, timeStep)	The periodic contact search, when due, and the central-difference update; see attemptIncrement().
beginStep(step, model, ...)	Start the step; see beginStep().
buildEquationSystem(model, step[, previous])	Build the equation system and everything sized by it.
buildMPCTransformation(model[, stepActions])	Collect the linear dependency records from all multi-point constraints of the model and assemble the master-slave condensation operator for the current equation system.
checkConvergence(R, ddU, dU, F, ...)	Check the convergence, individually for each field and for the scalar variables, based on the current flux residual and the correction, using the convergence criterion of the solver (see edelweissfe.solvers.base.convergencecriteria).
checkMPCDirichletConflicts(transformation, ...)	Raise if any Dirichlet boundary condition of the step prescribes a DOF that is a slave DOF of a multi-point constraint.
collectLoads(stepActions)	The loads among a step's actions: every action that is a LoadBase, whatever step action module it is defined in, so that a new kind of load needs no change to any solver.
computeElements(elements, U_np, dU, P, timeStep)	Loop over all elements, and evaluate them: their residual contributions, summed, are the internal force.
endStep(step, model)	Leave the model's state at the end time of the step, and report the step's performance timing.
energyOfTheModelsState()	The energy the model holds at the start of a cold step, which its energy balance starts from (see the module documentation, "Energy balance"): the kinetic energy of the model's velocity, read with the solution when the equation system was built, and the internal (strain) energy the materials publish for the state the elements hold.
formNetForceOfTheModelsState(dirichlets)	Form the net force a cold step starts from: that of the state the model holds, \(\boldsymbol{f}_\mathrm{net} = \boldsymbol{M} \boldsymbol{a}_0 + \boldsymbol{C} \boldsymbol{v}_0\), from its velocity and acceleration, read with the solution when the equation system was built, and the lumped operators; nothing is evaluated (see the module documentation, "The net force is solver state, and a step starts from the model's state").
getCriticalTimeStepForExplicitDynamics(model, U)	Compute the critical time step for explicit dynamics.
getRestartData()	The state this solver carries from one increment to the next, as it is; see checkpointedState.
internalInertialAndDampingForce(dirichlets)	The force this solver publishes as the node field result P of the increment just completed: the internal force with the inertial and damping forces, \(\boldsymbol{f}_\mathrm{int} + \boldsymbol{M} \boldsymbol{a} + \boldsymbol{C} \boldsymbol{v}\), without the loads and the constraints (see the module documentation, "The force it publishes").
isOutputIncrement(timeStep)	Only every output-frequency-th increment: an explicit run has millions.
kineticEnergy(V)	The kinetic energy \(\tfrac{1}{2} \sum m v^2\) of a velocity, with the unfolded lumped mass, over the degrees of freedom where it is an energy (ids_mechanicalEnergy): not on a first-order field, whose coefficient is a viscosity, nor on a field whose inertia is not mechanical.
linearSolve(A, b)	Solve the linear equation system.
locateConstrainedDofs(dirichlets)	Determine, up front, which global DOFs each Dirichlet BC constrains.
prepareIncrement(step, model, isRetry)	The topology update, when due, and the equation system; see prepareIncrement().
printResidualOutlierNodes(residualOutliers)	Print which nodes have the largest residuals.
reportTopologyChangeConservation(...)	Report what a topology change did to the quantities that ought to survive it.
schema	Option schema for this solver, per OptionSchemaProvider.
secondOrderMomentum(mass, V, model)	The linear momentum of the fields whose inertia is a mass, per spatial component; see linearMomentum().
setRestartData(data)	Restore the state getRestartData() returned, so that the resumed step continues exactly where the uninterrupted run would have.
solveIncrement(U_n, dU, V, netForce, ...)	Standard explicit update scheme to solve for an increment.
updateConstraintConnectivity(model)	Let the constraints whose connectivity is the outcome of a search re-run that search.
updateRigidBodies(model, timeStep)	Refresh the kinematics of all rigid bodies in the model after a converged increment.
updateTopologyAndConnectivity(model, step[, ...])	Run the topology update, then let every mesh-dependent consumer catch up on it.
validateModelCapabilities(model)	Refuse the model features this solver cannot integrate, on top of the base checks.
writeResultsToModel(U, dU, P, V, A)	Write the results of an increment to the model, through the model's layout of the equation system: the solution, its increment, the internal, inertial and damping force, the velocity and the acceleration to the node fields, and the solution to the scalar variables (see writeResultsToModel()); the incremental degrees of freedom are set to zero first, from where the next increment starts them (see zeroIncrementalCoefficients()).

schema

Option schema for this solver, per OptionSchemaProvider.

Attributes:

expectFirstOrderFields

expectSecondOrderFields

alias of NEDSchema

accelerationAtTheEndOfTheIncrement(dirichlets, timeIncrement)

The acceleration this solver publishes as the node field entry A of the increment just completed, the one at its end (see the module documentation, “The velocity and the acceleration it publishes”): with the velocity at that time, half an increment after the staggered one, it satisfies the equation of motion, \(\left( \boldsymbol{M} + \tfrac{1}{2} \Delta t \, \boldsymbol{C} \right) \boldsymbol{a} = \boldsymbol{f}_\mathrm{net} - \boldsymbol{C} \boldsymbol{v}_{n+1/2}\).

On a degree of freedom with a prescribed motion, whose acceleration the next increment’s prescribed velocity sets, it is zero. A multi-point-constraint slave, which integrates no equation of motion of its own, takes it from its masters.

Parameters:

• dirichlets – The Dirichlet boundary conditions of the step, located on the current equation system.

• timeIncrement (float) – The time increment just completed.

Returns:

The acceleration, the system’s vector.

Return type:

DofVector

acceptIncrement(step, model, timeStep)

Commit the increment to the model; see acceptIncrement().

Parameters:

• step – The step being solved.

• model (FEModel) – The model tree.

• timeStep (TimeStep) – The increment.

attemptIncrement(step, model, timeStep)

The periodic contact search, when due, and the central-difference update; see attemptIncrement().

Parameters:

• step – The step being solved.

• model (FEModel) – The model tree.

• timeStep (TimeStep) – The increment.

Raises:

StepFailed – If a material requests a cutback: an explicit increment is not retried smaller.

beginStep(step, model, fieldOutputController, outputmanagers)

Start the step; see beginStep(). The equation system is built by the first prepareIncrement().

Parameters:

• step – The step to be solved.

• model (FEModel) – The model tree.

• fieldOutputController (FieldOutputController) – The field output controller.

• outputmanagers (dict[str, OutputManagerBase]) – The output managers.

buildEquationSystem(model, step, previous=None)

Build the equation system and everything sized by it.

Called once before the increment loop, and again from inside it whenever a constraint reports that its DOF footprint changed – one method for both, so the path every model takes and the path only a contact model takes cannot drift apart. On that second path the DofManager is refreshed rather than rebuilt and the lumped operators are kept; see _operatorsReusable().

Parameters:

• model (FEModel) – The model tree.

• step – The step being solved; its actions are needed to check the multi-point constraints against the prescribed Dirichlet conditions.

• previous (ExplicitSystem) – The system being replaced, when this is a rebuild rather than the initial build. Its solution, velocity and net force are carried over verbatim rather than re-read from the model: they are the vectors the increment continues from. A rebuild triggered by a constraint’s connectivity leaves the mesh – hence the DOF layout – untouched, and that is checked rather than assumed: copying between two different layouts would mis-index every vector silently.

Returns:

The freshly built system.

Return type:

ExplicitSystem

computeElements(elements, U_np, dU, P, timeStep)

Loop over all elements, and evaluate them: their residual contributions, summed, are the internal force. Called by solveIncrement() once per increment.

Parameters:

• elements (list) – The finite elements, by number.

• U_np (DofVector) – The current solution vector.

• dU (DofVector) – The solution increment vector.

• P (DofVector) – The vector the elements’ residual contributions are added to; zero.

• timeStep (TimeStep) – The time step.

Returns:

• The internal force.

• The internal (strain) energy the elements publish.

Return type:

tuple[DofVector,float]

endStep(step, model)

Leave the model’s state at the end time of the step, and report the step’s performance timing.

The velocity the model holds lags half of the last increment behind the solution; it is advanced by that half increment with the acceleration (advanceVelocity()), so that the next step starts from the velocity at its start time (see the module documentation, “The velocity and the acceleration it publishes”). A step resumed at its very end does so as well: its checkpoint holds the staggered velocity and the last increment. This solver retries no increment, so its step never falls back to a checkpoint, which ends the abandoned attempt here too.

Parameters:

• step – The step that was solved.

• model (FEModel) – The model tree.

energyOfTheModelsState()

The energy the model holds at the start of a cold step, which its energy balance starts from (see the module documentation, “Energy balance”): the kinetic energy of the model’s velocity, read with the solution when the equation system was built, and the internal (strain) energy the materials publish for the state the elements hold. Nothing is evaluated.

Returns:

The kinetic and the internal energy.

Return type:

float

formNetForceOfTheModelsState(dirichlets)

Form the net force a cold step starts from: that of the state the model holds, \(\boldsymbol{f}_\mathrm{net} = \boldsymbol{M} \boldsymbol{a}_0 + \boldsymbol{C} \boldsymbol{v}_0\), from its velocity and acceleration, read with the solution when the equation system was built, and the lumped operators; nothing is evaluated (see the module documentation, “The net force is solver state, and a step starts from the model’s state”). Zero on a degree of freedom with a prescribed motion, which has no equation of motion, and on a multi-point-constraint slave, whose inertia and damping are folded onto its masters.

Parameters:

dirichlets – The Dirichlet boundary conditions of the step.

getCriticalTimeStepForExplicitDynamics(model, U)

Compute the critical time step for explicit dynamics.

Parameters:

• model (FEModel) – The model tree.

• U (DofVector)

Returns:

The critical time step for explicit dynamics.

Return type:

float

internalInertialAndDampingForce(dirichlets)

The force this solver publishes as the node field result P of the increment just completed: the internal force with the inertial and damping forces, \(\boldsymbol{f}_\mathrm{int} + \boldsymbol{M} \boldsymbol{a} + \boldsymbol{C} \boldsymbol{v}\), without the loads and the constraints (see the module documentation, “The force it publishes”).

On a degree of freedom with an equation of motion, the inertial and damping forces are the net force, since the next velocity update solves exactly \(\boldsymbol{M} \boldsymbol{a} + \boldsymbol{C} \boldsymbol{v} = \boldsymbol{f}_\mathrm{net}\) – so \(\boldsymbol{f}_\mathrm{int} + \boldsymbol{f}_\mathrm{net}\), which is the force of the loads and the constraints there. On a degree of freedom with a prescribed motion, whose acceleration the next increment’s prescribed velocity sets, it is the internal force alone. A multi-point-constraint slave has a zero net force (it is folded onto its masters), so it gets its internal force as well.

Parameters:

dirichlets – The Dirichlet boundary conditions of the step, located on the current equation system.

Returns:

The internal, inertial and damping force.

Return type:

DofVector

isOutputIncrement(timeStep)

Only every output-frequency-th increment: an explicit run has millions.

Parameters:

timeStep (TimeStep) – The accepted increment.

Returns:

True to write output.

Return type:

bool

kineticEnergy(V)

The kinetic energy \(\tfrac{1}{2} \sum m v^2\) of a velocity, with the unfolded lumped mass, over the degrees of freedom where it is an energy (ids_mechanicalEnergy): not on a first-order field, whose coefficient is a viscosity, nor on a field whose inertia is not mechanical.

Parameters:

V (DofVector) – The velocity vector.

Returns:

The kinetic energy.

Return type:

float

prepareIncrement(step, model, isRetry)

The topology update, when due, and the equation system; see prepareIncrement().

The topology update is due at the start of the step, and after every topology-check-frequency-th increment. At the start of the step it runs before anything sized by the equation system exists, so the mesh is final before the lumped mass, the multi-point-constraint condensation and the critical time step are derived from it. It runs here, before the next increment is proposed, because a refinement may lower the stable time increment that increment has to use.

Parameters:

• step – The step being solved.

• model (FEModel) – The model tree.

• isRetry (bool) – Never True: this solver does not retry increments.

reportTopologyChangeConservation(lumpedTotalsBefore, momentumBefore, kineticBefore, V, model)

Report what a topology change did to the quantities that ought to survive it.

Refinement interpolates the velocity onto new nodes and re-lumps every per-element scalar property, and the three invariants behave differently under that:

• Every field’s own lumped total is conserved exactly, by the same geometric identity regardless of what that field’s coefficient physically is: the children of a refined element tile it and carry the same value. Checked per FIELD, not per family, via ConservationCheck – summing across fields first, even ones that agree on units, would let a violation in a numerically small field hide inside a numerically large one. Violating any single field’s total raises.

• Linear momentum is conserved exactly for a spatially uniform velocity field, because the shape functions are a partition of unity and the child masses sum to the parent’s. For a general field the discrepancy is second order in the velocity gradient across the parent: discretisation error, not a defect. Reported, not enforced. Restricted to the fields whose inertia is a mass, per carriesLinearMomentum(): a viscosity or a non-mechanical inertia has no associated momentum, and a rotational inertia has an angular one that must not be added to a linear one. Across those fields the sum IS legitimate, because a real total system momentum is exactly the sum of its parts’ momenta.

• Kinetic energy is not conserved by interpolation plus re-lumping, and it is the most sensitive of the three, being quadratic in the interpolation error. Reported as a relative jump; more than roughly a percent is a reason to look at the transfer rather than to believe the physics. Restricted to the fields whose 0.5*m*v**2 is an energy – the mass-carrying ones and the rotational ones, a wider set than momentum admits – and reason, and summed across them for the same reason momentum is.

Parameters:

• lumpedTotalsBefore (dict[str, float]) – Every first- and second-order field’s own lumped total before the change, by field name; see _perFieldLumpedTotals().

• momentumBefore (ndarray) – Per-component momentum before the change.

• kineticBefore (float) – Kinetic energy before the change.

• V (DofVector) – The velocity vector of the rebuilt system.

• model (FEModel) – The model tree.

Raises:

RuntimeError – If any field’s own lumped total changed by more than its conservation tolerance.

secondOrderMomentum(mass, V, model)

The linear momentum of the fields whose inertia is a mass, per spatial component; see linearMomentum().

Parameters:

• mass (DofVector) – The lumped mass to weight with. Pass the unfolded mass: a multi-point-constraint slave carries real velocity, and folding its mass onto its masters would drop its momentum.

• V (DofVector) – The velocity vector.

• model (FEModel) – The model tree, for the fields’ spatial dimension.

Returns:

The momentum, one entry per spatial component.

Return type:

np.ndarray

solveIncrement(U_n, dU, V, netForce, internalForce, Minv, stepActions, model, timeStep, prevTimeStep)

Standard explicit update scheme to solve for an increment.

Parameters:

• U_n (DofVector) – The old solution vector, overwritten with the new one.

• dU (DofVector) – The solution increment, overwritten.

• V (DofVector) – The old velocity vector, overwritten with the new one.

• netForce (DofVector) – The net force of the last increment, which accelerates the nodes in this one, overwritten with the net force of this increment.

• internalForce (DofVector) – Overwritten with the elements’ sum of this increment, the internal force.

• Minv (DofVector) – The inverse lumped mass.

• stepActions (list) – The list of active step actions.

• model (FEModel) – The model tree.

• timeStep (TimeStep) – The time step.

• prevTimeStep (TimeStep) – The previous time step.

Returns:

A tuple containing

• the new solution vector

• the new velocity vector

• the net force of this increment

Return type:

tuple[DofVector,DofVector,DofVector]

updateConstraintConnectivity(model)

Let the constraints whose connectivity is the outcome of a search re-run that search.

Only the constraints in _dynamicConnectivityConstraints that are switched on are ticked (see updateConnectivity()), and the caller ticks them only every contact-update-frequency increments. Both matter: a node-to-surface search is O(slaves x facets) in Python, which is affordable once per increment of an implicit analysis – where it is amortised over a Newton loop and a linear solve – and not affordable tens of thousands of times. What makes throttling defensible rather than merely cheap is that an explicit time step is tiny: between two searches a node moves V * dT * frequency, orders of magnitude below a facet dimension. The caller reports the motion actually accumulated so that this can be checked against a given model instead of assumed.

Parameters:

model (FEModel) – The model tree.

Returns:

Whether any constraint’s DOF footprint changed, i.e. whether the equation system has to be rebuilt.

Return type:

bool

validateModelCapabilities(model)

Refuse the model features this solver cannot integrate, on top of the base checks.

Two beyond validateModelCapabilities():

A model modifier that would act after the analysis has started without configuring ``topology-check-frequency``. Modifiers that act only at simulation start are served by the initial topology update before the increment loop. Modifiers that act mid-run require topology-check-frequency to be set to a non-zero multiple of output-frequency; otherwise they are refused to prevent silently never adapting.

A constraint that introduces its own scalar variables (a Lagrange multiplier, an indirect-control unknown). Those DOFs carry no inertia, so their inverse lumped mass is zero and the explicit update leaves them untouched forever – the constraint would appear active and enforce nothing. Penalty-type constraints, which act through nodal forces alone, are supported: see assembleConstraintForces(). A constraint switched off from outside (switchedOn) is not refused: it takes no part in the step, and its scalar variables leave the equation system with it. One that is switched on is refused even while it is out of reach: it may come within reach during the step.

Parameters:

model (FEModel) – The model tree.

writeResultsToModel(U, dU, P, V, A)

Write the results of an increment to the model, through the model’s layout of the equation system: the solution, its increment, the internal, inertial and damping force, the velocity and the acceleration to the node fields, and the solution to the scalar variables (see writeResultsToModel()); the incremental degrees of freedom are set to zero first, from where the next increment starts them (see zeroIncrementalCoefficients()).

Written at every accepted increment, not only on output increments, for two reasons: an h-adaptivity event can fall on any increment and its interpolation reads these entries, and a restart checkpoint written from the node fields is the only way an explicit run can resume with its kinetic state intact. It is an O(nDof) copy.

Parameters:

• U (DofVector) – The solution vector.

• dU (DofVector) – The solution increment vector of the increment.

• P (DofVector | None) – The internal, inertial and damping force (see internalInertialAndDampingForce()); None to leave the model’s as it is.

• V (DofVector) – The velocity vector.

• A (DofVector | None) – The acceleration at the end of the increment (see accelerationAtTheEndOfTheIncrement()); None to leave the model’s as it is.

checkpointedState: dict | None = {'_conservationCheck': <class 'edelweissfe.solvers.base.conservationchecks.ConservationCheck'>, '_energyAtStepStart': <class 'float'>, '_externalWork': <class 'float'>, '_netForce': <class 'numpy.ndarray'>, '_warnedAboutMissingExternalWork': <class 'bool'>, '_warnedAboutMissingInternalEnergy': <class 'bool'>, 'prevTimeStep': <class 'edelweissfe.timesteppers.timestep.TimeStep'>}

The last completed increment and the net force, both read by the central-difference velocity update of the next increment, and the external work, an accumulator that cannot be reconstructed from the solution, with the energy at the step start it adds to. The velocity and the acceleration are the model’s (entries V and A of its node fields).

equationSystemLayout

What the model’s current equation system is made of (see buildEquationSystem()): the numbering of its degrees of freedom, and how the state and the results pass between the model and the vectors of the increments.

firstOrderFields

Fields integrated first order (forward Euler) and second order (central difference) in time, and their DOFs; derived in _classifyFieldsByScheme().

ids_mechanicalEnergy

Second-order DOFs whose 0.5*m*v**2 is a mechanical energy (mass-carrying and rotational); see phenomena.carriesKineticEnergy.

linearMomentumFields

Second-order fields whose inertia is a mass; the momentum diagnostic walks these.

nonMechanicalSecondOrderFields

Second-order fields whose inertia is neither a mass nor a rotational inertia; reported separately and kept out of both balances; see phenomena.inertiaKind.

prevTimeStep

The last completed increment; the velocity update reads 0.5 * (dT + dT_prev). None makes the first increment of a cold step the leapfrog’s half step; a resumed step uses the checkpointed one, since the checkpointed velocity already carries the half-step offset.

supportsMPC = True

Whether this solver supports master-slave condensation / multi-point constraints (e.g. surface ties). Subclasses supporting MPCs must set this to True.

supportsModelModifiers = True

Model modifiers act before the increment loop and every topology-check-frequency increments.

class edelweissfe.solvers.nonlinearexplicitdynamic.NEDSchema(courantNumber=0.8, outputFrequency=1000, topologyCheckFrequency=0, contactUpdateFrequency=100, expectSecondOrderFields=None, expectFirstOrderFields=None, reportPerformance=False, lumpedQuantityConservationTolerance=1e-06)

The options of the *solver datalines and of an >>options block routed to this solver, owned by this module and never mutated from outside it.

Mirrors NED.SolverSpecificOptions one-for-one; the plain self.options dict remains the actual source of truth consulted at runtime (see NISTSchema for why). The hyphenated option names are not valid Python identifiers, hence the optionName indirection.

Every option here is a property of the ANALYSIS – how far below the stability limit to run, how often to report, how often to let the mesh or a contact search change. Nothing here describes the model: which fields exist, which time derivative each carries, and what its coefficient means are all read from the model and from edelweissfe.config.phenomena, never asked of the deck. See NED._classifyFieldsByScheme().

The two expect-*-fields options are the one apparent exception, and are not one: they describe nothing and control nothing. They state what the deck’s author believes the derivation will produce, and NED._checkDerivedSchemeAgainstTheDeck() stops the run if it does not. Omitted – as in every deck that predates them – they assert nothing at all.

Attributes:

expectFirstOrderFields

expectSecondOrderFields

Parameters:

• courantNumber (float | None)

• outputFrequency (int | None)

• topologyCheckFrequency (int | None)

• contactUpdateFrequency (int | None)

• expectSecondOrderFields (list | None)

• expectFirstOrderFields (list | None)

• reportPerformance (bool)

• lumpedQuantityConservationTolerance (float | None)

NED

Option	Type / Default	Description
courant-number	float, default=0.8	The fraction of the critical time step actually used.
output-frequency	int, default=1000	The increment interval at which progress is logged and field outputs are finalized. Must be at least 1: it is a divisor, so 0 does not disable reporting.
topology-check-frequency	int, default=0	The increment interval at which model modifiers are offered a chance to change the mesh. 0 runs the topology update only once, before the increment loop. Must be a multiple of output-frequency, because a marker reads the last finalized field output.
contact-update-frequency	int, default=100	The increment interval at which constraints whose connectivity is the outcome of a search (contact) re-run that search. 0 disables periodic updates mid-run.
expect-second-order-fields	str, default=None	Fields the deck expects to be integrated SECOND order in time (central difference), as a comma-separated list. An ASSERTION, not a control: the scheme is always derived from the assembled inertia and damping, and this is compared against that derivation afterwards, never consulted to decide anything. A named field that turns out first order aborts the step. Optional; fields not named are not checked.
expect-first-order-fields	str, default=None	Fields the deck expects to be integrated FIRST order in time (forward Euler), as a comma-separated list. The counterpart of expect-second-order-fields and checked the same way. Optional; fields not named are not checked.
report-performance	bool, default=False	Print the cumulative performance table (totals since the start of the step) on every progress report. Off by default. An explicit run is millions of increments long, so the table printed at the end of the step is otherwise the only one anyone sees.
lumped-quantity-conservation-tolerance	float, default=1e-06	Relative-change tolerance for the per-topology-change conservation check on a row-sum-lumped quantity (mass, a first-order field’s viscosity, a second-order field’s non-mechanical inertia). Children of a refined element tile it and carry the same value, so this is conserved geometrically; the default is already generous relative to floating-point precision and exists to catch a genuinely wrong refinement or lumping, not to absorb ordinary quadrature noise. Raise this only when a specific, understood run trips it by a small margin (see edelweissfe.solvers.base.conservationchecks).

Rebuilding the equation system during a step

An explicit analysis with contact rebuilds its equation system while it runs: a contact search re-assigns which nodes a constraint couples, which changes that constraint’s degree-of-freedom footprint, and everything sized by the footprints has to follow.

Such a change is much narrower than it looks. It touches no node, no field, no scalar variable and no element, so the degree-of-freedom numbering and every element’s indices stay exactly as they were – and with them the lumped inertia, its inverse, the mass-proportional damping rate and the multi-point-constraint transformation, none of which is a function of the constraints’ footprints. The solver therefore re-locates the constraints alone (refreshConstraintIndices()) and keeps those operators, rather than constructing a new DofManager and assembling them again over every element.

That shortcut is only valid while everything the operators are a function of is unchanged, so it is taken only after checking the element set, the multi-point constraints, the node count and the scalar variables against what was recorded when they were last assembled, and only when no step action changes a material property mid-step. Anything that does not match falls back to a full build. A topology change – a mesh refinement – is not a connectivity change and never takes this path at all: it builds the system from scratch.

Work per increment outside the elements

An explicit increment is cheap, so the work around the element loop counts. Second-order degrees of freedom are updated by the central difference with mass-proportional damping, first-order ones by forward Euler. The damping rate and a second-order mask are formed once per equation system, so the update runs on whole vectors instead of indexed subsets, and gives the same result to the bit.

NEDParallel - Nonlinear Explicit Dynamic (parallel)

Parallel implementation of the NED solver.

class edelweissfe.solvers.nonlinearexplicitdynamicparallel.NEDParallel(jobInfo, journal, **kwargs)

Attributes:

convergenceTable

incrementStatus

mpcTransformation

Methods

accelerationAtTheEndOfTheIncrement(...)	The acceleration this solver publishes as the node field entry A of the increment just completed, the one at its end (see the module documentation, "The velocity and the acceleration it publishes"): with the velocity at that time, half an increment after the staggered one, it satisfies the equation of motion, \(\left( \boldsymbol{M} + \tfrac{1}{2} \Delta t \, \boldsymbol{C} \right) \boldsymbol{a} = \boldsymbol{f}_\mathrm{net} - \boldsymbol{C} \boldsymbol{v}_{n+1/2}\).
acceptIncrement(step, model, timeStep)	Commit the increment to the model; see acceptIncrement().
applyDirichletToResidual(timeStep, R, dirichlets)	Impose the Dirichlet BCs on the residual.
applyOptionsOverride(fieldValues)	Apply a partial override of this solver's own schema fields onto self.options.
applyStepActionsAtIncrementStart(model, ...)	Called when all step actions should be applied at the start of a step.
applyStepActionsAtStepEnd(model, stepActions)	Called when all step actions should finish a step.
applyStepActionsAtStepStart(model, step)	Apply every step action's step-start part (initial conditions, material initialization, prescribed fields, ...).
assembleConstraintForces(constraints, U_np, ...)	Evaluate every active constraint and add its residual contribution to the given vector; a constraint that is switched off (see active) is skipped.
assembleLoads(loads, U_np, dU, P, K, timeStep)	Assemble the residual contributions of all loads: each adds the negative of the external load it applies, and the tangent of that contribution, which is nonzero only for a load that follows the deformation; see The residual and the tangent.
assembleStiffnessCSR(K)	The assembled system matrix in CSR format, see assembledMatrix().
attemptIncrement(step, model, timeStep)	The periodic contact search, when due, and the central-difference update; see attemptIncrement().
beginStep(step, model, ...)	Start the step; see beginStep().
buildEquationSystem(model, step[, previous])	Build the equation system and everything sized by it.
buildMPCTransformation(model[, stepActions])	Collect the linear dependency records from all multi-point constraints of the model and assemble the master-slave condensation operator for the current equation system.
checkConvergence(R, ddU, dU, F, ...)	Check the convergence, individually for each field and for the scalar variables, based on the current flux residual and the correction, using the convergence criterion of the solver (see edelweissfe.solvers.base.convergencecriteria).
checkMPCDirichletConflicts(transformation, ...)	Raise if any Dirichlet boundary condition of the step prescribes a DOF that is a slave DOF of a multi-point constraint.
collectLoads(stepActions)	The loads among a step's actions: every action that is a LoadBase, whatever step action module it is defined in, so that a new kind of load needs no change to any solver.
endStep(step, model)	Leave the model's state at the end time of the step, and report the step's performance timing.
energyOfTheModelsState()	The energy the model holds at the start of a cold step, which its energy balance starts from (see the module documentation, "Energy balance"): the kinetic energy of the model's velocity, read with the solution when the equation system was built, and the internal (strain) energy the materials publish for the state the elements hold.
formNetForceOfTheModelsState(dirichlets)	Form the net force a cold step starts from: that of the state the model holds, \(\boldsymbol{f}_\mathrm{net} = \boldsymbol{M} \boldsymbol{a}_0 + \boldsymbol{C} \boldsymbol{v}_0\), from its velocity and acceleration, read with the solution when the equation system was built, and the lumped operators; nothing is evaluated (see the module documentation, "The net force is solver state, and a step starts from the model's state").
getCriticalTimeStepForExplicitDynamics(model, U)	Compute the critical time step for explicit dynamics.
getRestartData()	The state this solver carries from one increment to the next, as it is; see checkpointedState.
internalInertialAndDampingForce(dirichlets)	The force this solver publishes as the node field result P of the increment just completed: the internal force with the inertial and damping forces, \(\boldsymbol{f}_\mathrm{int} + \boldsymbol{M} \boldsymbol{a} + \boldsymbol{C} \boldsymbol{v}\), without the loads and the constraints (see the module documentation, "The force it publishes").
isOutputIncrement(timeStep)	Only every output-frequency-th increment: an explicit run has millions.
kineticEnergy(V)	The kinetic energy \(\tfrac{1}{2} \sum m v^2\) of a velocity, with the unfolded lumped mass, over the degrees of freedom where it is an energy (ids_mechanicalEnergy): not on a first-order field, whose coefficient is a viscosity, nor on a field whose inertia is not mechanical.
linearSolve(A, b)	Solve the linear equation system.
locateConstrainedDofs(dirichlets)	Determine, up front, which global DOFs each Dirichlet BC constrains.
prepareIncrement(step, model, isRetry)	The topology update, when due, and the equation system; see prepareIncrement().
printResidualOutlierNodes(residualOutliers)	Print which nodes have the largest residuals.
reportTopologyChangeConservation(...)	Report what a topology change did to the quantities that ought to survive it.
schema
secondOrderMomentum(mass, V, model)	The linear momentum of the fields whose inertia is a mass, per spatial component; see linearMomentum().
setRestartData(data)	Restore the state getRestartData() returned, so that the resumed step continues exactly where the uninterrupted run would have.
solveIncrement(U_n, dU, V, netForce, ...)	Standard explicit update scheme to solve for an increment.
updateConstraintConnectivity(model)	Let the constraints whose connectivity is the outcome of a search re-run that search.
updateRigidBodies(model, timeStep)	Refresh the kinematics of all rigid bodies in the model after a converged increment.
updateTopologyAndConnectivity(model, step[, ...])	Run the topology update, then let every mesh-dependent consumer catch up on it.
validateModelCapabilities(model)	Refuse the model features this solver cannot integrate, on top of the base checks.
writeResultsToModel(U, dU, P, V, A)	Write the results of an increment to the model, through the model's layout of the equation system: the solution, its increment, the internal, inertial and damping force, the velocity and the acceleration to the node fields, and the solution to the scalar variables (see writeResultsToModel()); the incremental degrees of freedom are set to zero first, from where the next increment starts them (see zeroIncrementalCoefficients()).

computeElements

beginStep(step, model, fieldOutputController, outputmanagers)

Start the step; see beginStep(). The equation system is built by the first prepareIncrement().

Parameters:

• step – The step to be solved.

• model – The model tree.

• fieldOutputController – The field output controller.

• outputmanagers – The output managers.

computeElements(elements, U_np, dU, P, timeStep)

Loop over all elements, and evaluate them: their residual contributions, summed, are the internal force. Called by solveIncrement() once per increment.

Parameters:

• elements (list) – The finite elements, by number.

• U_np (DofVector) – The current solution vector.

• dU (DofVector) – The solution increment vector.

• P (DofVector) – The vector the elements’ residual contributions are added to; zero.

• timeStep (TimeStep) – The time step.

Returns:

• The internal force.

• The internal (strain) energy the elements publish.

Return type:

tuple[DofVector,float]

NEDParallel

Option	Type / Default	Description
courant-number	float, default=0.8	The fraction of the critical time step actually used.
output-frequency	int, default=1000	The increment interval at which progress is logged and field outputs are finalized. Must be at least 1: it is a divisor, so 0 does not disable reporting.
topology-check-frequency	int, default=0	The increment interval at which model modifiers are offered a chance to change the mesh. 0 runs the topology update only once, before the increment loop. Must be a multiple of output-frequency, because a marker reads the last finalized field output.
contact-update-frequency	int, default=100	The increment interval at which constraints whose connectivity is the outcome of a search (contact) re-run that search. 0 disables periodic updates mid-run.
expect-second-order-fields	str, default=None	Fields the deck expects to be integrated SECOND order in time (central difference), as a comma-separated list. An ASSERTION, not a control: the scheme is always derived from the assembled inertia and damping, and this is compared against that derivation afterwards, never consulted to decide anything. A named field that turns out first order aborts the step. Optional; fields not named are not checked.
expect-first-order-fields	str, default=None	Fields the deck expects to be integrated FIRST order in time (forward Euler), as a comma-separated list. The counterpart of expect-second-order-fields and checked the same way. Optional; fields not named are not checked.
report-performance	bool, default=False	Print the cumulative performance table (totals since the start of the step) on every progress report. Off by default. An explicit run is millions of increments long, so the table printed at the end of the step is otherwise the only one anyone sees.
lumped-quantity-conservation-tolerance	float, default=1e-06	Relative-change tolerance for the per-topology-change conservation check on a row-sum-lumped quantity (mass, a first-order field’s viscosity, a second-order field’s non-mechanical inertia). Children of a refined element tile it and carry the same value, so this is conserved geometrically; the default is already generous relative to floating-point precision and exists to catch a genuinely wrong refinement or lumping, not to absorb ordinary quadrature noise. Raise this only when a specific, understood run trips it by a small margin (see edelweissfe.solvers.base.conservationchecks).