movement package

Submodules

movement.laplacian_smoothing module

class LaplacianSmoother(mesh, timestep, **kwargs)[source]

Bases: PrimeMover

Movement of a mesh is driven by a mesh velocity \(\mathbf{v}\), which is determined by solving a vector Laplace problem

\[\nabla^2_{\boldsymbol{\xi}}\mathbf{v} = \boldsymbol{0}, \quad \boldsymbol{\xi}\in\Omega,\]

under non-zero Dirichlet boundary conditions on a forced boundary section \(\partial\Omega_f\) and zero Dirichlet boundary conditions elsewhere:

\[\begin{split}\mathbf{v} = \left\{\begin{array}{rl} \mathbf{v}_D, & \boldsymbol{\xi}\in\partial\Omega_f\\ \boldsymbol{0}, & \boldsymbol{\xi}\in\partial\Omega\backslash\partial\Omega_f \end{array}\right.\end{split}\]

where the computational coordinates \(\boldsymbol{\xi} := \mathbf{x}(t_0)\) are the physical coordinates at the beginning of the simulation.

Parameters:
move(time, update_boundary_velocity=None, boundary_conditions=None, solver_parameters=None)[source]

Assemble and solve the Laplacian system and update the coordinates.

Parameters:

  • time (float) – the current time

  • update_boundary_velocity (Callable with a single argument of float type) – function that updates the boundary conditions at the current time

  • boundary_conditions (DirichletBC or list thereof) – Dirichlet boundary conditions to be enforced

  • solver_parameters (dict) – solver parameters for solving the Laplace equation

movement.math module

movement.monge_ampere module

Mesh movement based on solutions of equations of Monge-Ampère type.

MongeAmpereMover(mesh, monitor_function, method='relaxation', **kwargs)[source]

Factory function for generating Monge-Ampère mesh movers.

Movement of a mesh is determined by a monitor_function \(m\) and an equation of Monge-Ampère type,

\[m(\mathbf{x})\det(\mathbf{I} + \mathbf{H}(\phi)) = \theta,\]

for a convex scalar potential \(\phi=\phi(\boldsymbol{\xi})\), which is a function of the coordinates of the computational mesh. Here \(m=m(\mathbf{x})\) is a function of the coordinates of the physical mesh, \(\mathbf{I}\) is the identity matrix, \(\theta\) is a normalisation coefficient and \(\mathbf{H}(\phi)\) denotes the Hessian of \(\phi\) with respect to \(\boldsymbol{\xi}\).

Different implementations solve the Monge-Ampère equation in different ways. If the method argument is set to “relaxation” then it is solved in parabolised form in MongeAmpereMover_Relaxation. If the argument is set to “quasi_newton” then it is solved in its elliptic form using a quasi-Newton method in MongeAmpereMover_QuasiNewton. Descriptions of both methods may be found in [MCB18].

The physical mesh coordinates \(\mathbf{x}\) are updated according to

\[\mathbf{x} = \boldsymbol{\xi} + \nabla_{\boldsymbol{\xi}}\phi.\]
Parameters:

  • mesh (firedrake.mesh.MeshGeometry) – the physical mesh

  • monitor_function (Callable) – a Python function which takes a mesh as input

  • method (str) – choose from ‘relaxation’ and ‘quasi_newton’

  • phi_init (firedrake.function.Function) – initial guess for the scalar potential

  • H_init (firedrake.function.Function) – initial guess for the Hessian

  • maxiter (int) – maximum number of iterations for the solver

  • rtol (float) – relative tolerance for the residual

  • dtol (float) – divergence tolerance for the residual

  • pseudo_timestep (float) – pseudo-timestep (only relevant to relaxation method)

  • fixed_boundary_segments (list of str or int) – labels corresponding to boundary segments to be fixed with a zero Dirichlet condition. The ‘on_boundary’ label indicates the whole domain boundary

Returns:

the Monge-Ampère Mover object

Return type:

MongeAmpereMover_Relaxation or MongeAmpereMover_QuasiNewton

class MongeAmpereMover_QuasiNewton(mesh, monitor_function, phi_init=None, H_init=None, **kwargs)[source]

Bases: MongeAmpereMover_Base

The standard, elliptic form of the Monge-Ampère equation used for mesh movement is:

\[m(\mathbf{x})\det(\mathbf{I} + \mathbf{H}(\phi)) = \theta,\]

for a convex scalar potential \(\phi=\phi(\boldsymbol{\xi})\), which is a function of the coordinates of the computational mesh. Here \(m=m(\mathbf{x})\) is a user-provided monitor function, which is a function of the coordinates of the physical mesh. \(\mathbf{I}\) is the identity matrix, \(\theta\) is a normalisation coefficient and \(\mathbf{H}(\phi)\) denotes the Hessian of \(\phi\) with respect to \(\boldsymbol{\xi}\).

In this mesh mover, the elliptic Monge-Ampère equation is solved using a quasi-Newton method (see [MCB18] for details).

This approach typically takes fewer than ten iterations to converge, but each iteration is relatively expensive.

Parameters:

  • mesh (firedrake.mesh.MeshGeometry) – the physical mesh

  • monitor_function (Callable) – a Python function which takes a mesh as input

  • phi_init (firedrake.function.Function) – initial guess for the scalar potential

  • H_init (firedrake.function.Function) – initial guess for the Hessian

  • maxiter (int) – maximum number of iterations for the Quasi-Newton solver

  • rtol (float) – relative tolerance for the residual

  • dtol (float) – divergence tolerance for the residual

property equidistributor

Setup the equidistributor for the quasi-newton method.

The equidistributor solves

\[\int_{\Omega} \boldsymbol{\tau} \cdot \mathbf{H} \, \mathrm{d}x + \int_{\Omega} (\nabla \cdot \boldsymbol{\tau}) \cdot (\nabla \phi) \, \mathrm{d}x - \int_{\partial \Omega} (((\nabla \phi) \cdot \widehat{\mathbf{n}}) \cdot \boldsymbol{\tau}) \cdot \widehat{\mathbf{n}} \, \mathrm{d}s - \int_{\Omega} \psi (m \det(\mathbf{I} + \mathbf{H}) - \theta) \, \mathrm{d}x = 0, \quad \forall \boldsymbol{\tau} \in \mathbb{P}1^{d \times d}, \quad \forall \psi \in \mathbb{P}1,\]

for the potential \(\phi\) and its Hessian \(\mathbf{H}\), where \(d\) is the spatial dimension and \(\widehat{\mathbf{n}}\) is a normal vector to the boundary.

Returns:

the equidistributor

Return type:

NonlinearVariationalSolver

move()[source]

Run the quasi-Newton method to convergence and update the mesh.

Returns:

the iteration count

Return type:

int

class MongeAmpereMover_Relaxation(mesh, monitor_function, phi_init=None, H_init=None, **kwargs)[source]

Bases: MongeAmpereMover_Base

The standard, elliptic form of the Monge-Ampère equation used for mesh movement is:

\[m(\mathbf{x})\det(\mathbf{I} + \mathbf{H}(\phi)) = \theta,\]

for a convex scalar potential \(\phi=\phi(\boldsymbol{\xi})\), which is a function of the coordinates of the computational mesh. Here \(m=m(\mathbf{x})\) is a user-provided monitor function, which is a function of the coordinates of the physical mesh. \(\mathbf{I}\) is the identity matrix, \(\theta\) is a normalisation coefficient and \(\mathbf{H}(\phi)\) denotes the Hessian of \(\phi\) with respect to \(\boldsymbol{\xi}\).

In this mesh mover, the Monge-Ampère equation is instead solved in a parabolised form using a pseudo-time relaxation,

\[-\frac\partial{\partial\tau}\Delta\phi = m(\mathbf{x})\det(\mathbf{I} + \mathbf{H}(\phi)) - \theta,\]

where \(\tau\) is the pseudo-time variable. Forward Euler is used for the pseudo-time integration (see [MCB18] for details).

This approach typically takes tens or hundreds of iterations to converge, but each iteration is relatively cheap.

Parameters:

  • mesh (firedrake.mesh.MeshGeometry) – the physical mesh

  • monitor_function (Callable) – a Python function which takes a mesh as input

  • phi_init (firedrake.function.Function) – initial guess for the scalar potential

  • H_init (firedrake.function.Function) – initial guess for the Hessian

  • pseudo_timestep (float) – pseudo-timestep to use for the relaxation

  • maxiter (int) – maximum number of iterations for the relaxation

  • rtol (float) – relative tolerance for the residual

  • dtol (float) – divergence tolerance for the residual

  • fixed_boundary_segments (list of str or int) – labels corresponding to boundary segments to be fixed with a zero Dirichlet condition. The ‘on_boundary’ label indicates the whole domain boundary

property equidistributor

Setup the equidistributor for the relaxation method.

The equidistributor solves the following equation:

\[\int_{\Omega} \tau : \mathbf{H} \, \mathrm{d}x = -\int_{\Omega} (\nabla \cdot \tau) \cdot (\nabla \phi) \, \mathrm{d}x + \int_{\partial \Omega} ((\nabla \phi \cdot \widehat{\mathbf{n}}) \cdot \tau) \cdot \widehat{\mathbf{n}} \, \mathrm{d}s, \quad \forall \tau \in \mathbb{P}1^{d \times d}\]

for the Hessian \(\mathbf{H}\), where \(d\) is the spatial dimension and \(\widehat{\mathbf{n}}\) is a normal vector to the boundary.

Returns:

the equidistributor

Return type:

LinearVariationalSolver

move()[source]

Run the relaxation method to convergence and update the mesh.

Returns:

the iteration count

Return type:

int

property pseudotimestepper

Setup the pseudo-timestepper for the relaxation method.

Forward Euler is used for the pseudo-time integration (see [MCB18] for details). The pseudo-timestep may be set through the pseudo_timestep keyword argument to the constructor.

Returns:

the pseudo-timestepper

Return type:

LinearVariationalSolver

movement.mover module

class PrimeMover(mesh, monitor_function=None, raise_convergence_errors=True, tangling_check=None, quadrature_degree=None)[source]

Bases: ABC

Base class for all mesh movers.

Parameters:

  • mesh (firedrake.mesh.MeshGeometry) – the physical mesh

  • monitor_function (Callable) – a Python function which takes a mesh as input

  • raise_convergence_errors (bool) – convergence error handling behaviour: if True then ConvergenceErrors are raised, else warnings are raised and the program is allowed to continue

  • tangling_check (bool) – check whether the mesh has tangled elements (by default on in the 2D case and off otherwise)

  • quadrature_degree (int) – quadrature degree to be passed to Firedrakes measures

property coefficient_of_variation
Returns:

the coefficient of variation (σ/μ) of element volumes.

Return type:

float

abstract move()[source]

Move the mesh according to the method of choice.

to_computational_coordinates()[source]

Switch coordinates to correspond to the computational mesh mathcal{H}_C.

to_physical_coordinates()[source]

Switch coordinates to correspond to the physical mesh mathcal{H}_P.

property volume_ratio
Returns:

the ratio of the largest and smallest element volumes.

Return type:

float

movement.solver_parameters module

movement.spring module

SpringMover(*args, method='lineal', **kwargs)[source]

Movement of a mesh is determined by reinterpreting it as a structure of stiff beams and solving an associated discrete linear elasticity problem. (See [] for details.)

Parameters:
class SpringMover_Lineal(mesh, timestep, **kwargs)[source]

Bases: SpringMover_Base

Movement of a mesh is determined by reinterpreting it as a structure of stiff beams and solving an associated discrete linear elasticity problem.

We consider the ‘lineal’ case, as described in [].

Parameters:
move(time, update_boundary_displacement=None, boundary_conditions=None)[source]

Assemble and solve the lineal spring system and update the coordinates.

Parameters:

  • time (float) – the current time

  • update_boundary_displacement (Callable with a single argument of float type) – function that updates the boundary conditions at the current time

  • boundary_conditions (DirichletBC or list thereof) – Dirichlet boundary conditions to be enforced

class SpringMover_Torsional(*args, **kwargs)[source]

Bases: SpringMover_Lineal

Movement of a mesh is determined by reinterpreting it as a structure of stiff beams and solving an associated discrete linear elasticity problem.

We consider the ‘torsional’ case, as described in [].

Parameters:

movement.tangling module

class MeshTanglingChecker(mesh, raise_error=True)[source]

Bases: object

A class for tracking whether a mesh has tangled, i.e. whether any of its elements have become inverted.

Parameters:

  • mesh (firedrake.mesh.MeshGeometry) – the mesh to track if tangled

  • raise_error (bool) – if True, an error is raised if any element is tangled, otherwise a warning is raised

check()[source]

Check whether any element orientations have changed since the tangling checker was created.

property scaled_jacobian

Module contents