No ghost method for implicit boundary condition treatment of a scalar equation defined on cells. More...

Functions/Subroutines
subroutine	mod_add_cell_bc_noghost::add_cell_bc_noghost_rec (matrix, rhs, boundary_condition, equation_stencil, equation_ls_map)
	Set Robin boundary condition using ghost cells and polynomial extrapolation.

subroutine	mod_add_cell_bc_noghost::apply_extrapolation_on_matrix (matrix, rhs, ia, l, ghost_offsets, target_offsets, g, coefficients)
	Inject ghost cell contributions into the system matrix and RHS.

subroutine	mod_add_cell_bc_noghost_former::add_cell_bc_noghost_former (matrix, rhs, boundary_condition, equation_stencil, equation_ls_map)
	Set boundary condition to a cell-discretized linear system (without ghost cells)

Detailed Description

No ghost method for implicit boundary condition treatment of a scalar equation defined on cells.

Function/Subroutine Documentation

◆ add_cell_bc_noghost_former()

subroutine mod_add_cell_bc_noghost_former::add_cell_bc_noghost_former	(	double precision, dimension(:), intent(inout)	matrix,
		double precision, dimension(:), intent(inout)	rhs,
		type(t_boundary_condition), intent(in)	boundary_condition,
		type(t_cell_stencil), intent(in)	equation_stencil,
		type(t_ls_map), intent(in)	equation_ls_map )

Set boundary condition to a cell-discretized linear system (without ghost cells)

Deprecated: Use mod_add_cell_bc_noghost instead.

The Robin boundary condition is discretized:

\begin{align} g = \alpha \frac {\partial S_b} {\partial n} + \beta S_b \end{align}

The boundary condition is discretized with a centered scheme using exterior and interior nodes.

-----|·····  i: interior cell node,  e: exterior cell node,
  i  b  e    b: boundary face node.

The value of the unknwon on the exterior node can be expressed

\begin{align} S_e = A + B S_i \end{align}

with

\begin{align} A = \frac {g} {\alpha/\Delta + \beta/2} \end{align}

and

\begin{align} B = \frac {\alpha/\Delta - \beta/2} {\alpha/\Delta + \beta/2} \end{align}

Then, on the line of the matrix corresponding to the interior node, the matrix coefficient of the element corresponding to the interior node is set to zero, and matrix elements of the interior node are changed thanks to the relation find above. (as well as the RHS)

For instance, in 1D, on the left boundary:

\begin{align} a_0 S_i + a_{left} S_{i-1} + a_{right} S_{i+1} = Rhs_i \end{align}

is changed by:

\begin{align} (a_0 + Ba_{left}) S_i + a_{right} S_{i+1} = Rhs_i - a_{left} A \end{align}

Special care must be given is the stencil size is equal to 2. A second order extrapolation is done:

\begin{align} S_{i2}=2S_i - Se \end{align}

To be improved.

◆ add_cell_bc_noghost_rec()

subroutine mod_add_cell_bc_noghost::add_cell_bc_noghost_rec	(	double precision, dimension(:), intent(inout)	matrix,
		double precision, dimension(:), intent(inout)	rhs,
		type(t_boundary_condition), intent(in)	boundary_condition,
		type(t_cell_stencil), intent(in)	equation_stencil,
		type(t_ls_map), intent(in)	equation_ls_map )

private

Set Robin boundary condition using ghost cells and polynomial extrapolation.

Parameters

[in,out]	matrix	the system matrix
[in,out]	rhs	the system right-hand side
[in]	boundary_condition	the boundary condition
[in]	equation_stencil	the equation stencil
[in]	equation_ls_map	the linear system mapping

The objective of this routine is to reinject the matrix stencil's coefficients corresponding to exterior/ghost nodes back into the linear system, while ensuring the boundary condition. The routine handles stencil's extents of size 1 or 2. The method is based on polynomial reconstruction and average value computation. The order (linear/second, quadratic/third, cubic/fourth) of accuracy is piloted by the attribute boundary_condition%bc_scheme.

Step 1: Compute ghost cell average values/weights

For this step, we call mod_compute_ghost_weights_table::compute_ghost_weights_table for each necessary fictive ghost cell in order to compute the weights associated to each inner cell and the right-hand-side such that:

\begin{align} \bar{S}_{g_j} = \sum_{m=1}^{p} \omega_m^{(j)} \bar{S}_{i_m} + \omega_g^{(j)} g \end{align}

Step 2: Inject back ghost cell values/weights into the linear system

Once the ghost cell values/weights are computed, they are injected into the discrete linear system inner cells, while we nullify the coefficients associated to the ghost cells, modifying both the matrix and the right-hand side accordingly. See mod_add_cell_bc_noghost::apply_extrapolation_on_matrix routine for details.

If the matrix line is:

\begin{align} \gamma_{g_2} \bar{S}_{g_2} + \gamma_{g_1} \bar{S}_{g_1} + \sum_{m=1}^{N} \gamma_{i_m} \bar{S}_{i_m} = rhs \end{align}

then the injection consists in replacing \( \bar{S}_{g_1} \) and \( \bar{S}_{g_2} \) by the formulae gotten in Step 1, leading to:

\begin{align} \sum_{m=1}^{N} \left( \gamma_{i_m} + \gamma_{g_1} \omega_m^{(1)} + \gamma_{g_2} \omega_m^{(2)} \right) \bar{S}_{i_m} = rhs - \gamma_{g_1} \omega_g^{(1)} g - \gamma_{g_2} \omega_g^{(2)} g \end{align}

Note

The method is compatible with non-uniform grids.
This method is compatible with arbitrary high-order schemes. Higher order (more than fourth) schemes can be easily implemented.
The same strategy is used for filling ghost nodes for explicit schemes.

◆ apply_extrapolation_on_matrix()

subroutine mod_add_cell_bc_noghost::apply_extrapolation_on_matrix	(	double precision, dimension(:), intent(inout)	matrix,
		double precision, dimension(:), intent(inout)	rhs,
		integer, intent(in)	ia,
		integer, intent(in)	l,
		integer, dimension(:), intent(in)	ghost_offsets,
		integer, dimension(:), intent(in)	target_offsets,
		double precision, intent(in)	g,
		double precision, dimension(:,:), intent(in)	coefficients )

private

Inject ghost cell contributions into the system matrix and RHS.

Parameters

[in]	ia	the matrix line starting index
[in]	l	the rhs line index
[in]	ghost_offsets	the ghost cell offsets in the stencil
[in]	target_offsets	the target inner cell offsets in the stencil
[in]	g	the Robin's rhs boundary value
[in,out]	matrix	the system matrix
[in,out]	rhs	the system right-hand side
[in]	coefficients	the ghost cell coefficients table

This routine eliminates ghost unknowns from the discrete linear system by substituting their extrapolated expressions (from Robin BC extrapolation) into the equation rows where they appear. It is the post-processing step called after mod_compute_ghost_weights_table::compute_ghost_weights_table for the no ghost method.

Principle

Suppose a ghost cell average \(\bar{S}_g\) has been expressed as

\begin{align} \bar{S}_g = \sum_{k=1}^{n-1} a_k\,\bar{S}_{i_k} + a_g\, g , \end{align}

with coefficients stored in coefficients(ghost_id,:). Then any stencil entry of the form \( A_{row,g}\,\bar{S}_g \) in the system matrix is rewritten as:

\begin{align} A_{row,g}\,\bar{S}_g &= A_{row,g}\,\Big(\sum_{k=1}^{n-1} a_k\,\bar{S}_{i_k} + a_g\, g\Big) \\ &= \sum_{k=1}^{n-1}\,(a_k\,A_{row,g})\,\bar{S}_{i_k} \ +\ (a_g\,A_{row,g})\,g . \end{align}

The first term updates interior matrix entries,
The second term is added to the RHS (with sign convention depending on the boundary equation form).

Finally, the explicit ghost entry \(A_{row,g}\) in the matrix is set to zero.

Algorithm

For each ghost id:

Loop over all target ids and accumulate matrix(ia+target) += coeff(ghost,target) * matrix(ia+ghost)
Update RHS: rhs(l) -= c * coeff(ghost,n) * matrix(ia+ghost)
Set matrix(ia+ghost) = 0.

Symmetry (left/right boundaries)

The routine is independent on boundary side. By passing appropriate ghost_offsets and target_offsets, it applies equally to left/bottom/back or right/top/front boundaries.

Functions/Subroutines