pyamg.gallery#
Matrix gallery of model problems.
Functions#
poisson() : Poisson problem using Finite Differences
linear_elasticity() : Linear Elasticity using Finite Elements
stencil_grid() : General stencil generation from 1D, 2D, and 3D
diffusion_stencil_2d() : 2D rotated anisotropic FE/FD stencil
- pyamg.gallery.advection_2d(grid, theta=0.7853981633974483, l_bdry=1.0, b_bdry=1.0)[source]#
Generate matrix and RHS for upwind FD discretization of 2D advection.
(cos(theta),sin(theta)) dot grad(u) = 0,
with inflow boundaries on the left and bottom of the domain. Assume uniform grid spacing, dx=dy, even for grid[0] != grid[1].
- Parameters:
- gridtuple
(ny, nx) number of points in y and x
- thetafloat, optional
Rotation angle theta in radians defines direction of advection (cos(theta),sin(theta))
- l_bdryfloat, array
left boundary value. If float, then constant in-flow boundary value applied. If array, then length of array must be equal to ny = grid[0], and this array defines non-constant boundary value on the left.
- b_bdryfloat, array
bottom boundary value. If float, then constant in-flow boundary value applied. If array, then length of array must be equal to nx = grid[1], and this array defines non-constant boundary value on the bottom.
- Returns:
- Acsr_matrix
Defines 2D FD upwind discretization, with boundary
- rhsarray
Defines right-hand-side with boundary contributions
See also
Examples
>>> from numpy import pi >>> from pyamg.gallery import advection_2d >>> A, rhs = advection_2d( (4,4), theta=pi/4) >>> print(A.toarray().round(4)) [[ 1.4142 0. 0. -0.7071 0. 0. 0. 0. 0. ] [-0.7071 1.4142 0. 0. -0.7071 0. 0. 0. 0. ] [ 0. -0.7071 1.4142 0. 0. -0.7071 0. 0. 0. ] [ 0. 0. 0. 1.4142 0. 0. -0.7071 0. 0. ] [ 0. 0. 0. -0.7071 1.4142 0. 0. -0.7071 0. ] [ 0. 0. 0. 0. -0.7071 1.4142 0. 0. -0.7071] [ 0. 0. 0. 0. 0. 0. 1.4142 0. 0. ] [ 0. 0. 0. 0. 0. 0. -0.7071 1.4142 0. ] [ 0. 0. 0. 0. 0. 0. 0. -0.7071 1.4142]]
- pyamg.gallery.diffusion_stencil_2d(epsilon=1.0, theta=0.0, type='FE')[source]#
Rotated Anisotropic diffusion in 2d of the form.
-div Q A Q^T grad u
- Q = [cos(theta) -sin(theta)]
[sin(theta) cos(theta)]
- A = [1 0 ]
[0 eps ]
- Parameters:
- epsilonfloat, optional
Anisotropic diffusion coefficient: -div A grad u, where A = [1 0; 0 epsilon]. The default is isotropic, epsilon=1.0
- thetafloat, optional
Rotation angle theta in radians defines -div Q A Q^T grad, where Q = [cos(theta) -sin(theta); sin(theta) cos(theta)].
- type{‘FE’,’FD’}
Specifies the discretization as Q1 finite element (FE) or 2nd order finite difference (FD) The default is theta = 0.0
- Returns:
- stencilnumpy array
A 3x3 diffusion stencil
See also
Notes
Not all combinations are supported.
Examples
>>> import scipy as sp >>> from pyamg.gallery.diffusion import diffusion_stencil_2d >>> sten = diffusion_stencil_2d(epsilon=0.0001,theta=sp.pi/6,type='FD') >>> print(sten) [[-0.2164847 -0.750025 0.2164847] [-0.250075 2.0002 -0.250075 ] [ 0.2164847 -0.750025 -0.2164847]]
Consider a 2 x 4 grid ([x0, x1] x [y0, y1, y2, y3]). The first dimension of the stencil defines x.
>>> nx, ny = (2, 4) >>> sten = pyamg.gallery.diffusion_stencil_2d(epsilon=0.1, type='FD') >>> A = pyamg.gallery.stencil_grid(sten, (nx, ny)).toarray() >>> print(sten) [[-0. -1. 0. ] [-0.1 2.2 -0.1] [ 0. -1. -0. ]] >>> print(A) [[ 2.2 -0.1 0. 0. -1. 0. 0. 0. ] [-0.1 2.2 -0.1 0. 0. -1. 0. 0. ] [ 0. -0.1 2.2 -0.1 0. 0. -1. 0. ] [ 0. 0. -0.1 2.2 0. 0. 0. -1. ] [-1. 0. 0. 0. 2.2 -0.1 0. 0. ] [ 0. -1. 0. 0. -0.1 2.2 -0.1 0. ] [ 0. 0. -1. 0. 0. -0.1 2.2 -0.1] [ 0. 0. 0. -1. 0. 0. -0.1 2.2]]
- pyamg.gallery.gauge_laplacian(npts, spacing=1.0, beta=0.1)[source]#
Construct a Gauge Laplacian from Quantum Chromodynamics for regular 2D grids.
Note that this function is not written efficiently, but should be fine for N x N grids where N is in the low hundreds.
- Parameters:
- nptsint
number of pts in x and y directions
- spacingfloat
grid spacing between points
- betafloat
temperature Note that if beta=0, then we get the typical 5pt Laplacian stencil
- Returns:
- Acsr matrix
A is Hermitian positive definite for beta > 0.0 A is Symmetric semi-definite for beta = 0.0
References
[1]MacLachlan, S. and Oosterlee, C., “Algebraic Multigrid Solvers for Complex-Valued Matrices”, Vol. 30, SIAM J. Sci. Comp, 2008
Examples
>>> from pyamg.gallery import gauge_laplacian >>> A = gauge_laplacian(10)
- pyamg.gallery.linear_elasticity(grid, spacing=None, E=100000.0, nu=0.3, format=None)[source]#
Linear elasticity problem with Q1 finite elements on a regular rectangular grid.
- Parameters:
- gridtuple
length 2 tuple of grid sizes, e.g. (10, 10)
- spacingtuple
length 2 tuple of grid spacings, e.g. (1.0, 0.1)
- Efloat
Young’s modulus
- nufloat
Poisson’s ratio
- formatstring
Format of the returned sparse matrix (eg. ‘csr’, ‘bsr’, etc.)
- Returns:
- Acsr_matrix
FE Q1 stiffness matrix
- Barray
rigid body modes
See also
Notes
only 2d for now
References
[1]J. Alberty, C. Carstensen, S. A. Funken, and R. KloseDOI “Matlab implementation of the finite element method in elasticity” Computing, Volume 69, Issue 3 (November 2002) Pages: 239 - 263 http://www.math.hu-berlin.de/~cc/
Examples
>>> from pyamg.gallery import linear_elasticity >>> A, B = linear_elasticity((4, 4))
- pyamg.gallery.linear_elasticity_p1(vertices, elements, E=100000.0, nu=0.3, format=None)[source]#
P1 elements in 2 or 3 dimensions.
- Parameters:
- verticesarray_like
array of vertices of a triangle or tets
- elementsarray_like
array of vertex indices for tri or tet elements
- Efloat
Young’s modulus
- nufloat
Poisson’s ratio
- formatstring
‘csr’, ‘csc’, ‘coo’, ‘bsr’
- Returns:
- Acsr_matrix
FE Q1 stiffness matrix
Notes
works in both 2d and in 3d
References
[1]J. Alberty, C. Carstensen, S. A. Funken, and R. KloseDOI “Matlab implementation of the finite element method in elasticity” Computing, Volume 69, Issue 3 (November 2002) Pages: 239 - 263 http://www.math.hu-berlin.de/~cc/
Examples
>>> from pyamg.gallery import linear_elasticity_p1 >>> import numpy as np >>> E = np.array([[0, 1, 2],[1, 3, 2]]) >>> V = np.array([[0.0, 0.0],[1.0, 0.0],[0.0, 1.0],[1.0, 1.0]]) >>> A, B = linear_elasticity_p1(V, E)
- pyamg.gallery.load_example(name)[source]#
Load an example problem by name.
- Parameters:
- namestring (e.g. ‘airfoil’)
Name of the example to load
Notes
- Each example is stored in a dictionary with the following keys:
‘A’ : sparse matrix
‘B’ : near-nullspace candidates
‘vertices’ : dense array of nodal coordinates
‘elements’ : dense array of element indices
- Current example names are:
airfoil bar helmholtz_2D knot local_disc_galerkin_diffusion recirc_flow unit_cube unit_square
Examples
>>> from pyamg.gallery import load_example >>> ex = load_example('knot')
- pyamg.gallery.poisson(grid, dtype=<class 'float'>, format=None, type='FD')[source]#
Return a sparse matrix for the N-dimensional Poisson problem.
The matrix represents a finite Difference approximation to the Poisson problem on a regular n-dimensional grid with unit grid spacing and Dirichlet boundary conditions.
The last dimension is iterated over first: z, then y, then x. This should be used with np.mgrid() or np.ndenumerate.
- Parameters:
- gridtuple of integers
grid dimensions e.g. (100,100)
Notes
The matrix is symmetric and positive definite (SPD).
Examples
>>> from pyamg.gallery import poisson >>> # 4 nodes in one dimension >>> poisson((4,)).toarray() array([[ 2., -1., 0., 0.], [-1., 2., -1., 0.], [ 0., -1., 2., -1.], [ 0., 0., -1., 2.]])
>>> # rectangular two dimensional grid >>> poisson((2,3)).toarray() array([[ 4., -1., 0., -1., 0., 0.], [-1., 4., -1., 0., -1., 0.], [ 0., -1., 4., 0., 0., -1.], [-1., 0., 0., 4., -1., 0.], [ 0., -1., 0., -1., 4., -1.], [ 0., 0., -1., 0., -1., 4.]])
- pyamg.gallery.regular_triangle_mesh(nx, ny)[source]#
Construct a regular triangular mesh in the unit square.
- Parameters:
- nxint
Number of nodes in the x-direction
- nyint
Number of nodes in the y-direction
- Returns:
- Vertarray
nx*ny x 2 vertex list
- E2Varray
Nex x 3 element list
Examples
>>> from pyamg.gallery import regular_triangle_mesh >>> E2V,Vert = regular_triangle_mesh(3, 2)
- pyamg.gallery.sprand(m, n, density, format='csr')[source]#
Return a random sparse matrix.
- Parameters:
- m, nint
shape of the result
- densityfloat
target a matrix with nnz(A) = m*n*density, 0<=density<=1
- formatstring
sparse matrix format to return, e.g. ‘csr’, ‘coo’, etc.
- Returns:
- Asparse matrix
m x n sparse matrix
Examples
>>> from pyamg.gallery import sprand >>> A = sprand(5,5,3/5.0)
- pyamg.gallery.stencil_grid(S, grid, dtype=None, format=None)[source]#
Construct a sparse matrix form a local matrix stencil.
- Parameters:
- Sndarray
matrix stencil stored in N-d array
- gridtuple
tuple containing the N grid dimensions
- dtype
data type of the result
- formatstring
sparse matrix format to return, e.g. “csr”, “coo”, etc.
- Returns:
- Asparse matrix
Sparse matrix which represents the operator given by applying stencil S at each vertex of a regular grid with given dimensions.
Notes
The grid vertices are enumerated as arange(prod(grid)).reshape(grid). This implies that the last grid dimension cycles fastest, while the first dimension cycles slowest. For example, if grid=(2,3) then the grid vertices are ordered as (0,0), (0,1), (0,2), (1,0), (1,1), (1,2).
This coincides with the ordering used by the NumPy functions ndenumerate() and mgrid().
Examples
>>> from pyamg.gallery import stencil_grid >>> stencil = [-1,2,-1] # 1D Poisson stencil >>> grid = (5,) # 1D grid with 5 vertices >>> A = stencil_grid(stencil, grid, dtype=float, format='csr') >>> A.toarray() array([[ 2., -1., 0., 0., 0.], [-1., 2., -1., 0., 0.], [ 0., -1., 2., -1., 0.], [ 0., 0., -1., 2., -1.], [ 0., 0., 0., -1., 2.]])
>>> stencil = [[0,-1,0],[-1,4,-1],[0,-1,0]] # 2D Poisson stencil >>> grid = (3,3) # 2D grid with shape 3x3 >>> A = stencil_grid(stencil, grid, dtype=float, format='csr') >>> A.toarray() array([[ 4., -1., 0., -1., 0., 0., 0., 0., 0.], [-1., 4., -1., 0., -1., 0., 0., 0., 0.], [ 0., -1., 4., 0., 0., -1., 0., 0., 0.], [-1., 0., 0., 4., -1., 0., -1., 0., 0.], [ 0., -1., 0., -1., 4., -1., 0., -1., 0.], [ 0., 0., -1., 0., -1., 4., 0., 0., -1.], [ 0., 0., 0., -1., 0., 0., 4., -1., 0.], [ 0., 0., 0., 0., -1., 0., -1., 4., -1.], [ 0., 0., 0., 0., 0., -1., 0., -1., 4.]])