Utilities

EdelweissFE contains several classes and modules, which are not directly employed by the user, but very handy if fundamental new functionalities should be added.

Matrix conversion from COO sparse format to CSR format

Matrix assembly is a critical performance bottleneck in Finite Element simulations, particularly during non-linear implicit or explicit iterations where the system matrix structure (non-zero pattern) remains constant, but its numerical values are updated at every iteration.

EdelweissFE offers two generators to convert COO (Coordinate) format sparse matrices to CSR (Compressed Sparse Row) format efficiently without re-analysing the sparsity pattern:

  1. Legacy Generator (``csrgenerator``): a plain, single-threaded Cython implementation using a binary search for the COO-to-CSR mapping. Kept deliberately simple and short – not used by any solver, its purpose is teaching: reading it end to end is the fastest way to understand what the high-performance generator below does, before its parallelism, memory-layout tricks and the direct scatter-to-CSR alternative are layered on top.

  2. High-Performance Generator (``csrgeneratorv2``): A parallelized C++ engine with OpenMP support, thread-safe memory layouts, cache friendliness, and vectorized operations. This is what every solver actually uses.

See also

System matrix assembly compares this stage-then-gather route against direct scatter-to-CSR assembly, which avoids the staging array entirely – relevant whenever the number of stored contributions greatly exceeds the number of matrix entries, as it does for meshfree discretisations.

Legacy Generator (teaching reference)

Module edelweissfe.numerics.csrgenerator

class edelweissfe.numerics.csrgenerator.CSRGenerator

This cdef class generates Compressed Sparse Row Matrices from the COO format, and offers the possibility to update the matrix without reanalyzing the pattern (in contrast to SciPy).

Legacy, teaching-oriented generator: single-threaded, no memory-layout optimizations. Not used by any solver – see edelweissfe.numerics.csrgeneratorv2.CSRGenerator for the parallelized generator every solver actually uses, and doc/source/documentation/utils.rst for how the two relate.

Parameters:

systemMatrix – The system matrix in COO format, for the first initialization V can be a dummy (empty) vector.

Methods

updateCSR(V)

Get updated copies of the CSR matrix.
updateCSR(V)

Get updated copies of the CSR matrix.

Returns:

The system matrix in CSR format.

Return type:

csr_matrix

High-Performance Generator

Module edelweissfe.numerics.csrgeneratorv2

class edelweissfe.numerics.csrgeneratorv2.CSRGenerator

CSRGenerator class to create and manage a CSR matrix from COO format.

This class utilizes a C++ core for efficient conversion and updating of the CSR matrix.

Parameters:

  • systemMatrix (object) – An object containing COO format data with attributes I, J, and nDof.

  • patternOnly (bool) – Build only the CSR pattern (indptr/indices), not the gather map. Halves the sort array and, more importantly, removes the 32-bit limit on the pair count – gather_sources and assembly_ptr are int32 indices into the COO list, so a full build refuses more than INT32_MAX pairs. updateInPlace() and updateCSR() then raise, exactly as after releaseGatherMap(). Use it when the pattern is being borrowed by a DirectCSRAssembler and nothing will gather.

Attributes:
csrMatrix
gatherMapReleased

True once releaseGatherMap() has been called and this generator can no longer gather.

memoryBytes

Bytes held by the CSR pattern, the gather map and the CSR data array.

nnz

Number of stored entries in the CSR pattern.

Methods

addToPattern(entries)

Add entries to the CSR pattern that no COO pair contributes to; they are assembled as stored zeros.

releaseGatherMap()

Give back the gather map, keeping the CSR pattern -- and give up the ability to gather.

updateCSR(V)

Update the values of the CSR matrix and return an independent copy.

updateInPlace(V)

Update the values of the CSR matrix in-place based on the input vector V.
addToPattern(entries)

Add entries to the CSR pattern that no COO pair contributes to; they are assembled as stored zeros.

Makes room for values written into csrMatrix after the assembly, such as the multi-point-constraint condensation’s (see condensedPatternFillIn()).

csrMatrix is replaced by a new matrix with the extended pattern (and zero values). A matrix or pattern array obtained before keeps the old pattern and stays valid, but the generator no longer updates it.

Parameters:

entries – The entries to add (CSR, nDof x nDof), with sorted rows and none of them already in the pattern; only the pattern is read.

releaseGatherMap()

Give back the gather map, keeping the CSR pattern – and give up the ability to gather.

gather_sources is one int32 per COO pair, so the map is by far the largest thing this object holds (6.69 GiB at 43,350 DOF against 0.80 GiB for the pattern and the data array). A solver assembling directly into CSR borrows only the pattern, so holding the map would cancel most of what the direct path saves.

After this call updateInPlace() and updateCSR() raise. That is deliberate: the alternative is a nogil loop reading freed vectors, and a matrix that is quietly wrong is worse than one that is loudly unavailable. csrMatrix and the pattern stay valid and usable, and memoryBytes drops accordingly.

updateCSR(V)

Update the values of the CSR matrix and return an independent copy.

Use updateInPlace instead when the caller does not need to retain the matrix across subsequent assembly steps, to avoid the allocation cost of copying.

Parameters:

V (double[:]) – Input vector used to update the CSR matrix values.

Returns:

An independent copy of the updated CSR matrix.

Return type:

csr_matrix

updateInPlace(V)

Update the values of the CSR matrix in-place based on the input vector V.

Returns the internal CSR matrix directly (no copy). The caller must not retain the returned object across subsequent calls to updateInPlace or updateCSR, as the underlying data will be overwritten. It also must not be structurally mutated (eliminate_zeros(), prune(), resize(), a differently-shaped .data/.indices/.indptr reassignment, …) – the returned AliasedCSRMatrix raises rather than allowing this silently; see its docstring and GitHub issue #72.

Parameters:

V (double[:]) – Input vector used to update the CSR matrix values.

Returns:

A live view of the internal CSR matrix (not a copy).

Return type:

AliasedCSRMatrix

gatherMapReleased

True once releaseGatherMap() has been called and this generator can no longer gather.

memoryBytes

Bytes held by the CSR pattern, the gather map and the CSR data array.

The counterpart to DirectCSRAssembler.memoryBytes, on the same basis – what the object owns. The dominant term is the gather map’s gather_sources, one int32 per COO pair, which is what the direct path replaces with a uint16 offset per pair.

Neither figure counts the VIJ value array or the I/J index arrays: those belong to the DofManager, and I/J are needed by both paths.

nnz

Number of stored entries in the CSR pattern.

Under the hood, csrgeneratorv2 delegates sparse pattern building and assembly updates to a C++ engine class CSRCore (defined in _csrcore.h).

Algorithm & Optimization Details

  1. Row-based Partitioning: The rows are divided into partitions based on the number of threads (OpenMP max threads). By default, a heuristic of 4 partitions per thread is used to balance the load evenly across CPU cores.

  2. Parallel Histogram with Thread-Major Offset Layout: To calculate prefix sums and sort edges without atomic contention or cache conflicts: - A thread-local count of rows per partition is calculated. - The thread counts are reduced to global partition counts. - The write offset array write_offsets is transposed to a thread-major layout (indexed as [thread][partition]). Because each thread only writes to its own row in the offsets table, false sharing (cache line contention) is eliminated during the parallel scatter step under high thread counts.

  3. Parallel Scatter (Bucketing): COO elements (edges) are packed into 64-bit keys (High 32-bit: Row, Low 32-bit: Column) and scattered in parallel into partition buckets using the thread-major offsets.

  4. Parallel Sorting and Symbolic Counting: Each OpenMP thread takes ownership of a subset of partitions, sorts the packed edges using std::sort, and counts the unique non-zeros (nnz). Since partition boundaries map rows uniquely (the row partition assignment ensures no row overlaps partitions), thread-local writes to the indptr array are completely disjoint and synchronization-free.

  5. Vectorized Assembly Update: When updating the matrix via updateInPlace(V) or updateCSR(V), values corresponding to duplicated coordinate entries are summed. Because the lookup array gather_sources is contiguous in memory, the compiler is able to autovectorize the reduction loops using OpenMP SIMD directives:

    #pragma omp simd reduction(+:sum)
for (int32_t k = start; k < end; ++k) {
    sum += V_data[gather_sources[k]];
}

Memory Safety and Lifetime Management

  • Type Safety: The constructor of csrgeneratorv2.CSRGenerator explicitly casts the input coordinate arrays I and J to contiguous memoryviews of type int (representing 32-bit signed integers, np.intc). This ensures memory safety even if the coordinates originate from non-standard array types outside the dofmanager.

  • In-Place Matrix Updates: updateInPlace(V) modifies the internal SciPy csr_matrix.data buffer directly using raw C++ pointers, avoiding any Python object creation or memory allocations during iteration loops. Note that the returned matrix is a live view; its contents are mutated on subsequent update calls.

  • Lifetime of the C++ core: The pattern arrays (indptr, indices) of the generated SciPy csr_matrix point into the CSRCore. The core is owned by a small object, _CSRCoreOwner, held by the generator and by _CSRPattern, the NumPy base object of these arrays: the pattern stays valid as long as the generator or the matrix is referenced, whichever goes last; the gather map, which only the generator reads, goes with the generator. addToPattern moves the old pattern into its _CSRPattern, so a matrix obtained before keeps a valid pattern. Neither refers to the other, so no reference cycle keeps them alive until Python’s cyclic garbage collector runs: a replaced generator is freed right away.

Gathering efficiently of element results

Module edelweissfe.utils.elementresultcollector

class edelweissfe.utils.elementresultcollector.ElementResultCollector

A cdef class for collecting element results (by using the permanent results pointer (i.e., a numpy array) in large array of all elements and all quadrature points.

Collecting elemental results may be a performance critical part. This cdef class allows for the efficient gathering. A 3D array is assembled if multiple quadrature points are requested (shape [elements, quadraturePoints, resultVector]) or a 2D array for one quadrature point (shape [elements, resultVector]).

Method getCurrentResults() updates the assembly array and passes it back.

The caller is responsible to make a copy of it, if persistent results are needed!

Parameters:

  • elements – The list of elements for which the results should be collected.

  • quadraturePoints – The slice defining the desired quadrature points.

  • result – The name of the requested result.

Attributes:
resultsTable

Methods

getCurrentResults()

Update and get current results.

update()

Update all results.
getCurrentResults()

Update and get current results.

Returns:

The results array.

Return type:

np.ndarray

update()

Update all results.

Time stepping

A time stepper divides a simulation step into increments. It is a small state machine:

while not stepper.isFinished():
    timeStep = stepper.proposeTimeStep()      # the next increment; changes nothing
    try:
        solve the increment
    except failed:
        stepper.rejectTimeStep(cutbackFactor) # a smaller increment is proposed next
        continue
    stepper.acceptTimeStep(timeStep)          # the step advances by this increment
    write output (and restart checkpoints)

All the stepper knows is in its attributes, and after acceptTimeStep() they describe exactly the state the next increment starts from. A restart checkpoint is written after an increment is accepted, so it holds those attributes as they are, see checkpointedState.

class edelweissfe.timesteppers.base.timestepperbase.TimeStepperBase[source]

Base class for all time steppers.

It defines the interface which all solvers may rely on for controlling the incrementation of a simulation step.

Attributes:
checkpointedState

Methods

acceptTimeStep(timeStep)

Advance the step by timeStep, the increment just solved, and choose the size of the next.

changeIncrementSize(scaleFactor)

Modify the size of the next increment by a given scale factor within the bounds of the minimum and maximum increment size.

enforceTimeIncrement(timeIncrement)

Use timeIncrement for every following increment, instead of adapting it.

getRestartData()

The declared state, for a restart checkpoint.

isAtStepStart()

Whether no increment of the step was accepted yet.

isFinished()

Whether the step has reached its end.

numberOfIncrementsDone()

The number of increments accepted so far, counted against the maximum number of increments.

preventIncrementIncrease()

Keep the next increment from growing, e.g., in case of bad convergence.

proposeTimeStep()

The next increment, starting from the last accepted one.

rejectTimeStep(cutbackFactor)

Discard the proposed increment, which could not be solved; the next one proposed is smaller by cutbackFactor, within the bounds of the minimum and maximum increment size.

setRestartData(data)

Restore the declared state from a restart checkpoint.
abstractmethod acceptTimeStep(timeStep)[source]

Advance the step by timeStep, the increment just solved, and choose the size of the next.

Parameters:

timeStep (TimeStep) – The increment, as proposeTimeStep() returned it.

abstractmethod changeIncrementSize(scaleFactor)[source]

Modify the size of the next increment by a given scale factor within the bounds of the minimum and maximum increment size.

Parameters:

scaleFactor (float) – The factor for scaling based on the current increment.

enforceTimeIncrement(timeIncrement)[source]

Use timeIncrement for every following increment, instead of adapting it.

The caller is an explicit solver, whose stable time increment is a property of the mesh: set at the start of the step, and lowered when a refinement shrinks the smallest element.

Parameters:

timeIncrement (float) – The time increment to use.

Raises:

NotImplementedError – If this stepper cannot run on an enforced time increment.

getRestartData()[source]

The declared state, for a restart checkpoint.

Returns:

The state; see packState().

Return type:

dict[str, numpy.ndarray]

abstractmethod isAtStepStart()[source]

Whether no increment of the step was accepted yet. A step resumed from a checkpoint never is: a checkpoint is written after an accepted increment.

Returns:

True at the start of a step.

Return type:

bool

abstractmethod isFinished()[source]

Whether the step has reached its end.

Returns:

True if no increment is left.

Return type:

bool

abstractmethod numberOfIncrementsDone()[source]

The number of increments accepted so far, counted against the maximum number of increments.

Returns:

The number of increments.

Return type:

int

abstractmethod preventIncrementIncrease()[source]

Keep the next increment from growing, e.g., in case of bad convergence. Called before acceptTimeStep().

abstractmethod proposeTimeStep()[source]

The next increment, starting from the last accepted one. Changes nothing: proposing twice gives the same increment.

Returns:

The proposed increment.

Return type:

TimeStep

Raises:

ReachedMaxIncrements – If the step needs more increments than it is allowed.

abstractmethod rejectTimeStep(cutbackFactor)[source]

Discard the proposed increment, which could not be solved; the next one proposed is smaller by cutbackFactor, within the bounds of the minimum and maximum increment size.

Parameters:

cutbackFactor (float) – The factor scaling the increment size.

Raises:

ReachedMinIncrementSize – If the increment cannot be reduced any further: it is at or below the minimum size. The last increment of a step, the remainder of the step, may be below it.

setRestartData(data)[source]

Restore the declared state from a restart checkpoint.

Parameters:

data (dict[str, ndarray]) – What getRestartData() returned.

checkpointedState: dict | None = None

The attributes that make up the stepper’s progress within the step, by type; see checkpointedstate. Not its configuration, which a resumed run takes from its own input file.

identification = 'TimeStepper'

Overridden by every stepper; a default so base-class diagnostics can name their source.

Module edelweissfe.timesteppers.adaptivetimestepper

class edelweissfe.timesteppers.adaptivetimestepper.AdaptiveTimeStepper(currentTime, stepLength, startIncrement, maxIncrement, minIncrement, maxNumberIncrements, journal, increaseFactor=1.1, makeZeroIncrementFirst=True, goodIncrementsBeforeIncrease=3)[source]

Divides a step into increments whose size adapts to the convergence of the solver: it grows after good increments and is cut back after failed ones.

Parameters:

  • currentTime (float) – The start time of the step.

  • stepLength (float) – The total length of the step.

  • startIncrement (float) – The size of the start increment, relative to the step length.

  • maxIncrement (float) – The maximum size of an increment, relative to the step length.

  • minIncrement (float) – The minimum size of an increment, relative to the step length.

  • maxNumberIncrements (int) – The maximum number of allowed increments.

  • journal (Journal) – The journal instance for logging purposes.

  • increaseFactor (float) – The ratio to increase the increments in case of good convergence.

  • makeZeroIncrementFirst (bool) – If True, the first increment will be zero.

  • goodIncrementsBeforeIncrease (int) – The number of good increments that must already have passed before an accepted increment lets the next one grow: counted from the step start (the zero increment included), or from the last cutback, which itself counts as one. From then on, every further good increment lets the next one grow.

Methods

acceptTimeStep(timeStep)

Advance the step by timeStep, the increment just solved, and choose the size of the next.

changeIncrementSize(scaleFactor)

Modify the size of the next increment by a given scale factor within the bounds of the minimum and maximum increment size.

enforceTimeIncrement(timeIncrement)

Use timeIncrement for every following increment, instead of adapting it.

getRestartData()

The declared state, for a restart checkpoint.

isAtStepStart()

Whether no increment of the step was accepted yet.

isFinished()

Whether the step has reached its end.

numberOfIncrementsDone()

The number of increments accepted so far, counted against the maximum number of increments.

preventIncrementIncrease()

Keep the next increment from growing, e.g., in case of bad convergence.

proposeTimeStep()

The next increment, starting from the last accepted one.

rejectTimeStep(cutbackFactor)

Discard the proposed increment, which could not be solved; the next one proposed is smaller by cutbackFactor, within the bounds of the minimum and maximum increment size.

setRestartData(data)

Restore the declared state from a restart checkpoint.
acceptTimeStep(timeStep)[source]

Advance the step by timeStep, the increment just solved, and choose the size of the next.

Parameters:

timeStep (TimeStep) – The increment, as proposeTimeStep() returned it.

changeIncrementSize(scaleFactor)[source]

Modify the size of the next increment by a given scale factor within the bounds of the minimum and maximum increment size.

Parameters:

scaleFactor (float) – The factor for scaling based on the current increment.

isAtStepStart()[source]

Whether no increment of the step was accepted yet. A step resumed from a checkpoint never is: a checkpoint is written after an accepted increment.

Returns:

True at the start of a step.

Return type:

bool

isFinished()[source]

Whether the step has reached its end.

Returns:

True if no increment is left.

Return type:

bool

numberOfIncrementsDone()[source]

The number of increments accepted so far, counted against the maximum number of increments.

Returns:

The number of increments.

Return type:

int

preventIncrementIncrease()[source]

Keep the next increment from growing, e.g., in case of bad convergence. Called before acceptTimeStep().

proposeTimeStep()[source]

The next increment, starting from the last accepted one. Changes nothing: proposing twice gives the same increment.

Returns:

The proposed increment.

Return type:

TimeStep

Raises:

ReachedMaxIncrements – If the step needs more increments than it is allowed.

rejectTimeStep(cutbackFactor)[source]

Discard the proposed increment, which could not be solved; the next one proposed is smaller by cutbackFactor, within the bounds of the minimum and maximum increment size.

Parameters:

cutbackFactor (float) – The factor scaling the increment size.

Raises:

ReachedMinIncrementSize – If the increment cannot be reduced any further: it is at or below the minimum size. The last increment of a step, the remainder of the step, may be below it.

allowedToIncreasedNext

False if the solver asked to keep the next increment from growing.

checkpointedState: dict | None = {'allowedToIncreasedNext': <class 'bool'>, 'currentTime': <class 'float'>, 'finishedStepProgress': <class 'float'>, 'increment': <class 'float'>, 'incrementCounter': <class 'int'>, 'nPassedGoodIncrements': <class 'int'>}

The attributes that make up the stepper’s progress within the step, by type; see checkpointedstate. Not its configuration, which a resumed run takes from its own input file.

currentTime

The start time of the step.

finishedStepProgress

The accepted progress within the step, from 0 to 1.

identification = 'AdaptiveTimeStepper'

Overridden by every stepper; a default so base-class diagnostics can name their source.

increment

The size of the next increment, relative to the step length.

incrementCounter

The number of the next increment; the zero increment, if any, is number 0.

nPassedGoodIncrements

Counts the good increments; the increment grows once goodIncrementsBeforeIncrease have passed.

Module edelweissfe.timesteppers.simpletimestepper

class edelweissfe.timesteppers.simpletimestepper.SimpleTimeStepper(currentTime, stepLength, startIncrement, maxIncrement, minIncrement, maxNumberIncrements, journal)[source]

Divides a step into increments of a fixed size, or of a size enforced by the solver, after one zero increment.

The zero increment comes first, so that an explicit integrator can build its initial state (mass, internal force, contact) before taking a step. It is numbered 0; the increments after it are numbered from 1. A step resumed from a checkpoint does not repeat it: whether it was done is part of the stepper’s state.

Parameters:

  • currentTime (float) – The start time of the step.

  • stepLength (float) – The total length of the step.

  • startIncrement (float) – The size of the start increment, relative to the step length.

  • maxIncrement (float) – The maximum size of an increment, relative to the step length.

  • minIncrement (float) – The minimum size of an increment, relative to the step length.

  • maxNumberIncrements (float) – The maximum number of allowed increments.

  • journal (Journal) – The journal instance for logging purposes.

Methods

acceptTimeStep(timeStep)

Advance the step by timeStep, the increment just solved, and choose the size of the next.

changeIncrementSize(scaleFactor)

Modify the size of the next increment by a given scale factor within the bounds of the minimum and maximum increment size.

enforceTimeIncrement(timeIncrement)

Use timeIncrement for every following increment, instead of adapting it.

getRestartData()

The declared state, for a restart checkpoint.

isAtStepStart()

Whether no increment of the step was accepted yet.

isFinished()

Whether the step has reached its end.

numberOfIncrementsDone()

The number of increments accepted so far, counted against the maximum number of increments.

preventIncrementIncrease()

This time stepper never increases the increment size automatically, hence this is a no-op.

proposeTimeStep()

The next increment, starting from the last accepted one.

rejectTimeStep(cutbackFactor)

Discard the proposed increment, which could not be solved; the next one proposed is smaller by cutbackFactor, within the bounds of the minimum and maximum increment size.

setRestartData(data)

Restore the declared state from a restart checkpoint.
acceptTimeStep(timeStep)[source]

Advance the step by timeStep, the increment just solved, and choose the size of the next.

Parameters:

timeStep (TimeStep) – The increment, as proposeTimeStep() returned it.

changeIncrementSize(scaleFactor)[source]

Modify the size of the next increment by a given scale factor within the bounds of the minimum and maximum increment size.

Parameters:

scaleFactor (float) – The factor for scaling based on the current increment.

enforceTimeIncrement(timeIncrement)[source]

Use timeIncrement for every following increment, instead of adapting it.

The caller is an explicit solver, whose stable time increment is a property of the mesh: set at the start of the step, and lowered when a refinement shrinks the smallest element.

Parameters:

timeIncrement (float) – The time increment to use.

Raises:

NotImplementedError – If this stepper cannot run on an enforced time increment.

isAtStepStart()[source]

Whether no increment of the step was accepted yet. A step resumed from a checkpoint never is: a checkpoint is written after an accepted increment.

Returns:

True at the start of a step.

Return type:

bool

isFinished()[source]

Whether the step has reached its end.

Returns:

True if no increment is left.

Return type:

bool

numberOfIncrementsDone()[source]

The number of increments accepted so far, counted against the maximum number of increments.

Returns:

The number of increments.

Return type:

int

preventIncrementIncrease()[source]

This time stepper never increases the increment size automatically, hence this is a no-op.

proposeTimeStep()[source]

The next increment, starting from the last accepted one. Changes nothing: proposing twice gives the same increment.

Returns:

The proposed increment.

Return type:

TimeStep

Raises:

ReachedMaxIncrements – If the step needs more increments than it is allowed.

rejectTimeStep(cutbackFactor)[source]

Discard the proposed increment, which could not be solved; the next one proposed is smaller by cutbackFactor, within the bounds of the minimum and maximum increment size.

Parameters:

cutbackFactor (float) – The factor scaling the increment size.

Raises:

ReachedMinIncrementSize – If the increment cannot be reduced any further: it is at or below the minimum size. The last increment of a step, the remainder of the step, may be below it.

checkpointedState: dict | None = {'currentTime': <class 'float'>, 'enforcedTimeIncrement': <class 'float'>, 'finishedStepProgress': <class 'float'>, 'increment': <class 'float'>, 'totalIncrements': <class 'int'>, 'zeroIncrementDone': <class 'bool'>}

The attributes that make up the stepper’s progress within the step, by type; see checkpointedstate. Not its configuration, which a resumed run takes from its own input file.

currentTime

The start time of the step.

enforcedTimeIncrement

The time increment the solver enforces, or None.

finishedStepProgress

The accepted progress within the step, from 0 to 1.

identification = 'SimpleTimeStepper'

Overridden by every stepper; a default so base-class diagnostics can name their source.

increment

The size of the next increment, relative to the step length.

totalIncrements

The number of increments accepted after the zero increment.

zeroIncrementDone

Whether the zero increment was accepted.

Parent element faces and facet quadrature

The canonical face types behind the integrated contact formulation (The integrated formulation, step by step): the shape functions of a facet’s parent element face, and the quadrature rules used to integrate against them over a flat contact facet.

Module edelweissfe.utils.parentfacegeometry

edelweissfe.utils.parentfacegeometry.facetQuadratureRule(nFacetNodes, nPoints)[source]

The quadrature rule for one contact facet, in the facet’s own barycentric coordinates.

Parameters:

  • nFacetNodes (int) – 3 for a Tria3 facet, 2 for a Line2 facet.

  • nPoints (int) – The number of quadrature points requested.

Returns:

The barycentric coordinates, shape (nPoints, nFacetNodes), and the weights, which sum to one – so an integral over the facet is facetMeasure * sum(weights * f(points)).

Return type:

tuple[numpy.ndarray, numpy.ndarray]

edelweissfe.utils.parentfacegeometry.parentFaceShapeFunctions(faceType, parametricCoords)[source]

The parent face’s shape functions at one parametric point.

Parameters:

  • faceType (str) – One of the keys of PARENT_FACE_PARAMETRIC_COORDS.

  • parametricCoords (ndarray) – The evaluation point, shape (1,) for a line face or (2,) for a quad face.

Returns:

The shape function values, ordered like the canonical node ordering and summing to one. For quad8 and line3 these take negative values near the corners, which is the whole point of this module.

Return type:

numpy.ndarray

edelweissfe.utils.parentfacegeometry.PARENT_FACE_PARAMETRIC_COORDS = {'line2': array([[-1.],        [ 1.]]), 'line3': array([[-1.],        [ 1.],        [ 0.]]), 'quad4': array([[-1., -1.],        [ 1., -1.],        [ 1.,  1.],        [-1.,  1.]]), 'quad8': array([[-1., -1.],        [ 1., -1.],        [ 1.,  1.],        [-1.,  1.],        [ 0., -1.],        [ 1.,  0.],        [ 0.,  1.],        [-1.,  0.]])}

Parametric coordinates of each canonical parent-face node ordering, shape (nNodes, nDim - 1).