multigrid package

This is the pyro multigrid solver. THere are several versions.

MG implements a second-order discretization of a constant-coefficient Helmholtz equation:

(alpha - beta L) phi = f

variable_coeff_MG implements a variable-coefficient Poisson equation:

div { eta grad phi } = f

general_MG implements a more general elliptic equation:

alpha phi + div { beta grad phi } + gamma . grad phi = f

All use pure V-cycles to solve elliptic problems

Submodules

multigrid.MG module

The multigrid module provides a framework for solving elliptic problems. A multigrid object is just a list of grids, from the finest mesh down (by factors of two) to a single interior zone (each grid has the same number of guardcells).

The main multigrid class is setup to solve a constant-coefficient Helmholtz equation:

(alpha - beta L) phi = f

where L is the Laplacian and alpha and beta are constants. If alpha = 0 and beta = -1, then this is the Poisson equation.

We support Dirichlet or Neumann BCs, or a periodic domain.

The general usage is as follows:

a = multigrid.CellCenterMG2d(nx, ny, verbose=1, alpha=alpha, beta=beta)

this creates the multigrid object a, with a finest grid of nx by ny zones and the default boundary condition types. alpha and beta are the coefficients of the Helmholtz equation. Setting verbose = 1 causing debugging information to be output, so you can see the residual errors in each of the V-cycles.

Initialization is done as:

a.init_zeros()

this initializes the solution vector with zeros (this is not necessary if you just created the multigrid object, but it can be used to reset the solution between runs on the same object).

Next:

a.init_RHS(zeros((nx, ny), numpy.float64))

this initializes the RHS on the finest grid to 0 (Laplace’s equation). Any RHS can be set by passing through an array of (nx, ny) values here.

Then to solve, you just do:

a.solve(rtol = 1.e-10)

where rtol is the desired tolerance (residual norm / source norm)

to access the final solution, use the get_solution method:

v = a.get_solution()

For convenience, the grid information on the solution level is available as attributes to the class,

a.ilo, a.ihi, a.jlo, a.jhi are the indices bounding the interior of the solution array (i.e. excluding the ghost cells).

a.x and a.y are the coordinate arrays a.dx and a.dy are the grid spacings

class multigrid.MG.CellCenterMG2d(nx, ny, ng=1, xmin=0.0, xmax=1.0, ymin=0.0, ymax=1.0, xl_BC_type='dirichlet', xr_BC_type='dirichlet', yl_BC_type='dirichlet', yr_BC_type='dirichlet', xl_BC=None, xr_BC=None, yl_BC=None, yr_BC=None, alpha=0.0, beta=-1.0, nsmooth=10, nsmooth_bottom=50, verbose=0, aux_field=None, aux_bc=None, true_function=None, vis=0, vis_title='')[source]

Bases: object

The main multigrid class for cell-centered data.

We require that nx = ny be a power of 2 and dx = dy, for simplicity

get_solution(grid=None)[source]

Return the solution after doing the MG solve

If a grid object is passed in, then the solution is put on that grid – not the passed in grid must have the same dx and dy

Returns:out : ndarray
get_solution_gradient(grid=None)[source]

Return the gradient of the solution after doing the MG solve. The x- and y-components are returned in separate arrays.

If a grid object is passed in, then the gradient is computed on that grid. Note: the passed-in grid must have the same dx, dy

Returns:out : ndarray, ndarray
get_solution_object()[source]

Return the full solution data object at the finest resolution after doing the MG solve

Returns:out : CellCenterData2d object
grid_info(level, indent=0)[source]

Report simple grid information

init_RHS(data)[source]

Initialize the right hand side, f, of the Helmholtz equation (alpha - beta L) phi = f

Parameters:

data : ndarray

An array (of the same size as the finest MG level) with the values to initialize the solution to the elliptic problem.

init_solution(data)[source]

Initialize the solution to the elliptic problem by passing in a value for all defined zones

Parameters:

data : ndarray

An array (of the same size as the finest MG level) with the values to initialize the solution to the elliptic problem.

init_zeros()[source]

Set the initial solution to zero

smooth(level, nsmooth)[source]

Use red-black Gauss-Seidel iterations to smooth the solution at a given level. This is used at each stage of the V-cycle (up and down) in the MG solution, but it can also be called directly to solve the elliptic problem (although it will take many more iterations).

Parameters:

level : int

The level in the MG hierarchy to smooth the solution

nsmooth : int

The number of r-b Gauss-Seidel smoothing iterations to perform

solve(rtol=1e-11)[source]

The main driver for the multigrid solution of the Helmholtz equation. This controls the V-cycles, smoothing at each step of the way and uses simple smoothing at the coarsest level to perform the bottom solve.

Parameters:

rtol : float

The relative tolerance (residual norm / source norm) to solve to. Note that if the source norm is 0 (e.g. the righthand side of our equation is 0), then we just use the norm of the residual.

v_cycle(level)[source]

Perform a V-cycle for a single 2-level solve. This is applied recursively do V-cycle through the entire hierarchy.

multigrid.edge_coeffs module

class multigrid.edge_coeffs.EdgeCoeffs(g, eta, empty=False)[source]

Bases: object

a simple container class to hold edge-centered coefficients and restrict them to coarse levels

restrict()[source]

restrict the edge values to a coarser grid. Return a new EdgeCoeffs object

multigrid.general_MG module

This multigrid solver is build from multigrid/MG.py and implements a more general solver for an equation of the form:

alpha phi + div { beta grad phi } + gamma . grad phi = f

where alpha, beta, and gamma are defined on the same grid as phi. These should all come in as cell-centered quantities. The solver will put beta on edges. Note that gamma is a vector here, with x- and y-components.

A cell-centered discretization for phi is used throughout.

class multigrid.general_MG.GeneralMG2d(nx, ny, xmin=0.0, xmax=1.0, ymin=0.0, ymax=1.0, xl_BC_type='dirichlet', xr_BC_type='dirichlet', yl_BC_type='dirichlet', yr_BC_type='dirichlet', xl_BC=None, xr_BC=None, yl_BC=None, yr_BC=None, nsmooth=10, nsmooth_bottom=50, verbose=0, coeffs=None, true_function=None, vis=0, vis_title='')[source]

Bases: multigrid.MG.CellCenterMG2d

this is a multigrid solver that supports our general elliptic equation.

we need to accept a coefficient CellCenterData2d object with fields defined for alpha, beta, gamma_x, and gamma_y on the fine level.

We then restrict this data through the MG hierarchy (and average beta to the edges).

we need a new compute_residual() and smooth() function, that understands these coeffs.

smooth(level, nsmooth)[source]

Use red-black Gauss-Seidel iterations to smooth the solution at a given level. This is used at each stage of the V-cycle (up and down) in the MG solution, but it can also be called directly to solve the elliptic problem (although it will take many more iterations).

Parameters:

level : int

The level in the MG hierarchy to smooth the solution

nsmooth : int

The number of r-b Gauss-Seidel smoothing iterations to perform

multigrid.variable_coeff_MG module

This multigrid solver is build from multigrid/MG.py and implements a variable coefficient solver for an equation of the form:

div { eta grad phi } = f

where eta is defined on the same grid as phi.

A cell-centered discretization is used throughout.

class multigrid.variable_coeff_MG.VarCoeffCCMG2d(nx, ny, xmin=0.0, xmax=1.0, ymin=0.0, ymax=1.0, xl_BC_type='dirichlet', xr_BC_type='dirichlet', yl_BC_type='dirichlet', yr_BC_type='dirichlet', nsmooth=10, nsmooth_bottom=50, verbose=0, coeffs=None, coeffs_bc=None, true_function=None, vis=0, vis_title='')[source]

Bases: multigrid.MG.CellCenterMG2d

this is a multigrid solver that supports variable coefficients

we need to accept a coefficient array, coeffs, defined at each level. We can do this at the fine level and restrict it down the MG grids once.

we need a new compute_residual() and smooth() function, that understands coeffs.

smooth(level, nsmooth)[source]

Use red-black Gauss-Seidel iterations to smooth the solution at a given level. This is used at each stage of the V-cycle (up and down) in the MG solution, but it can also be called directly to solve the elliptic problem (although it will take many more iterations).

Parameters:

level : int

The level in the MG hierarchy to smooth the solution

nsmooth : int

The number of r-b Gauss-Seidel smoothing iterations to perform