Tools#

Optimization

newtonraphson([x0, fun, jac, solve, ...])

Find a root of a real function using the Newton-Raphson method.

tools.IterationState([iteration, fnorm, ...])

A class to keep track of the state of an iteration during evaluation.

tools.NewtonResult(x[, fun, jac, success, ...])

A data class which represents the result found by Newton's method.

Event Dispatcher

tools.EventDispatcher([plugins])

A class to dispatch events to plugins during evaluation.

tools.Context([job, step, substep])

A class to keep track of the context of a Job during evaluation.

Export of Results

save(region, field[, forces, gradient, ...])

Write field-data to a VTU file.

Convert Cell-Data to Point-Data

topoints(values, region[, average, mean])

Shift array of values located at quadrature points of cells to mesh-points.

project(values, region[, average, mean, dV, ...])

Project given values at quadrature-points to mesh-points.

tools.extrapolate(values, region[, average, ...])

Extrapolate (and optionally average) an array with values at quadrature points to mesh-points.

Reaction-Force and -Moment

tools.force(field, forces, boundary)

Evaluate the force vector sum on points of a boundary.

tools.moment(field, forces, boundary[, ...])

Evaluate the moment vector sum on points of a boundary at a given center point.

Detailed API Reference

felupe.newtonraphson(x0=None, fun=<function fun>, jac=<function jac>, solve=<function solve>, maxiter=16, update=<function update>, check=<function check>, args=(), kwargs=None, tol=np.float64(1.4901161193847656e-08), items=None, dof1=None, dof0=None, ext0=None, solver=<function spsolve>, verbose=None, callback=None, callback_kwargs=None, progress_bar=None, tqdm='tqdm', dispatcher=None)[ソース]#

Find a root of a real function using the Newton-Raphson method.

パラメータ:
  • x0 (felupe.FieldContainer, ndarray or None (optional)) -- Array or Field container with values of unknowns at a valid starting point (default is None).

  • fun (callable, optional) -- Callable which assembles the vector-valued objective function. Additional args and kwargs are passed. Function signature of a user-defined function has to be fun = lambda x, *args, **kwargs: f.

  • jac (callable, optional) -- Callable which assembles the matrix-valued Jacobian. Additional args and kwargs are passed. Function signature of a user-defined Jacobian has to be jac = lambda x, *args, **kwargs: K.

  • solve (callable, optional) -- Callable which prepares the linear equation system and solves it. If a keyword- argument from the list ["x", "dof1", "dof0", "ext0", "solver"] is found in the function-signature, then these arguments are passed to solve.

  • maxiter (int, optional) -- Maximum number of function iterations (default is 16).

  • update (callable, optional) -- Callable to update the unknowns. Function signature must be update = lambda x0, dx: x.

  • check (callable, optional) -- Callable to check the result with signature check = lambda xnorm, fnorm, success: dx, x, **kwargs.

  • tol (float, optional) -- Tolerance value to check if the function has converged (default is 1.490e-8).

  • items (list or None, optional) -- List with items which provide methods for assembly, e.g. like felupe.SolidBody (default is None).

  • dof1 (ndarray or None, optional) -- 1d-array of int with all active degrees of freedom (default is None).

  • dof0 (ndarray or None, optional) -- 1d-array of int with all prescribed degrees of freedom (default is None).

  • ext0 (ndarray or None, optional) -- Field values at mesh-points for the prescribed components of the unknowns based on dof0 (default is None).

  • solver (callable, optional) -- A sparse or dense solver (default is scipy.sparse.linalg.spsolve()). For a more performant alternative install PyPardiso and use pypardiso.spsolve().

  • verbose (bool or int or None, optional) -- Verbosity level to control how messages are printed during evaluation. If 1 or True and tqdm is installed, a progress bar is shown. If tqdm is missing or verbose is 2, more detailed text-based messages are printed. Default is None. If None, verbosity is set to True. If None and the environmental variable FELUPE_VERBOSE is set and its value is not true, then logging is turned off.

  • callback (callable or None, optional) -- An optional callback function with function signature callback = lambda dx, x, iteration, xnorm, fnorm, success: None, which is called after each completed iteration. Default is None.

  • progress_bar (tqdm or None, optional) -- Use an existing instance of a progress bar if verbose is True. If None and verbose is True, a new bar is created. Default is None.

  • tqdm (str, optional) -- If verbose is True, choose a backend for tqdm ("tqdm", "auto" or "notebook"). Default is "tqdm".

  • dispatcher (EventDispatcher or None, optional) -- An optional EventDispatcher to trigger events during evaluation. Default is None.

戻り値:

The result object.

戻り値の型:

felupe.tools.NewtonResult

メモ

Nonlinear equilibrium equations \(f(x)\) as a function of the unknowns \(x\) are solved by linearization of \(f\) at a valid starting point of given unknowns \(x_0\), see Eq. (1).

(1)#\[f(x_0) = 0\]

The linearization is given in Eq. (2)

(2)#\[f(x_0 + dx) \approx f(x_0) + K(x_0) \ dx \ (= 0)\]

with the Jacobian as in Eq. (3), evaluated at given unknowns \(x_0\),

(3)#\[K(x_0) = \frac{\partial f}{\partial x}(x_0)\]

and is rearranged to an equation system with left- and right-hand sides, see Eq. (4).

(4)#\[K(x_0) \ dx = -f(x_0)\]

After a solution is found, the unknowns are updated, see Eq. (5).

(5)#\[ \begin{align}\begin{aligned}dx &= \text{solve} \left( K(x_0), -f(x_0) \right) \nonumber\\x &= x_0 + dx\end{aligned}\end{align} \]

Repeated evaluations lead to an incrementally updated solution of \(x\). Herein, \(x_n\) refer to the initial unknowns whereas \(x\) are the updated unknowns, see Eq. (6).

注釈

The subscript \((\bullet)_{n+1}\) is dropped for easier readability.

(6)#\[ \begin{align}\begin{aligned}dx &= \text{solve} \left( K(x_n), -f(x_n) \right) \nonumber\\ x &= x_n + dx\end{aligned}\end{align} \]

Then, the nonlinear equilibrium equations are evaluated with the updated unknowns \(f(x)\). The procedure is repeated until convergence is reached.

サンプル

>>> import felupe as fem
>>>
>>> region = fem.RegionHexahedron(fem.Cube(n=6))
>>> field = fem.FieldContainer([fem.Field(region, dim=3)])
>>> boundaries = fem.dof.uniaxial(
...     field, move=0.2, clamped=True, return_loadcase=False
... )
>>> solid = fem.SolidBody(umat=fem.NeoHooke(mu=1.0, bulk=2.0), field=field)
>>> res = fem.newtonraphson(items=[solid], **loadcase)

Newton-Raphson solver
=====================

| # | norm(fun) |  norm(dx) |
|---|-----------|-----------|
| 1 | 7.553e-02 | 1.898e+00 |
| 2 | 1.310e-03 | 5.091e-02 |
| 3 | 3.086e-07 | 6.698e-04 |
| 4 | 2.255e-14 | 1.527e-07 |

Converged in 4 iterations ...

Newton's method had success

>>> res.success
True

and 4 iterations were needed to converge within the specified tolerance.

>>> res.iterations
4

The norm of the objective function for all active degrees of freedom is lower than 3e-15.

>>> np.linalg.norm(res.fun[loadcase["dof1"]])
2.7384964752762237e-15
class felupe.tools.NewtonResult(x, fun=None, jac=None, success=None, iterations=None, xnorms=None, fnorms=None)[ソース]#

A data class which represents the result found by Newton's method. All parameters are available as attributes.

パラメータ:
  • x (felupe.FieldContainer or ndarray) -- Array or Field container with values at a solution found by Newton's method.

  • fun (ndarray or None, optional) -- Values of objective function (default is None).

  • jac (ndarray or None, optional) -- Values of the Jacobian of the objective function (default is None).

  • success (bool or None, optional) -- A boolean flag which is True if the solution converged (default is None).

  • iterations (int or None, optional) -- Number of iterations until solution converged (default is None).

  • xnorms (array of float or None, optional) -- List with norms of the values of the solution (default is None).

  • fnorms (float or None, optional) -- List with norms of the objective function (default is None).

メモ

The objective function's norm is relative based on the function values on the prescribed degrees of freedom \((\bullet)_0\). A small number \(\varepsilon\) is added to avoid numeric instabilities.

\[\text{norm}(\boldsymbol{f}) = \frac{||\boldsymbol{f}_1||} {\varepsilon + ||\boldsymbol{f}_0||}\]
class felupe.tools.EventDispatcher(plugins=None)[ソース]#

A class to dispatch events to plugins during evaluation.

パラメータ:

plugins (list or None, optional) -- The list of plugins.

add_plugin(plugin)[ソース]#

Add a plugin to the dispatcher and reconfigure it.

configure(plugins)[ソース]#

Configure the dispatcher with the given plugins and return a dict with hooks as keys.

trigger(hook, context, state)[ソース]#

Trigger a hook with context and current state to all registered functions.

class felupe.tools.Context(job=None, step=None, substep=None)[ソース]#

A class to keep track of the context of a Job during evaluation.

パラメータ:
felupe.save(region, field, forces=None, gradient=None, filename='result.vtu', cell_data=None, point_data=None)[ソース]#

Write field-data to a VTU file.

パラメータ:
  • region (Region) -- The region to be saved.

  • field (FieldContainer) -- The field container to be saved.

  • forces (ndarray, optional) -- Array with reaction forces to be saved (default is None).

  • gradient (list of ndarray, optional) -- The result of umat.gradient() with the first Piola-Kirchhoff stress tensor as the first list item (default is None).

  • filename (str, optional) -- The filename for the results (default is "result.vtu").

  • cell_data (dict or None, optional) -- Additional dict with cell data (default is None).

  • point_data (dict or None, optional) -- Additional dict with point data (default is None).

サンプル

>>> import felupe as fem
>>>
>>> mesh = fem.Cube(n=6)
>>> region = fem.RegionHexahedron(mesh)
>>> field = fem.FieldContainer([fem.Field(region, dim=3)])
>>>
>>> boundaries = fem.dof.uniaxial(
...     field, clamped=True, move=0.3, return_loadcase=False
... )
>>>
>>> umat = fem.NeoHooke(mu=1)
>>> solid = fem.SolidBodyNearlyIncompressible(umat, field, bulk=5000)
>>> step = fem.Step(items=[solid], boundaries=boundaries)
>>> job = fem.CharacteristicCurve(steps=[step], boundary=boundaries["move"]).evaluate()
>>> fem.save(region, field, forces=job.res.fun, gradient=solid.results.stress)
felupe.topoints(values, region, average=True, mean=False)[ソース]#

Shift array of values located at quadrature points of cells to mesh-points.

パラメータ:
  • values (ndarray of shape (..., q, c)) -- Array with values located at the quadrature points q of cells c.

  • region (Region) -- A region used to shift the values to the mesh-points.

  • average (bool, optional) -- A flag to return values averaged at mesh-points if True and the mesh of the region is not already disconnected or to return values on a disconnected mesh with discontinuous values at cell transitions if False (default is True).

  • mean (bool, optional) -- A flag to shift the cell-means of the values to the mesh-points (default is False).

戻り値:

Field.values -- Array of values shifted to the mesh-points p.

戻り値の型:

ndarray of shape (p, ...)

メモ

If the number of quadrature-points per cell is greater than the number of points- per-cell, then the values are trimmed.

felupe.project(values, region, average=True, mean=False, dV=None, solver=<function spsolve>)[ソース]#

Project given values at quadrature-points to mesh-points.

パラメータ:
  • values (ndarray of shape (..., q, c)) -- Array with values located at the quadrature points q of cells c.

  • region (Region) -- A region used to project the values to the mesh-points.

  • average (bool, optional) -- A flag to return values averaged at mesh-points if True and the mesh of the region is not already disconnected or to return values on a disconnected mesh with discontinuous values at cell transitions if False (default is True).

  • mean (bool, optional) -- A flag to extrapolate the cell-means of the values to the mesh-points, see felupe.tools.extrapolate(). If True, dV and solver are ignored. Default is False.

  • dV (ndarray of shape (q, c) or None, optional) -- Differential volumes located at the quadrature points of the cells. If None, the differential volumes are taken from the region (default is None).

  • solver (callable, optional) -- A function for a sparse solver with signature x=solver(A, b). Default is scipy.sparse.linalg.spsolve().

戻り値:

Field.values -- Array of values projected to the mesh-points p.

戻り値の型:

ndarray of shape (p, ...)

メモ

The projection finds \(\hat{\boldsymbol{u}}\), located at mesh-points p, for values \(v\), evaluated at quadrature-points q, such that the variation of the minimization problem in Eq. (7) is solved.

(7)#\[ \begin{align}\begin{aligned}\Pi &= \int_V \frac{1}{2} (u - v) \cdot (u - v)\ dV \quad \rightarrow \quad \min\\\delta \Pi &= \int_V (u - v) \cdot \delta u\ dV = 0\end{aligned}\end{align} \]

This leads to assembled system-matrix and -vector(s) as shown in Eq. (8)

(8)#\[ \begin{align}\begin{aligned}\int_V v\ \cdot \delta u\ dV &\qquad \rightarrow \qquad \hat{\boldsymbol{A}}\\\int_V u\ \cdot \delta u\ dV &\qquad \rightarrow \qquad \hat{\boldsymbol{b}}\end{aligned}\end{align} \]

of an equation system to be solved, see Eq. (9). The right-arrows in Eq. (8) represent the assembly of the integral forms into the system matrix or vectors.

(9)#\[\hat{\boldsymbol{A}}\ \hat{\boldsymbol{u}} = \hat{\boldsymbol{b}}\]

The region, on which the values are projected to, may contain a mesh with a different (e.g. a disconnected) points- but same cells-array. The only requirement on the region is that it must use a compatible quadrature scheme. With the values at mesh-points at hand, new fields may be created to project between fields on different regions. If the array of differential volumes is not given, then the region must be created with Region(grad=True) to include the differential volumes.

警告

Triangular element formulations require regions with quadrature schemes with higher integration orders, the default choices are not sufficient for projection. For RegionTriangle and RegionTetra a second-order scheme is used.

import felupe as fem

mesh = fem.Rectangle(n=2).triangulate().add_midpoints_edges()
quadrature = fem.TriangleQuadrature(order=5)
fem.RegionQuadraticTriangle(mesh, quadrature=quadrature)

For RegionQuadraticTriangle and RegionQuadraticTetra use a fifth-order scheme.

mesh = fem.Cube(n=2).triangulate().add_midpoints_edges()
quadrature = fem.TetrahedronQuadrature(order=5)
fem.RegionQuadraticTetra(mesh, quadrature=quadrature)

サンプル

import felupe as fem

mesh = fem.Rectangle(n=2).triangulate()
region = fem.RegionTriangle(mesh)
field = fem.FieldAxisymmetric(region, dim=2)
values = field.extract()

region2 = region.copy()
region2.reload(quadrature=fem.TriangleQuadrature(order=2))
values_projected = project(values, region2)
felupe.tools.extrapolate(values, region, average=True, mean=False)[ソース]#

Extrapolate (and optionally average) an array with values at quadrature points to mesh-points.

パラメータ:
  • values (ndarray of shape (..., q, c)) -- Array with values located at the quadrature points q of cells c.

  • region (Region) -- A region used to extrapolate the values to the mesh-points.

  • average (bool, optional) -- A flag to return values averaged at mesh-points if True and the mesh of the region is not already disconnected or to return values on a disconnected mesh with discontinuous values at cell transitions if False (default is True).

  • mean (bool, optional) -- A flag to take the cell-means, averaged by the quadrature weights, of the values (default is False).

戻り値:

Field.values -- Array of values projected to the mesh-points p.

戻り値の型:

ndarray of shape (p, ...)

felupe.tools.force(field, forces, boundary)[ソース]#

Evaluate the force vector sum on points of a boundary.

felupe.tools.moment(field, forces, boundary, centerpoint=array([0., 0., 0.]))[ソース]#

Evaluate the moment vector sum on points of a boundary at a given center point.