Mesh and other FEM details

Base class

class phienics.mesh.Mesh(rs=1.0, r_min=0.0, r_max=1000000000.0, num_cells=500, too_fine=10000.0, Ntol=1e-08, linear_refine=0, linear_start=None, linear_stop=None, r_rm=None, A_rm=25.0, k_rm=20.0, adjust_transition=False)

Base class for mesh classes.

Arguments

rs

source radius (typically, 1.)

r_min

minimum radius (typically, 0.)

r_max

maximum radius box size (i.e. the simulation’s ‘infinity’). It must be large compared to the source radius and the Compton wavelength of the fields involved

num_cells

number of cells in the mesh

too_fine

this parameter is used to make sure points are not closer than machine precision allows to resolve. Points are required to be spaced at least \(\mathrm{too\_fine}\times\mathrm{DOLFIN\_EPS}\). See check_mesh()

linear_refine

after producing a non-uniform mesh, optionally slice cells in the [linear_start, linear_stop] interval this many times (use 0 to skip this step)

linear_start

refine linearly from here

linear_stop

refine linearly up to here

linear_cells

if linear_refine was chosen, how many cells were added by linear refinement

r_rm

optionally decluster points around this radius, using a transformation \(r = A_{\rm rm}/2 \arctan{\left( k_rm x \right)}\) (see section 3.4.3 of the accompanying paper )

A_rm

transformation parameter; only valid if r_rm is a valid radiusx

k_rm

transformation parameter; only valid if r_rm is a valid radius

adjust_transition

if True, radial coordinates are rescaled so that the minimum mesh spacing is placed exactly at \(r_s\) - this may be useful when working with a step or top-hat source. Otherwise, the spacing is determined by the input parameters alone. Default: False.

Inside a notebook, you can check the mesh parameters by using:

yourmesh.__dict__

Note

All distances are expressed in units of Mn as defined in UV.UVSolver or IR.IRSolver. For the (most sensible) choice of \(M_n={R_s}^{-1}\), distances are in units of the source radius.

Note

use utils.get_values() to output function values at mesh vertices, which takes care of this FEniCS bug: https://www.allanswered.com/post/jgvjw/fenics-plot-incorrectly-orders-coordinates-for-function-on-locally-refined-1d-mesh/.

apply_linear_refinement()

Refine the mesh linearly: slice all cells between a radius linear_start and a radius linear_stop in half, for a number of times specified by linear_refine.

build_mesh()

Build the mesh:

1. check that the input parameters for the mesh are valid (the function performing the test must be defined in the child mesh class)

2. set the mesh transformation. Depending on whether r_rm=None or r_rm=float, the transformation will be the baseline transformation (as defined in the child mesh class) or its updated form declustering mesh points at r_rm

3. create a starting linear mesh, in particular, compute its extrema by using the Newton method

4. obtain a nonlinear mesh by applying the transformation defined at (2)

5. apply additional linear refinement if required

6. check that the obtained mesh points are not closer than the user-set tolerance of \(\mathrm{too\_fine}\times\mathrm{DOLFIN\_EPS}\).

check_mesh()

Check that all points in the mesh are separated enough that machine precision can tell them apart. Instead of sheer machine precision, I employ a safer threshold too_fine \(\times\) DOLFIN_EPS, where DOLFIN_EPS is machine precision. The quantity too_fine is a mesh attribute.

transform_with_declustering()

Declusters points at a specified radius r_rm. See Sec. 3.4.3 of the accompanying paper for details.

In order to enable this feature in custom mesh classes, the baseline trasformation, first and second derivatives must be defined in the child mesh classes.

ArcTanExpMesh

class phienics.atexpmesh.ArcTanExpMesh(rs=1.0, r_min=0.0, r_max=1000000000.0, num_cells=500, too_fine=10000.0, Ntol=1e-08, linear_refine=0, linear_start=None, linear_stop=None, r_rm=None, A_rm=25.0, k_rm=20.0, adjust_transition=False, k=25.0, a=0.05, b=0.03)

This class implements the Arctan-exp mesh described in detail at Sec. 3.4.1 of the accompanying paper . This mesh applies the transformation:

\[r = T(x) = \frac{2}{\pi} r_{\rm s} \arctan{(kx)} \exp{(a x^3 + bx)}\]

to obtain a nonlinear radial mesh (\(r\) coordinate) starting from a uniform mesh (\(x\) coordinate).

Arguments

rs

source radius (typically, 1.)

r_min

minimum radius (typically, 0.)

r_max

maximum radius box size (i.e. the simulation’s ‘infinity’). It must be large compared to the source radius and the Compton wavelength of the fields involved

num_cells

number of cells in the mesh

too_fine

this parameter is used to make sure points are not closer than machine precision allows to resolve. Points are required to be spaced at least \(\mathrm{too\_fine}\times\mathrm{DOLFIN\_EPS}\). See check_mesh()

linear_refine

after producing a non-uniform mesh, optionally slice cells in the [linear_start, linear_stop] interval this many times (use 0 to skip this step)

linear_start

refine linearly from here

linear_stop

refine linearly up to here

linear_cells

if linear_refine was chosen, how many cells were added by linear refinement

r_rm

optionally remove points around this radius using a transformation \(r = A_{\rm rm}/2 \arctan{\left( k_rm x \right)}\) (see section 3.4.3 of the accompanying paper)

A_rm

transformation parameter; only valid if r_rm is a valid radiusx

k_rm

transformation parameter; only valid if r_rm is a valid radius

k

parameter \(k\) in the \(T(x)\) transformation

a

parameter \(a\) in the \(T(x)\) transformation

b

parameter \(b\) in the \(T(x)\) transformation

baseline_transform()

Defines the baseline trasformation

\[r = T(x) = \frac{2}{\pi} r_{\rm s} \arctan{(kx)} \exp{(a x^3 + bx)}\]

together with its first and second derivatives, and an approximation for the inverse transformation at small and large radii.

check_params()

Check that the input parameters \(a,b,k\) make sense: \(a\) and \(b\) must be positive or zero, but not simultaneously zero; \(k\) must be strictly positive.

ArcTanPowerLawMesh

class phienics.atplmesh.ArcTanPowerLawMesh(rs=1.0, r_min=0.0, r_max=1000000000.0, num_cells=2000, too_fine=10000.0, Ntol=1e-08, linear_refine=0, linear_start=None, linear_stop=None, r_rm=None, A_rm=25.0, k_rm=20.0, adjust_transition=False, k=25.0, gamma=6.5)

This class implements the Arctan-power law mesh described in detail at Sec. 3.4.2 of the accompanying paper . This mesh applies the transformation:

\[r = T(x) = \frac{2}{\pi} r_{\rm s} \arctan{(kx)} + x^{\gamma}\]

to obtain a nonlinear radial mesh (\(r\) coordinate) starting from a uniform mesh (\(x\) coordinate).

Arguments

rs

source radius (typically, 1.)

r_min

minimum radius (typically, 0.)

r_max

maximum radius box size (i.e. the simulation’s ‘infinity’). It must be large compared to the source radius and the Compton wavelength of the fields involved

num_cells

number of cells in the mesh

too_fine

this parameter is used to make sure points are not closer than machine precision allows to resolve. Points are required to be spaced at least \(\mathrm{too\_fine}\times\mathrm{DOLFIN\_EPS}\). See check_mesh()

linear_refine

after producing a non-uniform mesh, optionally slice cells in the [linear_start, linear_stop] interval this many times (use 0 to skip this step)

linear_start

refine linearly from here

linear_stop

refine linearly up to here

linear_cells

if linear_refine was chosen, how many cells were added by linear refinement

r_rm

optionally remove points around this radius using a transformation \(r = A_{\rm rm}/2 \arctan{\left( k_rm x \right)}\) (see section 3.4.3 of the accompanying paper )

A_rm

transformation parameter; only valid if r_rm is a valid radiusx

k_rm

transformation parameter; only valid if r_rm is a valid radius

k

parameter \(k\) in the \(T(x)\) transformation

gamma

parameter \(gamma\) in the \(T(x)\) transformation

baseline_transform()

Defines the baseline trasformation

\[r = T(x) = \frac{2}{\pi} r_{\rm s} \arctan{(kx)} + x^{\gamma}\]

together with its first and second derivatives, and an approximation for the inverse transformation at small and large radii.

check_params()

Check that the input parameters make sense: \(k\) and \(\gamma\) must satisfy \(k>0\) and \(\gamma \geq 1\).

Other FEM details

class phienics.fem.Fem(mesh, func_cont='CG', func_disc='DG', func_degree=4)

Base class with code-wide settings on the function spaces.

For further details (in particular, when continuous vs discontinuous function spaces are employed), see the accompanying paper .

Arguments

mesh

a mesh.Mesh object

func_cont

choice of interpolating function family for continuous function spaces. Default: ‘CG’, i.e. continuous Lagrange

func_disc

choice of interpolating function family for discontinuous function spaces. Default: ‘DG’, i.e. discontinuous Lagrange

func_degree

degree of interpolating polynomials