Model modifiers

ModelModifiers dynamically alter mesh topology, element/node sets, material states, or field allocations during an analysis step.

Unlike constraints, step actions or output managers – which act on a fixed mesh – a model modifier may change the mesh topology itself during an analysis: adding or removing nodes and elements, re-partitioning element/node sets and surfaces, and reallocating the solution fields. A modifier is declared with the *modelModifier keyword. At the start of every increment the solver runs all modifiers to a fixed point via update(), then lets mesh-dependent consumers catch up, then solves; when the topology changed, the equation system (DOF manager, sparsity pattern, solution vectors and any multi-point-constraint transformation) is rebuilt first. A modifier itself is written as two halves – plan(), which decides and may read solution state, and apply(), which carries the decision out and may not. See The topology pipeline: how a model changes itself for the full contract, why it is split that way, and what a new modifier must implement; this page covers the individual modifiers.

Topological containers have stable identity. NodeSet, ElementSet, EntityBasedSurface and NodeField are created once and mutated in place (replaceMembers/replaceData/resize) by a mesh-mutating modifier – they are never replaced with a new object under the same name. A component that simply caches one of these, e.g. self._nodes = model.nodeSets["mySet"] in __init__ and later iterates self._nodes, sees any later refinement automatically: no observer, no MeshDependent, no re-fetch. This is deliberate – it lets a student write a new constraint, step action or output without ever thinking about AMR.

A component that additionally pre-sizes a derived array or object to a container’s current size (e.g. a Dirichlet BC’s delta tiled to len(nSet), or a field output’s result collector pinning a snapshot element list) does need to notice a resize, since the array does not grow with the container on its own. Two mechanisms remain, narrowed to exactly these cases:

Lazy version check – each mutable container carries a _version counter, bumped on every in-place mutation. A component compares its own last-seen version against the container’s current one at its own per-increment entry point (e.g. StepActionBase._checkSetChanged, ConstraintBase._checkSetChanged) and recomputes the derived array only on a mismatch – see dirichlet, nodeforces and ElementFieldOutput for examples. This needs no registration and therefore has no observer lifecycle to leak.
Registered mesh dependent – for derived geometry that must be regenerated before the next equation-system rebuild (facet-based contact and tie; see below), a component registers itself via registerMeshDependent() and implements refresh(). Once per increment, after every modifier has settled, refreshMeshDependents() hands it the net change since it last looked – added/removed nodes and elements, the parent -> children map, the per-face child tiling, and which node/element sets or surfaces were touched (with touchesSurface/touchesNodeSet/touchesElementSet early-outs so a consumer can skip a change that doesn’t concern it).

The synchronous push notification that used to exist alongside this has been removed: with modifiers running to a fixed point, a per-mutation callback fires mid-pipeline and hands the consumer a state that no longer exists by the time the solve begins. See The topology pipeline: how a model changes itself.

`hAdaptivity` - Hanging-node h-adaptivity for HEX20

Module edelweissfe.modelmodifiers.adaptivity.hadaptivity

Dynamic h-adaptivity model modifier for HEX20 hanging-node AMR.

Dynamic adaptive \(h\)-refinement of 20-node serendipity hexahedra (GC3D20 / GC3D20R) in the small-strain, multifield (displacement + nonlocal damage) regime. Each increment the modifier evaluates one or more markers (see below), subdivides every marked element into splitFactor**3 children (default splitFactor=2, i.e. eight octree children; splitFactor=3 gives a 3x3x3 split into 27, and so on – honouring curved edges via the parent isoparametric map), enforces a one-level face-balance, transfers the converged nodal values (parent isoparametric interpolation) and quadrature-point history (via a pluggable state-transfer strategy, see edelweissfe.adaptivity.statetransfer) to the children, and couples the resulting non-conforming interface with an exact hanging-node multi-point constraint. Element/node sets, sections and element-based surfaces are propagated to the children so material assignment and surface loads stay consistent.

A child is assigned its section exactly like an element of the initial mesh: parameters set by >>materialParameterFromField (e.g. a random strength field) are evaluated at the child’s own center. The analytical field is a function of position only, so a refined region samples the same field realization more finely – just as a mesh that was refined from the start would.

The octree mirror only ever tracks the refineable 20-node elements: a model that also contains elements of a different kind (e.g. lower-order contact-facet elements bonded to the mesh) is left untouched by construction, since anything without exactly 20 nodes is skipped automatically. To restrict refinement explicitly (rather than relying on the node-count heuristic), or to refine only part of a purely 20-node mesh, set refineElSet (falls back to elSet if given) to the element set that should become the octree mirror. Marking (see below) is likewise restricted to the octree mirror’s own elements by default – there is nothing to gain from evaluating the marking expression on an element that can never be refined anyway, and most non-solid element types don’t expose most quadrature-point results in the first place.

The non-conforming 2:1 interface is coupled kinematically rather than by mortar: octree refinement is non-conforming but nested, the QUAD8 face-trace (and 3-node quadratic edge) spaces are invariant under the axis-aligned affine sub-maps of a uniform subdivision (of any factor, not only bisection), so pinning each hanging (slave) node to the coarse serendipity trace, \(u_s = \sum_a N_a(\xi_s)\, u_{m_a}\), is exact. The same field-independent weights apply to every field on the node (equal-order serendipity), so a single record per hanging node covers displacement and nonlocal damage alike. The constraint itself (*constraint, type=hangingnode) is documented under Constraints.

Vocabulary. A root element is an element of the mesh before refinement; every refined element lies inside exactly one. The reference cube of a root carries an integer reference lattice, fine enough that every node a subdivision can create sits on a lattice point, and every element covers an exact reference box of it. A root entity is a vertex, edge or face of the root mesh, named by the labels of its corner nodes. A node’s key is the lowest-dimensional root entity containing it plus its exact integer position there. A node hangs on an element if it lies on the element’s boundary without being one of its own nodes.

Topology is decided exactly, never from coordinates. Which nodes coincide, which node lies on which face or edge of an element, which elements share a face (2:1 balancing), which nodes hang – and with which weights – and which new nodes join a node set, are all decided topologically and in exact integer arithmetic, the way established AMR libraries do it (deal.II, p4est, MFEM):

every element knows its root element and its exact box in the root’s reference cube (from the subdivision order), so every node sits at an exact lattice point of its root;
the root mesh’s vertices, edges and faces are named by the labels of their corner nodes, and a point on them is measured in a frame those labels fix – two roots sharing a face therefore derive the identical name and coordinates for it in any relative orientation (RootEntityTable);
a hanging node is a node on the boundary of an active element that is not one of its own nodes; its master entity (edge or face) follows from its exact reference coordinates in that element, and its weights are the element’s own shape functions there, the exact coarse trace;
a new node joins a node set iff it lies, by its exact lattice position, on a parent face or edge whose nodes are all in the set.

No tolerance enters any of these decisions, so curved, warped and arbitrarily oriented unstructured meshes (e.g. from Cubit) behave exactly like axis-aligned ones. The requirements are a conforming root mesh and a refinement depth of at most 20 levels (the depth of the exact lattice). A collapsed root element – one that repeats a corner node, e.g. a degenerated wedge – or a non-manifold root mesh stops the analysis at setup with an error naming the element.

Node sets are inherited through faces and edges. A new node joins a node set iff it lies on a face or edge of its parent whose nodes are all in the set. Two consequences follow from node membership alone: a set that holds every node of a face – even if it was meant as the face’s perimeter only – gains the new nodes inside that face, and a set never gains the new nodes inside a parent element. Where a different membership is intended, define the boundary condition through an element set (or an element-based surface) instead, which refinement propagates to the children exactly.

2:1 balancing is face-based. After each refinement, any element with a face neighbour two or more levels finer is refined too, until every face-adjacent pair differs by at most one level. Elements that touch only along an edge or at a vertex are not balanced: next to a level-0 element, an element that shares only an edge with it may be at level 2 or deeper. Such an interface stays exactly conforming – the finer nodes on that edge hang on the coarsest element’s edge, through constraint chains that are resolved exactly and verified – so this affects mesh grading and local accuracy (an abrupt size jump along the edge), not correctness. With maxLevel=1 (the default) no two-level jump can occur at all. Edge or vertex balancing (as p4est’s edge/full balance) is not implemented; it would refine additional elements and change the results of multi-level runs.

Invariants, checked after every adaptation. A violation raises an error naming the offending nodes instead of producing a silently non-conforming mesh:

one node per point and one point per node (no split or merged nodes);
every node on the boundary of an active element is one of its nodes or a hanging-node slave (check_conformity());
every hanging node’s weights are verified in exact arithmetic: they vanish off the chosen edge or face and reproduce the constant, the node’s reference coordinates and every quadratic monomial;
the hanging-node constraint can never be overruled by another multi-point constraint claiming the same slave (it is registered first in model order; a dropped hanging-node record is an error).

class edelweissfe.adaptivity.rootentities.RootEntityTable(topology, latticeUnits)[source]

The vertices, edges and faces of the root mesh, keyed by root corner labels.

Parameters:

topology – The TopologyBase of the root elements.
latticeUnits (int) – \(M\): reference coordinate \(\xi\) of a root element is the lattice integer \((\xi + 1)\,M\).

Methods

`add_root`(rootEid, connectivity[, componentId])	Register a root element by its node labels (in the topology's local node order).
`key`(rootEid, point)	The `NodeKey` of a lattice point (three integers in [0, 2M]) of root `rootEid`.
`local_point`(rootEid, key)	The lattice point, in root element `rootEid`, of a key whose entity lies in that root's closure (see `roots_sharing()`).
`root_face_position`(rootEid, faceIndex, point)	The root face `faceIndex` (into `topology.faces`) of `rootEid` as a `RootFace`, and the position of `point` (a lattice point on it, possibly on its boundary) in the face's canonical frame -- identical from either root sharing the face.
`roots_sharing`(key)	The root elements whose closure contains the entity of `key`.
`roots_touching`(rootEid)	The root elements sharing at least one vertex (hence possibly an edge or a face) with `rootEid`, including itself.

add_root(rootEid, connectivity, componentId=0)[source]

Register a root element by its node labels (in the topology’s local node order).

Raises:

ValueError – For a collapsed root element (a corner label repeated), whose faces and edges cannot be named by their corners, and for a face owned by more than two roots (a non-manifold mesh).

Parameters:

rootEid (int)
componentId (int)

key(rootEid, point)[source]

The NodeKey of a lattice point (three integers in [0, 2M]) of root rootEid.

Parameters:: rootEid (int)
Return type:: NodeKey

local_point(rootEid, key)[source]

The lattice point, in root element rootEid, of a key whose entity lies in that root’s closure (see roots_sharing()).

Parameters:

rootEid (int)
key (NodeKey)

Return type:

tuple

root_face_position(rootEid, faceIndex, point)[source]

The root face faceIndex (into topology.faces) of rootEid as a RootFace, and the position of point (a lattice point on it, possibly on its boundary) in the face’s canonical frame – identical from either root sharing the face.

Parameters:

rootEid (int)
faceIndex (int)

Return type:

tuple

roots_sharing(key)[source]

The root elements whose closure contains the entity of key.

Parameters:: key (NodeKey)
Return type:: list

roots_touching(rootEid)[source]

The root elements sharing at least one vertex (hence possibly an edge or a face) with rootEid, including itself.

Parameters:: rootEid (int)
Return type:: set

class edelweissfe.adaptivity.refinement.AdaptiveMesh(splitFactor=2, topology=None, reserve_labels=None)[source]

Octree hierarchy of HEX20 elements: refinement, 2:1 balancing and hanging-node classification.

Every decision about how elements and nodes relate – which nodes coincide, which node lies on which face or edge, which elements share a face, which nodes hang and with which weights, which new nodes join a node set – is taken topologically and in exact integer arithmetic, never from physical coordinates: each element knows its root and its exact reference box, and the root mesh’s vertices, edges and faces are named by their corner labels (RootEntityTable). This holds on curved, warped and arbitrarily oriented unstructured meshes alike. The only requirement is a conforming root mesh without collapsed elements. check_conformity() verifies the result after every adaptation.

Parameters:

splitFactor (int) – n: a refined element is split into n parts per axis, i.e. n**3 children.
topology – The element TopologyBase (default: HEX20).
reserve_labels – The model’s node-label allocator, see NodeRegistry.

Methods

`add_root`(coords, labels[, componentId])	Add a root (level-0) element: its 20 node coordinates and node labels, in C3D20 order.
`balance_2to1`()	Refine coarser elements until no face-adjacent active pair differs by more than one level.
`check_conformity`(slaveLabels)	Raise unless the active mesh is conforming up to the given hanging-node slaves.
`classify_hanging`()	All hanging nodes of the active mesh, as `HangingNode`, sorted by slave label.
`hanging_mpc_records`()	The hanging-node constraints as `{slave: [(master, weight), ...]}`, every master an independent node: a chain -- a master that is itself hanging -- is substituted through, with the same flattening the multi-point-constraint transformation applies across all constraints.
`node_keys`(eid)	The `NodeKey` of each of the element's own nodes, in slot order.
`node_lattice_points`(eid)	The lattice points, in the element's root, of the element's own nodes (slot order).
`refine`(eid)	Split an active element into `splitFactor**3` children; deactivate it and keep element sets, surfaces and node sets consistent.

add_root(coords, labels, componentId=0)[source]

Add a root (level-0) element: its 20 node coordinates and node labels, in C3D20 order.

The labels are the model’s own, and each must have been registered with NodeRegistry.seed(). componentId is the connected body the element belongs to: node identities are kept per body, so two bodies sharing a flush interface are never welded together.

Raises:: ValueError – For an unregistered label, a collapsed root element, or a non-manifold root mesh.
Parameters:: componentId (int)
Return type:: int

balance_2to1()[source]

Refine coarser elements until no face-adjacent active pair differs by more than one level.

Face adjacency is exact and topological (see _facesOnEntities()), so balancing holds on curved and warped meshes as on flat ones. Returns the number of extra elements refined.

Only face neighbours are balanced. Elements touching only along an edge or at a vertex may differ by two or more levels; the interface stays conforming (the finer edge nodes hang on the coarsest element’s edge, through exactly resolved constraint chains), so this is a grading, not a correctness, limitation. Edge/vertex balancing is deliberately not implemented: it would refine additional elements and change the results of multi-level runs.

Return type:: int

check_conformity(slaveLabels)[source]

Raise unless the active mesh is conforming up to the given hanging-node slaves.

Exact, and independent of how the hanging nodes were found. Two conditions:

One node per point: a NodeKey maps to one label, and a label to one key. A split point would be two unconnected nodes; a merged one would weld distinct points.
Every node on the boundary of an active element is one of its own nodes or a hanging-node slave. Anything else is a free node on the element’s boundary: the interface is not conforming.

Parameters:: slaveLabels – The labels of all hanging-node slaves (e.g. the keys of hanging_mpc_records()).
Raises:: TopologyError – Naming the offending nodes and elements.

classify_hanging()[source]

All hanging nodes of the active mesh, as HangingNode, sorted by slave label.

A hanging node is a node on the boundary of an active element that is not one of its own nodes (_nodesOnElementBoundaries(): exact, exhaustive, no geometry). Its master entity is the lowest-dimensional one of that element containing it – an edge or a face – and, among several such elements, the coarsest, which guarantees continuity with the coarsest trace. Its weights are the master element’s own shape functions at the node’s exact reference point, which on that edge or face are non-zero only for the entity’s nodes: the exact coarse trace, on curved and warped elements alike. They are verified in exact arithmetic before use. Rational arithmetic is slow, but the weights depend only on the node’s point in the master’s box and the entity, and an octree has few distinct such pairs: each is computed and verified once.

Return type:: list

hanging_mpc_records()[source]

The hanging-node constraints as {slave: [(master, weight), ...]}, every master an independent node: a chain – a master that is itself hanging – is substituted through, with the same flattening the multi-point-constraint transformation applies across all constraints.

Return type:: dict

node_keys(eid)[source]

The NodeKey of each of the element’s own nodes, in slot order. Computed once: an element’s reference box never changes.

Return type:: list

node_lattice_points(eid)[source]

The lattice points, in the element’s root, of the element’s own nodes (slot order).

Return type:: list

refine(eid)[source]

Split an active element into splitFactor**3 children; deactivate it and keep element sets, surfaces and node sets consistent.

Children are returned in subdivision order, ix * n * n + iy * n + iz (see subdivision_children_param()).

Raises:: TopologyError – Beyond MAXIMUM_REFINEMENT_DEPTH.
Return type:: list

Refinement markers. Which elements are refined each increment is decided by one or more >>marker sub-keywords; the modifier refines the union of their marked sets. The available types are elementSet / nodeSet / surface (geometric, typically initialOnly for a fixed pre-refinement), fieldOutput (a boolean threshold on a perElement field output), and recoveryError (a Zienkiewicz–Zhu recovered-gradient error estimator on a nodal field, with Dörfler bulk marking – for gradient-enhanced damage). The theory of error-estimator marking, the averaging vs spr recovery, and why a reactive recoveryError marker is best paired with a predictive fieldOutput marker ahead of a propagating front are all covered under Adaptive refinement: error estimation and marking.

Options:
Option	Type / Default	Description
`moduleOptions`	dict, default={}	Internal
`elSet`	str, default=None	Fallback for ‘refineElSet’ if that is not given. Does not restrict marking itself – each ‘>>marker’ scopes its own eligible elements (a fieldOutput’s associated set, an elementSet/nodeSet/surface’s members).
`refineElSet`	str, default=None	Restrict the AMR octree mirror itself to this element set, e.g. the solid elements in a mesh that also contains contact-facet elements. Elements outside this set never become octree roots and are left untouched by refinement. Defaults to ‘elSet’ if given, otherwise to every 20-node (HEX20-family) element in the model.
`maxLevel`	int, default=1	Maximum refinement level.
`minMarkedElements`	int, default=1	Minimum number of eligible elements one evaluation of the markers must mark for a refinement pass to happen; with fewer, nothing is refined and no equation system is rebuilt. Marks are not carried over to later evaluations: they are a function of the current state only. Individual markers may cap their own marks (e.g. ‘maxRefinedFraction’) or expand them (e.g. ‘halo’) before this count is taken. Default 1 refines as soon as any element is marked.
`splitFactor`	int, default=2	Number of equal parts per axis a marked element is split into (2 = octree bisection into 8 children; 3 = 3x3x3 = 27 children, etc.). The hanging-node coupling stays exact for any factor.
`elementType`	str, default=None	Element type to instantiate for children (default: like parents).
`elementProvider`	str, default=’marmot’	Element provider.
`stateTransfer`	str, default=’nearestQp’	Quadrature-point state-transfer strategy for the whole state block: nearestQp\|projection\|virgin.
`stateTransferOverrides`	str, default=None	Per-state-variable overrides routing named variables to a different strategy, e.g. ‘strain:projection, stress:virgin’. Comma-separated ‘name:strategy’ pairs.
`marker`	tuple, default=()	AMR marker definition. At least one is required.

Keyword: marker
Option	Type / Default	Description
`type`	str (required)	Marker type, resolved through the ‘marker’ registry: fieldOutput, elementSet, nodeSet, surface, recoveryError. The type-specific options are defined by the selected marker’s own schema.
`initialOnly`	bool, default=False	Evaluate only once at simulation start

Example (dynamic refinement of a two-field GC3D20R cantilever): testfiles/marmot/AMR_DynamicRefinement/test.inp

** Dynamic h-adaptivity regression test -- exercises the full *modelModifier path end to end:
** marking -> octree refinement -> model mutation -> nodal warm-start -> quadrature-point state
** transfer -> hanging-node MPC -> equation-system rebuild -> solve on the refined mesh.
** Two coarse C3D20 elements form a cantilever; bending concentrates stress at the fixed end,
** so the stress marker refines exactly one element (2 -> 9 active, 13 hanging nodes), creating a
** 2:1 hanging interface that the MPC couples across the displacement field. (Multi-field
** coupling across a hanging interface is covered separately by AMR_TwoFieldGC3D20R.)
** Full integration, deliberately: a two-element mesh of reduced-integration C3D20R keeps a zero-energy
** mode the boundary conditions do not suppress, so its solution is not unique and any exact change to
** the linear system (e.g. how Dirichlet conditions are eliminated) would pick a different one.
** Pins the dynamic solution against U.ref. (A true adaptive-vs-uniform-fine PHYSICS benchmark --
** PLAN_AMR.md M3 -- needs a two-run comparison harness and remains future work.)
*material, name=linearelastic, id=mat
20000, 0.2

*section, name=sec, material=mat, type=solid
all

*modelModifier, type=hAdaptivity, name=amr
>>marker, type=fieldOutput, fieldOutput=stressForAMR, operator='>', threshold=0.6
maxLevel=1

*job, name=amr, domain=3d
*solver, solver=NIST, name=S
*fieldOutput
>>perNode, elSet=all, field=displacement, result=U, name=U
>>perElement, elSet=all, result=stress, quadraturePoint=0:27, name=stressForAMR, f(x)='np.abs(x)'

*step, solver=S
maxInc=0.5, minInc=1e-4, maxNumInc=50, maxIter=25,
>>dirichlet, name=x0, nSet=x0, field=displacement, 1=0.0, 2=0.0, 3=0.0
>>dirichlet, name=x1, nSet=x1, field=displacement, 3=0.0005

*Node
1, 0.0, 1.0, 1.0
2, 0.0, 0.0, 1.0
3, 0.0, 0.0, 0.0
4, 0.0, 1.0, 0.0
5, 1.0, 1.0, 1.0
6, 1.0, 0.0, 1.0
7, 1.0, 0.0, 0.0
8, 1.0, 1.0, 0.0
9, 0.0, 0.5, 1.0
10, 0.0, 0.0, 0.5
11, 0.0, 0.5, 0.0
12, 0.0, 1.0, 0.5
13, 1.0, 0.5, 1.0
14, 1.0, 0.0, 0.5
15, 1.0, 0.5, 0.0
16, 1.0, 1.0, 0.5
17, 0.5, 1.0, 1.0
18, 0.5, 0.0, 1.0
19, 0.5, 0.0, 0.0
20, 0.5, 1.0, 0.0
21, 2.0, 1.0, 1.0
22, 2.0, 0.0, 1.0
23, 2.0, 0.0, 0.0
24, 2.0, 1.0, 0.0
25, 2.0, 0.5, 1.0
26, 2.0, 0.0, 0.5
27, 2.0, 0.5, 0.0
28, 2.0, 1.0, 0.5
29, 1.5, 1.0, 1.0
30, 1.5, 0.0, 1.0
31, 1.5, 0.0, 0.0
32, 1.5, 1.0, 0.0
*Element, type=C3D20
1, 1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12, 13, 14, 15, 16, 17, 18, 19, 20
2, 5, 6, 7, 8, 21, 22, 23, 24, 13, 14, 15, 16, 25, 26, 27, 28, 29, 30, 31, 32
*elSet, elSet=all, generate
1, 2, 1
*nSet, nSet=x0
1, 2, 3, 4, 9, 10, 11, 12
*nSet, nSet=x1
21, 22, 23, 24, 25, 26, 27, 28
*nSet, nSet=yf
2, 3, 6, 7, 10, 14, 18, 19
22, 23, 26, 30, 31
*nSet, nSet=zf
3, 4, 7, 8, 11, 15, 19, 20
23, 24, 27, 31, 32

State-variable transfer strategies

Pluggable quadrature-point state-transfer strategies for adaptive mesh refinement.

The strategy for handing a refined parent element’s quadrature-point history to its children is selectable with the stateTransfer argument (default nearestQp):

nearestQp – NearestQuadraturePointCopy: each child quadrature point copies, verbatim, the nearest parent quadrature point (matched in the parent reference cube). Admissible by construction – the recommended default.
projection – PolynomialProjection: a tensor-product polynomial is fitted to the parent quadrature-point values by least squares and resampled at the children. Smooth across octants, but may produce an inadmissible internal state.
virgin – VirginState: children keep their freshly-initialised state; history is discarded (sound only when refining ahead of the process zone).

Different state variables generally require different treatment, and which ones may be projected, copied, or reset is entirely a property of the constitutive model – there is no universal choice. (For a hypoelastic model whose stress update depends only on the strain increment, for instance, a stored total strain may be irrelevant while the true history variables – stresses, back-stresses, plastic strains, damage, hardening – are what matter.) stateTransferOverrides therefore leaves the policy to the user: it routes named state variables to their own strategy while the stateTransfer default handles all the rest. The names are the material’s / element’s own state-variable names, located within a quadrature-point block via the element’s getStateVarSlice hook and dispatched by PerStateVarStateTransfer. For example, stateTransfer=nearestQp with stateTransferOverrides='<var1>:projection, <var2>:virgin' copies everything except <var1> (projected) and <var2> (reset to its initial value).

Example (per-state-variable transfer device; the routed variable is arbitrary, chosen only to exercise the mechanism): testfiles/marmot/AMR_DynamicRefinementProjection/test.inp

** Dynamic h-adaptivity exercising the PER-STATE-VARIABLE transfer device: one named state variable
** is routed to a different transfer strategy while all others are copied from the nearest parent
** quadrature point. This is a MECHANISM test -- the variable named below is arbitrary (chosen only
** because it exists on this element); WHICH state variables should be projected vs. copied vs. reset
** is entirely a property of the constitutive model and is left to the user. Exercises the
** getStateVarSlice element hook, the PerStateVarStateTransfer composite and one non-copy strategy
** (PolynomialProjection) end to end.
*material, name=linearelastic, id=mat
20000, 0.2

*section, name=sec, material=mat, type=solid
all

*modelModifier, type=hAdaptivity, name=amr
>>marker, type=fieldOutput, fieldOutput=stressForAMR, operator='>', threshold=0.6
maxLevel=1
stateTransferOverrides='strain:projection'

*job, name=amr, domain=3d
*solver, solver=NIST, name=S
*fieldOutput
>>perNode, elSet=all, field=displacement, result=U, name=U
>>perElement, elSet=all, result=stress, quadraturePoint=0:27, name=stressForAMR, f(x)='np.abs(x)'

*step, solver=S
maxInc=0.5, minInc=1e-4, maxNumInc=50, maxIter=25,
>>dirichlet, name=x0, nSet=x0, field=displacement, 1=0.0, 2=0.0, 3=0.0
>>dirichlet, name=x1, nSet=x1, field=displacement, 3=0.0005

*Node
1, 0.0, 1.0, 1.0
2, 0.0, 0.0, 1.0
3, 0.0, 0.0, 0.0
4, 0.0, 1.0, 0.0
5, 1.0, 1.0, 1.0
6, 1.0, 0.0, 1.0
7, 1.0, 0.0, 0.0
8, 1.0, 1.0, 0.0
9, 0.0, 0.5, 1.0
10, 0.0, 0.0, 0.5
11, 0.0, 0.5, 0.0
12, 0.0, 1.0, 0.5
13, 1.0, 0.5, 1.0
14, 1.0, 0.0, 0.5
15, 1.0, 0.5, 0.0
16, 1.0, 1.0, 0.5
17, 0.5, 1.0, 1.0
18, 0.5, 0.0, 1.0
19, 0.5, 0.0, 0.0
20, 0.5, 1.0, 0.0
21, 2.0, 1.0, 1.0
22, 2.0, 0.0, 1.0
23, 2.0, 0.0, 0.0
24, 2.0, 1.0, 0.0
25, 2.0, 0.5, 1.0
26, 2.0, 0.0, 0.5
27, 2.0, 0.5, 0.0
28, 2.0, 1.0, 0.5
29, 1.5, 1.0, 1.0
30, 1.5, 0.0, 1.0
31, 1.5, 0.0, 0.0
32, 1.5, 1.0, 0.0
*Element, type=C3D20
1, 1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12, 13, 14, 15, 16, 17, 18, 19, 20
2, 5, 6, 7, 8, 21, 22, 23, 24, 13, 14, 15, 16, 25, 26, 27, 28, 29, 30, 31, 32
*elSet, elSet=all, generate
1, 2, 1
*nSet, nSet=x0
1, 2, 3, 4, 9, 10, 11, 12
*nSet, nSet=x1
21, 22, 23, 24, 25, 26, 27, 28
*nSet, nSet=yf
2, 3, 6, 7, 10, 14, 18, 19
22, 23, 26, 30, 31
*nSet, nSet=zf
3, 4, 7, 8, 11, 15, 19, 20
23, 24, 27, 31, 32

Compatibility with facet-based contact and tie

Refining a solid whose surface feeds a facet-based contact or tie constraint works out of the box: the modifier keeps the relevant *surface definition in sync with the refined child faces, and the constraint – a MeshDependent – regenerates its facets from it. nodetodeformablesurfacepenalty notices via changesSince() at its own next connectivity update (a pull, since that tick already runs before the equation system is rebuilt); a tie has no such early tick of its own (its only hook is called from inside that rebuild, too late to safely swap in new facet elements), so it reconciles via the model’s push notification instead. Either way, no separate wiring is needed. A model with both solid elements to be refined and pre-existing contact-facet elements should restrict the octree mirror explicitly with refineElSet (see above), since a facet element is never itself 20-node but need not be excluded by name:

*modelModifier, type=hAdaptivity, name=amr
result=stress
expression='x > 300.0'
reducer=absmax
maxLevel=1
refineElSet=lower_all

Example (the master surface’s solid block is refined mid-run while contact is already engaged): testfiles/marmot/AMR_ContactRefinePatch/test.inp

** AMR + deformable-surface contact integration test (PLAN_AMR_CONTACT.md §5/§6 test 4): two
** GC3D20/C3D20 blocks in node-to-deformable-surface penalty contact; the lower (master) block is
** refined mid-run by hAdaptivity while contact is already engaged. Exercises the reconcile() path
** in nodeToDeformableSurfacePenalty: the master facets (built from surface 'lower_front') must be
** regenerated from the refined child faces and the constraint's cached facet list/reference
** coordinates rebound, without losing contact -- before PLAN_AMR_CONTACT.md's P2, the facets would
** silently go stale (still referencing the coarse element's now-interior corner nodes).
*material, name=linearelastic, id=mat
18000, 0.0

*section, name=sec, material=mat, type=solid
lower_all
upper_all

*modelGenerator, generator=boxGen, name=lower
nX      =1
nY      =1
nZ      =1
x0      =0
y0      =0
z0      =0
lX      =1
lY      =1
lZ      =1
elType  =C3D20

*modelGenerator, generator=boxGen, name=upper
nX      =1
nY      =1
nZ      =1
x0      =0
y0      =0
z0      =1.05
lX      =1
lY      =1
lZ      =1
elType  =C3D20

*modelGenerator, generator=surfaceElementGenerator, name=gen1
surface = upper_back
name    = slaveSurf

*modelGenerator, generator=surfaceElementGenerator, name=gen2
surface = lower_front
name    = masterSurf

*constraint, name=contact, type=nodeToDeformableSurfacePenalty
slaveSurface=slaveSurf_facets, masterSurface=masterSurf_facets, penalty=1e6, type=linear,
searchDistance=2.0, sliding=small

*modelModifier, type=hAdaptivity, name=amr
>>marker, type=fieldOutput, fieldOutput=stressForAMR, operator='>', threshold=300.0
maxLevel=1
refineElSet=lower_all

*job, name=amrcontact, domain=3d
*solver, solver=NIST, name=theSolver
*fieldOutput
>>perNode, elSet=lower_all, field=displacement, result=U, name=dispLower
>>perNode, elSet=upper_all, field=displacement, result=U, name=dispUpper
>>fromExpression, name=pN, nSet=slaveSurf_nodes, expression='model.constraints["contact"].getNormalPressures()'
>>fromExpression, name=pNmean, nSet=slaveSurf_nodes, expression='model.constraints["contact"].getNormalPressures()', f(x)='np.mean(x)'
>>perElement, elSet=lower_all, result=stress, quadraturePoint=0:27, name=stressForAMR, f(x)='np.abs(x)'

*step, solver=theSolver
maxInc=0.02, minInc=1e-6, maxNumInc=100, maxIter=25, stepLength=1
>>options, name=theSolver, extrapolation=off
>>dirichlet, name=fixLower, nSet=lower_back, field=displacement, 1=0.0, 2=0.0, 3=0.0
>>dirichlet, name=pushUpper, nSet=upper_front, field=displacement, 1=0.0, 2=0.0, 3=-0.15

Example (the master surface’s solid block is refined mid-run while already tied): testfiles/marmot/AMR_TieRefine/test.inp

** AMR + tie-constraint integration test (PLAN_AMR_CONTACT.md §5/§6 test 5): two C3D20 blocks tied
** together via DOF-elimination tie; the lower (master) block is refined mid-run by hAdaptivity
** while already tied. Exercises tie.reconcile(): the master facets (built from surface
** 'lower_front') must be regenerated from the refined child faces and every slave node
** re-projected onto its new (finer) closest facet -- before PLAN_AMR_CONTACT.md's P3, the tied
** records would silently go stale (still pointing at the coarse element's now-interior nodes).
*material, name=linearelastic, id=mat
18000, 0.0

*section, name=sec, material=mat, type=solid
lower_all
upper_all

*modelGenerator, generator=boxGen, name=lower
nX      =1
nY      =1
nZ      =1
x0      =0
y0      =0
z0      =0
lX      =1
lY      =1
lZ      =1
elType  =C3D20

*modelGenerator, generator=boxGen, name=upper
nX      =1
nY      =1
nZ      =1
x0      =0
y0      =0
z0      =1.05
lX      =1
lY      =1
lZ      =1
elType  =C3D20

*modelGenerator, generator=surfaceElementGenerator, name=gen1
surface = upper_back
name    = slaveSurf

*modelGenerator, generator=surfaceElementGenerator, name=gen2
surface = lower_front
name    = masterSurf

*constraint, name=tie, type=tie
slaveSurface=slaveSurf_facets, masterSurface=masterSurf_facets, adjust=False

*modelModifier, type=hAdaptivity, name=amr
>>marker, type=fieldOutput, fieldOutput=stressForAMR, operator='>', threshold=300.0
maxLevel=1
refineElSet=lower_all

*job, name=amrtie, domain=3d
*solver, solver=NIST, name=theSolver
*fieldOutput
>>perNode, elSet=lower_all, field=displacement, result=U, name=dispLower
>>perNode, elSet=upper_all, field=displacement, result=U, name=dispUpper
>>perElement, elSet=lower_all, result=stress, quadraturePoint=0:27, name=stressForAMR, f(x)='np.abs(x)'

*step, solver=theSolver
maxInc=0.02, minInc=1e-6, maxNumInc=100, maxIter=25, stepLength=1
>>options, name=theSolver, extrapolation=off
>>dirichlet, name=fixLower, nSet=lower_back, field=displacement, 1=0.0, 2=0.0, 3=0.0
>>dirichlet, name=pushUpper, nSet=upper_front, field=displacement, 1=0.0, 2=0.0, 3=-0.15

The rigid-body contact constraints (nodetorigidsurfacepenalty, nodetodiscreterigidbodypenalty) are likewise MeshDependent, but lighter still: their master geometry is rigid (an analytic plane, or a triangulated rigid body), so refinement only ever grows their watched slave nSet – no facet regeneration, no per-slave history, just a refreshed node list at the next updateConnectivity() tick.

Example (the block is refined mid-run while its face already rests against an analytic rigid wall): testfiles/marmot/AMR_RigidContactRefine/test.inp

** AMR + rigid-contact integration test (PLAN_AMR_CONTACT.md §5/P4 test): a single C3D20 block sits
** against an analytic rigid wall at its own left face (nodeToRigidSurfacePenalty, component=0,
** direction=-1, value=0.0), watching nSet 'gen_left'. The block is refined mid-run by hAdaptivity
** while the right face is being pushed left, compressing the block against the wall. Exercises
** Constraint.reconcile(): the refined left face gains new boundary nodes (face-interior +
** edge-midpoints of the octree split) that must be picked up into the watched node set and
** protected from penetrating the wall -- before PLAN_AMR_CONTACT.md's P4, they would silently stay
** unprotected (never added to the frozen slaveNodes list from construction).
*material, name=linearelastic, id=mat
30000, 0.15

*section, name=sec, material=mat, type=solid
all

*modelGenerator, generator=boxGen, name=gen
nX      =1
nY      =1
nZ      =1
lX      =1
lY      =1
lZ      =1
elType  =C3D20

*constraint, type=nodetorigidsurfacepenalty, name=wall
nSet=gen_left
field=displacement
component=0
penalty=1e6
value=0.0
direction=-1
type=linear

*modelModifier, type=hAdaptivity, name=amr
>>marker, type=fieldOutput, fieldOutput=stressForAMR, operator='>', threshold=300.0
maxLevel=1

*job, name=amrrigid, domain=3d
*solver, solver=NIST, name=theSolver
*fieldOutput
>>perNode, elSet=all, field=displacement, result=U, name=U
>>perElement, elSet=all, result=stress, quadraturePoint=0:27, name=stressForAMR, f(x)='np.abs(x)'

*step, solver=theSolver
maxInc=0.02, minInc=1e-6, maxNumInc=100, maxIter=25, stepLength=1
>>options, name=theSolver, extrapolation=off
>>dirichlet, name=fixBottomLeftFront, nSet=gen_bottomLeftFront, field=displacement, 2=0.0, 3=0.0
>>dirichlet, name=fixBottomLeftBack, nSet=gen_bottomLeftBack, field=displacement, 2=0.0
>>dirichlet, name=pushRight, nSet=gen_right, field=displacement, 1=-0.1

Transparency for ordinary (AMR-unaware) consumers

A constraint, step action or output that only caches a node/element set or a surface – e.g. self._nodes = model.nodeSets["mySet"] at construction – needs none of the above: it sees a later refinement automatically, because that container is mutated in place rather than replaced (see the note at the top of this page). testfiles/marmot/AMR_TransparencyProbe is the acceptance test for this: a minimal ConstraintBase (edelweissfe.constraints. amrtransparencyprobe) that does exactly this, registers no observer and implements no MeshDependent, yet raises if its cached node set is ever found not to have grown across a refinement it should have seen – guarding against a regression that reintroduces replacing a set instead of mutating it.

`surfaceSnap` - Analytical-surface node snapping for AMR-refined curved boundaries

Module edelweissfe.modelmodifiers.geometry.surfacesnap

Analytical-surface node snapping for AMR-refined curved boundaries (HEX20 only).

A coarse mesh represents a curved boundary (e.g. a borehole) as a polygon of flat facets. When AMR h-refines elements attached to that boundary, the new nodes it creates do NOT improve geometric fidelity toward the true curve – Hex20Topology.subdivide evaluates the PARENT element’s own (possibly already-flat) isoparametric map, so a straight parent edge only ever produces more straight sub-edges. This purely reactive model modifier reacts to an AMR modifier’s ModelChange and snaps exactly the newly created nodes on a tracked boundary onto a specified analytical surface (currently: a cylinder), leaving every pre-existing node untouched.

See the surfaceSnap section of Model modifiers for the full narrative (corner/midside modes, hanging-node and quality-safeguard behaviour, restart safety) and a runnable example.

Restart safety. Which element faces are “on the tracked surface”, and which candidate nodes are still waiting to be snapped, is decision-side state that must survive a checkpoint/resume exactly like hAdaptivity’s own tracked node sets do: by living in two genuine NodeSet instances (ModelModifier._wallSetName and ModelModifier._pendingSetName, both private to this modifier instance) that a restart rebuilds by REPLAYING this modifier’s own recorded SnapPlan history – the same mechanism, not a parallel one. This is why ModelModifier.apply() grows/shrinks those node sets independent of whether any node’s coordinates actually moved this round: a round can legitimately snap nothing (every candidate is a hanging-node collision, or the quality safeguard vetoed all of them) while still discovering real new wall-face membership – or resolving an earlier round’s pending one – that matters later. Touching either tracked set marks ModelChange non-empty (via changedNodeSets) even when movedNodes is empty, so that round is always recorded and always replays – no state is ever silently lost across a restart.

Retrying what a round could not snap. A candidate that could not be snapped this round – because it was an AMR hanging-node slave, or because snapping it (or a sibling on the same round) would have dropped an affected element’s quality below qualityDropThreshold – is not abandoned: it stays in the pending set and is reconsidered on every later round that changes anything, since either condition can resolve (a hanging slave stops being one once its coarse neighbour is itself refined; a quality veto can clear once local geometry changes).

A coarse HEX20 mesh represents a curved boundary – a borehole, a cylindrical bore, a fillet – as a polygon of flat facets. Hex20Topology.subdivide() (used by hAdaptivity, above) evaluates the parent element’s own isoparametric map when it creates new nodes, so a flat parent facet only ever produces more flat sub-facets: refinement makes the polygon finer-grained, never rounder. surfaceSnap is a purely reactive model modifier (initiatesTopologyChanges = False) that watches hAdaptivity’s own ModelChange and, for every newly-created node that lands on a tracked boundary, projects it exactly onto a specified analytical surface – currently a cylinder (originX/Y/Z + axisX/Y/Z + radius). Every pre-existing node is left untouched: the nodeSet option names a node set whose members must already lie exactly on the true surface (e.g. a borehole wall exported from an actual cylindrical CAD cut), and only faces entirely within that set – and their new AMR-created nodes – are ever moved.

Corners vs. midsides. A new corner node is always radially projected onto the cylinder, unconditionally. A new midside node has two modes, via midsideNodes:

straight (default) – recomputed as the mean of its two (already-snapped) corner endpoints, keeping the edge a straight chord. Faceted, but with correctly-placed vertices; cheap and safe.
curved – independently projected onto the cylinder too, giving a true curved quadratic edge. Better geometric fidelity per element, at the cost of a small extra deviation from a pure straight chord representation elsewhere in the model (see the theory note in the module docstring for the chord-vs-arc argument).

Hanging nodes and mesh quality are never compromised – and a skipped node is retried, not abandoned. A candidate node that is also an AMR hanging-node slave (Constraint) is skipped rather than snapped – moving it would break its own multi-point constraint – and a journal warning names it. Before committing a round, every affected element’s minimum corner Jacobian determinant is checked before and after; if any element would invert or drop below qualityDropThreshold (default 0.5, must be in (0, 1]) of its pre-snap value, the entire round’s snapping is skipped (nothing moves). Either way the skipped node stays in a modifier-private pending set and is reconsidered on every later round that changes anything: a hanging slave stops being one once its coarse neighbour is itself refined, and a quality veto can clear once local geometry changes.

Keeping AMR’s own mirror in sync. hAdaptivity keeps a private AdaptiveMesh mirror whose node coordinates are cached at the moment each element is created, and never re-read from the live model afterward. Without correcting that cache too, a LATER refinement of an element surfaceSnap already snapped would silently subdivide from the stale, pre-snap geometry. surfaceSnap calls HAdaptivityModelModifier.syncNodeCoordinates on every hAdaptivity instance in the model right after writing model.nodes, closing this gap.

Options:
Option	Type / Default	Description
`moduleOptions`	dict, default={}	Internal
`nodeSet`	str (required)	Node set naming the boundary to snap onto the analytical surface, e.g. a borehole wall exported by the mesh generator (its nodes must already lie exactly on the true surface – only nodes AMR creates later on faces fully within this set are ever moved).
`originX`	float, default=0.0	X coordinate of a point on the cylinder axis.
`originY`	float, default=0.0	Y coordinate of a point on the cylinder axis.
`originZ`	float, default=0.0	Z coordinate of a point on the cylinder axis.
`axisX`	float, default=0.0	X component of the cylinder axis direction.
`axisY`	float, default=1.0	Y component of the cylinder axis direction.
`axisZ`	float, default=0.0	Z component of the cylinder axis direction.
`radius`	float (required)	Cylinder radius.
`midsideNodes`	str, default=’straight’	‘straight’ (default): recompute a new midside node as the mean of its two (already-snapped) corner endpoints, keeping the edge a straight chord. ‘curved’: independently project the midside node onto the cylinder too, giving a true curved quadratic edge.
`qualityDropThreshold`	float, default=0.5	An affected element’s minimum corner Jacobian determinant must stay above this fraction of its pre-snap value (and stay positive), or the ENTIRE round’s snap for that element’s face is skipped this round (nothing is moved, but the boundary-face bookkeeping still advances, and the skipped nodes are retried on a later round). Must be in (0, 1].

Example (a coarse, geometrically exact-on-the-cylinder inner wall; refining it and snapping the new boundary nodes with midsideNodes=curved): testfiles/edelweiss-only/AMR_SurfaceSnap/test.inp

** Analytical-surface node snapping regression test -- exercises the full surfaceSnap
** *modelModifier path end to end alongside hAdaptivity: a coarse HEX20 sector mesh models a 90
** degree, 3-sector arc of a hollow cylindrical shell (inner radius 10, outer radius 16, height 4,
** y-axis). The inner ('wall') boundary is built geometrically EXACT on the true cylinder -- as a
** real Cubit export of an actual cylindrical surface would be -- but every other facet, and every
** node hAdaptivity's own refinement subsequently creates, is flat/faceted by construction
** (Hex20Topology.subdivide only ever interpolates the parent's own isoparametric map). Without
** surfaceSnap, refining an element touching 'wall' therefore leaves its new boundary nodes
** measurably off the true cylinder; with it, every new corner and every new midside
** (midsideNodes=curved here) is exactly reprojected. Pins the resulting elastic solution against
** U.ref.
*material, name=linearelastic, id=mat
30000.0, 0.2

*section, name=sec, material=mat, type=solid
ring

*modelModifier, type=hAdaptivity, name=amr
>>marker, type=nodeSet, nSet=wall, initialOnly=True
maxLevel=1
splitFactor=2

*modelModifier, type=surfaceSnap, name=snap
nodeSet=wall
originX=0.0, originY=0.0, originZ=0.0
axisX=0.0, axisY=1.0, axisZ=0.0
radius=10.0
midsideNodes=curved

*job, name=amrSurfaceSnap, domain=3d
*solver, solver=NIST, name=S
*fieldOutput
>>perNode, elSet=ring, field=displacement, result=U, name=U

*step, solver=S
maxInc=1.0, minInc=1.0, maxNumInc=1, maxIter=25, stepLength=1
>>dirichlet, name=fixBase, nSet=base, field=displacement, 1=0.0, 2=0.0, 3=0.0
>>dirichlet, name=stretchTop, nSet=top, field=displacement, 2=0.01

*Node
1, 10.0000000000, 0.0000000000, 0.0000000000
2, 9.6592582629, 0.0000000000, 2.5881904510
3, 8.6602540378, 0.0000000000, 5.0000000000
4, 7.0710678119, 0.0000000000, 7.0710678119
5, 5.0000000000, 0.0000000000, 8.6602540378
6, 2.5881904510, 0.0000000000, 9.6592582629
7, 0.0000000000, 0.0000000000, 10.0000000000
8, 10.0000000000, 2.0000000000, 0.0000000000
9, 8.6602540378, 2.0000000000, 5.0000000000
10, 5.0000000000, 2.0000000000, 8.6602540378
11, 0.0000000000, 2.0000000000, 10.0000000000
12, 10.0000000000, 4.0000000000, 0.0000000000
13, 9.6592582629, 4.0000000000, 2.5881904510
14, 8.6602540378, 4.0000000000, 5.0000000000
15, 7.0710678119, 4.0000000000, 7.0710678119
16, 5.0000000000, 4.0000000000, 8.6602540378
17, 2.5881904510, 4.0000000000, 9.6592582629
18, 0.0000000000, 4.0000000000, 10.0000000000
19, 13.0000000000, 0.0000000000, 0.0000000000
20, 11.2583302492, 0.0000000000, 6.5000000000
21, 6.5000000000, 0.0000000000, 11.2583302492
22, 0.0000000000, 0.0000000000, 13.0000000000
23, 13.0000000000, 4.0000000000, 0.0000000000
24, 11.2583302492, 4.0000000000, 6.5000000000
25, 6.5000000000, 4.0000000000, 11.2583302492
26, 0.0000000000, 4.0000000000, 13.0000000000
27, 16.0000000000, 0.0000000000, 0.0000000000
28, 14.9282032303, 0.0000000000, 4.0000000000
29, 13.8564064606, 0.0000000000, 8.0000000000
30, 10.9282032303, 0.0000000000, 10.9282032303
31, 8.0000000000, 0.0000000000, 13.8564064606
32, 4.0000000000, 0.0000000000, 14.9282032303
33, 0.0000000000, 0.0000000000, 16.0000000000
34, 16.0000000000, 2.0000000000, 0.0000000000
35, 13.8564064606, 2.0000000000, 8.0000000000
36, 8.0000000000, 2.0000000000, 13.8564064606
37, 0.0000000000, 2.0000000000, 16.0000000000
38, 16.0000000000, 4.0000000000, 0.0000000000
39, 14.9282032303, 4.0000000000, 4.0000000000
40, 13.8564064606, 4.0000000000, 8.0000000000
41, 10.9282032303, 4.0000000000, 10.9282032303
42, 8.0000000000, 4.0000000000, 13.8564064606
43, 4.0000000000, 4.0000000000, 14.9282032303
44, 0.0000000000, 4.0000000000, 16.0000000000
*Element, type=C3D20
1, 14, 3, 1, 12, 40, 29, 27, 38, 9, 2, 8, 13, 35, 28, 34, 39, 24, 20, 19, 23
2, 16, 5, 3, 14, 42, 31, 29, 40, 10, 4, 9, 15, 36, 30, 35, 41, 25, 21, 20, 24
3, 18, 7, 5, 16, 44, 33, 31, 42, 11, 6, 10, 17, 37, 32, 36, 43, 26, 22, 21, 25
*elSet, elSet=ring, generate
1, 3, 1
*nSet, nSet=wall
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
*nSet, nSet=base
1
2
3
4
5
6
7
19
20
21
22
27
28
29
30
31
32
33
*nSet, nSet=top
12
13
14
15
16
17
18
23
24
25
26
38
39
40
41
42
43
44

Restart safety. Which element faces lie on the tracked surface, and which candidates are still pending a retry, is decision-side state that must survive a checkpoint/resume exactly like hAdaptivity’s own tracked node sets do. Rather than private Python state, both live in genuine, modifier-private NodeSet instances that a restart rebuilds by replaying this modifier’s own recorded decisions – the same mechanism described below for hAdaptivity, not a parallel one. Because a round can legitimately discover new boundary-face membership (or resolve an earlier round’s pending node) while moving zero node coordinates, ModelChange carries a dedicated movedNodes field and ModelChangeType a GEOMETRY_CHANGE kind, so that a coordinate-only (or even a zero-coordinate, node-set-only) mutation is never mistaken for an empty, unrecorded one.

Restart / checkpointing

A model modifier that mutates topology, like hAdaptivity, cannot rely on the plain reconstruct-then-overwrite scheme every other checkpointed component uses (see *restart): a refined mesh’s new elements and nodes aren’t in the .inp file to rebuild from. Instead, the checkpoint records only the decisions that drove each past change – and the resumed run replays them through the modifier’s own apply(), the very same code the live run executed. The marker evaluation that produced a decision is never re-run.

A modifier does not implement its own restart. It implements encodePlan() and decodePlan() so that its decision survives a checkpoint; FEModel records every applied decision in history and replays it. An earlier design had each modifier serializing its own history and implementing its own replay, which is precisely how a resumed run came to rebuild a differently-numbered mesh – two implementations of one mutation always drift. See The topology pipeline: how a model changes itself.

Implementing your own model modifiers

Subclass from the model-modifier base class in module edelweissfe.modelmodifiers.base.modelmodifierbase

Abstract base class for all ModelModifier entities in EdelweissFE.

class edelweissfe.modelmodifiers.base.modelmodifierbase.ModelModifierBase(name, model, journal, **kwargs)[source]

Abstract base class for entities that dynamically mutate the FEModel topology, mesh, or state variables during analysis steps.

Attributes:

actsOnlyAtSimulationStart: Whether every topology change this modifier can plan happens on the very first topology update of the analysis.
checkpointedState
name: Name of the model modifier.
schema

Parameters:

name (str)
model (FEModel)
journal (Journal)

Methods

`apply`(model, plan)	Carry out `plan`, mutating the model.
`declaredDomain`(model)	The element numbers this modifier claims exclusive authority over.
`decodePlan`(data)	Inverse of `encodePlan()`.
`encodePlan`(plan)	Serialize a plan for the topology history, as a flat mapping of name to `numpy.ndarray` (or scalar).
`getRestartData`()	The declared state; see `checkpointedState`.
`plan`(model, change, step)	Decide what, if anything, this modifier wants to change about the model -- without changing it.
`setRestartData`(data)	Restore the state `getRestartData()` returned.

abstractmethod apply(model, plan)[source]

Carry out plan, mutating the model.

Must not read solution state. A pure function of (model, plan) – everything the decision depended on belongs in the plan. This is the single code path shared by a live run and a restart replay, and the reason a replayed mesh is numbered identically: there is no second implementation to drift apart from this one.

Runs inside an open topology window (see changes()), so it may create and delete elements – through reserveElementNumbers() and createElement(), never by writing model.elements directly.

Parameters:

model (FEModel) – The FEModel object, mutated in place.
plan – A plan previously returned by plan() (live) or restored from the recorded topology history (replay).

Returns:

The changeset this mutation produced.

Return type:

ModelChange

declaredDomain(model)[source]

The element numbers this modifier claims exclusive authority over.

Two modifiers that both own an element will fight over it: each mutates model.elements directly, so one deleting or subdividing an element the other still tracks leaves the second holding a stale reference, which later corrupts element-set membership and can surface as a node that is simultaneously Dirichlet-prescribed and a multi-point-constraint slave.

checkModelModifierDomains() compares these pairwise once, at the end of setup, and refuses the model rather than letting the conflict appear deep in the solve loop. The default claims nothing, which is correct for a modifier that only ever adds entities.

Parameters:: model (FEModel)
Return type:: set

abstractmethod decodePlan(data)[source]

Inverse of encodePlan().

Parameters:: data (dict)

abstractmethod encodePlan(plan)[source]

Serialize a plan for the topology history, as a flat mapping of name to numpy.ndarray (or scalar).

Everything apply() needs must be here: on a restart, the plan is rebuilt from this and nothing else – plan() is not called.

Return type:: dict

getRestartData()[source]

The declared state; see checkpointedState.

Returns:: A flat mapping of array name to array.
Return type:: dict[str, numpy.ndarray]

abstractmethod plan(model, change, step)[source]

Decide what, if anything, this modifier wants to change about the model – without changing it.

This half is allowed to read solution state: markers, node fields, the current time, the step. It must return a serializable description of the decision, which apply() can carry out on its own, or None if there is nothing to do.

A restart replay never calls this. It reads the recorded plan and calls apply() directly, which is what makes a resumed run reproduce the original topology exactly instead of re-deciding from state that may not be identical.

Parameters:

model (FEModel) – The FEModel object. Read only – do not mutate it here.
change – The net ModelChange since this modifier last planned within the current topology update, or None on the first round. Return ``None`` when it does not touch this modifier’s domain – that is what lets the pipeline reach a fixed point instead of looping (see update()).
step – The current step.

Returns:

A plan for apply(), or None.

Return type:

object | None

setRestartData(data)[source]

Restore the state getRestartData() returned.

Parameters:: data (dict[str, ndarray]) – The mapping of arrays.

property actsOnlyAtSimulationStart: bool

Whether every topology change this modifier can plan happens on the very first topology update of the analysis.

True means a solver that runs the topology update exactly once, before its increment loop, is enough for this modifier; False means the update has to be re-run as the solution evolves. The distinction exists because those are genuinely different capabilities: re-running the update mid-analysis invalidates everything derived from the mesh, which an implicit solver rebuilds from its node fields but an explicit one cannot, its velocity state being solver-local. See validateModelCapabilities().

A purely reactive modifier (initiatesTopologyChanges == False) cannot act on its own at all, so it never needs a later update and answers True.

Returns:: Whether this modifier is fully served by a single topology update at the start.
Return type:: bool

checkpointedState: dict | None = None: The state this modifier carries from one decision to the next, by attribute name and type (see checkpointedstate); not the mesh, which a restart replays through apply(). Every modifier declares it, an empty mapping if it carries none.

initiatesTopologyChanges = True: Whether this modifier can start a topology change by itself. A purely reactive modifier sets this to False, so it does not disqualify a solver without topology updates (see NonlinearSolverBase.validateModelCapabilities).

property name: str: Name of the model modifier.

Model modifiers

`hAdaptivity` - Hanging-node h-adaptivity for HEX20

Batching refinement across increments

State-variable transfer strategies

Predictor after a refinement

Re-equilibration after a refinement

Refinement under a dynamic solver

Compatibility with facet-based contact and tie

Transparency for ordinary (AMR-unaware) consumers

`surfaceSnap` - Analytical-surface node snapping for AMR-refined curved boundaries

Restart / checkpointing

Implementing your own model modifiers

Model modifiers

hAdaptivity - Hanging-node h-adaptivity for HEX20

Batching refinement across increments

State-variable transfer strategies

Predictor after a refinement

Re-equilibration after a refinement

Refinement under a dynamic solver

Compatibility with facet-based contact and tie

Transparency for ordinary (AMR-unaware) consumers

surfaceSnap - Analytical-surface node snapping for AMR-refined curved boundaries

Restart / checkpointing

Implementing your own model modifiers

`hAdaptivity` - Hanging-node h-adaptivity for HEX20

`surfaceSnap` - Analytical-surface node snapping for AMR-refined curved boundaries