Grid & AbstractGrid
Grid
Ferrite.generate_grid
— Functiongenerate_grid(celltype::Cell, nel::NTuple, [left::Vec, right::Vec)
Return a Grid
for a rectangle in 1, 2 or 3 dimensions. celltype
defined the type of cells, e.g. Triangle
or Hexahedron
. nel
is a tuple of the number of elements in each direction. left
and right
are optional endpoints of the domain. Defaults to -1 and 1 in all directions.
Ferrite.Node
— TypeNode{dim, T}
A Node
is a point in space.
Fields
x::Vec{dim,T}
: stores the coordinates
Ferrite.CellIndex
— TypeA CellIndex
wraps an Int and corresponds to a cell with that number in the mesh
Ferrite.VertexIndex
— TypeA VertexIndex
wraps an (Int, Int) and defines a local vertex by pointing to a (cell, vert).
Ferrite.EdgeIndex
— TypeA EdgeIndex
wraps an (Int, Int) and defines a local edge by pointing to a (cell, edge).
Ferrite.FaceIndex
— TypeA FaceIndex
wraps an (Int, Int) and defines a local face by pointing to a (cell, face).
Ferrite.FacetIndex
— TypeA FacetIndex
wraps an (Int, Int) and defines a local facet by pointing to a (cell, facet).
Ferrite.Grid
— TypeGrid{dim, C<:AbstractCell, T<:Real} <: AbstractGrid}
A Grid
is a collection of Ferrite.AbstractCell
s and Ferrite.Node
s which covers the computational domain. Helper structures for applying boundary conditions or define subdomains are gathered in cellsets
, nodesets
, facetsets
, and vertexsets
.
Fields
cells::Vector{C}
: stores all cells of the gridnodes::Vector{Node{dim,T}}
: stores thedim
dimensional nodes of the gridcellsets::Dict{String, OrderedSet{Int}}
: maps aString
key to anOrderedSet
of cell idsnodesets::Dict{String, OrderedSet{Int}}
: maps aString
key to anOrderedSet
of global node idsfacetsets::Dict{String, OrderedSet{FacetIndex}}
: maps aString
to anOrderedSet
ofFacetIndex
vertexsets::Dict{String, OrderedSet{VertexIndex}}
: maps aString
key to anOrderedSet
ofVertexIndex
Utility Functions
Ferrite.getcells
— Functiongetcells(grid::AbstractGrid)
getcells(grid::AbstractGrid, v::Union{Int,Vector{Int}}
getcells(grid::AbstractGrid, setname::String)
Returns either all cells::Collection{C<:AbstractCell}
of a <:AbstractGrid
or a subset based on an Int
, Vector{Int}
or String
. Whereas the last option tries to call a cellset
of the grid
. Collection
can be any indexable type, for Grid
it is Vector{C<:AbstractCell}
.
Ferrite.getncells
— FunctionReturns the number of cells in the <:AbstractGrid
.
Ferrite.getnodes
— Functiongetnodes(grid::AbstractGrid)
getnodes(grid::AbstractGrid, v::Union{Int,Vector{Int}}
getnodes(grid::AbstractGrid, setname::String)
Returns either all nodes::Collection{N}
of a <:AbstractGrid
or a subset based on an Int
, Vector{Int}
or String
. The last option tries to call a nodeset
of the <:AbstractGrid
. Collection{N}
refers to some indexable collection where each element corresponds to a Node.
Ferrite.getnnodes
— FunctionReturns the number of nodes in the grid.
Ferrite.nnodes_per_cell
— FunctionReturns the number of nodes of the i
-th cell.
Ferrite.getcellset
— Functiongetcellset(grid::AbstractGrid, setname::String)
Returns all cells as cellid in the set with name setname
.
Ferrite.getnodeset
— Functiongetnodeset(grid::AbstractGrid, setname::String)
Returns all nodes as nodeid in the set with name setname
.
Ferrite.getfacetset
— Functiongetfacetset(grid::AbstractGrid, setname::String)
Returns all faces as FacetIndex
in the set with name setname
.
Ferrite.getvertexset
— Functiongetvertexset(grid::AbstractGrid, setname::String)
Returns all vertices as VertexIndex
in the set with name setname
.
Ferrite.transform_coordinates!
— Functiontransform_coordinates!(grid::Abstractgrid, f::Function)
Transform the coordinates of all nodes of the grid
based on some transformation function f(x)
.
Ferrite.getcoordinates
— Functiongetcoordinates(grid::AbstractGrid, idx::Union{Int,CellIndex})
getcoordinates(cache::CellCache)
Get a vector with the coordinates of the cell corresponding to idx
or cache
Ferrite.getcoordinates!
— Functiongetcoordinates!(x::Vector{<:Vec}, grid::AbstractGrid, idx::Union{Int,CellIndex})
getcoordinates!(x::Vector{<:Vec}, grid::AbstractGrid, cell::AbstractCell)
Mutate x
to the coordinates of the cell corresponding to idx
or cell
.
Ferrite.geometric_interpolation
— Methodgeometric_interpolation(::AbstractCell)::ScalarInterpolation
geometric_interpolation(::Type{<:AbstractCell})::ScalarInterpolation
Each AbstractCell
type has a unique geometric interpolation describing its geometry. This function returns that interpolation, which is always a scalar interpolation.
Ferrite.get_node_coordinate
— Functionget_node_coordinate(::Node)
Get the value of the node coordinate.
get_node_coordinate(grid::AbstractGrid, n::Int)
Return the coordinate of the n
th node in grid
Ferrite.getspatialdim
— MethodFerrite.getspatialdim(grid::AbstractGrid)
Get the spatial dimension of the grid, corresponding to the vector dimension of the grid's coordinates.
Ferrite.getrefdim
— MethodFerrite.getrefdim(cell::AbstractCell)
Ferrite.getrefdim(::Type{<:AbstractCell})
Get the reference dimension of the cell, i.e. the dimension of the cell's reference shape.
Topology
Ferrite.ExclusiveTopology
— TypeExclusiveTopology(grid::AbstractGrid)
The experimental feature ExclusiveTopology
saves topological (connectivity/neighborhood) data of the grid. Only the highest dimensional neighborhood is saved. I.e., if something is connected by a face and an edge, only the face neighborhood is saved. The lower dimensional neighborhood is recomputed when calling getneighborhood if needed.
Fields
vertex_to_cell::AbstractArray{AbstractVector{Int}, 1}
: global vertex id to all cells containing the vertexcell_neighbor::AbstractArray{AbstractVector{Int}, 1}
: cellid to all connected cellsface_neighbor::AbstractArray{AbstractVector{FaceIndex}, 2}
:face_neighbor[cellid, local_face_id]
-> neighboring facesedge_neighbor::AbstractArray{AbstractVector{EdgeIndex}, 2}
:edge_neighbor[cellid, local_edge_id]
-> neighboring edgesvertex_neighbor::AbstractArray{AbstractVector{VertexIndex}, 2}
:vertex_neighbor[cellid, local_vertex_id]
-> neighboring verticesface_skeleton::Union{Vector{FaceIndex}, Nothing}
: List of unique faces in the grid given asFaceIndex
edge_skeleton::Union{Vector{EdgeIndex}, Nothing}
: List of unique edges in the grid given asEdgeIndex
vertex_skeleton::Union{Vector{VertexIndex}, Nothing}
: List of unique vertices in the grid given asVertexIndex
The implementation only works with conforming grids, i.e. grids without "hanging nodes". Non-conforming grids will give unexpected results. Grids with embedded cells (different reference dimension compared to the spatial dimension) are not supported, and will error on construction.
Ferrite.getneighborhood
— Functiongetneighborhood(topology, grid::AbstractGrid, cellidx::CellIndex, include_self=false)
getneighborhood(topology, grid::AbstractGrid, faceidx::FaceIndex, include_self=false)
getneighborhood(topology, grid::AbstractGrid, vertexidx::VertexIndex, include_self=false)
getneighborhood(topology, grid::AbstractGrid, edgeidx::EdgeIndex, include_self=false)
Returns all connected entities of the same type as defined by the respective topology. If include_self
is true, the given entity is included in the returned list as well.
Ferrite.facetskeleton
— Functionfacetskeleton(top::ExclusiveTopology, grid::AbstractGrid)
Materializes the skeleton from the neighborhood
information by returning an iterable over the unique facets in the grid, described by FacetIndex
.
Ferrite.vertex_star_stencils
— Functionvertex_star_stencils(top::ExclusiveTopology, grid::Grid) -> AbstractVector{AbstractVector{VertexIndex}}
Computes the stencils induced by the edge connectivity of the vertices.
Ferrite.getstencil
— Functiongetstencil(top::ArrayOfVectorViews{VertexIndex, 1}, grid::AbstractGrid, vertex_idx::VertexIndex) -> AbstractVector{VertexIndex}
Get an iterateable over the stencil members for a given local entity.
Grid Sets Utility
Ferrite.addcellset!
— Functionaddcellset!(grid::AbstractGrid, name::String, cellid::AbstractVecOrSet{Int})
addcellset!(grid::AbstractGrid, name::String, f::function; all::Bool=true)
Adds a cellset to the grid with key name
. Cellsets are typically used to define subdomains of the problem, e.g. two materials in the computational domain. The DofHandler
can construct different fields which live not on the whole domain, but rather on a cellset. all=true
implies that f(x)
must return true
for all nodal coordinates x
in the cell if the cell should be added to the set, otherwise it suffices that f(x)
returns true
for one node.
addcellset!(grid, "left", Set((1,3))) #add cells with id 1 and 3 to cellset left
addcellset!(grid, "right", x -> norm(x[1]) < 2.0 ) #add cell to cellset right, if x[1] of each cell's node is smaller than 2.0
Ferrite.addfacetset!
— Functionaddfacetset!(grid::AbstractGrid, name::String, faceid::AbstractVecOrSet{FacetIndex})
addfacetset!(grid::AbstractGrid, name::String, f::Function; all::Bool=true)
Adds a facetset to the grid with key name
. A facetset maps a String
key to a OrderedSet
of tuples corresponding to (global_cell_id, local_facet_id)
. Facetsets can be used to initialize Dirichlet
boundary conditions for the ConstraintHandler
. all=true
implies that f(x)
must return true
for all nodal coordinates x
on the facet if the facet should be added to the set, otherwise it suffices that f(x)
returns true
for one node.
addfacetset!(grid, "right", Set((FacetIndex(2,2), FacetIndex(4,2)))) #see grid manual example for reference
addfacetset!(grid, "clamped", x -> norm(x[1]) ≈ 0.0) #see incompressible elasticity example for reference
Ferrite.addboundaryfacetset!
— Functionaddboundaryfacetset!(grid::AbstractGrid, topology::ExclusiveTopology, name::String, f::Function; all::Bool=true)
Adds a boundary facetset to the grid with key name
. A facetset maps a String
key to a OrderedSet
of tuples corresponding to (global_cell_id, local_facet_id)
. Facetsets are used to initialize Dirichlet
structs, that are needed to specify the boundary for the ConstraintHandler
. all=true
implies that f(x)
must return true
for all nodal coordinates x
on the facet if the facet should be added to the set, otherwise it suffices that f(x)
returns true
for one node.
Ferrite.addvertexset!
— Functionaddvertexset!(grid::AbstractGrid, name::String, faceid::AbstractVecOrSet{FaceIndex})
addvertexset!(grid::AbstractGrid, name::String, f::Function)
Adds a vertexset to the grid with key name
. A vertexset maps a String
key to a OrderedSet
of tuples corresponding to (global_cell_id, local_vertex_id)
. Vertexsets can be used to initialize Dirichlet
boundary conditions for the ConstraintHandler
.
addvertexset!(grid, "right", Set((VertexIndex(2,2), VertexIndex(4,2))))
addvertexset!(grid, "clamped", x -> norm(x[1]) ≈ 0.0)
Ferrite.addboundaryvertexset!
— Functionaddboundaryvertexset!(grid::AbstractGrid, topology::ExclusiveTopology, name::String, f::Function; all::Bool=true)
Adds a boundary vertexset to the grid with key name
. A vertexset maps a String
key to an OrderedSet
of tuples corresponding to (global_cell_id, local_vertex_id)
. all=true
implies that f(x)
must return true
for all nodal coordinates x
on the face if the face should be added to the set, otherwise it suffices that f(x)
returns true
for one node.
Ferrite.addnodeset!
— Functionaddnodeset!(grid::AbstractGrid, name::String, nodeid::AbstractVecOrSet{Int})
addnodeset!(grid::AbstractGrid, name::String, f::Function)
Adds a nodeset::OrderedSet{Int}
to the grid
's nodesets
with key name
. Has the same interface as addcellset
. However, instead of mapping a cell id to the String
key, a set of node ids is returned.
Multithreaded Assembly
Ferrite.create_coloring
— Functioncreate_coloring(g::Grid, cellset=1:getncells(g); alg::ColoringAlgorithm)
Create a coloring of the cells in grid g
such that no neighboring cells have the same color. If only a subset of cells should be colored, the cells to color can be specified by cellset
.
Returns a vector of vectors with cell indexes, e.g.:
ret = [
[1, 3, 5, 10, ...], # cells for color 1
[2, 4, 6, 12, ...], # cells for color 2
]
Two different algorithms are available, specified with the alg
keyword argument:
alg = ColoringAlgorithm.WorkStream
(default): Three step algorithm from Turcksin et al. [11], albeit with a greedy coloring in the second step. Generally results in more colors thanColoringAlgorithm.Greedy
, however the cells are more equally distributed among the colors.alg = ColoringAlgorithm.Greedy
: greedy algorithm that works well for structured quadrilateral grids such as e.g. quadrilateral grids fromgenerate_grid
.
The resulting colors can be visualized using Ferrite.write_cell_colors
.
In a previous version of Ferrite this function returned a dictionary mapping cell ID to color numbers as the first argument. If you need this mapping you can create it using the following construct:
colors = create_coloring(...)
cell_colormap = Dict{Int,Int}(
cellid => color for (color, cellids) in enumerate(final_colors) for cellid in cellids
)
References
- [11] Turcksin et al. ACM Trans. Math. Softw. 43 (2016).