Index
StartUpDG.BoundaryTagPlotter
StartUpDG.CurvedMesh
StartUpDG.CutCellMesh
StartUpDG.MeshData
StartUpDG.MeshData
StartUpDG.MeshData
StartUpDG.MeshImportOptions
StartUpDG.MeshPlotter
StartUpDG.MultidimensionalQuadrature
StartUpDG.MultipleRefElemData
StartUpDG.NamedArrayPartition
StartUpDG.NonConformingMesh
StartUpDG.PhysicalFrame
StartUpDG.Polynomial
StartUpDG.RefElemData
StartUpDG.RefElemData
StartUpDG.RefElemData
StartUpDG.RefElemData
StartUpDG.RefElemData
StartUpDG.RefElemData
StartUpDG.RefElemData
StartUpDG.RefElemData
StartUpDG.TensorProductGaussCollocation
StartUpDG.TensorProductQuadrature
StartUpDG.VertexMappedMesh
StartUpDG.VertexMeshPlotter
NodesAndModes.equi_nodes
StartUpDG.MeshData_to_vtk
StartUpDG.MeshData_to_vtk
StartUpDG.SUD_to_vtk_order
StartUpDG.boundary_face_centroids
StartUpDG.build_node_maps
StartUpDG.caratheodory_pruning_qr
StartUpDG.ck45
StartUpDG.connect_mesh
StartUpDG.connect_mesh
StartUpDG.construct_cut_surface_quadrature
StartUpDG.construct_cut_volume_quadrature
StartUpDG.construct_physical_frame_elements
StartUpDG.estimate_h
StartUpDG.findline
StartUpDG.get_boundary_face_labels
StartUpDG.get_node_boundary_tags
StartUpDG.get_num_elements
StartUpDG.hex_vtk_order
StartUpDG.hybridized_SBP_operators
StartUpDG.hybridized_SBP_operators
StartUpDG.inverse_trace_constant
StartUpDG.make_periodic
StartUpDG.n_verts_between
StartUpDG.quad_vtk_order
StartUpDG.read_Gmsh_2D
StartUpDG.read_Gmsh_2D_v2
StartUpDG.read_Gmsh_2D_v4
StartUpDG.read_Gmsh_2D_v4
StartUpDG.refine
StartUpDG.remap_element_grouping
StartUpDG.sort_by_axis
StartUpDG.sparse_low_order_SBP_operators
StartUpDG.subcell_limiting_operators
StartUpDG.subtriangulated_cutcell_quadrature
StartUpDG.tag_boundary_faces
StartUpDG.tag_boundary_faces
StartUpDG.triangle_vtk_order
StartUpDG.triangulateIO_to_VXYEToV
StartUpDG.type_to_vtk
StartUpDG.type_to_vtk
StartUpDG.type_to_vtk
StartUpDG.type_to_vtk
StartUpDG.uniform_mesh
StartUpDG.vtk_order
StartUpDG.vtk_order
StartUpDG.vtk_order
StartUpDG.vtk_order
StartUpDG.wedge_vtk_order
Triangulate.triangulate
Functions
StartUpDG.BoundaryTagPlotter
— TypeBoundaryTagPlotter(triout::TriangulateIO)
Plot recipe to visualize boundary tags by color. Usage: plot(BoundaryTagPlotter(triout))
StartUpDG.CurvedMesh
— Typestruct CurvedMesh{T}
Mesh type indicating that the mesh has been curved. Stores the original mesh type as a field.
Fields
originalmeshtype :: T
StartUpDG.CutCellMesh
— TypeCutCellMesh
is used in the MeshData
field mesh_type
for cut cell meshes.
The field physical_frame_elements
is a container with shifting/scaling information for each element. We evaluate the physical basis over each element by applying a shifting and scaling of the physical coordinates. The resulting shifted/scaled coordinates then fall into the reference element and can be used to evaluate a reference element basis.
The field cut_face_nodes
is a container whose elements are indices of face nodes for a cut element. In other words, md.xf.cut[cut_face_nodes[1]]
returns the face nodes of the first element.
We assume all cut elements have the same number of volume quadrature points (which is at least the dimension of a degree 2N polynomial space).
The field objects
contains a tuple of the objects used to define the cut region.
The field cut_cell_operators
contains optionally precomputed operators (mass, differntiation, face interpolation, and lifting operators).
The field cut_cell_data
contains additional data from PathIntersections.
StartUpDG.MeshData
— TypeMeshData(rd::RefElemData, objects,
vx::AbstractVector, vy::AbstractVector,
quadrature_type=Subtriangulation();
quad_rule_face=get_1d_quadrature(rd),
precompute_operators=false)
Constructor for MeshData utilizing moment fitting. Does not guarantee positive quadrature weights, and is slower due to the use of adaptive sampling to construct
!!! Warning: this may be deprecated or removed in future versions.
StartUpDG.MeshData
— Typestruct MeshData{Dim, Tv, Ti}
MeshData: contains info for a high order piecewise polynomial discretization on an unstructured mesh.
Example:
N, K1D = 3, 2
rd = RefElemData(Tri(), N)
VXY, EToV = uniform_mesh(Tri(), K1D)
md = MeshData(VXY, EToV, rd)
(; x, y ) = md
StartUpDG.MeshData
— MethodMeshData(VXYZ, EToV, rd::RefElemData)
MeshData((VXYZ, EToV), rd::RefElemData)
Returns a MeshData struct with high order DG mesh information from the unstructured mesh information (VXYZ..., EToV).
MeshData(rd::RefElemData, md::MeshData, xyz...)
Given new nodal positions xyz...
(e.g., from mesh curving), recomputes geometric terms and outputs a new MeshData struct. Only fields modified are the coordinate-dependent terms xyz
, xyzf
, xyzq
, rstxyzJ
, J
, nxyzJ
, Jf
.
StartUpDG.MeshImportOptions
— TypeMeshImportOptions
This struct allows the user to opt for supported features when importing a Gmsh 4.1 .msh file.
Support
- grouping::Bool | On import would you like to include physical group assignements of 2D elements?
- remap_group_name::Bool | On import would you like to maintain or remap physical group ID? Remap results in groupIds in the range 1:number_group_ids.
StartUpDG.MeshPlotter
— TypeMeshPlotter(rd::RefElemData, md::RefElemData)
Plot recipe to plot a (possibly curved) quadrilateral or triangular mesh. Usage: plot(MeshPlotter(...))
StartUpDG.MultidimensionalQuadrature
— TypeMultidimensionalQuadrature
A type parameter for Polynomial
indicating that the quadrature has no specific structure. Example usage:
# these are both equivalent
approximation_type = Polynomial{MultidimensionalQuadrature}()
approximation_type = Polynomial(MultidimensionalQuadrature())
StartUpDG.MultipleRefElemData
— Typestruct MultipleRefElemData{T <: NamedTuple}
data::T
end
Holds multiple RefElemData
objects in data
where typeof(data) <: NamedTuple
.
Individual RefElemData
can be accessed via getproperty
, for example rds.Tri
.
StartUpDG.NamedArrayPartition
— TypeNamedArrayPartition(; kwargs...)
NamedArrayPartition(x::NamedTuple)
Similar to an ArrayPartition
but the individual arrays can be accessed via the constructor-specified names. However, unlike ArrayPartition
, each individual array must have the same element type.
StartUpDG.NonConformingMesh
— TypeThis is an experimental feature and may change without warning in future releases.
This is a proof of concept implementation of a non-conforming mesh in StartUpDG.jl. The intended usage is as follows:
rd = RefElemData(Quad(), N=7)
md = MeshData(NonConformingQuadMeshExample(), rd)
(; x, y) = md
u = @. sin(pi * x) * sin(pi * y)
# interpolate to faces
num_total_faces = num_faces(rd.element_type) * md.num_elements
uf = reshape(rd.Vf * u, :, num_total_faces)
# interpolate faces to mortars (`uf` denotes mortar faces for `NonConformingMesh` types)
(; nonconforming_faces, mortar_interpolation_matrix) = md.mesh_type
u_mortar = reshape(mortar_interpolation_matrix * uf[:, nonconforming_faces], :, 2 * length(nonconforming_faces))
# construct interior (uM = u⁻ "minus") values and exterior (uP = u⁺ "plus") values
uM = hcat(uf, u_mortar) # uM = both element faces and mortar faces
uP = uM[md.mapP]
The mortar_projection_matrix
similarly maps values from 2 mortar faces back to values on the original non-conforming face. These can be used to create DG solvers on non-conforming meshes.
StartUpDG.PhysicalFrame
— Type`PhysicalFrame{NDIMS} <: AbstractElemShape{NDIMS}`
PhysicalFrame
element type. Uses a total degree N approximation space, but is computed with a tensor product Legendre basis as opposed to a triangular PKDO basis. Stores fields shifting
and scaling
to shift/scale physical coordinates so that they are on the reference element.
PhysicalFrame()
PhysicalFrame(x, y)
PhysicalFrame(x, y, vx, vy): stores coordinates `vx, vy` of background Cartesian cell
Constructors for a PhysicalFrame object (optionally uses arrays of points x
, y
on a cut element).
StartUpDG.Polynomial
— TypePolynomial{T}
Represents polynomial approximation types (as opposed to finite differences). By default, Polynomial()
constructs a Polynomial{StartUpDG.DefaultPolynomialType}
. Specifying a type parameters allows for dispatch on additional structure within a polynomial approximation (e.g., collocation, tensor product quadrature, etc).
StartUpDG.RefElemData
— Typestruct RefElemData
RefElemData: contains info (interpolation points, volume/face quadrature, operators) for a high order nodal basis on a given reference element.
Example:
N = 3
rd = RefElemData(Tri(), N)
(; r, s ) = rd
StartUpDG.RefElemData
— MethodRefElemData(elem; N, kwargs...)
RefElemData(elem, approx_type; N, kwargs...)
Keyword argument constructor for RefElemData (to "label" N
via rd = RefElemData(Line(), N=3)
)
StartUpDG.RefElemData
— MethodRefElemData(elem::Line, approximation_type, N;
quad_rule_vol = quad_nodes(elem, N+1))
RefElemData(elem, approximation_type, N;
quad_rule_vol = quad_nodes(elem, N),
quad_rule_face = quad_nodes(face_type(elem), N))
Constructor for RefElemData
for different element types.
StartUpDG.RefElemData
— MethodRefElemData(elementType::Line, approxType::SBP, N)
RefElemData(elementType::Quad, approxType::SBP, N)
RefElemData(elementType::Hex, approxType::SBP, N)
RefElemData(elementType::Tri, approxType::SBP, N)
SBP reference element data for Quad()
, Hex()
, and Tri()
elements.
For Line()
, Quad()
, and Hex()
, approxType
is SBP{TensorProductLobatto}
.
For Tri()
, approxType
can be SBP{Kubatko{LobattoFaceNodes}}
, SBP{Kubatko{LegendreFaceNodes}}
, or SBP{Hicken}
.
StartUpDG.RefElemData
— MethodRefElemData(elem::Pyr,
approximation_type::Polynomial, N;
quad_rule_vol=quad_nodes(elem, N),
quad_rule_face_quad=quad_nodes(Quad(), N),
quad_rule_face_tri=quad_nodes(Tri(), N),
quad_rule_face=(quad_rule_face_quad, quad_rule_face_tri),
Nplot=10)
Builds operators for pyramids.
StartUpDG.RefElemData
— MethodRefElemData(elem::Union{Line, Quad, Hex}, approximation_type::Polynomial{Gauss}, N)
Builds a rd::RefElemData
with (N+1)-point Gauss quadrature in each dimension.
StartUpDG.RefElemData
— MethodRefElemData(elem::Union{Tri, Tet, Pyr}, approx_type::Polynomial{<:TensorProductQuadrature}, N; kwargs...)
RefElemData(elem::Union{Wedge},
approx_type::Polynomial{<:TensorProductQuadrature}, N;
quad_rule_tri = stroud_quad_nodes(Tri(), 2 * N),
quad_rule_line = gauss_quad(0, 0, N),
kwargs...)
Uses collapsed coordinate volume quadrature. Should be called via
RefElemData(Tri(), Polynomial(TensorProductQuadrature()), N)
StartUpDG.RefElemData
— MethodRefElemData(elem::Wedge, approximation_type::Polynomial, N;
quad_rule_vol=quad_nodes(elem, N),
quad_rule_face_quad=quad_nodes(Quad(), N),
quad_rule_face_tri=quad_nodes(Tri(), N),
quad_rule_face=(quad_rule_face_quad, quad_rule_face_tri),
Nplot=10)
Builds operators for prisms/wedges
StartUpDG.TensorProductGaussCollocation
— TypeTensorProductGaussCollocation
Polynomial{TensorProductGaussCollocation} type indicates a tensor product
(N+1)-point Gauss quadrature on tensor product elements.
StartUpDG.TensorProductQuadrature
— TypeTensorProductQuadrature{T}
A type parameter to Polynomial
indicating that the quadrature has a tensor product structure. Example usage:
# these are both equivalent
approximation_type = Polynomial{TensorProductQuadrature}(gauss_quad(0, 0, 1))
approximation_type = Polynomial(TensorProductQuadrature(gauss_quad(0, 0, 1)))
StartUpDG.VertexMappedMesh
— Typestruct VertexMappedMesh
The default MeshData
mesh type, represents a mesh which is defined purely by vertex locations and element-to-vertex connectivities. For example, these include affine triangular meshes or bilinear quadrilateral or trilinear hexahedral meshes.
Fields
element_type :: TE <: AbstractElemShape
VXYZ :: TV
EToV :: TEV
StartUpDG.VertexMeshPlotter
— TypeVertexMeshPlotter((VX, VY), EToV, fv)
VertexMeshPlotter(triout::TriangulateIO)
Plot recipe to plot a quadrilateral or triangular mesh. Usage: plot(VertexMeshPlotter(...))
NodesAndModes.equi_nodes
— Methodfunction NodesAndModes.equi_nodes(elem::PhysicalFrame, curve, N)
Returns back Np(N)
equally spaced nodes on the background quadrilateral corresponding to elem
, with points inside of curve
removed.
StartUpDG.MeshData_to_vtk
— FunctionMeshDatatovtk(md, rd, dim, data, dataname, datatype, filename, writedata = false, equidist_nodes = true)
Translate the given mesh into a vtk-file. md
holds a MeshData
object rd
holds a reference element data/RefElemData
of a TensorProductWedge data
holds an array of arrays (of size num_nodes
by num_elements
) with plotting data dataname
is an array of strings with name of the associated data write_data
, flag if data should be written or not (e.g., if data is not written, only the mesh will be saved as output) equi_dist_nodes
flag if points should be interpolated to equidstant nodes
StartUpDG.MeshData_to_vtk
— MethodMeshData_to_vtk(md, rd, dim, data, dataname, datatype, filename, write_data = false, equi_dist_nodes = true)
Translate the given mesh into a vtk-file. md
holds a MeshData
object rd
holds a reference element data/RefElemData
object. data
holds an array of arrays (of size num_nodes
by num_elements
) with plotting data dataname
is an array of strings with name of the associated data write_data
, flag if data should be written or not (e.g., if data is not written, only the mesh will be saved as output) equi_dist_nodes
flag if points should be interpolated to equidstant nodes
StartUpDG.SUD_to_vtk_order
— MethodSUD_to_vtk_order(rd::RefElemData, dim)
Compute the permutation of the nodes between StartUpDG and VTK
StartUpDG.boundary_face_centroids
— Methodfunction boundary_face_centroids(md)
Returns face centroids and boundary_face_ids
on the boundaries of the domain given by md::MeshData.
StartUpDG.build_node_maps
— Methodbuild_node_maps(FToF, Xf)
Intialize the connectivity table along all edges and boundary node tables of all elements. mapM
- map minus (interior). mapP
- map plus (exterior).
Xf = (xf, yf, zf)
and FToF
is size (Nfaces * K)
and FToF[face]
= face neighbor
mapM
, mapP
are size Nfp
x (Nfaces*K)
Examples
julia> mapM, mapP, mapB = build_node_maps(FToF, (xf, yf))
StartUpDG.caratheodory_pruning_qr
— Methodcaratheodory_pruning_qr(V, w_in)
This performs Caratheodory pruning using a naive QR-based algorithm. Returns (w, inds), where inds
denotes sub-selected indices for a reduced quadrature rule, and w
is a vector of reduced positive weights.
The original Matlab code this was based on was authored by Akil Narayan.
StartUpDG.ck45
— Methodck45()
Returns coefficients rka,rkb,rkc for the 4th order 5-stage low storage Carpenter/Kennedy Runge Kutta method. Coefficients evolve the residual, solution, and local time, e.g.,
Example
for i in eachindex(rk4a, rk4b)
@. res = rk4a[i] * res + dt * rhs # i = RK stage
@. u += rk4b[i] * res
end
StartUpDG.connect_mesh
— Methodconnect_mesh(EToV,fv)
Inputs:
EToV
is anum_elements
byNv
matrix whose rows identify theNv
vertices
which make up one of the num_elements
elements.
fv
(an array of arrays containing unordered indices of face vertices).
Output: FToF
, an length(fv)
by num_elements
index array containing face-to-face connectivity.
StartUpDG.connect_mesh
— Methodconnect_mesh(rd, face_centroids, region_flags, cutcells; tol = 1e2 * eps())
Connects faces of a cut mesh to each other, returns FToF
such that face f
is connected to FToF[f]
.
Inputs:
- rd::RefElemData
- facecentroids = (facecentroidsx, facecentroidsy), where `facecentroids_x/y` are vectors of coordinates of face centroids
region_flags
,cutcells
are return arguments fromPathIntersections.define_regions
The keyword argument tol
is the tolerance for matches between face centroids.
StartUpDG.construct_cut_surface_quadrature
— Functionconstruct_cut_surface_quadrature(N, cutcells, quad_rule_1D = gauss_quad(0, 0, N))
Constructs cut surface quadrature using a degree N
geometric mapping and a reference quadrature rule quad_rule_1D
. Returns xf, yf, nxJ, nyJ, wf
which are vectors, and face_node_indices
, which is a Vector{Vector{Int}}
of global face node indices (which index into xf.cut
, yf.cut
, etc) for each face of each cut element.
On boundaries of cut cells, the surface quadrature is taken to be exact for degree N^2 + (N-1)
polynomials. This ensures satisfaction of a weak GSBP property.
StartUpDG.construct_cut_volume_quadrature
— Methodconstruct_cut_volume_quadrature(N, cutcells; target_degree = 2 * N - 1)
Constructs volume quadrature using subtriangulations of cut cells and Caratheodory pruning. The resulting quadrature is exact for polynomials of degree target_degree
.
We set target_degree
to 2 * N - 1
by default, which is sufficient to ensure that ∫du/dx * v is integrated exactly so that integration by parts holds under the generated cut cell quadrature.
StartUpDG.construct_physical_frame_elements
— Methodconstruct_physical_frame_elements(region_flags, cutcells)
Computes physical frame shifting and scaling parameters from the vertices of cut cells and the background cell location.
StartUpDG.estimate_h
— Methodestimate_h(rd::RefElemData, md::MeshData)
estimate_h(e, rd::RefElemData, md::MeshData) # e = element index
Estimates the mesh size via min sizeofdomain * |J|/|Jf|, since |J| = O(hᵈ) and |Jf| = O(hᵈ⁻¹).
StartUpDG.findline
— Methodfindline(word::String, lines)
Outputs the line number of word
in lines
.
It is assumed that the word exists at least once in the file.
StartUpDG.get_boundary_face_labels
— Methodfunction get_boundary_face_labels(triout::TriangulateIO, md::MeshData{2})
Find Triangle segment labels of boundary faces. Returns two arguments:
boundary_face_tags
: tags of faces on the boundaryboundary_faces
: list of faces on the boundary of the domain
StartUpDG.get_node_boundary_tags
— Methodfunction get_node_boundary_tags(triout::TriangulateIO,md::MeshData{2},rd::RefElemData{2,Tri})
Computes node_tags
= Nfp
x Nfaces * num_elements
array where each entry is a Triangulate.jl tag number.
StartUpDG.get_num_elements
— Functionreturns the number of elements in a .msh file of a specified dimension
Notes: Gmsh includes elements in a .msh file of multiple dimensions. We want a count of how many
2D elements are in our file. This corisponds to the number of elements in our tri mesh.
StartUpDG.hex_vtk_order
— Functionhex_vtk_order(corner_verts, order, dim, skip = false)
Compute the coordinates of a VTKLAGRANGEHEXAHEDRON of a hex of order order
defined by the coordinates of the vertices given in corner_verts
. dim
is the dimension of the coordinates given. If skip
is set to true, the coordinates of the vertex- and edge-points aren't computed, which can be used to compute points of a VTK_LAGRANGE_HEXHEDRON
Inspired by: https://github.com/ju-kreber/paraview-scripts/blob/master/node_ordering.py
VTK node numbering of a hexagon:
8+------+7
/| /|
/ | / |
5+------+6 |
z | 4+–-|–+3 | y | / | / |/ |/ |/ 0–> x 1+–––+2
StartUpDG.hybridized_SBP_operators
— Methodfunction hybridized_SBP_operators(rd::RefElemData{DIMS})
Constructs hybridized SBP operators given a RefElemData
. Returns operators Qrsth..., VhP, Ph
.
StartUpDG.hybridized_SBP_operators
— Methodhybridized_SBP_operators(md::MeshData{2, <:CutCellMesh})
This constructs hybridized SBP operators using the approach taken in Chan (2019), "Skew-Symmetric Entropy Stable Modal Discontinuous Galerkin Formulations". https://doi.org/10.1007/s10915-019-01026-w
This function returns hybridized_operators::Vector{Tuple{<:Matrix, <:Matrix}}
and project_and_interp_operators, projection_operators, interpolation_operators
, which are all Vector{<:Matrix}
, where each entry corresponds to a cut element.
StartUpDG.inverse_trace_constant
— Methodfunction inverse_trace_constant(rd::RefElemData)
Returns the degree-dependent constant in the inverse trace equality over the reference element (as reported in "GPU-accelerated dG methods on hybrid meshes" by Chan, Wang, Modave, Remacle, Warburton 2016).
Can be used to estimate dependence of maximum stable timestep on degree of approximation.
StartUpDG.make_periodic
— Methodmake_periodic(md::MeshData{Dim}, is_periodic...) where {Dim}
make_periodic(md::MeshData{Dim}, is_periodic = ntuple(x -> true, Dim)) where {Dim}
make_periodic(md::MeshData, is_periodic = true)
Returns new MeshData such that the node maps mapP
and face maps FToF
are now periodic. Here, is_periodic
is a tuple of Bool
indicating whether or not to impose periodic BCs in the x
,y
, or z
coordinate.
StartUpDG.n_verts_between
— Methodn_verts_between(n, from, to, dim)
Compute the coordinates of n equally distributed points between the points given by from
and to
. dim
is the dimension of from
and to
. Inspired by: https://github.com/ju-kreber/paraview-scripts/blob/master/node_ordering.py
StartUpDG.quad_vtk_order
— Functionquad_vtk_order(corner_verts, order, dim, skip = false)
Compute the coordinates of a VTKLAGRANGEQUAD of a quad of order order
defined by the coordinates of the vertices given in corner_verts
. dim
is the dimension of the coordinates given. If skip
is set to true, the coordinates of the vertex- and edge-points aren't computed, which can be used to compute points of a VTK_LAGRANGE_QUAD
Inspired by: https://github.com/ju-kreber/paraview-scripts/blob/master/node_ordering.py
StartUpDG.read_Gmsh_2D
— Methodread_Gmsh_2D(filename, args...)
Reads a 2D triangular Gmsh file. Mesh formats 2.2 and 4.1 supported. Returns (VX, VY), EToV.
Examples
VXY, EToV = read_Gmsh_2D("eulerSquareCylinder2D.msh") # v2.2 file format
VXY, EToV = read_Gmsh_2D("test/testset_Gmsh_meshes/periodicity_mesh_v4.msh") # v4.1 file format
# if MeshImportOptions.grouping=true, then a third variable `grouping` is returned
VXY, EToV, grouping = read_Gmsh_2D("test/testset_Gmsh_meshes/periodicity_mesh_v4.msh", MeshImportOptions(true, false))
VXY, EToV, grouping = read_Gmsh_2D("test/testset_Gmsh_meshes/periodicity_mesh_v4.msh", true) # same as above
See also
https://gmsh.info/doc/texinfo/gmsh.html#MSH-file-format
https://gmsh.info/doc/texinfo/gmsh.html#MSH-file-format-version-2-0028Legacy0029
StartUpDG.read_Gmsh_2D_v2
— Methodread_Gmsh_2D_v2(filename)
Reads triangular GMSH 2D file format 2.2 0 8. Returns (VX, VY), EToV.
Examples
VXY, EToV = read_Gmsh_2D_v2("eulerSquareCylinder2D.msh")
https://gmsh.info/doc/texinfo/gmsh.html#MSH-file-format-version-2-0028Legacy0029
StartUpDG.read_Gmsh_2D_v4
— FunctionFor brevity when grouping is the only supported feature.
example: VXY, EToV, grouping = read_Gmsh_2D_v4("file.msh",true)
example: VXY, EToV = read_Gmsh_2D_v4("file.msh",false)
StartUpDG.read_Gmsh_2D_v4
— Methodfunction read_Gmsh_2D_v4(filename, options)
reads triangular GMSH 2D .msh files.
Output
This depends on if grouping is opted for or not
- returns: (VX, VY), EToV
- returns: (VX, VY), EToV, grouping
Supported formats and features:
- version 4.1 'physical group support 'remap group ids
grouping application
When modeling the wave equation you might want wave speeds to vary across your domain. By assigning physical groups in Gmsh we can maintain such groupings upon importing the .msh file. Each imported element will be a member of a phyical group.
VXY, EToV = read_Gmsh_2D_v4("eulerSquareCylinder2D.msh")
VXY, EToV = read_Gmsh_2D_v4("eulerSquareCylinder2D.msh",false)
VXY, EToV, grouping = read_Gmsh_2D_v4("eulerSquareCylinder2D.msh", true)
option = MeshImportOption(true)
VXY, EToV, grouping = read_Gmsh_2D_v4("eulerSquareCylinder2D.msh", option)
https://gmsh.info/doc/texinfo/gmsh.html#MSH-file-format
Notes: the version 4 format has a more detailed block data format this leads to more complicated parser.
StartUpDG.refine
— Functionfunction refine(triout, h, href = h/2)
Refinement of a previous mesh given the current mesh size h
. Preserves boundary/volume tags.
StartUpDG.remap_element_grouping
— Methodremapelementgrouping!(eg::Vector{Int}) GMSH uses integers for naming conventions. This function remaps the Gmsh ids to a list of ids 1:numGroups. This just cleans up a little after Gmsh
Example output
remapelementgrouping([16,16,17,17]) -> [1,1,2,2]
StartUpDG.sort_by_axis
— Methodsort_by_axis(corner_verts)
Given the points 'corner_verts' sort them in a lexicographical order and return the permutated points.
StartUpDG.sparse_low_order_SBP_operators
— Methodfunction sparse_low_order_SBP_operators(rd; factor=1.01)
Constructs sparse low order SBP operators given a RefElemData
. Returns operators Qrst..., E ≈ Vf * Pq
that satisfy a generalized summation-by-parts (GSBP) property:
`Q_i + Q_i^T = E' * B_i * E`
factor
is a scaling which determines how close a node must be to another node to be considered a neighbor.
StartUpDG.subcell_limiting_operators
— MethodΔrst, Rrst = subcell_limiting_operators(rd::RefElemData)
Returns tuples of subcell limiting operators Drst = (Δr, Δs, ...) and R = (Rr, Rs, ...) such that for r where sum(r) = 0, sum(D * Diagonal(θ) * R * r) = 0 for any choice of θ. These operators are useful for conservative subcell limiting techniques (see https://doi.org/10.1016/j.compfluid.2022.105627 for an example of such an approach on tensor product elements).
Sparse SBP operators used in an intermediate step when buidling these subcell limiting operators; by default, these operators are constructed using sparse_low_order_SBP_operators
. To construct subcell limiting operators for a general SBP operator, one can use the following:
Δ, R = subcell_limiting_operators(Q::AbstractMatrix; tol = 100 * eps())
StartUpDG.subtriangulated_cutcell_quadrature
— Functionsubtriangulated_cutcell_quadrature(cutcell, rd_tri::RefElemData,
r1D = gauss_lobatto_quad(0,0,rd_tri.N))
Constructs a quadrature from subtriangulations. The degree of both the quadrature rule and isoparametric mapping are specified by rd_tri
.
The optional argument r1D
specifies the nodal points used for constructing a curved mapping via interpolatory warp and blend.
StartUpDG.tag_boundary_faces
— Methodfunction tag_boundary_faces(md, boundary_name::Symbol = :entire_boundary)
function tag_boundary_faces(md, boundary_list::Dict{Symbol, <:Function})
When called without arguments, just returns Dict(:entire_boundary => boundary_faces)
`.
Example usage:
julia> rd = RefElemData(Tri(), N=1)
julia> md = MeshData(uniform_mesh(Tri(), 2)..., rd)
julia> on_bottom_boundary(x, y, tol = 1e-13) = abs(y+1) < tol
julia> on_top_boundary(x, y, tol = 1e-13) = abs(y-1) < tol
julia> tag_boundary_faces(Dict(:bottom => on_bottom_boundary,
:top => on_top_boundary), md)
StartUpDG.tag_boundary_faces
— Methodfunction tag_boundary_faces(triout::TriangulateIO,
rd::RefElemData{2,Tri}, md::MeshData{2},
boundary_list::Union{NamedTuple,Dict{Symbol,Int}})
Here, boundary_list
is a Dict
(or NamedTuple
) whose values are the boundary tags for a TriangulateIO
mesh format. The output is a Dict
or NamedTuple
with keys given by boundary_list
and values
equal to vectors of faces on that given boundary.
Example usage:
julia> using Triangulate, StartUpDG
julia> triout = scramjet()
julia> rd = RefElemData(Tri(),N=1)
julia> md = MeshData(triangulateIO_to_VXYEToV(triout)...,rd)
julia> tag_boundary_faces(triout,rd,md, Dict(:wall=>1, :inflow=>2, :outflow=>3))
StartUpDG.triangle_vtk_order
— Functiontriangle_vtk_order(corner_verts, order, dim, skip = false)
Compute the coordinates of a VTK_LAGRANGE_TRIANGLE
of a triangle or order order
defined by the coordinates of the vertices given in corner_verts
. dim
is the dimension of the coordinates given. If skip
is set to true, the coordinates of the vertex- and edge-points aren't computed, which can be used to compute points of a VTK_LAGRANGE_TRIANGLE
Inspired by: https://github.com/ju-kreber/paraview-scripts/blob/master/node_ordering.py
StartUpDG.triangulateIO_to_VXYEToV
— Methodfunction triangulateIO_to_VXYEToV(triout::TriangulateIO)
Computes VX
,VY
,EToV
from a TriangulateIO
object.
StartUpDG.type_to_vtk
— Methodtype_to_vtk(elem::Hex)
return the VTK-type
StartUpDG.type_to_vtk
— Methodtype_to_vtk(elem::Quad)
return the VTK-type
StartUpDG.type_to_vtk
— Methodtype_to_vtk(elem::Tri)
return the VTK-type
StartUpDG.type_to_vtk
— Methodtype_to_vtk(elem::Wedge)
return the VTK-type
StartUpDG.uniform_mesh
— Methoduniform_mesh(elem::Line,Kx)
uniform_mesh(elem::Tri,Kx,Ky)
uniform_mesh(elem::Quad,Kx,Ky)
uniform_mesh(elem::Hex,Kx,Ky,Kz)
uniform_mesh(elem, K)
Uniform Kx
(by Ky
by Kz
) mesh on $[-1,1]^d$, where d
is the spatial dimension. Returns (VX,VY,VZ)
, EToV
. When only one K
is specified, it assumes a uniform mesh with K
elements in each coordinate direction.
K
can also be specified using a keyword argument K1D
, e.g., uniform_mesh(elem; K1D = 16)
.
StartUpDG.vtk_order
— Methodvtk_order(elem::Hex, order)
Construct all node-points of a VTKLAGRANGEHEXAHEDRON of order order
. The corner-nodes are given by the reference hexahedron used by StartUpDG in the order defined by vtk.
StartUpDG.vtk_order
— Methodvtk_order(elem::Quad, order)
Construct all node-points of a VTKLAGRANGEQUAD of order order
. The corner-nodes are given by the reference quadrilateral used by StartUpDG in the order defined by vtk
StartUpDG.vtk_order
— Methodvtk_order(elem::Tri, order)
Construct all node-points of a VTK_LAGRANGE_TRIANGLE
of order order
. The corner-nodes are given by the reference-triangle used by StartUpDG in the order defined by vtk
StartUpDG.vtk_order
— Methodvtk_order(elem::Wedge, order)
Construct all node-points of a VTKLAGRANGEWEDGE of order order
. The corner-nodes are given by the reference-wedge used by StartUpDG
StartUpDG.wedge_vtk_order
— Methodwedge_vtk_order(corner_verts, order, dim)
Compute the coordinates of a VTKLAGRANGEWEDGE of order order
defined by the coordinates of the vertices given in corner_verts
. dim
is the dimension of the coordinates given. Inspired by: https://github.com/ju-kreber/paraview-scripts/blob/master/node_ordering.py
Triangulate.triangulate
— Functionfunction Triangulate.triangulate(triin::TriangulateIO, maxarea, minangle=20)
Convenience routine to avoid writing out @sprintf
each time. Returns a TriangulateIO
object.