CombineDofs
CombineDofs
provides a mechanism to couple degrees of freedom (DOFs) in a finite element system. This is especially useful for enforcing periodic boundary conditions, multi-point constraints, or other situations where DOFs from different locations or unknowns should be treated as identical in the assembled system.
ExtendableFEM.CombineDofs
— Typefunction CombineDofs(uX, uY, dofsX, dofsY, factors; kwargs...)
When assembled, the dofsX of the unknown uX will be coupled with the dofsY of uY, e.g., for periodic boundary conditions.
Keyword arguments:
name
: name for operator used in printouts. Default: ''CombineDofs''penalty
: penalty for fixed degrees of freedom. Default: 1.0e30verbosity
: verbosity level. Default: 0
ExtendableFEM.CombineDofs
— MethodCombineDofs(
uX,
uY,
coupling_matrix::AbstractMatrix;
kwargs...
) -> CombineDofs{_A, CT, Vector{Int64}} where {_A, CT<:(AbstractMatrix)}
Creates an operator that enforces linear constraints between degrees of freedom (dofs) of two unknowns, e.g. for periodic boundary conditions.
Arguments
uX
: The unknown (or block index) whose dofs are to be constrained (the "from" side).uY
: The unknown (or block index) whose dofs are used as targets (the "to" side).coupling_matrix
: A precomputed sparse matrix encoding the coupling between dofs (as fromget_periodic_coupling_matrix
).
Keyword Arguments
name
: name for operator used in printouts. Default: ''CombineDofs''penalty
: penalty for fixed degrees of freedom. Default: 1.0e30verbosity
: verbosity level. Default: 0
Returns
A CombineDofs
operator that, when assembled, enforces the specified linear relationships between dofs by adding penalty terms to the system matrix.
Periodic Boundary Conditions
To set up periodic boundary conditions, you often need to determine which DOFs should be coupled. The following utility functions can help:
ExtendableFEM.get_periodic_coupling_info
— Functionfunction get_periodic_coupling_info(FES, xgrid, b1, b2, is_opposite::Function; factor_vectordofs = "auto")
computes the dofs that have to be coupled for periodic boundary conditions on the given xgrid for boundary regions b1, b2. The isopposite function evaluates if two provided face midpoints are on opposite sides to each other (the mesh xgrid should be appropriate). For vector-valued FETypes the user can provide factorvectordofs to incorporate a sign change if needed. This is automatically done for all Hdiv-conforming elements and (for the normal-weighted face bubbles of) the Bernardi-Raugel element H1BR.
ExtendableFEM.get_periodic_coupling_matrix
— Functionget_periodic_coupling_matrix(
FES::FESpace,
b_from,
b_to,
give_opposite!::Function;
mask = :auto,
sparsity_tol = 1.0e-12
)
Compute a coupling information for each dof on one boundary as a linear combination of dofs on another boundary
Input:
- FES: FE space to be coupled (on its dofgrid)
- bfrom: boundary region(s) of the grid which dofs should be replaced in terms of dofs on bto
- bto: boundary region(s) of the grid with dofs to replace the dofs in bfrom
- give_opposite! Function in (y,x)
- mask: (optional) vector of masking components
- sparsity_tol: threshold for treating an interpolated value as zero
giveopposite!(y,x) has to be defined in a way that for each x ∈ bfrom the resulting y is in the opposite boundary. For each x in the grid, the resulting y has to be in the grid, too: incorporate some mirroring of the coordinates. Example: If bfrom is at x[1] = 0 and the opposite boundary is at y[1] = 1, then giveopposite!(y,x) = y .= [ 1-x[1], x[2] ]
The return value is a (𝑛 × 𝑛) sparse matrix 𝐴 (𝑛 is the total number of dofs) containing the periodic coupling information. The relation ship between the degrees of freedome is dofᵢ = ∑ⱼ Aⱼᵢ ⋅ dofⱼ. It is guaranteed that i) Aᵢⱼ=0 if dofᵢ is 𝑛𝑜𝑡 on the boundary b_from. ii) Aᵢⱼ=0 if the opposite of dofᵢ is not in the same grid cell as dofⱼ.
Example: Periodic Boundary Coupling (extract from Example212)
Suppose you want to enforce periodicity between the left and right boundaries of a 2D domain. You can use get_periodic_coupling_matrix
to find the corresponding DOFs, and then use CombineDofs
to couple them in the system. The following code is adapted from Example212_PeriodicBoundary2D.jl:
# Define a function to map points on the left to the right boundary
function give_opposite!(y, x)
y .= x
y[1] = width - x[1]
return nothing
end
# Compute the coupling matrix
coupling_matrix = get_periodic_coupling_matrix(FES, reg_left, reg_right, give_opposite!)
# Assign the CombineDofs operator to the problem description
assign_operator!(PD, CombineDofs(u, u, coupling_matrix))
See also Example212 and Example312 for more advanced use cases and details on periodic boundary conditions.