helpers

helpers

Small helper functions used throughout the library

Functions:

Name	Description
mesh_from_file	Reads a mesh from a file to be used in dolfinx.
rescale_mesh	copies and rescales/shifts a mesh
rescale_mesh_in_place	rescales/shifts a mesh in place.
create_periodic_boundary_conditions	creates periodic boundary conditions on a 2D box mesh.

hommx.helpers.PoissonFEM

solves the poisson equation

__init__

__init__(msh: mesh.Mesh, A: Callable[[ufl.SpatialCoordinate], ufl.Form], f: ufl.form)

Initializes the solver, meshes and boundary condtions.

Parameters:

Name	Type	Description	Default
msh	mesh.Mesh	The macro mesh, on which we want to solve the oscillatory Poisson equation. This mesh can live on MPI.COMM_WORLD and the cell problems are automatically solved on each process in parallel.	required
A	Callable[[ufl.SpatialCoordinate], ufl.Form]	The coefficient, it should be callable like: A(x)(y), where \(x\) is a spatial coordinate on the macro mesh (the cell center \(c_T\)) and \(y\) is a ufl.SpatialCoordinate on the microscopic mesh, that is passed to dolfinx to solve the cell problem. A needs to be 1-periodic in y, at least for the theory to work.	required
f	ufl.form	The right hand side of the Poisson problem.	required
msh_micro	mesh.Mesh	The microscopic mesh, this needs to be the unit-square. Further it needs to live on MPI.COMM_SELF, since every process owns a whole copy of the microscopic mesh. If any other communicator but MPI.COMM_SELF is used the results will most likely be rubish.	required
eps		\(\varepsilon\), the microscopic scaling. Note that this needs to be small enough, so that the cells live entirely within their corresponding element. If this is not the case, results may be rubish.	required
petsc_options_global_solve	optional	PETSc solver options for the global solver, see	required
petsc_options_cell_problem	optional	PETSc solver options for the global solver, see	required
petsc_options_prefix	optional	options prefix used for PETSc options. Defaults to "hommx_PoissonHMM".	required

solve

solve() -> fem.Function

Assemble the LHS, RHS and solve the problem

This method assembles the HMM stiffness matrix, so depending on the problem it might run for some time.

plot_solution

plot_solution(u: fem.Function | None = None)

Simple plot of the solution using pyvista.

Solve needs to be run before calling this. On parallel methods each process only plots the local part.

hommx.helpers.solve_diffusion_1d

solve_diffusion_1d(epsilon: float, nx: int, A_callable: callable) -> np.ndarray[tuple[int], np.dtype[float]]

Solves diffusion equation with multiscale coefficient. Solves the diffusion equation:

\[ \mathrm{div}(A \\nabla u) = f \]

Where \(A\) is of the shape:

\[ A(x) = 1 + 0.5*\\sin(x\\frac{2*\pi}{\\varepsilon}) \]

Parameters:

Name	Type	Description	Default
epsilon	float	\(\\varepsilon\)	required
nx	int	discretization of the FEM solution	required
A_callable	callable	callable, such that A_callable(x) returns a proper ufl.Form where x are the Spatial Coordinates	required

hommx.helpers.solve_diffusion_2d

solve_diffusion_2d(epsilon: float, nx: int, A_callable: callable) -> np.ndarray[tuple[int], np.dtype[float]]

Solves diffusion equation with multiscale coefficient. Solves the diffusion equation:

\[ \mathrm{div}(A \\nabla u) = f \]

Where \(A\) is of the shape:

\[ A(x) = 1 + 0.5*\\sin(x\\frac{2*\pi}{\\varepsilon}) \]

Parameters:

Name	Type	Description	Default
epsilon	float	\(\\varepsilon\)	required
nx	int	discretization of the FEM solution	required
A_callable	callable	callable, such that A_callable(x) returns a proper ufl.Form where x are the Spatial Coordinates	required

hommx.helpers.mesh_from_file

mesh_from_file(filename)

hommx.helpers.mesh_from_delaunay

mesh_from_delaunay(points: np.ndarray[tuple[int, int], np.dtype[float]], triangles: np.ndarray[tuple[int, int], np.dtype[int]]) -> mesh.Mesh

Creates mesh from list of points and triangles.

Parameters:

Name	Type	Description	Default
points	np.ndarray[tuple[int, int], np.dtype[float]]	(N_points, 2) shape array points on the plane.	required
triangles	np.ndarray[tuple[int, int], np.dtype[int]]	(N_triangles, 3) shape array containing indices into the points array where each row represents one triangle.	required

Notes

Setting up the unit square could look something like this:

from scipy.spatial import Delaunay
import numpy as np

x = np.linspace(0, 1, 10)
y = np.linspace(0, 1, 10)
X, Y = np.meshgrid(x, y)
points - np.stack([X, Y], axis=-1).reshape(-1, 2)
triangles = Delaunay(points).simplices

hommx.helpers.rescale_mesh

rescale_mesh(msh: mesh.Mesh, scale: float = 1, shift: np.ndarray[tuple[int], np.dtype[float]] = np.array([0, 0, 0])) -> mesh.Mesh

Creates a rescaled and shifted copy of the mesh.

Parameters:

Name	Type	Description	Default
msh	mesh.Mesh	Mesh	required
scale	float	constant scalign factor	1
shift	np.ndarray[tuple[int], np.dtype[float]]	(3,) shape array containing the shift.	np.array([0, 0, 0])

For example consider the unit-square \([0,1]^2\). This is remapped to \(\text{scale}*[0, 1]^2 + \text{shift} = [\text{shift}[0], \text{scale}+\text{shift}[1]]^2\). This is used to map between the Unit cell \(Y=[0,1]^d\) and \(Y_\varepsilon\) in the HMM.

hommx.helpers.rescale_mesh_in_place

rescale_mesh_in_place(msh: mesh.Mesh, scale: float = 1, shift: np.ndarray[tuple[int], np.dtype[float]] = np.array([0, 0, 0])) -> mesh.Mesh

Rescales and shifts the mesh in place.

Parameters:

Name	Type	Description	Default
msh	mesh.Mesh	Mesh	required
scale	float	constant scalign factor	1
shift	np.ndarray[tuple[int], np.dtype[float]]	(3,) shape array containing the shift.	np.array([0, 0, 0])

For example consider the unit-square \([0,1]^2\). This is remapped to \(\text{scale}*[0, 1]^2 + \text{shift} = [\text{shift}[0], \text{scale}+\text{shift}[1]]^2\). This is used to map between the Unit cell \(Y=[0,1]^d\) and \(Y_\varepsilon\) in the HMM.

hommx.helpers.plot_fem_function

plot_fem_function(V: fem.FunctionSpace, u: fem.Function)