#!/usr/bin/env python3
# -*- coding: utf-8 -*-
# ---------------------------------------------------------------------
#
# _____ _ _ _ _____ _____
# | ____|__| | ___| |_ _____(_)___ ___| ___| ____|
# | _| / _` |/ _ \ \ \ /\ / / _ \ / __/ __| |_ | _|
# | |__| (_| | __/ |\ V V / __/ \__ \__ \ _| | |___
# |_____\__,_|\___|_| \_/\_/ \___|_|___/___/_| |_____|
#
#
# Unit of Strength of Materials and Structural Analysis
# University of Innsbruck,
# 2017 - today
#
# Matthias Neuner matthias.neuner@uibk.ac.at
#
# This file is part of EdelweissFE.
#
# This library is free software; you can redistribute it and/or
# modify it under the terms of the GNU Lesser General Public
# License as published by the Free Software Foundation; either
# version 2.1 of the License, or (at your option) any later version.
#
# The full text of the license can be found in the file LICENSE.md at
# the top level directory of EdelweissFE.
# ---------------------------------------------------------------------
# Created on Fri Jan 27 19:53:45 2017
# @author: Matthias Neuner
import itertools
import textwrap
from operator import attrgetter
import h5py
import numpy as np
from edelweissfe.config.phenomena import getFieldSize, phenomena
from edelweissfe.fields.nodefield import NodeField
from edelweissfe.journal.journal import Journal
from edelweissfe.models.equationsystemlayout import EquationSystemLayout
from edelweissfe.models.modelchange import TopologyRecord
from edelweissfe.models.topologypipeline import TopologyPipeline
from edelweissfe.numerics.dofmanager import DofManager
from edelweissfe.numerics.parallelizationutilities import (
getNumberOfWorkerThreads,
mapChunksInParallel,
)
from edelweissfe.utils.checkpoint import readRestartDataInto, writeRestartDataOf
from edelweissfe.utils.exceptions import (
MeshChangedSinceCheckpoint,
RestartError,
TopologyError,
)
from edelweissfe.variables.fieldvariable import FieldVariable
from edelweissfe.variables.scalarvariable import ScalarVariable
[docs]
class FEModel:
"""This is is a standard finite element model tree.
It takes care of the correct number of variables,
for nodes and scalar degrees of freedem, and it manages the fields.
Parameters
----------
dimension
The dimension of the model.
"""
identification = "FEModel"
def __init__(self, dimension: int):
self.time = 0.0 #: Current time of the model.
self.nodes = {} #: Nodes in the model.
self.elements = {} #: Elements in the model.
self.nodeSets = {} #: NodeSets in the model.
self.nodeFields = {} #: NodeFields in the model.
self.elementSets = {} #: ElementSets in the model.
self.sections = {} #: Sections in the model.
self.surfaces = {} #: Surface definitions in the model.
self.constraints = {} #: Constraints in the model.
self.constraintSets = {} #: ConstraintsSets in the model.
self.multiPointConstraints = {} #: Multi-point (DOF-elimination) constraints in the model.
self.modelModifiers = {} #: Model modifiers (dynamic topology / mesh mutation entities) in the model.
self.contactFacetRecipes = {} #: facet elSet name -> (surfaceName, prefix, triangulation, nodalWeights).
self.materials = {} #: Materials in the model.
self.analyticalFields = {} #: AnalyticalFields in the model.
self.scalarVariables = {} #: ScalarVariables in the model.
self.additionalParameters = {} #: Additional information.
self.rigidBodies = {} #: RigidBodies in the model.
self.elementProperties = [] #: Element properties.
self.domainSize = dimension #: Spatial dimension of the model
self.fieldOutputController = None #: Set once by the driver; lets in-model entities (e.g. AMR markers) look up a named *fieldOutput by value, not just by declaration.
#: How the mesh may change during a run; see :class:`~edelweissfe.models.topologypipeline.TopologyPipeline`.
self.topology = TopologyPipeline(self)
def __copy__(self):
"""A shallow copy of the model, with a topology pipeline that acts on the copy.
A shallow copy shares every attribute with the original, which would include the pipeline --
and a pipeline changes the model it was created for. The copy therefore gets its own pipeline
bound to it (see :meth:`~edelweissfe.models.topologypipeline.TopologyPipeline.boundTo`), which
is what a reduced model needs, e.g. EdelweissMeshfree's active sub-model.
"""
clone = object.__new__(type(self))
clone.__dict__.update(self.__dict__)
clone.topology = self.topology.boundTo(clone)
return clone
[docs]
def createNode(self, node):
"""Add a freshly created node to the model.
Parameters
----------
node
The node, already carrying a label obtained from :meth:`~edelweissfe.models.topologypipeline.TopologyPipeline.reserveNodeNumbers`.
"""
if not self.topology.isOpen:
raise TopologyError(
"node {:} was created outside a topology change: only model modifiers may create "
"or delete nodes, inside TopologyPipeline.changes()".format(node.label)
)
if node.label in self.nodes:
raise TopologyError(
"node label {:} is already taken -- node labels are reserved via "
"TopologyPipeline.reserveNodeNumbers() and never recycled".format(node.label)
)
self.nodes[node.label] = node
[docs]
def createElement(self, element):
"""Add a freshly created element to the model.
Parameters
----------
element
The element, already carrying a number obtained from :meth:`~edelweissfe.models.topologypipeline.TopologyPipeline.reserveElementNumbers`.
"""
if not self.topology.isOpen:
raise TopologyError(
"element {:} was created outside a topology change: only model modifiers may create "
"or delete elements, inside TopologyPipeline.changes()".format(element.elNumber)
)
if element.elNumber in self.elements:
raise TopologyError(
"element number {:} is already taken -- element numbers are reserved via "
"TopologyPipeline.reserveElementNumbers() and never recycled".format(element.elNumber)
)
self.elements[element.elNumber] = element
[docs]
def removeElement(self, elNumber: int):
"""Remove an element from the model. Its number is retired, never reissued.
Parameters
----------
elNumber
The number of the element to remove.
"""
if not self.topology.isOpen:
raise TopologyError(
"element {:} was deleted outside a topology change: only model modifiers may create "
"or delete elements, inside TopologyPipeline.changes()".format(elNumber)
)
del self.elements[elNumber]
def _populateNodeFieldVariablesFromElements(
self,
):
"""Creates FieldVariables on Nodes depending on the all
elements.
"""
for element in self.elements.values():
for node, elementNodeFields in zip(element.nodes, element.fields):
for field in elementNodeFields:
if field not in node.fields:
node.fields[field] = FieldVariable(node, field)
def _populateNodeFieldVariablesFromConstraints(
self,
):
"""Creates FieldVariables on Nodes depending on the all
constraints.
"""
for constraint in itertools.chain(self.constraints.values(), self.multiPointConstraints.values()):
for node, nodeFields in zip(constraint.nodes, constraint.fieldsOnNodes):
for field in nodeFields:
if field not in node.fields:
node.fields[field] = FieldVariable(node, field)
def _createNodeFieldsFromNodes(self, nodes: list, nodeSets: list) -> dict[str, NodeField]:
"""Bundle nodal FieldVariables together in contiguous NodeFields.
Parameters
----------
nodes
The list of Nodes from which the NodeFields should be created.
nodeSets
The list of NodeSets, which should be considered in the index map of the NodeFields.
Returns
-------
dict[str,NodeField]
The dictionary containing the NodeField instances for every active field."""
domainSize = self.domainSize
theNodeFields = dict()
for field in phenomena.keys():
theNodeField = NodeField(field, getFieldSize(field, domainSize), nodes)
if theNodeField.nodes:
theNodeFields[field] = theNodeField
return theNodeFields
def _linkFieldVariableObjects(self, nodes):
"""Link NodeFields to individual FieldVariable objects.
Parameters
----------
nodes
Nodes to be linked
Inside a deferring :meth:`~edelweissfe.models.topologypipeline.TopologyPipeline.changes` window this only notes the request; the link is
made once, for these nodes as they are then, when the window closes.
"""
if self.topology.isDeferringFieldBookkeeping:
self.topology.deferFieldVariableLink(nodes)
return
for node in nodes:
for field, fieldVariable in node.fields.items():
nodeField = self.nodeFields[field]
idx = nodeField._indicesOfNodesInArray[node]
fieldVariable.values = self.nodeFields[field]["U"][idx, :]
return
def _requestAdditionalScalarVariable(self, name: str):
"""Create a new scalar variables
Parameters
----------
name
The name of the variable.
Returns
-------
ScalarVariable
The instance of the variable.
"""
if name in self.scalarVariables:
raise Exception("ScalarVariable with name {:} already exists!".format(name))
self.scalarVariables[name] = ScalarVariable(name)
return self.scalarVariables[name]
def _createAndAssignScalarVariableForConstraints(self, journal: Journal):
"""Create ScalarVariables for constraints.
Parameters
----------
journal
The journal intance.
"""
# we may have additional scalar degrees of freedom, not associated with any node (e.g, lagrangian multipliers of constraints)
for constraintName, constraint in self.constraints.items():
nAdditionalScalarVariables = constraint.getNumberOfAdditionalNeededScalarVariables()
if nAdditionalScalarVariables > 0:
journal.message(
"Constraint {:} requests {:} additional ScalarVariables".format(
constraintName, nAdditionalScalarVariables
),
self.identification,
2,
)
scalarVariables = [
self._requestAdditionalScalarVariable("{:}_{:}".format(constraintName, i))
for i in range(nAdditionalScalarVariables)
]
constraint.assignAdditionalScalarVariables(scalarVariables)
def _prepareVariablesAndFields(self, journal):
"""Prepare all variables and fields for a simulation.
Parameters
----------
journal
The journal instance.
"""
journal.message(
"Activating fields on nodes from Elements and Constraints",
self.identification,
)
self._populateNodeFieldVariablesFromElements()
self._populateNodeFieldVariablesFromConstraints()
journal.message("Bundling fields on nodes to NodeFields", self.identification)
self.nodeFields = self._createNodeFieldsFromNodes(self.nodeSets["all"], self.nodeFields.values())
journal.message("Assembling ScalarVariables", self.identification)
self.scalarVariables = dict()
self._createAndAssignScalarVariableForConstraints(journal)
def _resizeNodeFieldsForNodes(self, journal: Journal):
"""Resize the existing NodeFields in place for the current ``self.nodeSets["all"]``,
instead of rebuilding :attr:`nodeFields` from scratch as :meth:`_prepareVariablesAndFields`
does. Used by mesh mutators (e.g. AMR's ``hadaptivity._materialize``) so that NodeField
identity -- and hence any :class:`~edelweissfe.fields.nodefield.NodeFieldSubset` or
reference a consumer cached -- survives a topology change. ScalarVariables (e.g. Lagrange
multipliers of constraints) are rebuilt, but values are preserved by name for constraints
that still exist, so their converged state survives the refinement too.
Everything here is recomputed from the current elements, constraints and ``nodeSets["all"]``
-- the node order of the resized fields is the model's node creation order, which is why one
call after several mutations yields the same layout as one call per mutation. Inside a
deferring :meth:`~edelweissfe.models.topologypipeline.TopologyPipeline.changes` window the call is therefore only noted, and made once when
the window closes.
Parameters
----------
journal
The journal instance.
"""
if self.topology.isDeferringFieldBookkeeping:
self.topology.deferNodeFieldResize(journal)
return
journal.message(
"Activating fields on nodes from Elements and Constraints",
self.identification,
)
self._populateNodeFieldVariablesFromElements()
self._populateNodeFieldVariablesFromConstraints()
journal.message("Resizing NodeFields", self.identification)
nodes = self.nodeSets["all"]
for nodeField in self.nodeFields.values():
nodeField.resize(nodes)
# a phenomenon activated for the first time (e.g. by a newly materialized constraint) has
# no NodeField yet -- create it exactly as _prepareVariablesAndFields would
for field in phenomena.keys():
if field not in self.nodeFields:
newNodeField = NodeField(field, getFieldSize(field, self.domainSize), nodes)
if newNodeField.nodes:
self.nodeFields[field] = newNodeField
journal.message("Assembling ScalarVariables", self.identification)
# Rebuilding cold-starts every value: restore the converged Lagrange multipliers of unchanged
# constraints by name, so that only new ones start at zero.
previousScalarVariableValues = {name: v.value for name, v in self.scalarVariables.items()}
self.scalarVariables = dict()
self._createAndAssignScalarVariableForConstraints(journal)
for name, v in self.scalarVariables.items():
if name in previousScalarVariableValues:
v.value = previousScalarVariableValues[name]
def _prepareElements(self, journal: Journal):
"""Prepare elements for a simulation.
In detail, sections are assigned.
Parameters
----------
journal
The journal instance.
"""
for section in self.sections.values():
section.assignSectionPropertiesToModel(self)
for elementProperty in self.elementProperties:
elementProperty.assignElementPropertiesToModel(self)
# check if all elements are assigned a material
materialAssigned = np.fromiter(map(attrgetter("hasMaterial"), self.elements.values()), dtype=bool)
if not materialAssigned.all():
elementIds = np.array([str(elId) for elId in self.elements.keys()])[np.logical_not(materialAssigned)]
raise Exception(f"No material was assigned to element(s) with id(s) {', '.join(elementIds)}.")
[docs]
def prepareYourself(self, journal: Journal):
"""Prepare the model for a simulation.
Creates the variables, bundles the fields,
and initializes elements.
Parameters
----------
journal
The journal instance.
"""
self.topology.adoptSetupElementNumbers()
self.topology.ensureSurfaceFacetModifier(journal)
self.topology.checkModelModifierDomains()
self._prepareVariablesAndFields(journal)
self._prepareElements(journal)
self.topology.recordSetupFingerprint()
[docs]
def activeConstraints(self) -> list:
"""The constraints that take part in the equation system: all constraints, in their order, that are active,
i.e. switched on and within reach (see
:attr:`~edelweissfe.constraints.base.constraintbase.ConstraintBase.active`).
Returns
-------
list[ConstraintBase]
The active constraints.
"""
return [constraint for constraint in self.constraints.values() if constraint.active]
[docs]
def updateConnectivity(self) -> bool:
"""Bring the coupling of the degrees of freedom up to date with the converged state, before an increment: which
nodes the entities of the equation system couple. Called by the solvers' topology update at the start of every
increment and of every retry after a cutback, before the increment is proposed (see
:meth:`~edelweissfe.solvers.base.nonlinearsolverbase.NonlinearSolverBase.updateTopologyAndConnectivity`). If it
reports a change, an implicit solver rebuilds its equation system (:meth:`buildEquationSystem`).
In a finite element model, the elements keep their nodes, and only a constraint may change the nodes it
couples: a contact constraint searches for the surface it touches. A meshfree model of EdelweissMeshfree moves
its material points through a background grid, and its particles carry their kernel functions along: it
reassigns them here first, then lets its constraints follow.
Only the constraints that are switched on are located anew
(:attr:`~edelweissfe.constraints.base.constraintbase.ConstraintBase.switchedOn`), those out of reach included,
so that they can find themselves within reach again. A constraint switched off from outside rests until it is
switched on again, and is located anew at the start of the next increment then.
A change of the coupling is not a change of the mesh: it neither moves nodes nor creates them, so a solver
rebuilds its equation system, but does not re-equilibrate (``equilibrateAfterModelChange``).
Returns
-------
bool
True if a coupling changed, so that the equation system has to be rebuilt.
"""
# A list, not a generator in any(): every constraint must update, not just those up to the first change.
return any(
[constraint.updateConnectivity(self) for constraint in self.constraints.values() if constraint.switchedOn]
)
[docs]
def buildEquationSystem(self, withSystemMatrix: bool = True) -> EquationSystemLayout:
"""Lay out the global equation system of the model as it is now: number its degrees of freedom, and say which
elements contribute to it, see :class:`~edelweissfe.models.equationsystemlayout.EquationSystemLayout`. Called
by a solver whenever it (re)builds its equation system: at the start of a step, and after the topology update
reported a change of the mesh or of the coupling (:meth:`updateConnectivity`).
The equation system of a finite element model holds the node fields, the elements, and the active constraints
with their scalar variables (e.g. Lagrange multipliers). A constraint that is not active (see
:attr:`~edelweissfe.constraints.base.constraintbase.ConstraintBase.active`) is left out together with its
scalar variables, which would otherwise be unknowns without an equation.
Parameters
----------
withSystemMatrix
False for an explicit solver, which assembles no system matrix: then neither the pattern of the system
matrix is built, nor the map from the degrees of freedom back to their nodes, by which the Newton iterations
report the nodes with the largest residuals.
Returns
-------
EquationSystemLayout
The layout of the equation system.
"""
constraints = self.activeConstraints()
dofManager = DofManager(
self.nodeFields.values(),
[scalarVariable for constraint in constraints for scalarVariable in constraint.scalarVariables],
self.elements.values(),
constraints,
self.nodeSets.values(),
initializeVIJPattern=withSystemMatrix,
determiningIndexToHostObjectMapping=withSystemMatrix,
)
return EquationSystemLayout(self, dofManager, list(self.elements.values()))
[docs]
def advanceToTime(self, time: float):
"""Accept the current state of the model and sub instances, and
set the new time.
Of the constraints, only the active ones accept their state (see :meth:`activeConstraints`): a constraint that
took no part in the increment has no new state, and keeps its history as it was, e.g. the multipliers and the
frictional history of a contact switched off for a step.
Parameters
----------
time
The new time.
"""
self.time = time
self._acceptElementStates()
# Serial: negligible next to the element loop.
for constraint in self.activeConstraints():
constraint.acceptLastState()
for mpc in self.multiPointConstraints.values():
mpc.acceptLastState()
[docs]
def advanceVelocity(self, timeIncrement: float):
"""Advance the velocity of the model by its acceleration over a time increment, :math:`\\boldsymbol{V}
\\leftarrow \\boldsymbol{V} + \\Delta t \\, \\boldsymbol{A}`, in every node field, on every node.
An explicit solver integrates the velocity staggered half an increment behind the solution within its step,
and publishes the acceleration at the end of every increment. At the end of its step, it advances the velocity
by the half increment it lags behind, so that the step leaves the velocity at its end time, as every solver
does (see the documentation page *Node field results*). Where the acceleration is zero, the velocity stays as
it is. The model does this itself, without an equation system: a step resumed at its very end has built none.
Parameters
----------
timeIncrement
The time increment.
"""
for nodeField in self.nodeFields.values():
velocity = nodeField["V"]
velocity += timeIncrement * nodeField["A"]
def _acceptElementStates(self):
"""Let every element accept its computed state, across the available threads.
An element's ``acceptLastState`` touches only that element's own state buffers -- for a
Marmot element it is the single ``self._stateVars[:] = self._stateVarsTemp`` copy -- so the
elements are independent of one another and this loop parallelizes without any coordination.
It is worth parallelizing because of how the explicit solver uses it, not because of how the
implicit ones do. An implicit analysis calls this once per *converged* increment, amortised
over a Newton loop and a linear solve, where 19 ms disappears into the noise. Explicit
dynamics calls it on every one of millions of increments: on the anchor pry-out model it was
15.6 % of the whole step, the single largest cost outside the element kernels themselves, and
every millisecond of it was serial Python holding 31 of 32 threads idle.
"""
elements = self.elements
def acceptChunk(chunk):
for element in chunk:
element.acceptLastState()
if len(elements) < getNumberOfWorkerThreads(): # too few to be worth handing to the pool
acceptChunk(elements.values())
else:
mapChunksInParallel(elements.values(), acceptChunk)
[docs]
def writeRestart(self, restartFile: h5py.File):
"""Write the current (converged) state of the model to a restart checkpoint.
Does not serialize model topology (nodes/elements/sets/sections/materials) -- only state
that a rebuild from the original ``.inp`` file cannot reproduce: node field values, element
quadrature-point history, scalar variables (e.g. Lagrange multipliers), stateful
constraints' internal history (e.g. frictional contact), and the additional parameters
(e.g. the load parameter of the arc length method).
Parameters
----------
restartFile
An open, writable :class:`h5py.File` (or group) to write the checkpoint into.
"""
f = restartFile
f.attrs["time"] = self.time
nodeFieldsGroup = f.create_group("nodeFields")
for nf in self.nodeFields.values():
nodeFieldGroup = nodeFieldsGroup.create_group(nf.name)
for entryName in NodeField.ENTRIES:
nodeFieldGroup.create_dataset(entryName, data=nf[entryName])
scalarVariablesGroup = f.create_group("scalarVariables")
for name, scalarVariable in self.scalarVariables.items():
scalarVariablesGroup.attrs[name] = scalarVariable.value
additionalParametersGroup = f.create_group("additionalParameters")
for name, value in self.additionalParameters.items():
additionalParametersGroup.attrs[name] = value
writeRestartDataOf(f.create_group("modelModifiers"), self.modelModifiers)
elementsGroup = f.create_group("elements")
for elNumber, element in self.elements.items():
elementsGroup.create_dataset(str(elNumber), data=element.getStateVars())
writeRestartDataOf(f.create_group("constraints"), self.constraints)
# The ordered record of every applied model-modifier decision (e.g. AMR refinements). A resumed run
# replays these (see readRestart), so that a single code path mutates topology, live or replayed.
historyGroup = f.create_group("topologyHistory")
historyGroup.attrs["setupFingerprint"] = self.topology.setupFingerprint()
historyGroup.attrs["count"] = len(self.topology.history)
for index, record in enumerate(self.topology.history):
recordGroup = historyGroup.create_group("{:06d}".format(index))
recordGroup.attrs["modifier"] = record.modifier
recordGroup.attrs["roundNumber"] = record.roundNumber
recordGroup.attrs["time"] = record.time
recordGroup.attrs["fingerprint"] = record.fingerprint
for entryName, entryValues in record.plan.items():
recordGroup.create_dataset(entryName, data=entryValues)
[docs]
def readRestart(self, restartFile: h5py.File, journal: Journal = None):
"""Read the state of the model from a restart checkpoint written by :meth:`writeRestart`.
Either the model has been rebuilt from the original ``.inp`` file and :meth:`prepareYourself`
called -- a resumed run -- or it is the very model the checkpoint was written from, going back
to it within the run (see :meth:`~edelweissfe.steps.base.stepbase.StepBase.fallBackToCheckpoint`).
Of the recorded topology decisions, those the model has not applied yet are replayed: all of
them for a rebuilt model, none for a model going back. A model whose mesh has changed since
the checkpoint cannot go back to it, since a mesh change cannot be undone.
Parameters
----------
restartFile
An open, readable :class:`h5py.File` (or group) to read the checkpoint from.
journal
Optional Journal for the progress messages of the topology replay.
Raises
------
RestartError
If the model is not the one the checkpoint was written from, before anything is restored, or
if the checkpoint does not hold every entry of a node field (see
:attr:`~edelweissfe.fields.nodefield.NodeField.ENTRIES`).
MeshChangedSinceCheckpoint
If the model's mesh has changed since the checkpoint was written, before anything is restored.
"""
f = restartFile
# Checked before anything is restored, so that a refused checkpoint leaves the model as it is.
historyGroup = f["topologyHistory"]
if self.topology.setupFingerprint() != historyGroup.attrs["setupFingerprint"]:
raise RestartError(
"the model rebuilt from the input file is not the one the checkpoint was written from: "
"its mesh (nodes, elements, coordinates, numbering) differs"
)
recordGroups = [historyGroup["{:06d}".format(index)] for index in range(int(historyGroup.attrs["count"]))]
appliedRecords = self.topology.history
if any(
(applied.modifier, applied.roundNumber, applied.fingerprint)
!= (str(group.attrs["modifier"]), int(group.attrs["roundNumber"]), str(group.attrs["fingerprint"]))
for applied, group in zip(appliedRecords, recordGroups)
):
raise RestartError("the topology decisions of the model differ from those the checkpoint records")
if len(appliedRecords) > len(recordGroups):
raise MeshChangedSinceCheckpoint(appliedRecords[len(recordGroups)])
for fieldName, storedField in f["nodeFields"].items():
if set(storedField.keys()) != set(NodeField.ENTRIES):
raise RestartError(
"the checkpoint holds the entries {:} of node field {:}, not {:}".format(
", ".join(sorted(storedField.keys())), fieldName, ", ".join(NodeField.ENTRIES)
)
)
self.time = f.attrs["time"]
# Must run before every restore below: the replay materializes elements/nodes (e.g. AMR-refined
# children) which those restores address by label.
records = [
TopologyRecord(
modifier=str(recordGroup.attrs["modifier"]),
roundNumber=int(recordGroup.attrs["roundNumber"]),
time=float(recordGroup.attrs["time"]),
plan={entryName: values[:] for entryName, values in recordGroup.items()},
fingerprint=str(recordGroup.attrs["fingerprint"]),
)
for recordGroup in recordGroups[len(appliedRecords) :]
]
self.topology.replayHistory(records, journal)
readRestartDataInto(f["modelModifiers"], self.modelModifiers)
# Bring the mesh-dependent consumers (ties, contact surfaces) up to date with the replayed mesh, as
# they were when the checkpoint was written.
self.topology.refreshMeshDependents()
for nf in self.nodeFields.values():
storedField = f["nodeFields"].get(nf.name)
if storedField is None:
raise RestartError("the checkpoint holds no node field {:}".format(nf.name))
for entryName in NodeField.ENTRIES:
storedField[entryName].read_direct(nf[entryName])
for name, scalarVariable in self.scalarVariables.items():
scalarVariable.value = f["scalarVariables"].attrs[name]
self.additionalParameters = dict(f["additionalParameters"].attrs)
# Sound only because the replay reproduces the original numbering exactly (verified by fingerprint):
# a missing element means a mismatched model and must be reported.
for elementKey, stateVars in f["elements"].items():
elNumber = int(elementKey)
element = self.elements.get(elNumber)
if element is None:
raise RestartError(
"the checkpoint holds state for element {:}, which does not exist after the "
"topology replay -- the replayed model does not match the one checkpointed".format(elNumber)
)
element.setStateVars(stateVars[:])
# A rigid body moves its surface nodes into the current configuration, a pure function of
# its reference point's restored DOFs: put them where the uninterrupted run had them.
for rigidBody in self.rigidBodies.values():
rigidBody.updateKinematics()
readRestartDataInto(f["constraints"], self.constraints)
def printPrettyModelSummary(model: FEModel, journal: Journal):
identification = "PrettyModelSummary"
def wrapList(theList):
for line in textwrap.wrap("[" + ", ".join(theList) + "]"):
journal.message(
" {:<20} ".format(line),
identification,
0,
)
journal.message(
"Finite element model with spatial dimension {:} has".format(model.domainSize),
identification,
0,
)
journal.message(
" {:<20}{:<15} ".format("nodes:", len(model.nodes)),
identification,
0,
)
journal.message(
" {:<20}{:<15} ".format("node sets:", len(model.nodeSets)),
identification,
0,
)
wrapList(model.nodeSets.keys())
journal.message(
" {:<20}{:<15} ".format("node fields:", len(model.nodeFields)),
identification,
0,
)
wrapList(["{:} ({:})".format(k, len(v.nodes)) for k, v in model.nodeFields.items()])
journal.message(
" {:<20}{:<15}".format("elements: ", len(model.elements)),
identification,
0,
)
journal.message(
" {:<20}{:<15}".format("element sets: ", len(model.elementSets)),
identification,
0,
)
wrapList(model.elementSets.keys())
if model.constraints:
journal.message(
" {:<20}{:<15}".format("constraints: ", len(model.constraints)),
identification,
0,
)
if model.multiPointConstraints:
journal.message(
" {:<20}{:<15}".format("multi-point constraints: ", len(model.multiPointConstraints)),
identification,
0,
)
if model.scalarVariables:
journal.message(
" {:<20}{:<15}".format("scalar variables: ", len(model.scalarVariables)),
identification,
0,
)
journal.message(
" {:<20}{:<15}".format("materials: ", len(model.materials)),
identification,
0,
)
wrapList(model.materials.keys())
if model.analyticalFields:
journal.message(
" {:<20}{:<15}".format("analytical fields: ", len(model.analyticalFields)),
identification,
0,
)
wrapList(model.analyticalFields.keys())