lapy.Solver¶

class lapy.Solver(geometry, lump=False, aniso=None, aniso_smooth=10, use_cholmod=False)[source]¶

FEM solver for Laplace Eigenvalue and Poisson Equation.

Inputs can be geometry classes which have vertices and elements. Currently TriaMesh and TetMesh are implemented. FEM matrices (stiffness (or A) and mass matrix (or B)) are computed during the construction. After that the Eigenvalue solver (lapy.Solver.eigs) or Poisson Solver (lapy.Solver.poisson) can be called.

Parameters:

geometryTriaMesh | TetMesh: Mesh geometry.
lumpbool: If True, lump the mass matrix (diagonal).
anisofloat | tuple of shape (2,): Anisotropy for curvature based anisotopic Laplace. If a tuple (a_min, a_max), differentially affects the minimum and maximum curvature directions. e.g. ``(0, 50) will set scaling to 1 into the minimum curvature direction, even if the maximum curvature is large in those regions ( i.e. isotropic in regions with large maximum curvature and minimum curvature close to 0, i.e. a concave cylinder).
aniso_smoothint | None: Number of smoothing iterations for curvature computation on vertices.
use_cholmodbool, default: False: If True, attempts to use the Cholesky decomposition for improved execution speed. Requires the scikit-sparse library. If it can not be found, an error will be thrown. If False, will use slower LU decomposition.

Methods

`eigs`([k])	Compute the linear finite-element method Laplace-Beltrami spectrum.
`fem_tria_mass`(tria[, lump])	Compute the sparse symmetric mass matrix of the Laplace-Beltrami operator for a given triangle mesh.
`poisson`([h, dtup, ntup])	Solver for the Poisson equation with boundary conditions.

Notes

The class has a static member to create the mass matrix of TriaMesh for external function that do not need stiffness.

eigs(k=10)[source]¶

Compute the linear finite-element method Laplace-Beltrami spectrum.

Parameters:

kint, default=10: The number of eigenvalues and eigenvectors desired. k must be smaller than N (number of vertices). It is not possible to compute all eigenvectors of a matrix.

Returns:

eigenvaluesnp.ndarray: Array of k eigenvalues, shape (k,). For closed meshes or Neumann boundary condition, 0 will be the first eigenvalue (with constant eigenvector).
eigenvectorsnp.ndarray: Array representing the k eigenvectors, shape (n_vertices, k). The column eigenvectors[:, i] is the eigenvector corresponding to eigenvalues[i].

static fem_tria_mass(tria, lump=False)[source]¶

Compute the sparse symmetric mass matrix of the Laplace-Beltrami operator for a given triangle mesh.

The sparse symmetric matrix is computed for a given triangle mesh using the linear finite element method (assuming a closed mesh or Neumann boundary condition). This function is faster than the constructor above and can be used when only a mass matrix is needed.

Parameters:

triaTriaMesh: Triangle mesh.
lumpbool, default=False: If True, B should be lumped (diagonal).

Returns:

scipy.sparse.csc_matrix: Sparse symmetric positive definite mass matrix B, shape (n, n).

Notes

This only returns the mass matrix B of the Eigenvalue problem: A x = lambda B x. The area of the surface mesh can be obtained via B.sum().

poisson(h=0.0, dtup=(), ntup=())[source]¶

Solver for the Poisson equation with boundary conditions.

This solver is based on the A and B Laplace matrices where A x = B h and A is a sparse symmetric positive semi-definite matrix of shape (n, n) and B is a sparse symmetric positive definite matrix of shape (n, n).

Parameters:

hfloat or np.ndarray, default=0.0: Right hand side, can be constant or array with vertex values. The default 0.0 corresponds to Laplace equation A x = 0.
dtuptuple, default=(): Dirichlet boundary condition as a tuple containing the index and data arrays of same length. The default, an empty tuple, corresponds to no Dirichlet condition.
ntuptuple, default=(): Neumann boundary condition as a tuple containing the index and data arrays of same length. The default, an empty tuple, corresponds to Neumann on all boundaries.

Returns:

np.ndarray: Array with vertex values of the solution, shape (n_vertices,).

Raises:

ValueError: If input matrices are not square or have different dimensions. If h is not scalar or array matching matrix dimensions. If dtup or ntup don’t contain exactly 2 arrays. If dtup indices are not unique. If dtup or ntup arrays have mismatched lengths. If matrix format is not CSC or CSR.

Notes

A and B are obtained via computeAB for either triangle or tetrahedral mesh.