Boundary conditions
Ferrite.ConstraintHandler
Ferrite.Dirichlet
Ferrite.PeriodicDirichlet
Ferrite.RHSData
Ferrite.add!
Ferrite.apply!
Ferrite.apply_analytical!
Ferrite.apply_assemble!
Ferrite.apply_local!
Ferrite.apply_rhs!
Ferrite.apply_zero!
Ferrite.close!
Ferrite.collect_periodic_facets
Ferrite.collect_periodic_facets!
Ferrite.get_rhs_data
Ferrite.update!
Ferrite.ConstraintHandler
— TypeConstraintHandler([T=Float64], dh::AbstractDofHandler)
A collection of constraints associated with the dof handler dh
. T
is the numeric type for stored values.
Ferrite.Dirichlet
— TypeDirichlet(u::Symbol, ∂Ω::AbstractVecOrSet, f::Function, components=nothing)
Create a Dirichlet boundary condition on u
on the ∂Ω
part of the boundary. f
is a function of the form f(x)
or f(x, t)
where x
is the spatial coordinate and t
is the current time, and returns the prescribed value. components
specify the components of u
that are prescribed by this condition. By default all components of u
are prescribed.
The set, ∂Ω
, can be an AbstractSet
or AbstractVector
with elements of type FacetIndex
, FaceIndex
, EdgeIndex
, VertexIndex
, or Int
. For most cases, the element type is FacetIndex
, as shown below. To constrain a single point, using VertexIndex
is recommended, but it is also possible to constrain a specific nodes by giving the node numbers via Int
elements. To constrain e.g. an edge in 3d EdgeIndex
elements can be given.
For example, here we create a Dirichlet condition for the :u
field, on the facetset called ∂Ω
and the value given by the sin
function:
Examples
# Obtain the facetset from the grid
∂Ω = getfacetset(grid, "boundary-1")
# Prescribe scalar field :s on ∂Ω to sin(t)
dbc = Dirichlet(:s, ∂Ω, (x, t) -> sin(t))
# Prescribe all components of vector field :v on ∂Ω to 0
dbc = Dirichlet(:v, ∂Ω, x -> 0 * x)
# Prescribe component 2 and 3 of vector field :v on ∂Ω to [sin(t), cos(t)]
dbc = Dirichlet(:v, ∂Ω, (x, t) -> [sin(t), cos(t)], [2, 3])
Dirichlet
boundary conditions are added to a ConstraintHandler
which applies the condition via apply!
and/or apply_zero!
.
Ferrite.PeriodicDirichlet
— TypePeriodicDirichlet(u::Symbol, facet_mapping, components=nothing)
PeriodicDirichlet(u::Symbol, facet_mapping, R::AbstractMatrix, components=nothing)
PeriodicDirichlet(u::Symbol, facet_mapping, f::Function, components=nothing)
Create a periodic Dirichlet boundary condition for the field u
on the facet-pairs given in facet_mapping
. The mapping can be computed with collect_periodic_facets
. The constraint ensures that degrees-of-freedom on the mirror facet are constrained to the corresponding degrees-of-freedom on the image facet. components
specify the components of u
that are prescribed by this condition. By default all components of u
are prescribed.
If the mapping is not aligned with the coordinate axis (e.g. rotated) a rotation matrix R
should be passed to the constructor. This matrix rotates dofs on the mirror facet to the image facet. Note that this is only applicable for vector-valued problems.
To construct an inhomogeneous periodic constraint it is possible to pass a function f
. Note that this is currently only supported when the periodicity is aligned with the coordinate axes.
See the manual section on Periodic boundary conditions for more information.
Ferrite.collect_periodic_facets
— Functioncollect_periodic_facets(grid::Grid, mset, iset, transform::Union{Function,Nothing}=nothing; tol=1e-12)
Match all mirror facets in mset
with a corresponding image facet in iset
. Return a dictionary which maps each mirror facet to a image facet. The result can then be passed to PeriodicDirichlet
.
mset
and iset
can be given as a String
(an existing facet set in the grid) or as a AbstractSet{FacetIndex}
directly.
By default this function looks for a matching facet in the directions of the coordinate system. For other types of periodicities the transform
function can be used. The transform
function is applied on the coordinates of the image facet, and is expected to transform the coordinates to the matching locations in the mirror set.
The keyword tol
specifies the tolerance (i.e. distance and deviation in facet-normals) between a image-facet and mirror-facet, for them to be considered matched.
See also: collect_periodic_facets!
, PeriodicDirichlet
.
collect_periodic_facets(grid::Grid, all_facets::Union{AbstractSet{FacetIndex},String,Nothing}=nothing; tol=1e-12)
Split all facets in all_facets
into image and mirror sets. For each matching pair, the facet located further along the vector (1, 1, 1)
becomes the image facet.
If no set is given, all facets on the outer boundary of the grid (i.e. all facets that do not have a neighbor) is used.
See also: collect_periodic_facets!
, PeriodicDirichlet
.
Ferrite.collect_periodic_facets!
— Functioncollect_periodic_facets!(facet_map::Vector{PeriodicFacetPair}, grid::Grid, mset, iset, transform::Union{Function,Nothing}; tol=1e-12)
Same as collect_periodic_facets
but adds all matches to the existing facet_map
.
Ferrite.add!
— Functionadd!(sdh::SubDofHandler, name::Symbol, ip::Interpolation)
Add a field called name
approximated by ip
to the SubDofHandler sdh
.
add!(dh::DofHandler, name::Symbol, ip::Interpolation)
Add a field called name
approximated by ip
to the DofHandler dh
.
The field is added to all cells of the underlying grid, use SubDofHandler
s if the grid contains multiple cell types, or to add the field to subset of all the cells.
add!(ch::ConstraintHandler, ac::AffineConstraint)
Add the AffineConstraint
to the ConstraintHandler
.
add!(ch::ConstraintHandler, dbc::Dirichlet)
Add a Dirichlet
boundary condition to the ConstraintHandler
.
add!(proj::L2Projector, set::AbstractVecOrSet{Int}, ip::Interpolation;
qr_rhs, [qr_lhs])
Add an interpolation ip
on the cells in set
to the L2Projector
proj
.
qr_rhs
sets the quadrature rule used to later integrate the right-hand-side of the projection equation, when callingproject
. It should match the quadrature points used when creating the quadrature-point variables to project.- The optional
qr_lhs
sets the quadrature rule used to integrate the left-hand-side of the projection equation, and defaults to a quadrature rule that integrates the mass-matrix exactly for the given interpolationip
.
Ferrite.close!
— Functionclose!(dh::AbstractDofHandler)
Closes dh
and creates degrees of freedom for each cell.
close!(ch::ConstraintHandler)
Close and finalize the ConstraintHandler
.
close!(proj::L2Projector)
Close proj
which assembles and calculates the left-hand-side of the projection equation, before doing a Cholesky factorization of the mass-matrix.
Ferrite.update!
— Functionupdate!(ch::ConstraintHandler, time::Real=0.0)
Update time-dependent inhomogeneities for the new time. This calls f(x)
or f(x, t)
when applicable, where f
is the function(s) corresponding to the constraints in the handler, to compute the inhomogeneities.
Note that this is called implicitly in close!(::ConstraintHandler)
.
Ferrite.apply!
— Functionapply!(K::SparseMatrixCSC, rhs::AbstractVector, ch::ConstraintHandler)
Adjust the matrix K
and right hand side rhs
to account for the Dirichlet boundary conditions specified in ch
such that K \ rhs
gives the expected solution.
apply!(K, rhs, ch)
essentially calculates
rhs[free] = rhs[free] - K[constrained, constrained] * a[constrained]
where a[constrained]
are the inhomogeneities. Consequently, the sign of rhs
matters (in contrast with apply_zero!
).
apply!(v::AbstractVector, ch::ConstraintHandler)
Apply Dirichlet boundary conditions and affine constraints, specified in ch
, to the solution vector v
.
Examples
K, f = assemble_system(...) # Assemble system
apply!(K, f, ch) # Adjust K and f to account for boundary conditions
u = K \ f # Solve the system, u should be "approximately correct"
apply!(u, ch) # Explicitly make sure bcs are correct
The last operation is not strictly necessary since the boundary conditions should already be fulfilled after apply!(K, f, ch)
. However, solvers of linear systems are not exact, and thus apply!(u, ch)
can be used to make sure the boundary conditions are fulfilled exactly.
Ferrite.apply_zero!
— Functionapply_zero!(K::SparseMatrixCSC, rhs::AbstractVector, ch::ConstraintHandler)
Adjust the matrix K
and the right hand side rhs
to account for prescribed Dirichlet boundary conditions and affine constraints such that du = K \ rhs
gives the expected result (e.g. du
zero for all prescribed degrees of freedom).
apply_zero!(v::AbstractVector, ch::ConstraintHandler)
Zero-out values in v
corresponding to prescribed degrees of freedom and update values prescribed by affine constraints, such that if a
fulfills the constraints, a ± v
also will.
These methods are typically used in e.g. a Newton solver where the increment, du
, should be prescribed to zero even for non-homogeneouos boundary conditions.
See also: apply!
.
Examples
u = un + Δu # Current guess
K, g = assemble_system(...) # Assemble residual and tangent for current guess
apply_zero!(K, g, ch) # Adjust tangent and residual to take prescribed values into account
ΔΔu = K \ g # Compute the (negative) increment, prescribed values are "approximately" zero
apply_zero!(ΔΔu, ch) # Make sure values are exactly zero
Δu .-= ΔΔu # Update current guess
The last call to apply_zero!
is only strictly necessary for affine constraints. However, even if the Dirichlet boundary conditions should be fulfilled after apply!(K, g, ch)
, solvers of linear systems are not exact. apply!(ΔΔu, ch)
can be used to make sure the values for the prescribed degrees of freedom are fulfilled exactly.
Ferrite.apply_local!
— Functionapply_local!(
local_matrix::AbstractMatrix, local_vector::AbstractVector,
global_dofs::AbstractVector, ch::ConstraintHandler;
apply_zero::Bool = false
)
Similar to apply!
but perform condensation of constrained degrees-of-freedom locally in local_matrix
and local_vector
before they are to be assembled into the global system.
When the keyword argument apply_zero
is true
all inhomogeneities are set to 0
(cf. apply!
vs apply_zero!
).
This method can only be used if all constraints are "local", i.e. no constraint couples with dofs outside of the element dofs (global_dofs
) since condensation of such constraints requires writing to entries in the global matrix/vector. For such a case, apply_assemble!
can be used instead.
Note that this method is destructive since it, by definition, modifies local_matrix
and local_vector
.
Ferrite.apply_assemble!
— Functionapply_assemble!(
assembler::AbstractAssembler, ch::ConstraintHandler,
global_dofs::AbstractVector{Int},
local_matrix::AbstractMatrix, local_vector::AbstractVector;
apply_zero::Bool = false
)
Assemble local_matrix
and local_vector
into the global system in assembler
by first doing constraint condensation using apply_local!
.
This is similar to using apply_local!
followed by assemble!
with the advantage that non-local constraints can be handled, since this method can write to entries of the global matrix and vector outside of the indices in global_dofs
.
When the keyword argument apply_zero
is true
all inhomogeneities are set to 0
(cf. apply!
vs apply_zero!
).
Note that this method is destructive since it modifies local_matrix
and local_vector
.
Ferrite.get_rhs_data
— Functionget_rhs_data(ch::ConstraintHandler, A::SparseMatrixCSC) -> RHSData
Returns the needed RHSData
for apply_rhs!
.
This must be used when the same stiffness matrix is reused for multiple steps, for example when timestepping, with different non-homogeneouos Dirichlet boundary conditions.
Ferrite.apply_rhs!
— Functionapply_rhs!(data::RHSData, f::AbstractVector, ch::ConstraintHandler, applyzero::Bool=false)
Applies the boundary condition to the right-hand-side vector without modifying the stiffness matrix.
See also: get_rhs_data
.
Ferrite.RHSData
— TypeRHSData
Stores the constrained columns and mean of the diagonal of stiffness matrix A
.
Initial conditions
Ferrite.apply_analytical!
— Functionapply_analytical!(
a::AbstractVector, dh::AbstractDofHandler, fieldname::Symbol,
f::Function, cellset=1:getncells(get_grid(dh)))
Apply a solution f(x)
by modifying the values in the degree of freedom vector a
pertaining to the field fieldname
for all cells in cellset
. The function f(x)
are given the spatial coordinate of the degree of freedom. For scalar fields, f(x)::Number
, and for vector fields with dimension dim
, f(x)::Vec{dim}
.
This function can be used to apply initial conditions for time dependent problems.
This function only works for standard nodal finite element interpolations when the function value at the (algebraic) node is equal to the corresponding degree of freedom value. This holds for e.g. Lagrange and Serendipity interpolations, including sub- and superparametric elements.