Index

Functions

StartUpDG.BoundaryTagPlotterType
BoundaryTagPlotter(triout::TriangulateIO)

Plot recipe to visualize boundary tags by color. Usage: plot(BoundaryTagPlotter(triout))

source
StartUpDG.CurvedMeshType
struct CurvedMesh{T}

Mesh type indicating that the mesh has been curved. Stores the original mesh type as a field.

Fields

originalmeshtype :: T

source
StartUpDG.CutCellMeshType

CutCellMesh 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.

source
StartUpDG.MeshDataType
MeshData(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.

source
StartUpDG.MeshDataType
struct 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
source
StartUpDG.MeshDataMethod
MeshData(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.

source
StartUpDG.MeshImportOptionsType
MeshImportOptions

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.
source
StartUpDG.MeshPlotterType
MeshPlotter(rd::RefElemData, md::RefElemData)

Plot recipe to plot a (possibly curved) quadrilateral or triangular mesh. Usage: plot(MeshPlotter(...))

source
StartUpDG.MultidimensionalQuadratureType
MultidimensionalQuadrature

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())
source
StartUpDG.MultipleRefElemDataType
struct 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.

source
StartUpDG.NamedArrayPartitionType
NamedArrayPartition(; 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.

source
StartUpDG.NonConformingMeshType
Experimental implementation

This 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.

source
StartUpDG.PhysicalFrameType
`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).

source
StartUpDG.PolynomialType
Polynomial{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).

source
StartUpDG.RefElemDataType
struct 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
source
StartUpDG.RefElemDataMethod
RefElemData(elem; N, kwargs...)
RefElemData(elem, approx_type; N, kwargs...)

Keyword argument constructor for RefElemData (to "label" N via rd = RefElemData(Line(), N=3))

source
StartUpDG.RefElemDataMethod
RefElemData(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.

source
StartUpDG.RefElemDataMethod
RefElemData(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}.

source
StartUpDG.RefElemDataMethod
RefElemData(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.

source
StartUpDG.RefElemDataMethod
RefElemData(elem::Union{Line, Quad, Hex}, approximation_type::Polynomial{Gauss}, N)

Builds a rd::RefElemData with (N+1)-point Gauss quadrature in each dimension.

source
StartUpDG.RefElemDataMethod
RefElemData(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)
source
StartUpDG.RefElemDataMethod
RefElemData(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

source
StartUpDG.TensorProductQuadratureType
TensorProductQuadrature{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)))
source
StartUpDG.VertexMappedMeshType
struct 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

source
StartUpDG.VertexMeshPlotterType
VertexMeshPlotter((VX, VY), EToV, fv)
VertexMeshPlotter(triout::TriangulateIO)

Plot recipe to plot a quadrilateral or triangular mesh. Usage: plot(VertexMeshPlotter(...))

source
NodesAndModes.equi_nodesMethod
function 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.

source
StartUpDG.MeshData_to_vtkFunction

MeshDatatovtk(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

source
StartUpDG.MeshData_to_vtkMethod
MeshData_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

source
StartUpDG.build_node_mapsMethod
build_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))
source
StartUpDG.caratheodory_pruning_qrMethod
caratheodory_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.

source
StartUpDG.ck45Method
ck45()

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
source
StartUpDG.connect_meshMethod
connect_mesh(EToV,fv)

Inputs:

  • EToV is a num_elements by Nv matrix whose rows identify the Nv 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.

source
StartUpDG.connect_meshMethod
connect_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 from PathIntersections.define_regions

The keyword argument tol is the tolerance for matches between face centroids.

source
StartUpDG.construct_cut_surface_quadratureFunction
construct_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.

source
StartUpDG.construct_cut_volume_quadratureMethod
construct_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.

source
StartUpDG.estimate_hMethod
estimate_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ᵈ⁻¹).

source
StartUpDG.findlineMethod

findline(word::String, lines)

Outputs the line number of word in lines.

It is assumed that the word exists at least once in the file.

source
StartUpDG.get_boundary_face_labelsMethod
function 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 boundary
  • boundary_faces: list of faces on the boundary of the domain
source
StartUpDG.get_node_boundary_tagsMethod
function 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.

source
StartUpDG.get_num_elementsFunction

returns 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.

source
StartUpDG.hex_vtk_orderFunction
hex_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

source
StartUpDG.hybridized_SBP_operatorsMethod
function hybridized_SBP_operators(rd::RefElemData{DIMS})

Constructs hybridized SBP operators given a RefElemData. Returns operators Qrsth..., VhP, Ph.

source
StartUpDG.hybridized_SBP_operatorsMethod
hybridized_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.

source
StartUpDG.make_periodicMethod
make_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.

source
StartUpDG.n_verts_betweenMethod
n_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

source
StartUpDG.quad_vtk_orderFunction
quad_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

source
StartUpDG.read_Gmsh_2DMethod
read_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

source
StartUpDG.read_Gmsh_2D_v2Method
read_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

source
StartUpDG.read_Gmsh_2D_v4Function

For 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)
source
StartUpDG.read_Gmsh_2D_v4Method
function 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.

source
StartUpDG.refineFunction
function refine(triout, h, href = h/2)

Refinement of a previous mesh given the current mesh size h. Preserves boundary/volume tags.

source
StartUpDG.remap_element_groupingMethod

remapelementgrouping!(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]

source
StartUpDG.sort_by_axisMethod
sort_by_axis(corner_verts)

Given the points 'corner_verts' sort them in a lexicographical order and return the permutated points.

source
StartUpDG.sparse_low_order_SBP_operatorsMethod
function 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.

source
StartUpDG.subcell_limiting_operatorsMethod
Δ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())
source
StartUpDG.subtriangulated_cutcell_quadratureFunction
subtriangulated_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.

source
StartUpDG.tag_boundary_facesMethod
function 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)
source
StartUpDG.tag_boundary_facesMethod
function 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))
source
StartUpDG.triangle_vtk_orderFunction
triangle_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

source
StartUpDG.uniform_meshMethod
uniform_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).

source
StartUpDG.vtk_orderMethod
vtk_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.

source
StartUpDG.vtk_orderMethod
vtk_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

source
StartUpDG.vtk_orderMethod
vtk_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

source
StartUpDG.vtk_orderMethod
vtk_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

source
StartUpDG.wedge_vtk_orderMethod
wedge_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

source
Triangulate.triangulateFunction
function Triangulate.triangulate(triin::TriangulateIO, maxarea, minangle=20)

Convenience routine to avoid writing out @sprintf each time. Returns a TriangulateIO object.

source