Preconditioner: BDDC
LFAToolkit.BDDC
— TypeBDDC(operator, injectiontype)
BDDC preconditioner for finite element operators
Arguments:
operator::Operator
: finite element operator to preconditioninjectiontype::BDDCInjectionType
: type of injection into subassembled space to use
Returns:
- BDDC preconditioner object
LFAToolkit.getprimalnodes
— Methodgetprimalnodes(bddc)
Compute or retrieve the primal nodes for a BDDC preconditioner
Returns:
- vector of primal nodes for BDDC preconditioner
Example:
# setup
mesh = Mesh2D(1.0, 1.0);
p = 4
diffusion = GalleryOperator("diffusion", p, p, mesh);
bddc = LumpedBDDC(diffusion);
# verify primal nodes
@assert bddc.primalnodes == [1, p, p^2 - p + 1, p^2]
# output
LFAToolkit.getsubassemblednodes
— Methodgetsubassemblednodes(bddc)
Compute or retrieve the subassembled nodes for a BDDC preconditioner
Returns:
- vector of subassembled nodes for BDDC preconditioner
Example:
# setup
mesh = Mesh2D(1.0, 1.0);
p = 4
diffusion = GalleryOperator("diffusion", p, p, mesh);
bddc = LumpedBDDC(diffusion);
# verify subassembled nodes
@assert bddc.subassemblednodes == setdiff(1:p^2, [1, p, p^2 - p + 1, p^2])
# output
LFAToolkit.getinterfacenodes
— Methodgetinterfacenodes(bddc)
Compute or retrieve the interface nodes for a BDDC preconditioner
Returns:
- vector of interface nodes for BDDC preconditioner
Example:
# setup
mesh = Mesh2D(1.0, 1.0);
p = 4
diffusion = GalleryOperator("diffusion", p, p, mesh);
bddc = LumpedBDDC(diffusion);
# get true interface nodes
trueinterfacenodes = [1:p..., p^2-p+1:p^2...]
for i = 1:p-2
push!(trueinterfacenodes, i * p + 1)
push!(trueinterfacenodes, (i + 1) * p)
end
trueinterfacenodes = sort(trueinterfacenodes)
trueinterfacenodes = setdiff(trueinterfacenodes, bddc.primalnodes)
# verify interface nodes
@assert bddc.interfacenodes == trueinterfacenodes
# output
LFAToolkit.getinteriornodes
— Methodgetinteriornodes(bddc)
Compute or retrieve the interior nodes for a BDDC preconditioner
Returns:
- vector of interior nodes for BDDC preconditioner
Example:
# setup
mesh = Mesh2D(1.0, 1.0);
p = 4
diffusion = GalleryOperator("diffusion", p, p, mesh);
bddc = LumpedBDDC(diffusion);
# get true interface nodes
trueinterfacenodes = [1:p..., p^2-p+1:p^2...]
for i = 1:p-2
push!(trueinterfacenodes, i * p + 1)
push!(trueinterfacenodes, (i + 1) * p)
end
trueinterfacenodes = sort(trueinterfacenodes)
trueinteriornodes = setdiff(bddc.subassemblednodes, trueinterfacenodes)
# verify
@assert bddc.interiornodes == trueinteriornodes
# output
LFAToolkit.getprimalmodes
— Functiongetprimalmodes(bddc)
Compute or retrieve the primal modes for a BDDC preconditioner
Returns:
- vector of primal modes for BDDC preconditioner
Example:
# setup
mesh = Mesh2D(1.0, 1.0);
p = 4
diffusion = GalleryOperator("diffusion", p, p, mesh);
bddc = LumpedBDDC(diffusion);
# verify primal nodes
@assert bddc.primalmodes == [1]
# output
LFAToolkit.getsubassembledmodes
— Functiongetsubassembledmodes(bddc)
Compute or retrieve the subassembled modes for a BDDC preconditioner
Returns:
- vector of subassembled modes for BDDC preconditioner
Example:
# setup
mesh = Mesh2D(1.0, 1.0);
p = 4
diffusion = GalleryOperator("diffusion", p, p, mesh);
bddc = LumpedBDDC(diffusion);
# verify subassembled modes
@assert bddc.subassembledmodes == setdiff(1:(p-1)^2, [1])
# output
LFAToolkit.getinterfacemodes
— Functiongetinterfacemodes(bddc)
Compute or retrieve the interface modes for a BDDC preconditioner
Returns:
- vector of interface modes for BDDC preconditioner
Example:
# setup
mesh = Mesh2D(1.0, 1.0);
p = 4
diffusion = GalleryOperator("diffusion", p, p, mesh);
bddc = LumpedBDDC(diffusion);
# get true interface modes
trueinterfacemodes = [1:p-1...]
for i = 1:p-2
push!(trueinterfacemodes, i * (p - 1) + 1)
end
trueinterfacemodes = sort(trueinterfacemodes)
trueinterfacemodes = setdiff(trueinterfacemodes, bddc.primalmodes)
# verify interface modes
@assert bddc.interfacemodes == trueinterfacemodes
# output
LFAToolkit.getinteriormodes
— Functiongetinteriormodes(bddc)
Compute or retrieve the interior modes for a BDDC preconditioner
Returns:
- vector of interior modes for BDDC preconditioner
Example:
# setup
mesh = Mesh2D(1.0, 1.0);
p = 4
diffusion = GalleryOperator("diffusion", p, p, mesh);
bddc = LumpedBDDC(diffusion);
# verify interior modes
@assert bddc.interiormodes == setdiff(bddc.subassembledmodes, bddc.interfacemodes)
# output
LFAToolkit.getsubassembledinverse
— Functiongetsubassembledinverse(bddc)
Compute or retrieve the solver for subdomain for a BDDC preconditioner
Returns:
- solver for subdomain
Example:
# setup
mesh = Mesh2D(1.0, 1.0);
p = 4
diffusion = GalleryOperator("diffusion", p, p, mesh);
bddc = LumpedBDDC(diffusion);
# verify subassembled inverse
@assert bddc.subassembledinverse ≈
Matrix(diffusion.elementmatrix[bddc.subassemblednodes, bddc.subassemblednodes])^-1
# output
LFAToolkit.getinteriorinverse
— Functiongetinteriorinverse(bddc)
Compute or retrieve the solver for subdomain interior for a BDDC preconditioner
Returns:
- solver for subdomain interior
Example:
# setup
mesh = Mesh2D(1.0, 1.0);
p = 4
diffusion = GalleryOperator("diffusion", p, p, mesh);
bddc = LumpedBDDC(diffusion);
# verify interior inverse
@assert bddc.interiorinverse ≈
Matrix(diffusion.elementmatrix[bddc.interiornodes, bddc.interiornodes])^-1
# output
LFAToolkit.getschur
— Functiongetschur(bddc)
Compute or retrieve the Schur complement matrix for a BDDC preconditioner
Returns:
- matrix for Schur complement
Example:
# setup
mesh = Mesh2D(1.0, 1.0);
p = 4
diffusion = GalleryOperator("diffusion", p, p, mesh);
bddc = LumpedBDDC(diffusion);
# verify Schur complement
@assert bddc.schur ≈
diffusion.elementmatrix[bddc.primalnodes, bddc.primalnodes] -
diffusion.elementmatrix[bddc.primalnodes, bddc.subassemblednodes] *
Matrix(diffusion.elementmatrix[bddc.subassemblednodes, bddc.subassemblednodes])^-1 *
diffusion.elementmatrix[bddc.subassemblednodes, bddc.primalnodes]
# output
LFAToolkit.getmixedmultiplicity
— Functiongetmixedmultiplicity(bddc)
Compute or retrieve the diagonal matrix of mixed interface node and primal mode multiplicity for the BDDC preconditioner
Returns:
- matrix of mixed multiplicity for the BDDC preconditioner
LFAToolkit.getJDT
— FunctiongetJDT(bddc)
Compute or retrieve the matrix mapping the jump over subdomain interface modes
Returns:
- matrix mapping the jump over subdomain interfacemodes
LFAToolkit.getprimalrowmodemap
— Functiongetprimalrowmodemap(bddc)
Compute or retrieve the matrix mapping the rows of the primal BDDC matrix to the primal symbol matrix
Returns:
- matrix mapping rows of primal BDDC matrix to primal symbol matrix
Example:
# setup
mesh = Mesh2D(1.0, 1.0);
diffusion = GalleryOperator("diffusion", 3, 3, mesh);
bddc = LumpedBDDC(diffusion);
# verify primal row mode map
@assert bddc.primalrowmodemap ≈ [1 1 1 1]
# output
LFAToolkit.getprimalcolumnmodemap
— Functiongetprimalcolumnmodemap(bddc)
Compute or retrieve the matrix mapping the columns of the primal BDDC matrix to the primal symbol matrix
Returns:
- matrix mapping columns of primal BDDC matrix to primal symbol matrix
Example:
# setup
mesh = Mesh2D(1.0, 1.0);
diffusion = GalleryOperator("diffusion", 3, 3, mesh);
bddc = LumpedBDDC(diffusion);
# verify primal column mode map
@assert bddc.primalcolumnmodemap ≈ [1; 1; 1; 1]
# output
LFAToolkit.getsubassembledrowmodemap
— Functiongetsubassembledrowmodemap(bddc)
Compute or retrieve the matrix mapping the rows of the subassembled BDDC matrix to the subassembled symbol matrix
Returns:
- matrix mapping rows of subassembled BDDC matrix to subassembled symbol matrix
Example:
# setup
mesh = Mesh2D(1.0, 1.0);
diffusion = GalleryOperator("diffusion", 3, 3, mesh);
bddc = LumpedBDDC(diffusion);
# verify subassembled row mode map
@assert bddc.subassembledrowmodemap ≈
bddc.operator.rowmodemap[bddc.subassembledmodes, bddc.subassemblednodes]
# output
LFAToolkit.getsubassembledcolumnmodemap
— Functiongetsubassembledcolumnmodemap(bddc)
Compute or retrieve the matrix mapping the columns of the subassembled BDDC matrix to the subassembled symbol matrix
Returns:
- matrix mapping columns of subassembled BDDC matrix to subassembled symbol matrix
Example:
# setup
mesh = Mesh2D(1.0, 1.0);
diffusion = GalleryOperator("diffusion", 3, 3, mesh);
bddc = LumpedBDDC(diffusion);
# verify subassembled column mode map
@assert bddc.subassembledcolumnmodemap ≈
bddc.operator.columnmodemap[bddc.subassemblednodes, bddc.subassembledmodes]
# output
LFAToolkit.getmixedrowmodemap
— Functiongetmixedrowmodemap(bddc)
Compute or retrieve the matrix mapping the rows of the mixed BDDC matrix to the symbol matrix
Returns:
- matrix mapping rows of mixed BDDC matrix to symbol matrix
LFAToolkit.getmixedcolumnmodemap
— Functiongetmixedcolumnmodemap(bddc)
Compute or retrieve the matrix mapping the columns of the mixed BDDC matrix to the symbol matrix
Returns:
- matrix mapping columns of mixed BDDC matrix to symbol matrix
LFAToolkit.getmodepermutation
— Functiongetmodepermutation(bddc)
Compute or retrieve the matrix permuting multi-component modes to standard ordering
Returns:
- matrix mapping BDDC mode ordering to standard ordering
LFAToolkit.computesymbolsinjection
— Methodcomputesymbolsinjection(bddc)
Compute or retrieve the injection operator for the BDDC symbol matrix
Returns:
- matrix providing the injection operator of BDDC symbol matrix