Gridap.ReferenceFEs
Gridap.ReferenceFEs
— ModuleThe exported names are
BDM
BDMRefFE
Bezier
BezierRefFE
CDConformity
CONT
Conformity
ContraVariantPiolaMap
CurlConformity
DISC
DivConforming
DivConformity
Dof
ExtrusionPolytope
GeneralPolytope
GenericLagrangianRefFE
GenericQuadrature
GenericRefFE
GradConformity
H1Conformity
HEX
HEX8
HEX_AXIS
INVALID_PERM
L2Conformity
Lagrangian
LagrangianDofBasis
LagrangianRefFE
ModalC0
ModalC0RefFE
MomentBasedDofBasis
Nedelec
NedelecRefFE
PYRAMID
Polygon
Polyhedron
Polytope
QUAD
QUAD4
Quadrature
QuadratureName
RaviartThomas
RaviartThomasRefFE
ReferenceFE
ReferenceFEName
SEG2
SEGMENT
SerendipityRefFE
TET
TET4
TET_AXIS
TRI
TRI3
VERTEX
VERTEX1
WEDGE
bdm
bezier
compress_cell_data
compute_cell_to_modalC0_reffe
compute_face_orders
compute_lagrangian_reffaces
compute_monomial_basis
compute_nodes
compute_own_nodes
compute_own_nodes_permutations
compute_shapefuns
duffy
evaluate!
evaluate_dof
expand_cell_data
get_bounding_box
get_coordinates
get_dimrange
get_dimranges
get_dof_basis
get_dof_to_comp
get_dof_to_node
get_edge_tangent
get_extrusion
get_face_coordinates
get_face_dimranges
get_face_dofs
get_face_moments
get_face_nodes
get_face_nodes_dofs
get_face_own_dofs
get_face_own_dofs_permutations
get_face_own_nodes
get_face_own_nodes_permutations
get_face_type
get_face_vertex_permutations
get_face_vertices
get_facedims
get_faces
get_facet_normal
get_facet_orientations
get_name
get_node_and_comp_to_dof
get_node_coordinates
get_nodes
get_offset
get_offsets
get_order
get_orders
get_own_dofs_permutations
get_own_nodes_permutations
get_polytope
get_prebasis
get_reffaces
get_shapefuns
get_vertex_coordinates
get_vertex_node
get_vertex_permutations
get_weights
is_P
is_Q
is_S
is_first_order
is_n_cube
is_simplex
lagrangian
linear_combination
modalC0
nedelec
num_cell_dims
num_dims
num_dofs
num_edges
num_faces
num_facets
num_nodes
num_point_dims
num_points
num_vertices
raviart_thomas
return_cache
return_type
simplexify
strang
tensor_product
test_dof
test_dof_array
test_lagrangian_reference_fe
test_polytope
test_quadrature
test_reference_fe
Gridap.ReferenceFEs.HEX
— Constantconst HEX = Polytope(HEX_AXIS,HEX_AXIS,HEX_AXIS)
Gridap.ReferenceFEs.HEX8
— Constantconst HEX8 = LagrangianRefFE(Float64,HEX,1)
Gridap.ReferenceFEs.HEX_AXIS
— ConstantConstant to be used in order to indicate a hex-like extrusion axis.
Gridap.ReferenceFEs.INVALID_PERM
— ConstantConstant of type Int
used to signal that a permutation is not valid.
Gridap.ReferenceFEs.PYRAMID
— Constantconst PYRAMID = Polytope(HEX_AXIS,HEX_AXIS,TET_AXIS)
Gridap.ReferenceFEs.QUAD
— Constantconst QUAD = Polytope(HEX_AXIS,HEX_AXIS)
Gridap.ReferenceFEs.QUAD4
— Constantconst QUAD4 = LagrangianRefFE(Float64,QUAD,1)
Gridap.ReferenceFEs.SEG2
— Constantconst SEG2 = LagrangianRefFE(Float64,SEGMENT,1)
Gridap.ReferenceFEs.SEGMENT
— Constantconst SEGMENT = Polytope(HEX_AXIS)
Gridap.ReferenceFEs.TET
— Constantconst TET = Polytope(TET_AXIS,TET_AXIS,TET_AXIS)
Gridap.ReferenceFEs.TET4
— Constantconst TET4 = LagrangianRefFE(Float64,TET,1)
Gridap.ReferenceFEs.TET_AXIS
— ConstantConstant to be used in order to indicate a tet-like extrusion axis.
Gridap.ReferenceFEs.TRI
— Constantconst TRI = Polytope(TET_AXIS,TET_AXIS)
Gridap.ReferenceFEs.TRI3
— Constantconst TRI3 = LagrangianRefFE(Float64,TRI,1)
Gridap.ReferenceFEs.VERTEX
— Constantconst VERTEX = Polytope()
Gridap.ReferenceFEs.VERTEX1
— Constantconst VERTEX1 = LagrangianRefFE(Float64,VERTEX,1)
Gridap.ReferenceFEs.WEDGE
— Constantconst WEDGE = Polytope(TET_AXIS,TET_AXIS,HEX_AXIS)
Gridap.ReferenceFEs.Conformity
— MethodConformity(reffe::ReferenceFE) -> Conformity
Gridap.ReferenceFEs.ExtrusionPolytope
— Typestruct ExtrusionPolytope{D} <: Polytope{D}
extrusion::Point{D,Int}
# + private fields
end
Concrete type for polytopes that can be represented with an "extrusion" tuple. The underlying extrusion is available in the field extrusion
. Instances of this type can be obtained with the constructors
Gridap.ReferenceFEs.ExtrusionPolytope
— MethodExtrusionPolytope(extrusion::Int...)
Generates an ExtrusionPolytope
from the tuple extrusion
. The values in extrusion
are either equal to the constant HEX_AXIS
or the constant TET_AXIS
.
Examples
Creating a quadrilateral, a triangle, and a wedge
using Gridap.ReferenceFEs
quad = ExtrusionPolytope(HEX_AXIS,HEX_AXIS)
tri = ExtrusionPolytope(TET_AXIS,TET_AXIS)
wedge = ExtrusionPolytope(TET_AXIS,TET_AXIS,HEX_AXIS)
println(quad == QUAD)
println(tri == TRI)
println(wedge == WEDGE)
# output
true
true
true
Gridap.ReferenceFEs.GeneralPolytope
— Typestruct GeneralPolytope{D,Dp,Tp} <: Polytope{D}
The GeneralPolytope
is definded defined by a set of vertices and a rototation system (a planar oriented graph). This polytopal representation can represent any polytope in 2 and 3 dimensions.
In 2 dimensions (Polygon
), the representation of the polygon is a closed polyline.
In 3 dimensions (Polyhedron
), the rotation system generates the connectivities, each facet is a closed cycle of the graph. This construction allows complex geometrical operations, e.g., intersecting polytopes by halfspaces. See also,
K. Sugihara, "A robust and consistent algorithm for intersecting convex polyhedra", Comput. Graph. Forum 13 (3) (1994) 45–54, doi: 10.1111/1467-8659.1330045
D. Powell, T. Abel, "An exact general remeshing scheme applied to physically conservative voxelization", J. Comput. Phys. 297 (Sept. 2015) 340–356, doi: [10.1016/j.jcp.2015.05.022](https://doi.org/10.1016/j.jcp.2015.05.022.
S. Badia, P. A. Martorell, F. Verdugo. "Geometrical discretisations for unfitted finite elements on explicit boundary representations", J.Comput. Phys. 460 (2022): 111162. doi: 10.1016/j.jcp.2022.111162
Gridap.ReferenceFEs.GeneralPolytope
— MethodGeneralPolytope{D}(vertices,graph;kwargs...)
Constructor of a GeneralPolytope
that generates a polytope of D dimensions with the given vertices
and graph
of connectivities.
Gridap.ReferenceFEs.GenericLagrangianRefFE
— Typestruct GenericLagrangianRefFE{C,D} <: LagrangianRefFE{D} reffe::GenericRefFE{C,D} face_nodes::Vector{Vector{Int}} end
Gridap.ReferenceFEs.GenericQuadrature
— Typestruct GenericQuadrature{D,T,C <: AbstractVector{Point{D,T}},W <: AbstractVector{T}} <: Quadrature{D,T}
coordinates::Vector{Point{D,T}}
weights::Vector{T}
name::String
end
Gridap.ReferenceFEs.GenericRefFE
— Typestruct GenericRefFE{T,D} <: ReferenceFE{D}
# + private fields
end
This type is a materialization of the ReferenceFE
interface. That is, it is a struct
that stores the values of all abstract methods in the ReferenceFE
interface. This type is useful to build reference FEs from the underlying ingredients without the need to create a new type.
Note that some fields in this struct
are type unstable deliberately in order to simplify the type signature. Don't access them in computationally expensive functions, instead extract the required fields before and pass them to the computationally expensive function.
Gridap.ReferenceFEs.GenericRefFE
— Method GenericRefFE{T}(
ndofs::Int,
polytope::Polytope{D},
prebasis::AbstractVector{<:Field},
dofs::AbstractVector{<:Dof},
conformity::Conformity,
metadata,
face_dofs::Vector{Vector{Int}},
shapefuns::AbstractVector{<:Field}=compute_shapefuns(dofs,prebasis)) where {T,D}
Constructs a GenericRefFE
object with the provided data.
Gridap.ReferenceFEs.GenericRefFE
— Method GenericRefFE{T}(
ndofs::Int,
polytope::Polytope{D},
prebasis::AbstractVector{<:Field},
predofs::AbstractVector{<:Dof},
conformity::Conformity,
metadata,
face_dofs::Vector{Vector{Int}},
shapefuns::AbstractVector{<:Field},
dofs::AbstractVector{<:Dof}=compute_dofs(predofs,shapefuns)) where {T,D}
Constructs a GenericRefFE
object with the provided data.
Gridap.ReferenceFEs.LagrangianDofBasis
— Typestruct LagrangianDofBasis{P,V} <: AbstractArray{<:Dof}
nodes::Vector{P}
dof_to_node::Vector{Int}
dof_to_comp::Vector{Int}
node_and_comp_to_dof::Vector{V}
end
Type that implements a Lagrangian dof basis.
Fields:
nodes::Vector{P}
vector of points (P<:Point
) storing the nodal coordinatesnode_and_comp_to_dof::Vector{V}
vector such thatnode_and_comp_to_dof[node][comp]
returns the dof associated with nodenode
and the componentcomp
in the typeV
.dof_to_node::Vector{Int}
vector of integers such thatdof_to_node[dof]
returns the node id associated with dof iddof
.dof_to_comp::Vector{Int}
vector of integers such thatdof_to_comp[dof]
returns the component id associated with dof iddof
.
Gridap.ReferenceFEs.LagrangianDofBasis
— MethodLagrangianDofBasis(::Type{T},nodes::Vector{<:Point}) where T
Creates a LagrangianDofBasis
for fields of value type T
associated with the vector of nodal coordinates nodes
.
Gridap.ReferenceFEs.LagrangianRefFE
— Typeabstract type LagrangianRefFE{D} <: ReferenceFE{D}
Abstract type representing a Lagrangian reference FE. Lagrangian in the sense that get_dof_basis
returns an instance of LagrangianDofBasis
. The interface for this type is defined with the methods of ReferenceFE
plus the following ones
Gridap.ReferenceFEs.LagrangianRefFE
— MethodLagrangianRefFE(::Type{T},p::Polytope,orders) where T
LagrangianRefFE(::Type{T},p::Polytope,order::Int) where T
Builds a LagrangianRefFE
object on top of the given polytope. T
is the type of the value of the approximation space (e.g., T=Float64
for scalar-valued problems, T=VectorValue{N,Float64}
for vector-valued problems with N
components). The arguments order
or orders
are for the polynomial order of the resulting space, which allows isotropic or anisotropic orders respectively (provided that the cell topology allows the given anisotropic order). The argument orders
should be an indexable collection of D
integers (e.g., a tuple or a vector), being D
the number of space dimensions.
In order to be able to use this function, the type of the provided polytope p
has to implement the following additional methods. They have been implemented for ExtrusionPolytope
in the library. They need to be implemented for new polytope types in order to build Lagangian reference elements on top of them.
compute_monomial_basis(::Type{T},p::Polytope,orders) where T
compute_own_nodes(p::Polytope,orders)
compute_face_orders(p::Polytope,face::Polytope,iface::Int,orders)
The following methods are also used in the construction of the LagrangianRefFE
object. A default implementation of them is available in terms of the three previous methods. However, the user can also implement them for new polytope types increasing customization possibilities.
Gridap.ReferenceFEs.LinearCombinationDofVector
— Typestruct LinearCombinationDofVector{T} <: AbstractVector{Dof}
change_of_basis::Matrix{T}
dof_basis::AbstractVector{<:Dof}
end
Type that implements a dof basis (a) as the linear combination of a dof basis (b). The dofs are first evaluated at dof basis (b) (field dof_basis
) and the dof values are next mapped to dof basis (a) applying a change of basis (field change_of_basis
).
Fields:
change_of_basis::Matrix{T}
the matrix of the change from dof basis (b) to (a)dof_basis::AbstractVector{<:Dof}
A type representing dof basis (b)
Gridap.ReferenceFEs.Polygon
— TypePolygon = GeneralPolytope{2}
A polygon is a GeneralPolytope
in 2 dimensions.
Gridap.ReferenceFEs.Polyhedron
— TypePolyhedron = GeneralPolytope{3}
A polyhedron is a GeneralPolytope
in 3 dimensions.
Gridap.ReferenceFEs.Polytope
— Typeabstract type Polytope{D} <: GridapType
Abstract type representing a polytope (i.e., a polyhedron in arbitrary dimensions). D
is the environment dimension (typically, 0, 1, 2, or 3). This type parameter is needed since there are functions in the Polytope
interface that return containers with Point{D}
objects. We adopt the usual nomenclature for polytope-related objects. All objects in a polytope (from vertices to the polytope itself) are called n-faces or simply faces. The notation n-faces is used only when it is needed to refer to the object dimension n. Otherwise we simply use face. In addition, we say
- vertex (pl. vertices): for 0-faces
- edge: for 1-faces
- facet: for (
D-1
)-faces
The Polytope
interface is defined by overloading the following functions
get_faces(p::Polytope)
get_dimranges(p::Polytope)
Polytope{N}(p::Polytope,faceid::Integer) where N
get_vertex_coordinates(p::Polytope)
(==)(a::Polytope{D},b::Polytope{D}) where D
And optionally these ones:
get_edge_tangent(p::Polytope)
get_facet_normal(p::Polytope)
get_facet_orientations(p::Polytope)
get_vertex_permutations(p::Polytope)
is_n_cube(p::Polytope)
is_simplex(p::Polytope)
simplexify(p::Polytope)
The interface can be tested with the function
Gridap.ReferenceFEs.Polytope
— MethodPolytope(extrusion::Int...)
Equivalent to ExtrusionPolytope(extrusion...)
Gridap.ReferenceFEs.Polytope
— MethodPolytope{N}(p::Polytope,faceid::Integer) where N
Returns a Polytope{N}
object representing the "reference" polytope of the N
-face with id faceid
. The value faceid
refers to the numeration restricted to the dimension N
(it starts with 1 for the first N
-face).
Gridap.ReferenceFEs.Quadrature
— Typeabstract type Quadrature{D,T} <: GridapType end
-get_coordinates(q::Quadrature)
-get_weights(q::Quadrature)
-test_quadrature
Gridap.ReferenceFEs.Quadrature
— MethodQuadrature(polytope::Polytope{D},degree) where D
Gridap.ReferenceFEs.ReferenceFE
— Typeabstract type ReferenceFE{D} <: GridapType
Abstract type representing a Reference finite element. D
is the underlying coordinate space dimension. We follow the Ciarlet definition. A reference finite element is defined by a polytope (cell topology), a basis of an interpolation space of top of this polytope (denoted here as the prebasis), and a basis of the dual of this space (i.e. the degrees of freedom). From this information one can compute the shape functions (i.e, the canonical basis of w.r.t. the degrees of freedom) with a simple change of basis. In addition, we also encode in this type information about how the interpolation space in a reference finite element is "glued" with neighbors in order to build conforming cell-wise spaces.
The ReferenceFE
interface is defined by overloading these methods:
num_dofs(reffe::ReferenceFE)
get_polytope(reffe::ReferenceFE)
get_prebasis(reffe::ReferenceFE)
get_dof_basis(reffe::ReferenceFE)
Conformity(reffe::ReferenceFE)
get_face_own_dofs(reffe::ReferenceFE,conf::Conformity)
get_face_own_dofs_permutations(reffe::ReferenceFE,conf::Conformity)
get_face_dofs(reffe::ReferenceFE)
The interface is tested with
Gridap.ReferenceFEs.ReferenceFE
— MethodReferenceFE{N}(reffe::GenericLagrangianRefFE{GradConformity},iface::Integer) where N
Base.:==
— Method(==)(a::LagrangianRefFE{D}, b::LagrangianRefFE{D}) where D
Base.:==
— Method(==)(a::Polytope{D},b::Polytope{D}) where D
Returns true
if the polytopes a
and b
are equivalent. Otherwise, it returns false
. Note that the operator ==
returns false
by default for polytopes of different dimensions. Thus, this function has to be overloaded only for the case of polytopes a
and b
of same dimension.
Base.isopen
— Methodisopen(p::GeneralPolytope) -> Bool
In return whether the polytope is watter tight or not.
Gridap.Polynomials.get_order
— Methodget_order(reffe::LagrangianRefFE)
Gridap.Polynomials.get_orders
— Methodget_orders(reffe::LagrangianRefFE)
Gridap.ReferenceFEs.BDMRefFE
— MethodBDMRefFE(::Type{et},p::Polytope,order::Integer) where et
The order
argument has the following meaning: the divergence of the functions in this basis is in the P space of degree order-1
.
Gridap.ReferenceFEs.ModalC0RefFE
— MethodModalC0RefFE(::Type{T},p::Polytope{D},orders) where {T,D}
Returns an instance of GenericRefFE{ModalC0}
representing a ReferenceFE with Modal C0-continuous shape functions (multivariate scalar-valued, vector-valued, or tensor-valued, iso- or aniso-tropic).
For more details about the shape functions, see Section 1.1.5 in
Ern, A., & Guermond, J. L. (2013). Theory and practice of finite elements (Vol. 159). Springer Science & Business Media.
and references therein.
The constructor is only implemented for for n-cubes and the minimum order in any dimension must be greater than one. The DoFs are numbered by n-faces in the same way as with CLagrangianRefFEs.
Gridap.ReferenceFEs.NedelecRefFE
— MethodNedelecRefFE(::Type{et},p::Polytope,order::Integer) where et
The order
argument has the following meaning: the curl of the functions in this basis is in the Q space of degree order
.
Gridap.ReferenceFEs.RaviartThomasRefFE
— MethodRaviartThomasRefFE(::Type{et},p::Polytope,order::Integer) where et
The order
argument has the following meaning: the divergence of the functions in this basis is in the Q space of degree order
.
Gridap.ReferenceFEs.SerendipityRefFE
— MethodSerendipityRefFE(::Type{T},p::Polytope,order::Int) where T
SerendipityRefFE(::Type{T},p::Polytope,orders::Tuple) where T
Returns an instance of LagrangianRefFE
, whose underlying approximation space is the serendipity space of order order
. Implemented for order from 1 to 4. The type of the polytope p
has to implement all the queries detailed in the constructor LagrangianRefFE(::Type{T},p::Polytope{D},orders) where {T,D}
.
Examples
using Gridap.ReferenceFEs
order = 2
reffe = SerendipityRefFE(Float64,QUAD,order)
println( num_dofs(reffe) )
# output
8
Gridap.ReferenceFEs.check_polytope_graph
— Methodcheck_polytope_graph(p::GeneralPolytope) -> Bool
It checks whether the graph is well-constructed. The graph must be oriented and planar.
Gridap.ReferenceFEs.compute_dofs
— Methodcompute_dofs(predofs,shapefuns)
Helper function used to compute the dof basis associated with the dof basis predofs
and the basis shapefuns
.
It is equivalent to
change = inv(evaluate(predofs,shapefuns))
linear_combination(change,predofs) # i.e. transpose(change)*predofs
Gridap.ReferenceFEs.compute_face_orders
— Methodcompute_face_orders(p::Polytope,face::Polytope,iface::Int,orders)
Returns a vector or a tuple with the order per direction at the face face
of the polytope p
when restricting the order per direction orders
to this face. iface
is the face id of face
in the numeration restricted to the face dimension.
Gridap.ReferenceFEs.compute_lagrangian_reffaces
— Methodcompute_lagrangian_reffaces(::Type{T},p::Polytope,orders) where T
Returns a tuple of length D
being the number of space dimensions. The entry d+1
of this tuple contains a vector of LagrangianRefFE
one for each face of dimension d
on the boundary of the polytope.
Gridap.ReferenceFEs.compute_monomial_basis
— Methodcompute_monomial_basis(::Type{T},p::Polytope,orders) where T -> MonomialBasis
Returns the monomial basis of value type T
and order per direction described by orders
on top of the polytope p
.
Gridap.ReferenceFEs.compute_nodes
— Methodcompute_nodes(p::Polytope,orders)
When called
node_coords, face_own_nodes = compute_nodes(p,orders)
Returns node_coords
, the nodal coordinates of all the Lagrangian nodes associated with the order per direction orders
, and face_own_nodes
, being a vector of vectors indicating which nodes are owned by each of the faces of the polytope p
.
Gridap.ReferenceFEs.compute_own_nodes
— Methodcompute_own_nodes(p::Polytope{D},orders) where D -> Vector{Point{D,Float64}}
Returns the coordinates of the nodes owned by the interior of the polytope associated with a Lagrangian space with the order per direction described by orders
.
Gridap.ReferenceFEs.compute_own_nodes_permutations
— Methodcompute_own_nodes_permutations(
p::Polytope, own_nodes_coordinates) -> Vector{Vector{Int}}
Returns a vector of vectors with the permutations of the nodes owned by the interior of the polytope.
Gridap.ReferenceFEs.compute_shapefuns
— Methodcompute_shapefuns(dofs,prebasis)
Helper function used to compute the shape function basis associated with the dof basis dofs
and the basis prebasis
.
It is equivalent to
change = inv(evaluate(dofs,prebasis))
linear_combination(change,prebasis) # i.e. transpose(change)*prebasis
Gridap.ReferenceFEs.expand_cell_data
— MethodGridap.ReferenceFEs.get_bounding_box
— Methodget_bounding_box(p::Polytope{D}) where D
Gridap.ReferenceFEs.get_coordinates
— Methodget_coordinates(q::Quadrature)
Gridap.ReferenceFEs.get_dimrange
— Methodget_dimrange(p::Polytope,d::Integer)
Equivalent to
get_dimranges(p)[d+1]
Gridap.ReferenceFEs.get_dimranges
— Methodget_dimranges(p::Polytope) -> Vector{UnitRange{Int}}
Given a polytope p
it returns a vector of ranges. The entry d+1
in this vector contains the range of face ids for the faces of dimension d
.
Examples
using Gridap.ReferenceFEs
ranges = get_dimranges(SEGMENT)
println(ranges)
# output
UnitRange{Int}[1:2, 3:3]
Face ids for the vertices in the segment range from 1 to 2 (2 vertices), the face ids for edges in the segment range from 3 to 3 (only one edge with id 3).
Gridap.ReferenceFEs.get_dof_basis
— Methodget_dof_basis(reffe::ReferenceFE) -> Dof
Returns the underlying dof basis encoded in a Dof
object.
Gridap.ReferenceFEs.get_dof_to_comp
— Methodget_dof_to_comp(reffe::LagrangianRefFE)
Gridap.ReferenceFEs.get_dof_to_node
— Methodget_dof_to_node(reffe::LagrangianRefFE)
Gridap.ReferenceFEs.get_edge_tangent
— Methodget_edge_tangent(p::Polytope) -> Vector{VectorValue{D,Float64}}
Given a polytope p
, returns a vector of VectorValue
objects representing the unit tangent vectors to the polytope edges.
Gridap.ReferenceFEs.get_extrusion
— Methodget_extrusion(p::ExtrusionPolytope)
Equivalent to p.extrusion
.
Gridap.ReferenceFEs.get_face_coordinates
— Methodget_face_coordinates(p::Polytope)
get_face_coordinates(p::Polytope,d::Integer)
Gridap.ReferenceFEs.get_face_dimranges
— Methodget_face_dimranges(p::Polytope,d::Integer)
Gridap.ReferenceFEs.get_face_dofs
— Methodget_face_dofs(reffe::ReferenceFE,d::Integer)
Gridap.ReferenceFEs.get_face_dofs
— Methodget_face_dofs(reffe::ReferenceFE) -> Vector{Vector{Int}}
Returns a vector of vector that, for each face, stores the dofids in the closure of the face.
Gridap.ReferenceFEs.get_face_nodes
— Methodget_face_nodes(reffe::LagrangianRefFE,d::Integer)
Gridap.ReferenceFEs.get_face_nodes
— Methodget_face_nodes(reffe::LagrangianRefFE)
Gridap.ReferenceFEs.get_face_own_dofs
— Methodget_face_own_dofs(reffe::ReferenceFE,conf::Conformity,d::Integer)
Gridap.ReferenceFEs.get_face_own_dofs
— Methodget_face_own_dofs(reffe::ReferenceFE,conf::Conformity) -> Vector{Vector{Int}}
Gridap.ReferenceFEs.get_face_own_dofs
— Methodget_face_own_dofs(reffe::ReferenceFE,d::Integer)
Gridap.ReferenceFEs.get_face_own_dofs
— Methodget_face_own_dofs(reffe::ReferenceFE) -> Vector{Vector{Int}}
Gridap.ReferenceFEs.get_face_own_dofs_permutations
— Methodget_face_own_dofs_permutations(reffe::ReferenceFE,conf::Conformity,d::Integer)
Gridap.ReferenceFEs.get_face_own_dofs_permutations
— Methodget_face_own_dofs_permutations(reffe::ReferenceFE,conf::Conformity) -> Vector{Vector{Vector{Int}}}
Gridap.ReferenceFEs.get_face_own_dofs_permutations
— Methodget_face_own_dofs_permutations(reffe::ReferenceFE,d::Integer)
Gridap.ReferenceFEs.get_face_own_dofs_permutations
— Methodget_face_own_dofs_permutations(reffe::ReferenceFE) -> Vector{Vector{Vector{Int}}}
Gridap.ReferenceFEs.get_face_own_nodes
— Methodget_face_own_nodes(reffe::LagrangianRefFE,conf::Conformity,d::Integer)
Gridap.ReferenceFEs.get_face_own_nodes
— Methodget_face_own_nodes(reffe::LagrangianRefFE,conf::Conformity)
Gridap.ReferenceFEs.get_face_own_nodes
— Methodget_face_own_nodes(reffe::LagrangianRefFE,d::Integer)
Gridap.ReferenceFEs.get_face_own_nodes
— Methodget_face_own_nodes(reffe::LagrangianRefFE)
Gridap.ReferenceFEs.get_face_own_nodes_permutations
— Methodget_face_own_nodes_permutations(reffe::LagrangianRefFE,conf::Conformity,d::Integer)
Gridap.ReferenceFEs.get_face_own_nodes_permutations
— Methodget_face_own_nodes_permutations(reffe::LagrangianRefFE,conf::Conformity)
Gridap.ReferenceFEs.get_face_own_nodes_permutations
— Methodget_face_own_nodes_permutations(reffe::LagrangianRefFE,d::Integer)
Gridap.ReferenceFEs.get_face_own_nodes_permutations
— Methodget_face_own_nodes_permutations(reffe::LagrangianRefFE)
Gridap.ReferenceFEs.get_face_type
— Methodget_face_type(reffe::GenericLagrangianRefFE{GradConformity}, d::Integer) -> Vector{Int}
Gridap.ReferenceFEs.get_face_type
— Methodget_face_type(p::Polytope,d::Integer) -> Vector{Int}
Return a vector of integers denoting, for each face of dimension d
, an index to the vector get_reffaces(Polytope{d},p)
Examples
Get the unique polytopes for the facets of a wedge and identify of which type each face is.
using Gridap.ReferenceFEs
reffaces = get_reffaces(Polytope{2},WEDGE)
face_types = get_face_type(WEDGE,2)
println(reffaces)
println(face_types)
# output
Gridap.ReferenceFEs.ExtrusionPolytope{2}[TRI, QUAD]
[1, 1, 2, 2, 2]
The three first facets are of type 1
, i.e, QUAD
, and the last ones of type 2
, i.e., TRI
.
Gridap.ReferenceFEs.get_face_vertex_permutations
— Methodget_face_vertex_permutations(p::Polytope)
get_face_vertex_permutations(p::Polytope,d::Integer)
Gridap.ReferenceFEs.get_face_vertices
— Methodget_face_vertices(p::Polytope) -> Vector{Vector{Int}}
get_face_vertices(p::Polytope,dim::Integer) -> Vector{Vector{Int}}
Gridap.ReferenceFEs.get_facedims
— Methodget_facedims(p::Polytope) -> Vector{Int}
Given a polytope p
, returns a vector indicating the dimension of each face in the polytope
Examples
using Gridap.ReferenceFEs
dims = get_facedims(SEGMENT)
println(dims)
# output
[0, 0, 1]
The first two faces in the segment (the two vertices) have dimension 0 and the third face (the segment itself) has dimension 1
Gridap.ReferenceFEs.get_faces
— Methodget_faces(p::Polytope,dimfrom::Integer,dimto::Integer) -> Vector{Vector{Int}}
For dimfrom >= dimto
returns a vector that for each face of dimension dimfrom
stores a vector of the ids of faces of dimension dimto
on its boundary.
For dimfrom < dimto
returns a vector that for each face of dimfrom
stores a vector of the face ids of faces of dimension dimto
that touch it.
The numerations used in this function are the ones restricted to each dimension.
using Gridap.ReferenceFEs
edge_to_vertices = get_faces(QUAD,1,0)
println(edge_to_vertices)
vertex_to_edges_around = get_faces(QUAD,0,1)
println(vertex_to_edges_around)
# output
Array{Int,1}[[1, 2], [3, 4], [1, 3], [2, 4]]
Array{Int,1}[[1, 3], [1, 4], [2, 3], [2, 4]]
Gridap.ReferenceFEs.get_faces
— Methodget_faces(p::Polytope) -> Vector{Vector{Int}}
Given a polytope p
the function returns a vector of vectors defining the incidence relation of the faces in the polytope.
Each face in the polytope receives a unique integer id. The id 1 is assigned to the first 0-face. Consecutive increasing ids are assigned to the other 0-faces, then to 1-faces, and so on. The polytope itself receives the largest id which coincides with num_faces(p)
. For a face id iface
, get_faces(p)[iface]
is a vector of face ids, corresponding to the faces that are incident with the face labeled with iface
. That is, faces that are either on its boundary or the face itself. In this vector of incident face ids, faces are ordered by dimension, starting with 0-faces. Within each dimension, the labels are ordered in a consistent way with the polyope object for the face iface
itself.
Examples
using Gridap.ReferenceFEs
faces = get_faces(SEGMENT)
println(faces)
# output
Array{Int,1}[[1], [2], [1, 2, 3]]
The constant SEGMENT
is bound to a predefined instance of polytope that represents a segment. The face labels associated with a segment are [1,2,3]
, being 1
and 2
for the vertices and 3
for the segment itself. In this case, this function returns the vector of vectors [[1],[2],[1,2,3]]
meaning that vertex 1
is incident with vertex 1
(idem for vertex 2), and that the segment (id 3
) is incident with the vertices 1
and 2
and the segment itself.
Gridap.ReferenceFEs.get_facet_normal
— Methodget_facet_normal(p::Polytope) -> Vector{VectorValue{D,Float64}}
Given a polytope p
, returns a vector of VectorValue
objects representing the unit outward normal vectors to the polytope facets.
Gridap.ReferenceFEs.get_facet_orientations
— Methodget_facet_orientations(p::Polytope) -> Vector{Int}
Given a polytope p
returns a vector of integers of length num_facets(p)
. Facets, whose vertices are ordered consistently with the outwards normal vector, receive value 1
in this vector. Otherwise, facets receive value -1
.
Gridap.ReferenceFEs.get_graph
— Methodget_graph(p::GeneralPolytope) -> Vector{Vector{Int32}}
It returns the edge-vertex graph of the polytope p
.
Gridap.ReferenceFEs.get_metadata
— Methodget_metadata(p::GeneralPolytope)
It return the metadata stored in the polytope p
.
Gridap.ReferenceFEs.get_name
— Methodget_name(q::Quadrature)
Gridap.ReferenceFEs.get_node_and_comp_to_dof
— Methodget_node_and_comp_to_dof(reffe::LagrangianRefFE)
Gridap.ReferenceFEs.get_node_coordinates
— Methodget_node_coordinates(reffe::LagrangianRefFE)
Gridap.ReferenceFEs.get_offset
— Methodget_offset(p::Polytope,d::Integer)
Equivalent to get_offsets(p)[d+1]
.
Gridap.ReferenceFEs.get_offsets
— Methodget_offsets(p::Polytope) -> Vector{Int}
Given a polytope p
, it returns a vector of integers. The position in the d+1
entry in this vector is the offset that transforms a face id in the global numeration in the polytope to the numeration restricted to faces to dimension d
.
Examples
using Gridap.ReferenceFEs
offsets = get_offsets(SEGMENT)
println(offsets)
# output
[0, 2]
Gridap.ReferenceFEs.get_own_dofs_permutations
— Methodget_own_dofs_permutations(reffe::ReferenceFE,conf::Conformity)
Gridap.ReferenceFEs.get_own_dofs_permutations
— Methodget_own_dofs_permutations(reffe::ReferenceFE)
Gridap.ReferenceFEs.get_own_nodes_permutations
— Methodget_own_nodes_permutations(reffe::LagrangianRefFE,conf::Conformity)
Gridap.ReferenceFEs.get_own_nodes_permutations
— Methodget_own_nodes_permutations(reffe::LagrangianRefFE)
Gridap.ReferenceFEs.get_polytope
— Methodget_polytope(reffe::ReferenceFE) -> Polytope
Returns the underlying polytope object.
Gridap.ReferenceFEs.get_prebasis
— Methodget_prebasis(reffe::ReferenceFE) -> Field
Returns the underlying prebasis encoded as a Field
object.
Gridap.ReferenceFEs.get_reffaces
— Methodget_reffaces(::Type{Polytope{d}},p::Polytope) where d -> Vector{Polytope{d}}
Get a vector of the unique polytopes for the faces of dimension d
.
Examples
Get the unique polytopes for the facets of a wedge.
using Gridap.ReferenceFEs
reffaces = get_reffaces(Polytope{2},WEDGE)
println(reffaces)
# output
Gridap.ReferenceFEs.ExtrusionPolytope{2}[TRI, QUAD]
Gridap.ReferenceFEs.get_reffaces
— Methodget_reffaces(
::Type{ReferenceFE{d}},
reffe::GenericLagrangianRefFE{GradConformity}) where d -> Vector{GenericLagrangianRefFE{GradConformity,M,d}}
Gridap.ReferenceFEs.get_shapefuns
— Methodget_shapefuns(reffe::ReferenceFE) -> Field
Returns the basis of shape functions (i.e. the canonical basis) associated with the reference FE. The result is encoded as a Field
object.
Gridap.ReferenceFEs.get_vertex_coordinates
— Methodget_vertex_coordinates(p::Polytope) -> Vector{Point{D,Float64}}
Given a polytope p
return a vector of points representing containing the coordinates of the vertices.
Gridap.ReferenceFEs.get_vertex_node
— Methodget_vertex_node(reffe::LagrangianRefFE,conf::Conformity) -> Vector{Int}
Gridap.ReferenceFEs.get_vertex_node
— Methodget_vertex_node(reffe::LagrangianRefFE) -> Vector{Int}
Gridap.ReferenceFEs.get_vertex_permutations
— Methodget_vertex_permutations(p::Polytope) -> Vector{Vector{Int}}
Given a polytope p
, returns a vector of vectors containing all admissible permutations of the polytope vertices. An admissible permutation is one such that, if the vertices of the polytope are re-labeled according to this permutation, the resulting polytope preserves the shape of the original one.
Examples
using Gridap.ReferenceFEs
perms = get_vertex_permutations(SEGMENT)
println(perms)
# output
Array{Int,1}[[1, 2], [2, 1]]
The first admissible permutation for a segment is [1,2]
,i.e., the identity. The second one is [2,1]
, i.e., the first vertex is relabeled as 2
and the second vertex is relabeled as 1
.
Gridap.ReferenceFEs.get_weights
— Methodget_weights(q::Quadrature)
Gridap.ReferenceFEs.is_P
— Methodis_P(reffe::GenericLagrangianRefFE{GradConformity})
Gridap.ReferenceFEs.is_Q
— Methodis_Q(reffe::GenericLagrangianRefFE{GradConformity})
Gridap.ReferenceFEs.is_S
— Methodis_S(reffe::GenericLagrangianRefFE{GradConformity})
Gridap.ReferenceFEs.is_first_order
— Methodis_first_order(reffe::GenericLagrangianRefFE{GradConformity}) -> Bool
Gridap.ReferenceFEs.is_n_cube
— Methodis_n_cube(p::Polytope) -> Bool
Gridap.ReferenceFEs.is_simplex
— Methodis_simplex(p::Polytope) -> Bool
Gridap.ReferenceFEs.isactive
— Methodisactive(p::GeneralPolytope,vertex::Integer) -> Bool
It returns whether a vertex is connected to any other vertex.
Gridap.ReferenceFEs.next_corner!
— MethodSweep through the dim
-dimensional hypercube recursively, collecting all simplices.
We represent vertices as bit patterns. In dim
dimensions, the lowermost dim
bits are either zero or one. Interpreted as integer, this labels the vertices of the hypercube from the origin ("bottom") 0
to the diagonally opposite vertex ("top") 2^dim-1
.
Each simplex contains both the bottom vertex 0
as well as the top vertex 2^dim-1
. Its other vertices trace a path from the bottom to the top. The algorithm below finds all possible paths.
simplices
is the accumulator where the simplices are collected.vertices
is the current set of vertices as we sweep from the origin to the diagonally opposide vertex.corner
is the current corner.
Gridap.ReferenceFEs.num_cell_dims
— Methodnum_cell_dims(::Type{<:ReferenceFE{D}}) where D
num_cell_dims(reffe::ReferenceFE{D}) where D
Returns D
.
Gridap.ReferenceFEs.num_dims
— Methodnum_dims(::Type{<:Polytope{D}}) where D
num_dims(p::Polytope{D}) where D
Returns D
.
Gridap.ReferenceFEs.num_dims
— Methodnum_dims(::Type{<:ReferenceFE{D}}) where D
num_dims(reffe::ReferenceFE{D}) where D
Returns D
.
Gridap.ReferenceFEs.num_dims
— Methodnum_dims(::Quadrature{D}) where D where D
num_dims(::Type{<:Quadrature{D}}) where D
Gridap.ReferenceFEs.num_dofs
— Methodnum_dofs(reffe::ReferenceFE) -> Int
Returns the number of DOFs.
Gridap.ReferenceFEs.num_edges
— Methodnum_edges(p::Polytope)
Returns the number of edges in the polytope p
.
Gridap.ReferenceFEs.num_edges
— Methodnum_edges(reffe::ReferenceFE)
Gridap.ReferenceFEs.num_faces
— Methodnum_faces(p::Polytope,dim::Integer)
Returns the number of faces of dimension dim
in polytope p
.
Gridap.ReferenceFEs.num_faces
— Methodnum_faces(p::Polytope)
Returns the total number of faces in polytope p
(from vertices to the polytope itself).
Gridap.ReferenceFEs.num_faces
— Methodnum_faces(reffe::ReferenceFE)
num_faces(reffe::ReferenceFE,d::Integer)
Gridap.ReferenceFEs.num_facets
— Methodnum_facets(p::Polytope)
Returns the number of facets in the polytope p
.
Gridap.ReferenceFEs.num_facets
— Methodnum_facets(reffe::ReferenceFE)
Gridap.ReferenceFEs.num_nodes
— Methodnum_nodes(reffe::LagrangianRefFE)
Gridap.ReferenceFEs.num_point_dims
— Methodnum_point_dims(::Type{<:ReferenceFE{D}}) where D
num_point_dims(reffe::ReferenceFE{D}) where D
Returns D
.
Gridap.ReferenceFEs.num_point_dims
— Methodnum_point_dims(::Quadrature{D}) where D
num_point_dims(::Type{<:Quadrature{D}}) where D
Gridap.ReferenceFEs.num_points
— Methodnum_points(q::Quadrature)
Gridap.ReferenceFEs.num_vertices
— Methodnum_vertices(p::Polytope)
Returns the number of vertices in the polytope p
.
Gridap.ReferenceFEs.num_vertices
— Methodnum_vertices(reffe::ReferenceFE)
Gridap.ReferenceFEs.simplexify
— Methodsimplexify(p::Polytope) -> Tuple{Vector{Vector{Int}},Polytope}
Gridap.ReferenceFEs.simplexify_interior
— Methodsimplexify_interior(p::Polyhedron)
simplex_interior
computes a simplex partition of the volume inside the Polyhedron p
. It returns a vector of coordinates and an array of connectivitties.
Gridap.ReferenceFEs.simplexify_surface
— Methodsimplexify_surface(p::Polyhedron)
simplex_surface
computes a simplex partition of the surface bounding the Polyhedron p
. It returns a vector of coordinates and an array of connectivitties.
Gridap.ReferenceFEs.test_dof
— Methodtest_dof(dof,field,v;cmp::Function=(==))
Test that the Dof
interface is properly implemented for object dof
. It also checks if the object dof
when evaluated at the field field
returns the same value as v
. Comparison is made with the comp
function.
Gridap.ReferenceFEs.test_lagrangian_reference_fe
— Methodtest_lagrangian_reference_fe(reffe::LagrangianRefFE)
Gridap.ReferenceFEs.test_polytope
— Methodtest_polytope(p::Polytope{D}; optional::Bool=false) where D
Function that stresses out the functions in the Polytope
interface. It tests whether the function in the polytope interface are defined for the given object, and whether they return objects of the expected type. With optional=false
(the default), only the mandatory functions are checked. With optional=true
, the optional functions are also tested.
Gridap.ReferenceFEs.test_quadrature
— Methodtest_quadrature(q::Quadrature{D,T}) where {D,T}
Gridap.ReferenceFEs.test_reference_fe
— Methodtest_reference_fe(reffe::ReferenceFE{D}) where D
Test if the methods in the ReferenceFE
interface are defined for the object reffe
.