Grid & AbstractGrid

Grid

Ferrite.generate_gridFunction
generate_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.

source
Ferrite.NodeType
Node{dim, T}

A Node is a point in space.

Fields

  • x::Vec{dim,T}: stores the coordinates
source
Ferrite.VertexIndexType

A VertexIndex wraps an (Int, Int) and defines a local vertex by pointing to a (cell, vert).

source
Ferrite.EdgeIndexType

A EdgeIndex wraps an (Int, Int) and defines a local edge by pointing to a (cell, edge).

source
Ferrite.FaceIndexType

A FaceIndex wraps an (Int, Int) and defines a local face by pointing to a (cell, face).

source
Ferrite.FacetIndexType

A FacetIndex wraps an (Int, Int) and defines a local facet by pointing to a (cell, facet).

source
Ferrite.GridType
Grid{dim, C<:AbstractCell, T<:Real} <: AbstractGrid}

A Grid is a collection of Ferrite.AbstractCells and Ferrite.Nodes 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 grid
  • nodes::Vector{Node{dim,T}}: stores the dim dimensional nodes of the grid
  • cellsets::Dict{String, OrderedSet{Int}}: maps a String key to an OrderedSet of cell ids
  • nodesets::Dict{String, OrderedSet{Int}}: maps a String key to an OrderedSet of global node ids
  • facetsets::Dict{String, OrderedSet{FacetIndex}}: maps a String to an OrderedSet of FacetIndex
  • vertexsets::Dict{String, OrderedSet{VertexIndex}}: maps a String key to an OrderedSet of VertexIndex
source

Utility Functions

Ferrite.getcellsFunction
getcells(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}.

source
Ferrite.getnodesFunction
getnodes(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.

source
Ferrite.getcellsetFunction
getcellset(grid::AbstractGrid, setname::String)

Returns all cells as cellid in the set with name setname.

source
Ferrite.getnodesetFunction
getnodeset(grid::AbstractGrid, setname::String)

Returns all nodes as nodeid in the set with name setname.

source
Ferrite.getfacetsetFunction
getfacetset(grid::AbstractGrid, setname::String)

Returns all faces as FacetIndex in the set with name setname.

source
Ferrite.getvertexsetFunction
getvertexset(grid::AbstractGrid, setname::String)

Returns all vertices as VertexIndex in the set with name setname.

source
Ferrite.transform_coordinates!Function
transform_coordinates!(grid::Abstractgrid, f::Function)

Transform the coordinates of all nodes of the grid based on some transformation function f(x).

source
Ferrite.getcoordinatesFunction
getcoordinates(grid::AbstractGrid, idx::Union{Int,CellIndex})
getcoordinates(cache::CellCache)

Get a vector with the coordinates of the cell corresponding to idx or cache

source
Ferrite.getcoordinates!Function
getcoordinates!(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.

source
Ferrite.geometric_interpolationMethod
geometric_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.

source
Ferrite.get_node_coordinateFunction
get_node_coordinate(::Node)

Get the value of the node coordinate.

source
get_node_coordinate(grid::AbstractGrid, n::Int)

Return the coordinate of the nth node in grid

source
Ferrite.getspatialdimMethod
Ferrite.getspatialdim(grid::AbstractGrid)

Get the spatial dimension of the grid, corresponding to the vector dimension of the grid's coordinates.

source
Ferrite.getrefdimMethod
Ferrite.getrefdim(cell::AbstractCell)
Ferrite.getrefdim(::Type{<:AbstractCell})

Get the reference dimension of the cell, i.e. the dimension of the cell's reference shape.

source

Topology

Ferrite.ExclusiveTopologyType
ExclusiveTopology(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 vertex
  • cell_neighbor::AbstractArray{AbstractVector{Int}, 1}: cellid to all connected cells
  • face_neighbor::AbstractArray{AbstractVector{FaceIndex}, 2}: face_neighbor[cellid, local_face_id] -> neighboring faces
  • edge_neighbor::AbstractArray{AbstractVector{EdgeIndex}, 2}: edge_neighbor[cellid, local_edge_id] -> neighboring edges
  • vertex_neighbor::AbstractArray{AbstractVector{VertexIndex}, 2}: vertex_neighbor[cellid, local_vertex_id] -> neighboring vertices
  • face_skeleton::Union{Vector{FaceIndex}, Nothing}: List of unique faces in the grid given as FaceIndex
  • edge_skeleton::Union{Vector{EdgeIndex}, Nothing}: List of unique edges in the grid given as EdgeIndex
  • vertex_skeleton::Union{Vector{VertexIndex}, Nothing}: List of unique vertices in the grid given as VertexIndex
Limitations

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.

source
Ferrite.getneighborhoodFunction
getneighborhood(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.

source
Ferrite.facetskeletonFunction
facetskeleton(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.

source
Ferrite.vertex_star_stencilsFunction
vertex_star_stencils(top::ExclusiveTopology, grid::Grid) -> AbstractVector{AbstractVector{VertexIndex}}

Computes the stencils induced by the edge connectivity of the vertices.

source
Ferrite.getstencilFunction
getstencil(top::ArrayOfVectorViews{VertexIndex, 1}, grid::AbstractGrid, vertex_idx::VertexIndex) -> AbstractVector{VertexIndex}

Get an iterateable over the stencil members for a given local entity.

source

Grid Sets Utility

Ferrite.addcellset!Function
addcellset!(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
source
Ferrite.addfacetset!Function
addfacetset!(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
source
Ferrite.addboundaryfacetset!Function

addboundaryfacetset!(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.

source
Ferrite.addvertexset!Function
addvertexset!(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)
source
Ferrite.addboundaryvertexset!Function

addboundaryvertexset!(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.

source
Ferrite.addnodeset!Function
addnodeset!(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.

source

Multithreaded Assembly

Ferrite.create_coloringFunction
create_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 than ColoringAlgorithm.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 from generate_grid.

The resulting colors can be visualized using Ferrite.write_cell_colors.

Cell to color mapping

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).
source