ExtendableFEMBase.AbstractFiniteElement
— Typeabstract type AbstractFiniteElement
reconstruction identity operator: evaluates a reconstructed version of the finite element function.
FEreconst specifies the reconstruction space and reconstruction algorithm if it is defined for the finite element that it is applied to.
ExtendableFEMBase.AbstractFunctionOperator
— Typeabstract type AbstractFunctionOperator
root type for FunctionOperators.
ExtendableFEMBase.AccumulatingVector
— Typestruct AccumulatingVector{T} <: AbstractArray{T, 2}
block of an FEVector that carries coefficients for an associated FESpace and can be assigned as an AbstractArray (getindex, setindex, size, length)
ExtendableFEMBase.Curl2D
— Typeabstract type Curl2D <: ??
evaluates the curl of some two-dimensional vector field, i.e. Curl2D((u1,u2)) = du2/dx1 - du1/dx2
ExtendableFEMBase.Curl3D
— Typeabstract type Curl3D <: ??
evaluates the curl of some three-dimensional vector field, i.e. Curl3D(u) = abla imes u
ExtendableFEMBase.CurlScalar
— Typeabstract type CurlScalar <: ??
evaluates the curl of some scalar function in 2D, i.e. the rotated gradient.
ExtendableFEMBase.Divergence
— Typeabstract type Divergence <: ??
evaluates the divergence of the finite element function.
ExtendableFEMBase.DofMap
— Typeabstract type DofMap <: ExtendableGrids.AbstractGridAdjacency
Dofmaps are stored as an ExtendableGrids.AbstractGridAdjacency in the finite element space and collect information with respect to different AssemblyTypes. They are generated automatically on demand and the dofmaps associated to each subtype can be accessed via FESpace[DofMap].
ExtendableFEMBase.FEEvaluator
— Methodfunction FEEvaluator(FE::FESpace, operator::AbstractFunctionOperator, qrule::QuadratureRule; T = Float64, AT = ON_CELLS, L2G = nothing)
Constructs a FEEvaluator that handles evaluations of finite element basis function evaluation for the given FESpace, operator at the quadrature points of the given QuadratureRule. It has an update! function to update the evaluation upon entry to a new cell. Evaluations can be accessed via FEEvaluator.cvals[j,k,i] where i is the quadrature point id, k is the local dof number and j is the component.
Note that matrix-valued operators evaluations, e.g. for Gradient, are given as a long vector (in component-wise order).
ExtendableFEMBase.FEMatrix
— Typestruct FEMatrix{TvM, TiM, TvG, TiG, nbrow, nbcol, nbtotal} <: SparseArrays.AbstractSparseArray{TvM, TiM, 2}
an AbstractMatrix (e.g. an ExtendableSparseMatrix) with an additional layer of several FEMatrixBlock subdivisions each carrying coefficients for their associated pair of FESpaces
ExtendableFEMBase.FEMatrix
— MethodFEMatrix{TvM,TiM}(FESX, FESY; name = "auto")
Creates FEMatrix with one rectangular block (FESX,FESY) if FESX and FESY are single FESpaces, or a rectangular block matrix with blocks corresponding to the entries of the FESpace vectors FESX and FESY. Optionally a name for the matrix can be given.
ExtendableFEMBase.FEMatrix
— MethodFEMatrix{TvM,TiM}(name::String, FES::FESpace{TvG,TiG,FETypeX,APTX}) where {TvG,TiG,FETypeX,APTX}
Creates FEMatrix with one square block (FES,FES).
ExtendableFEMBase.FEMatrix
— MethodFEMatrix{TvM,TiM}(FESX, FESY; name = "auto")
Creates an FEMatrix with blocks coressponding to the ndofs of FESX (rows) and FESY (columns).
ExtendableFEMBase.FEMatrixBlock
— Typestruct FEMatrixBlock{TvM, TiM, TvG, TiG, FETypeX, FETypeY, APTX, APTY} <: AbstractArray{TvM, 2}
block of an FEMatrix that carries coefficients for an associated pair of FESpaces and can be assigned as an two-dimensional AbstractArray (getindex, setindex, size)
ExtendableFEMBase.FESpace
— Typestruct FESpace{Tv, Ti, FEType<:AbstractFiniteElement,AT<:AssemblyType}
name::String # full name of finite element space (used in messages)
broken::Bool # if true, broken dofmaps are generated
ndofs::Int # total number of dofs
coffset::Int # offset for component dofs
xgrid::ExtendableGrid[Tv,Ti} # link to xgrid
dofmaps::Dict{Type{<:AbstractGridComponent},Any} # backpack with dofmaps
end
A struct that has a finite element type as parameter and carries dofmaps (CellDofs, FaceDofs, BFaceDofs) plus additional grid information and access to arrays holding coefficients if needed.
ExtendableFEMBase.FESpace
— Methodfunction FESpace{FEType<:AbstractFiniteElement,AT<:AssemblyType}(
xgrid::ExtendableGrid{Tv,Ti};
name = "",
broken::Bool = false)
Constructor for FESpace of the given FEType, AT = ONCELLS/ONFACES/ONEDGES generates a finite elements space on the cells/faces/edges of the provided xgrid (if omitted ONCELLS is used as default). The broken switch allows to generate a broken finite element space (that is piecewise H1/Hdiv/HCurl). If no name is provided it is generated automatically from FEType. If no AT is provided, the space is generated ON_CELLS.
ExtendableFEMBase.FEVector
— Typestruct FEVector{T, Tv, Ti}
a plain array but with an additional layer of several FEVectorBlock subdivisions each carrying coefficients for their associated FESpace. The j-th block can be accessed by getindex(::FEVector, j) or by getindex(::FEVector, tag) if tags are associated. The full vector can be accessed via FEVector.entries
ExtendableFEMBase.FEVector
— MethodFEVector{T}(FES; name = nothing, tags = nothing, kwargs...) where T <: Real
Creates FEVector that has one block if FES is a single FESpace, and a blockwise FEVector if FES is a vector of FESpaces. Optionally a name for the vector (as a String) or each of the blocks (as a vector of Strings), or tags (as an Array{Any}) for the blocks can be specified.
ExtendableFEMBase.FEVectorBlock
— Typestruct FEVectorBlock{T, Tv, Ti, FEType, APT} <: AbstractArray{T, 1}
block of an FEVector that carries coefficients for an associated FESpace and can be assigned as an AbstractArray (getindex, setindex, size, length)
ExtendableFEMBase.Gradient
— Typeabstract type Gradient <: ??
evaluates the gradient of the finite element function.
ExtendableFEMBase.H1BR
— Typeabstract type H1BR{edim} <: AbstractH1FiniteElementWithCoefficients where {edim<:Int}
vector-valued (ncomponents = edim) Bernardi–Raugel element (first-order polynomials + normal-weighted face bubbles)
allowed ElementGeometries:
- Triangle2D (piecewise linear + normal-weighted face bubbles)
- Quadrilateral2D (Q1 space + normal-weighted face bubbles)
- Tetrahedron3D (piecewise linear + normal-weighted face bubbles)
ExtendableFEMBase.H1BUBBLE
— Typeabstract type H1BUBBLE{ncomponents} <: AbstractH1FiniteElement where {ncomponents<:Int}
Piecewise bubbles (=zero at boundary)
allowed element geometries:
- Edge1D (one quadratic bubble)
- Triangle2D (one cubic bubble)
- Quadrilateral2D (one quartic bubble)
- Tetrahedron3D (one cubic bubble)
ExtendableFEMBase.H1CR
— Typeabstract type H1CR{ncomponents} <: AbstractH1FiniteElement where {ncomponents<:Int}
Crouzeix-Raviart element (only continuous at face centers).
allowed ElementGeometries:
- Triangle2D (piecewise linear, similar to P1)
- Quadrilateral2D (similar to Q1 space)
- Tetrahedron3D (piecewise linear, similar to P1)
ExtendableFEMBase.H1MINI
— Typeabstract type H1MINI{ncomponents,edim} <: AbstractH1FiniteElement where {ncomponents<:Int,edim<:Int}
Mini finite element.
allowed element geometries:
- Triangle2D (linear polynomials + cubic cell bubble)
- Quadrilateral2D (Q1 space + quartic cell bubble)
- Tetrahedron3D (linear polynomials + cubic cell bubble)
ExtendableFEMBase.H1P1
— Typeabstract type H1P1{ncomponents} <: AbstractH1FiniteElement where {ncomponents<:Int}
Continuous piecewise first-order linear polynomials.
allowed ElementGeometries:
- Edge1D
- Triangle2D
- Tetrahedron3D
ExtendableFEMBase.H1P1TEB
— Typeabstract type H1P1TEB{edim} <: AbstractH1FiniteElementWithCoefficients where {edim<:Int}
vector-valued (ncomponents = edim) element that uses P1 functions + tangential-weighted edge bubbles as suggested by [Diening, L., Storn, J. & Tscherpel, T., "Fortin operator for the Taylor–Hood element", Numer. Math. 150, 671–689 (2022)]
(is inf-sup stable for Stokes if paired with continuous P1 pressure space, less degrees of freedom than MINI)
allowed ElementGeometries:
- Triangle2D
- Tetrahedron3D
ExtendableFEMBase.H1P2
— Typeabstract type H1P2{ncomponents,edim} <: AbstractH1FiniteElement where {ncomponents<:Int,edim<:Int}
Continuous piecewise second-order polynomials.
allowed ElementGeometries:
- Edge1D
- Triangle2D
- Tetrahedron3D
ExtendableFEMBase.H1P2B
— Typeabstract type H1P2B{ncomponents,edim} <: AbstractH1FiniteElement where {ncomponents<:Int,edim<:Int}
Continuous piecewise second-order polynomials.
allowed ElementGeometries:
- Triangle2D
ExtendableFEMBase.H1P3
— Typeabstract type H1P3{ncomponents,edim} <: AbstractH1FiniteElement where {ncomponents<:Int,edim<:Int}
Continuous piecewise third-order polynomials.
allowed ElementGeometries:
- Edge1D
- Triangle2D
- Tetrahedron3D
ExtendableFEMBase.H1Pk
— Typeabstract type H1PK{ncomponents,edim,order} <: AbstractH1FiniteElement where {ncomponents<:Int,edim<:Int,order<:Int}
Continuous piecewise polynomials of arbitrary order >= 1 with ncomponents components in edim space dimensions.
allowed ElementGeometries:
- Edge1D
- Triangle2D
ExtendableFEMBase.H1Q1
— Typeabstract type Q1P1{ncomponents} <: AbstractH1FiniteElement where {ncomponents<:Int}
Continuous piecewise first-order polynomials on simplices and quads, can be used for mixed geometries.
allowed ElementGeometries:
- Edge1D (P1 space)
- Triangle2D (P1 space)
- Quadrilateral2D (Q1 space)
- Tetrahedron3D (P1 space)
- Hexahedron3D (Q1 space)
ExtendableFEMBase.H1Q2
— Typeabstract type H1Q2{ncomponents,edim} <: AbstractH1FiniteElement where {ncomponents<:Int,edim<:Int}
Continuous piecewise second-order polynomials on simplices and quads. Can be used with mixed geometries (in 2D).
allowed ElementGeometries:
- Edge1D (P2 space)
- Triangle2D (P2 space)
- Quadrilateral2D (Q2 space with cell bubble)
- Tetrahedron3D (P2 space)
ExtendableFEMBase.HCURLN0
— Typeabstract type HCURLN0{edim} <: AbstractHcurlFiniteElement where {edim<:Int}
Hcurl-conforming vector-valued (ncomponents = edim) lowest-order Nedelec space of first kind.
allowed ElementGeometries:
- Triangle2D
- Quadrilateral2D
- Tetrahedron3D
ExtendableFEMBase.HCURLN1
— Typeabstract type HCURLN1{edim} <: AbstractHcurlFiniteElement where {edim<:Int}
Hcurl-conforming vector-valued (ncomponents = edim) Nedelec space of first kind and order 1.
allowed ElementGeometries:
- Triangle2D
ExtendableFEMBase.HDIVBDM1
— Typeabstract type HDIVBDM1{edim} <: AbstractHdivFiniteElement where {edim<:Int}
Hdiv-conforming vector-valued (ncomponents = edim) lowest-order Brezzi-Douglas-Marini space
allowed ElementGeometries:
- Triangle2D
- Quadrilateral2D
- Tetrahedron3D
ExtendableFEMBase.HDIVBDM2
— Typeabstract type HDIVBDM2{edim} <: AbstractHdivFiniteElement where {edim<:Int}
Hdiv-conforming vector-valued (ncomponents = edim) Brezzi-Douglas-Marini space of order 2
allowed ElementGeometries:
- Triangle2D
ExtendableFEMBase.HDIVRT0
— Typeabstract type HDIVRT0{edim} <: AbstractHdivFiniteElement where {edim<:Int}
Hdiv-conforming vector-valued (ncomponents = edim) lowest-order Raviart-Thomas space.
allowed ElementGeometries:
- Triangle2D
- Quadrilateral2D
- Tetrahedron3D
- Hexahedron3D
ExtendableFEMBase.HDIVRT1
— Typeabstract type HDIVRT1{edim} <: AbstractHdivFiniteElement where {edim<:Int}
Hdiv-conforming vector-valued (ncomponents = edim) Raviart-Thomas space of order 1.
allowed ElementGeometries:
- Triangle2D
- Tetrahedron3D
ExtendableFEMBase.Hessian
— Typeabstract type Hessian <: ??
evaluates the full Hessian of some (possibly vector-valued) finite element function.
ExtendableFEMBase.Identity
— Typeabstract type Identity <: ??
identity operator: evaluates finite element function.
ExtendableFEMBase.IdentityComponent
— Typeabstract type IdentityComponent{c} <: ??
identity operator: evaluates only the c-th component of the finite element function.
ExtendableFEMBase.L2P0
— Typeabstract type L2P0{ncomponents} <: AbstractH1FiniteElement where {ncomponents<:Int}
Piecewise constant polynomials on cells.
allowed ElementGeometries:
- any
ExtendableFEMBase.L2P1
— Typeabstract type L2P1{ncomponents} <: AbstractH1FiniteElement where {ncomponents<:Int}
Discontinuous piecewise first-order linear polynomials.
allowed ElementGeometries:
- any
ExtendableFEMBase.Laplacian
— Typeabstract type Laplacian <: ??
evaluates the Laplacian of some (possibly vector-valued) finite element function.
ExtendableFEMBase.NormalFlux
— Typeabstract type NormalFlux <: ??
evaluates the normal-flux of the finite element function.
only available on FACES/BFACES and currently only for H1 and Hdiv elements
ExtendableFEMBase.OperatorPair
— Typeabstract type OperatorPair{<:StandardFunctionOperator,<:StandardFunctionOperator} <: StandardFunctionOperator
allows to evaluate two operators in place of one, e.g. OperatorPair{Identity,Gradient}.
ExtendableFEMBase.OperatorTriple
— Typeabstract type OperatorTriple{<:StandardFunctionOperator,<:StandardFunctionOperator} <: StandardFunctionOperator
allows to evaluate three operators in place of one, e.g. OperatorTriple{Identity,Gradient,Hessian}.
ExtendableFEMBase.PointEvaluator
— Typefunction Pointevaluator(
[kernel!::Function],
oa_args::Array{<:Tuple{<:Any, DataType},1};
kwargs...)
Generates a PointEvaluator that can evaluate the specified operator evaluations at arbitrary points. If no kernel function is given, the arguments are given directly. If a kernel is provided, the arguments are postprocessed accordingly and the kernel has be conform to the interface
kernel!(result, eval_args, qpinfo)
where qpinfo allows to access information at the current evaluation point. Additionally the length of the result needs to be specified via the kwargs.
Evaluation can be triggered via the evaluate function after an initialize! call.
Operator evaluations are tuples that pair a tag (to identify an unknown or the position in the vector) with a FunctionOperator.
Keyword arguments:
resultdim: dimension of result field (default = length of operators). Default: 0
params: array of parameters that should be made available in qpinfo argument of kernel function. Default: nothing
name: name for operator used in printouts. Default: ''PointEvaluator''
verbosity: verbosity level. Default: 0
ExtendableFEMBase.QuadratureRule
— Typeabstract type QuadratureRule{T<:Real, ET<:ExtendableGrids.AbstractElementGeometry}
A struct that contains the name of the quadrature rule, the reference points and the weights for the parameter-determined element geometry.
ExtendableFEMBase.QuadratureRule
— Methodfunction QuadratureRule{T,ET}(order::Int) where {T<:Real, ET <: AbstractElementGeometry0D}
Constructs 0D quadrature rule of specified order (always point evaluation).
ExtendableFEMBase.QuadratureRule
— Methodfunction QuadratureRule{T,ET}(order::Int) where {T<:Real, ET <: AbstractElementGeometry1D}
Constructs 1D quadrature rule of specified order.
ExtendableFEMBase.QuadratureRule
— Methodfunction QuadratureRule{T,ET}(order::Int) where {T<:Real, ET <: Parallelepiped3D}
Constructs quadrature rule on Parallelepiped3D of specified order.
ExtendableFEMBase.QuadratureRule
— Methodfunction QuadratureRule{T,ET}(order::Int) where {T<:Real, ET <: Parallelogram2D}
Constructs quadrature rule on Parallelogram2D of specified order.
ExtendableFEMBase.QuadratureRule
— Methodfunction QuadratureRule{T,ET}(order::Int) where {T<:Real, ET <: Tetrahedron3D}
Constructs quadrature rule on Tetrahedron3D of specified order.
ExtendableFEMBase.QuadratureRule
— Methodfunction QuadratureRule{T,ET}(order::Int) where {T<:Real, ET <: Triangle2D}
Constructs quadrature rule on Triangle2D of specified order.
ExtendableFEMBase.Reconstruct
— Typeabstract type Reconstruct{FETypeR, O} <: ??
reconstruction operator: evaluates a reconstructed version of the finite element function.
FETypeR specifies the reconstruction space (needs to be defined for the finite element that it is applied to). O specifies the StandardFunctionOperator that shall be evaluated.
ExtendableFEMBase.SegmentIntegrator
— Methodfunction SegmentIntegrator(
EG::ElementGeometry,
[kernel!::Function],
oa_args::Array{<:Tuple{<:Any, DataType},1};
kwargs...)
Generates an SegmentIntegrator that can intgrate over segments of the specified geometry EG. To do so, it evaluates, at each quadrature point, the specified operator evaluations, postprocesses them with the kernel function (if provided) and accumulates the results with the quadrature weights. If no kernel is given, the arguments are integrated directly. If a kernel is provided it has be conform to the interface
kernel!(result, eval_args, qpinfo)
where qpinfo allows to access information at the current quadrature point. Additionally the length of the result needs to be specified via the kwargs.
Evaluation can be triggered via the integrate_segment! function after an initialize!
Operator evaluations are tuples that pair a tag (to identify an unknown or the position in the vector) with a FunctionOperator.
Keyword arguments:
factor: factor that should be multiplied during assembly. Default: 1
resultdim: dimension of result field (default = length of arguments). Default: 0
matrix_mode: integrator integrates basis functions of FEspace seperately to assembly a matrix that maps solution to segment integrations (requires that kernel is linear). Default: false
name: name for operator used in printouts. Default: ''SegmentIntegrator''
params: array of parameters that should be made available in qpinfo argument of kernel function. Default: nothing
quadorder: quadrature order. Default: ''auto''
bonus_quadorder: additional quadrature order added to quadorder. Default: 0
entrytolerance: threshold to add entry to sparse matrix (only in matrixmode). Default: 0
verbosity: verbosity level. Default: 0
ExtendableFEMBase.SymmetricGradient
— Typeabstract type SymmetricGradient{offdiagval} <: ??
evaluates the symmetric part of the gradient of the finite element function and returns its Voigt compression (off-diagonals on position [j,k] and [k,j] are added together and weighted with offdiagval).
ExtendableFEMBase.SymmetricHessian
— Typeabstract type SymmetricHessian{offdiagval} <: ??
evaluates only the symmetric part of the Hessian of some (possibly vector-valued) finite element function. A concatenation of Voigt compressed Hessians for each component of the finite element functions is returned. The weighting of the mixed derivatives can be steered with offdiagval (e.g. √2 or 1 depending on the use case).
ExtendableFEMBase.TangentFlux
— Typeabstract type TangentFlux <: ??
evaluates the tangent-flux of the finite element function.
ExtendableFEMBase.TangentialGradient
— Typeabstract type TangentialGradient <: ??
evaluates the gradient of the tangential part of some vector-valued finite element function.
Base.append!
— Methodappend!(
FEF::FEVector{T},
FES::FESpace{Tv, Ti, FEType, APT};
name,
tag
) -> Int64
Overloaded append
function for FEVector
that adds a FEVectorBlock at the end.
Base.eltype
— Methodeltype(
_::FESpace{Tv, Ti, FEType<:AbstractFiniteElement, APT}
) -> Type{FEType} where FEType<:AbstractFiniteElement
Custom eltype
function for FESpace
returns the finite element type parameter of the finite element space.
Base.eltype
— Methodeltype(
_::QuadratureRule{T<:Real, ET<:ExtendableGrids.AbstractElementGeometry}
) -> Any
Custom eltype
function for QuadratureRule{T,ET}
.
Base.fill!
— Methodfill!(b::FEVectorBlock, value)
Overloaded fill
function for FEVectorBlock
(only fills the block, not the complete FEVector).
Base.fill!
— Methodfill!(B::FEMatrixBlock{Tv, Ti}, value)
Custom fill
function for FEMatrixBlock
(only fills the already present nzval in the block, not the complete FEMatrix).
Base.get!
— Methodget!(FES::FESpace, DM::Type{<:DofMap}) -> Any
To be called by getindex. This triggers lazy creation of non-existing dofmaps
Base.getindex
— MethodBase.getindex(FES::FESpace,DM::Type{<:DofMap})
Generic method for obtaining dofmap. This method is mutating in the sense that non-existing dofmaps are created on demand. Due to the fact that components are stored as Any the return value triggers type instability.
Base.length
— Methodlength(FEB::FEVectorBlock) -> Int64
Custom length
function for FEVectorBlock
that gives the coressponding number of degrees of freedoms of the associated FESpace
Base.length
— Methodlength(FEF::FEVector) -> Int64
Custom length
function for FEVector
that gives the number of defined FEMatrixBlocks in it
Base.length
— Methodlength(
_::FEMatrix{TvM, TiM, TvG, TiG, nbrow, nbcol, nbtotal}
) -> Any
Custom length
function for FEMatrix
that gives the total number of defined FEMatrixBlocks in it
Base.setindex!
— Methodsetindex!(FES::FESpace, v, DM::Type{<:DofMap}) -> Any
Set new dofmap
Base.show
— Methodshow(io::IO, FEF::FEVector)
Custom show
function for FEVector
that prints some information on its blocks.
Base.show
— Methodshow(io::IO, Q::QuadratureRule)
Custom show
function for QuadratureRule{T,ET}
that prints some information.
Base.show
— Methodshow(
io::IO,
FES::FESpace{Tv, Ti, FEType<:AbstractFiniteElement, APT}
)
Custom show
function for FESpace
that prints some information and all available dofmaps.
Base.show
— Methodshow(
io::IO,
FEM::FEMatrix{TvM, TiM, TvG, TiG, nbrow, nbcol, nbtotal}
)
Custom show
function for FEMatrix
that prints some information on its blocks.
Base.size
— Methodsize(FEB::FEMatrixBlock) -> Tuple{Int64, Int64}
Custom size
function for FEMatrixBlock
that gives a tuple with the size of the block (that coressponds to the number of degrees of freedoms in X and Y)
Base.size
— Methodsize(
_::FEMatrix{TvM, TiM, TvG, TiG, nbrow, nbcol, nbtotal}
) -> Tuple{Any, Any}
Custom size
function for FEMatrix
that gives a tuple with the number of rows and columns of the FEBlock overlay
Base.view
— Methodreturns a view of the part of the full FEVector that coressponds to the block.
ExtendableFEMBase.FESpaces
— MethodFESpaces(FEV::FEVector{T, Tv, Ti}) -> Any
Returns the vector of FEspaces for the blocks of the given FEVector.
ExtendableFEMBase.add!
— Methodadd!(A::FEMatrix{Tv, Ti}, B::FEMatrix{Tv, Ti}; kwargs...)
Adds FEMatrix/ExtendableSparseMatrix/CSCMatrix B to FEMatrix A.
ExtendableFEMBase.addblock!
— Methodaddblock!(
a::FEVectorBlock,
b::AbstractVector;
factor,
offset
)
Adds Array b to FEVectorBlock a.
ExtendableFEMBase.addblock!
— Methodaddblock!(a::FEVectorBlock, b::FEVectorBlock; factor)
Adds FEVectorBlock b to FEVectorBlock a.
ExtendableFEMBase.addblock!
— Methodaddblock!(
A::FEMatrixBlock{Tv, Ti},
B::FEMatrixBlock{Tv, Ti};
factor,
transpose
)
Adds FEMatrixBlock B to FEMatrixBlock A.
ExtendableFEMBase.addblock!
— Methodaddblock!(
A::FEMatrixBlock{Tv},
B::ExtendableSparse.ExtendableSparseMatrix{Tv, Ti<:Integer};
factor,
transpose
)
Adds ExtendableSparseMatrix B to FEMatrixBlock A.
ExtendableFEMBase.addblock!
— Methodaddblock!(
A::FEMatrixBlock{Tv},
cscmat::SparseArrays.SparseMatrixCSC{Tv, Ti<:Integer};
factor,
transpose
)
Adds SparseMatrixCSC B to FEMatrixBlock A.
ExtendableFEMBase.addblock_matmul!
— Methodaddblock_matmul!(
a::AbstractArray{Tv, 1},
B::FEMatrixBlock{Tv, Ti},
b::AbstractArray{Tv, 1};
factor,
transposed
)
Adds matrix-vector product B times b to FEVectorBlock a.
ExtendableFEMBase.addblock_matmul!
— Methodaddblock_matmul!(
A::FEMatrixBlock{Tv},
cscmatB::SparseArrays.SparseMatrixCSC{Tv, Ti},
cscmatC::SparseArrays.SparseMatrixCSC{Tv, Ti};
factor,
transposed
)
Adds matrix-matrix product B times C to FEMatrixBlock A.
ExtendableFEMBase.addblock_matmul!
— Methodaddblock_matmul!(
a::FEVectorBlock{Tv},
B::ExtendableSparse.ExtendableSparseMatrix{Tv, Ti<:Integer},
b::FEVectorBlock{Tv};
factor
)
Adds matrix-vector product B times b to FEVectorBlock a.
ExtendableFEMBase.addblock_matmul!
— Methodaddblock_matmul!(
a::FEVectorBlock{Tv},
B::FEMatrixBlock{Tv, Ti},
b::FEVectorBlock{Tv};
factor,
transposed
)
Adds matrix-vector product B times b (or B' times b if transposed = true) to FEVectorBlock a.
ExtendableFEMBase.assemblytype
— Methodassemblytype(
_::FESpace{Tv, Ti, FEType<:AbstractFiniteElement, APT}
) -> Any
returns the assembly type parameter of the finite element space, i.e. on which entities of the grid the finite element is defined.
ExtendableFEMBase.displace_mesh!
— Methodfunction displace_mesh!(xgrid::ExtendableGrid, source::FEVectorBlock; magnify = 1)
Moves all nodes of the grid by adding the displacement field in source (expects a vector-valued finite element) times a magnify value.
ExtendableFEMBase.eval_febe!
— Function eval_febe!(result, FEBE::FEBasisEvaluator, j::Int, i::Int, offset::Int = 0, factor = 1)
Evaluates the linear combination of the basisfunction with given coefficients at the i-th quadrature point and writes the (possibly vector-valued) evaluation into result (beginning at offset and with the specified factor).
ExtendableFEMBase.eval_febe!
— Function eval_febe!(result, FEBE::FEBasisEvaluator, j::Int, i::Int, offset::Int = 0, factor = 1)
Evaluate the j-th basis function of the FEBasisEvaluator at the i-th quadrature point and writes the (possibly vector-valued) evaluation into result (beginning at offset and with the specified factor).
ExtendableFEMBase.eval_func
— Methodfunction eval_func(PE::PointEvaluator)
Yields the function (result, x) -> evaluate!(result,PE,x).
ExtendableFEMBase.eval_func_bary
— Methodfunction eval_func_bary(PE::PointEvaluator)
Yields the function (result, xref, item) -> evaluate_bary!(result,PE,xref,item).
ExtendableFEMBase.evaluate!
— Methodfunction evaluate!(
result,
PE::PointEvaluator,
x
)
Evaluates the PointEvaluator at the specified coordinates x. (To do so it internally calls CellFinder to find the cell and the barycentric coordinates of x and calls evaluate_bary!.)
ExtendableFEMBase.evaluate_bary!
— Methodfunction evaluate_bary!(
result,
PE::PointEvaluator,
xref,
item
)
Evaluates the PointEvaluator at the specified reference coordinates in the cell with the specified item number.
ExtendableFEMBase.initialize!
— Methodfunction initialize!(
O::PointEvaluator,
sol;
time = 0,
kwargs...)
Initializes the PointEvaluator for the specified solution.
ExtendableFEMBase.initialize!
— Methodfunction initialize!(
O::SegmentIntegrator{T, UT},
sol;
time = 0,
kwargs...)
Initializes the SegmentIntegrator for the specified solution.
ExtendableFEMBase.integrate!
— Methodintegrate!(
integral4items::AbstractArray{T},
grid::ExtendableGrids.ExtendableGrid{Tv, Ti},
AT::Type{<:ExtendableGrids.AssemblyType},
integrand;
offset,
bonus_quadorder,
quadorder,
time,
items,
force_quadrature_rule,
kwargs...
)
Integration that writes result on every item into integral4items.
ExtendableFEMBase.integrate
— Methodintegrate(
grid::ExtendableGrids.ExtendableGrid,
AT::Type{<:ExtendableGrids.AssemblyType},
integrand!,
resultdim::Int64;
T,
kwargs...
) -> Union{Float64, Vector{Float64}}
Integration that returns total integral.
ExtendableFEMBase.integrate_segment!
— Methodfunction integrate_segment!(
result::Array{T,1},
SI::SegmentIntegrator,
w::Array{Array{T,1},1},
b::Array{Array{T,1},1},
item
) where {T}
Integrate a segment with world coordinates w and barycentric coordinates b in the cell with the given item number.
ExtendableFEMBase.lazy_interpolate!
— Methodfunction lazy_interpolate!(
target::FEVectorBlock{T1,Tv,Ti},
source::FEVectorBlock{T2,Tv,Ti};
operator = Identity,
postprocess = nothing,
xtrafo = nothing,
items = [],
not_in_domain_value = 1e30,
eps = 1e-13,
use_cellparents::Bool = false) where {T1,T2,Tv,Ti}
Interpolates (operator-evaluations of) the given finite element function into the finite element space assigned to the target FEVectorBlock. (Currently not the most efficient way as it is based on the PointEvaluation pattern and cell search. If CellParents are available in the grid components of the target grid, these parent cell information can be used to improve the search. To activate this put 'use_cellparents' = true). By some given kernel function that is conforming to the interface
kernel!(result, input, qpinfo)
the operator evaluation (=input) can be further postprocessed. The qpinfo argument allows to access information at the current quadrature point.
Note: discontinuous quantities at vertices of the target grid will be evaluted in the first found cell of the source grid. No averaging is performed. With eps the tolerances of the cell search via ExtendableGrids.CellFinder can be steered.
ExtendableFEMBase.ldrdmatmul
— Methodldrdmatmul(
a1::AbstractArray{Tv, 1},
a2::AbstractArray{Tv, 1},
B::ExtendableSparse.ExtendableSparseMatrix{Tv, Ti<:Integer},
b1::AbstractArray{Tv, 1},
b2::AbstractArray{Tv, 1};
factor
) -> Any
Computes vector'-matrix-vector product (a1-a2)'B(b1-b2).
ExtendableFEMBase.lrmatmul
— Methodlrmatmul(
a::AbstractArray{Tv, 1},
B::ExtendableSparse.ExtendableSparseMatrix{Tv, Ti<:Integer},
b::AbstractArray{Tv, 1};
factor
) -> Any
Computes vector'-matrix-vector product a'Bb.
ExtendableFEMBase.nbcols
— Methodnbcols(
_::FEMatrix{TvM, TiM, TvG, TiG, nbrow, nbcol, nbtotal}
) -> Any
Gives the number of FEMatrixBlocks in each row.
ExtendableFEMBase.nbrows
— Methodnbrows(
_::FEMatrix{TvM, TiM, TvG, TiG, nbrow, nbcol, nbtotal}
) -> Any
Gives the number of FEMatrixBlocks in each column.
ExtendableFEMBase.nodevalues!
— Functionfunction nodevalues!(
target::AbstractArray{<:Real,2},
source::FEVectorBlock,
operator::Type{<:AbstractFunctionOperator} = Identity;
regions::Array{Int,1} = [0],
abs::Bool = false,
factor = 1,
cellwise = false, # return cellwise nodevalues ncells x nnodes_on_cell
target_offset::Int = 0, # start to write into target after offset
zero_target::Bool = true, # target vector is zeroed
continuous::Bool = false)
Evaluates the finite element function with the coefficient vector source and the specified FunctionOperator at all the nodes of the (specified regions of the) grid and writes the values into target. Discontinuous (continuous = false) quantities are evaluated in all neighbouring cells of each node and then averaged. Continuous (continuous = true) quantities are only evaluated once at each node.
ExtendableFEMBase.nodevalues!
— Methodfunction nodevalues!(
target::AbstractArray{<:Real,2},
source::AbstractArray{T,1},
FE::FESpace{Tv,Ti,FEType,AT},
operator::Type{<:AbstractFunctionOperator} = Identity;
regions::Array{Int,1} = [0],
abs::Bool = false,
factor = 1,
target_offset::Int = 0, # start to write into target after offset
zero_target::Bool = true, # target vector is zeroed
continuous::Bool = false)
Evaluates the finite element function with the coefficient vector source (interpreted as a coefficient vector for the FESpace FE) and the specified FunctionOperator at all the nodes of the (specified regions of the) grid and writes the values into target. Discontinuous (continuous = false) quantities are evaluated in all neighbouring cells of each node and then averaged. Continuous (continuous = true) quantities are only evaluated once at each node.
ExtendableFEMBase.nodevalues
— Methodfunction nodevalues(
source::FEVectorBlock,
operator::Type{<:AbstractFunctionOperator} = Identity;
regions::Array{Int,1} = [0],
abs::Bool = false,
factor = 1,
nodes = [],
cellwise = false, # return cellwise nodevalues ncells x nnodes_on_cell (only if nodes == [])
target_offset::Int = 0, # start to write into target after offset
zero_target::Bool = true, # target vector is zeroed
continuous::Bool = false)
Evaluates the finite element function with the coefficient vector source and the specified FunctionOperator at the specified list of nodes of the grid (default = all nodes) and writes the values in that order into target. Nodes that are not part of the specified regions (default = all regions) are set to zero. Discontinuous (continuous = false) quantities are evaluated in all neighbouring cells of each node and then averaged. Continuous (continuous = true) quantities are only evaluated once at each node.
ExtendableFEMBase.nodevalues_subset!
— Methodfunction nodevalues_subset!(
target::AbstractArray{<:Real,2},
source::AbstractArray{T,1},
FE::FESpace{Tv,Ti,FEType,AT},
operator::Type{<:AbstractFunctionOperator} = Identity;
regions::Array{Int,1} = [0],
abs::Bool = false,
factor = 1,
nodes = [],
target_offset::Int = 0, # start to write into target after offset
zero_target::Bool = true, # target vector is zeroed
continuous::Bool = false)
Evaluates the finite element function with the coefficient vector source (interpreted as a coefficient vector for the FESpace FE) and the specified FunctionOperator at the specified list of nodes of the grid and writes the values in that order into target. Node values for nodes that are not part of the specified regions (default = all regions) are set to zero. Discontinuous (continuous = false) quantities are evaluated in all neighbouring cells of each node and then averaged. Continuous (continuous = true) quantities are only evaluated once at each node.
ExtendableFEMBase.nodevalues_view
— Methodfunction nodevalues_view(
source::FEVectorBlock,
operator::Type{<:AbstractFunctionOperator} = Identity)
Returns a vector of views of the nodal values of the source block (currently works for unbroken H1-conforming elements) that directly accesses the coefficients.
ExtendableFEMBase.ref_integrate!
— Methodref_integrate!(
integral::AbstractArray,
EG::Type{<:ExtendableGrids.AbstractElementGeometry},
order::Int64,
integrand::Function
)
Integration for reference basis functions on reference domains (merely for testing stuff).
Note: area of reference geometry is not multiplied
ExtendableFEMBase.unicode_gridplot
— Methodfunction unicode_gridplot(
xgrid::ExtendableGrid;
title = "gridplot",
resolution = (40,20),
color = (200,200,200),
bface_color = (255,0,0),
CanvasType = BrailleCanvas,
plot_based = ON_CELLS, # or ON_FACES/ON_EDGES
kwargs...
Plots the grid on a UnicodePlots canvas (default: BrailleCanvas) by drawing all edges in the triangulation.
ExtendableFEMBase.unicode_scalarplot
— Methodfunction unicode_scalarplot(
u::FEVectorBlock;
components = 1:get_ncomponents(u),
abs = false,
resolution = (30,30),
colormap = :viridis,
title = u.name,
kwargs...)
Plots components of the finite element function in the FEVectorBlock u by using a lazy_interpolate! onto a coarse uniform mesh and UnicodePlots.jl lineplot or heatmap for 1D or 2D, respectively.
In 1D all components all plotted in the same lineplot, while in 2D all components are plotted in a separate heatmap.
If abs = true, only the absolute value over the components is plotted.
ExtendableFEMBase.update_basis!
— Methodfunction update_basis!(FEBE::FEEvaluator, item::Integer)
Sets FEBE.citem[] = item and updates the basis.
ExtendableFEMBase.update_basis!
— Methodfunction update_basis!(FEBE::SingleFEEvaluator)
Updates the basis for the current item FEBE.citem[].
ExtendableGrids.interpolate!
— Methodfunction ExtendableGrids.interpolate!(target::FEVectorBlock,
source::Function;
items = [],
bonus_quadorder = 0,
time = 0,
kwargs...)
Interpolates the given source function into the finite element space assigned to the target FEVectorBlock.
The source functions should adhere to the interface
source!(result, qpinfo)
The qpinfo argument communicates vast information of the current quadrature/evaluation point.
The bonus_quadorder argument can be used to steer the quadrature order of integrals that needs to be computed for the interpolation (the default quadrature order corressponds to the polynomial order of the finite element).
ExtendableGrids.interpolate!
— Methodfunction ExtendableGrids.interpolate!(target::FEVectorBlock,
AT::Type{<:AssemblyType},
source!::Function;
items = [],
bonus_quadorder = 0,
time = 0,
kwargs...)
Interpolates the given source into the finite elements space assigned to the target FEVectorBlock with the specified AssemblyType (usualy ON_CELLS).
The source functions should adhere to the interface
source!(result, qpinfo)
The qpinfo argument communicates vast information of the current quadrature/evaluation point.
The bonus_quadorder argument can be used to steer the quadrature order of integrals that needs to be computed for the interpolation (the default quadrature order corressponds to the polynomial order of the finite element).
LinearAlgebra.dot
— Methoddot(a::FEVectorBlock{T}, b::FEVectorBlock{T}) -> Any
Scalar product between two FEVEctorBlocks