Finite Element Operator

LFAToolkit.getdimensionFunction
getdimension(operator)

Retrieve the dimension of an operator

Returns:

  • dimension of the operator

Example

# setup
mesh = Mesh2D(1.0, 1.0);
mass = GalleryOperator("mass", 4, 4, mesh);

# verify
@assert mass.dimension == 2

# output
source
LFAToolkit.getelementmatrixFunction
getelementmatrix(operator)

Compute or retrieve the element matrix of operator for computing the symbol

Returns:

  • assembled element matrix

Mass matrix example:

# setup
mesh = Mesh2D(1.0, 1.0);
mass = GalleryOperator("mass", 4, 4, mesh);

# element matrix computation
elementmatrix = mass.elementmatrix;

# verify
u = ones(4 * 4);
v = elementmatrix * u;

total = sum(v);
@assert total ≈ 1.0

# output

Diffusion matrix example:

# setup
mesh = Mesh2D(1.0, 1.0);
diffusion = GalleryOperator("diffusion", 4, 4, mesh);

# element matrix computation
elementmatrix = diffusion.elementmatrix;

# verify
u = ones(4 * 4);
v = elementmatrix * u;

total = sum(v);
@assert abs(total) < 1e-14

# output
source
LFAToolkit.getmultiplicityFunction
getmultiplicity(operator)

Compute or retrieve the vector of node multiplicity for the operator

Returns:

  • vector of node multiplicity for the operator

Example

for dimension = 1:3
    # setup
    mesh = []
    if dimension == 1
        mesh = Mesh1D(1.0)
    elseif dimension == 2
        mesh = Mesh2D(1.0, 1.0)
    elseif dimension == 3
        mesh = Mesh3D(1.0, 1.0, 1.0)
    end
    diffusion = GalleryOperator("diffusion", 5, 5, mesh)

    # compute multiplicity
    multiplicity = diffusion.multiplicity

    # verify
    multiplicity1D = [2.0 1.0 1.0 1.0 2.0]'
    if dimension == 1
        @assert multiplicity ≈ multiplicity1D
    elseif dimension == 2
        @assert multiplicity ≈ kron(multiplicity1D, multiplicity1D)
    elseif dimension == 3
        @assert multiplicity ≈ kron(multiplicity1D, multiplicity1D, multiplicity1D)
    end
end

# output
source
LFAToolkit.getrowmodemapMethod
getrowmodemap(operator)

Compute or retrieve the matrix mapping the rows of the element matrix to the symbol matrix

Returns:

  • matrix mapping rows of element matrix to symbol matrix

Example:

# setup
mesh = Mesh1D(1.0);
mass = GalleryOperator("mass", 4, 4, mesh);

# verify row mode map
@assert mass.rowmodemap ≈ [1 0 0 1; 0 1 0 0; 0 0 1 0]

# output
source
LFAToolkit.getcolumnmodemapMethod
getcolumnmodemap(operator)

Compute or retrieve the matrix mapping the columns of the element matrix to the symbol matrix

Returns:

  • matrix mapping columns of element matrix to symbol matrix

Example:

# setup
mesh = Mesh1D(1.0);
mass = GalleryOperator("mass", 4, 4, mesh);

# verify column mode map
@assert mass.columnmodemap ≈ [1 0 0; 0 1 0; 0 0 1; 1 0 0]

# output
source
LFAToolkit.getnodecoordinatedifferencesMethod
getnodecoordinatedifferences(operator)

Compute or retrieve the array of differences in coordinates between nodes

Returns:

  • array of differences in coordinates between nodes

Example:

# setup
mesh = Mesh1D(1.0);
mass = GalleryOperator("mass", 4, 4, mesh);

# compute node coordinate differences
nodedifferences = mass.nodecoordinatedifferences;

# verify
truenodes = [-1, -√(1 / 5), √(1 / 5), 1];
truenodedifferences = [(truenodes[j] - truenodes[i]) / 2.0 for i = 1:4, j = 1:4];
@assert nodedifferences ≈ truenodedifferences

# output
source
LFAToolkit.massoperatorFunction
massoperator(basis, mesh)

Convenience constructor for mass operator

Weak form:

  • $\int v u$

Arguments:

  • basis::AbstractBasis: basis for all operator fields to use
  • mesh::Mesh: mesh for operator

Returns:

  • mass matrix operator with basis on mesh

Example:

# mass operator
mesh = Mesh2D(1.0, 1.0);
basis = TensorH1LagrangeBasis(3, 4, 1, mesh.dimension);
mass = LFAToolkit.massoperator(basis, mesh);

# verify
println(mass)

# output

finite element operator:
2d mesh:
    dx: 1.0
    dy: 1.0

2 inputs:
operator field:
  tensor product basis:
    numbernodes1d: 3
    numberquadraturepoints1d: 4
    numbercomponents: 1
    dimension: 2
  evaluation mode:
    interpolation
operator field:
  tensor product basis:
    numbernodes1d: 3
    numberquadraturepoints1d: 4
    numbercomponents: 1
    dimension: 2
  evaluation mode:
    quadratureweights

1 output:
operator field:
  tensor product basis:
    numbernodes1d: 3
    numberquadraturepoints1d: 4
    numbercomponents: 1
    dimension: 2
  evaluation mode:
    interpolation
source
LFAToolkit.diffusionoperatorFunction
diffusionoperator(basis, mesh)

Convenience constructor for diffusion operator

Weak form:

  • $\int \nabla v \nabla u$

Arguments:

  • basis::AbstractBasis: basis for all operator fields to use
  • mesh::Mesh: mesh for operator

Returns:

  • diffusion operator with basis on mesh

Example:

# diffusion operator
mesh = Mesh2D(1.0, 1.0);
basis = TensorH1LagrangeBasis(3, 4, 1, mesh.dimension);
diffusion = LFAToolkit.diffusionoperator(basis, mesh);

# verify
println(diffusion)

# output

finite element operator:
2d mesh:
    dx: 1.0
    dy: 1.0

2 inputs:
operator field:
  tensor product basis:
    numbernodes1d: 3
    numberquadraturepoints1d: 4
    numbercomponents: 1
    dimension: 2
  evaluation mode:
    gradient
operator field:
  tensor product basis:
    numbernodes1d: 3
    numberquadraturepoints1d: 4
    numbercomponents: 1
    dimension: 2
  evaluation mode:
    quadratureweights

1 output:
operator field:
  tensor product basis:
    numbernodes1d: 3
    numberquadraturepoints1d: 4
    numbercomponents: 1
    dimension: 2
  evaluation mode:
    gradient
source
LFAToolkit.advectionoperatorFunction
advectionoperator(basis, mesh; parameters)

Convenience constructor for advection operator

Weak form:

  • $\int \nabla v u$

Arguments:

  • basis::AbstractBasis: basis for all operator fields to use
  • mesh::Mesh: mesh for operator

Keyword Arguments:

  • parameters::NamedTuple = ([wind = [1., 1.],): named tuple of model parameters, defines wind speed

Returns:

  • advection operator with basis on mesh

Example:

# advection operator
mesh = Mesh2D(1.0, 1.0);
mapping = hale_trefethen_strip_transformation(1.4);
basis = TensorH1LagrangeBasis(3, 4, 1, mesh.dimension, mapping = mapping);
parameters = (wind = [1.0, 1.0],);
advection = LFAToolkit.advectionoperator(basis, mesh; parameters = parameters);

# verify
println(advection)

# output

finite element operator:
2d mesh:
    dx: 1.0
    dy: 1.0

2 inputs:
operator field:
  tensor product basis:
    numbernodes1d: 3
    numberquadraturepoints1d: 4
    numbercomponents: 1
    dimension: 2
  evaluation mode:
    interpolation
operator field:
  tensor product basis:
    numbernodes1d: 3
    numberquadraturepoints1d: 4
    numbercomponents: 1
    dimension: 2
  evaluation mode:
    quadratureweights

1 output:
operator field:
  tensor product basis:
    numbernodes1d: 3
    numberquadraturepoints1d: 4
    numbercomponents: 1
    dimension: 2
  evaluation mode:
    gradient
source
LFAToolkit.supgadvectionoperatorFunction
supgadvectionoperator(basis, mesh; parameters)

Convenience constructor for SUPG advection operator

Weak form: Right hand side

  • $\int wind \nabla v u - wind wind τ \nabla v \nabla u$

Arguments:

  • basis::AbstractBasis: basis for all operator fields to use
  • mesh::Mesh: mesh for operator

Keyword Arguments:

  • parameters::NamedTuple = ([wind = [1., 1.], τ = 1.0): named tuple of model parameters, defines wind speed and SUPG scaling

Returns:

  • SUPG advection operator with basis on mesh

Example:

# supg advection operator
mesh = Mesh2D(1.0, 1.0);
basis = TensorH1LagrangeBasis(3, 4, 1, mesh.dimension);
parameters = (wind = [1.0, 1.0], τ = 1.0);
supgadvection = LFAToolkit.supgadvectionoperator(basis, mesh; parameters = parameters);

# verify
println(supgadvection)

# output

finite element operator:
2d mesh:
    dx: 1.0
    dy: 1.0

2 inputs:
operator field:
  tensor product basis:
    numbernodes1d: 3
    numberquadraturepoints1d: 4
    numbercomponents: 1
    dimension: 2
  evaluation modes:
    interpolation
    gradient
operator field:
  tensor product basis:
    numbernodes1d: 3
    numberquadraturepoints1d: 4
    numbercomponents: 1
    dimension: 2
  evaluation mode:
    quadratureweights

1 output:
operator field:
  tensor product basis:
    numbernodes1d: 3
    numberquadraturepoints1d: 4
    numbercomponents: 1
    dimension: 2
  evaluation mode:
    gradient
source
LFAToolkit.supgmassoperatorFunction
supgmassoperator(basis, mesh; parameters)

Convenience constructor for SUPG mass matrix operator

Weak form: Left hand side

  • $\int v u_t + wind τ u_t \nabla v$

Arguments:

  • basis::AbstractBasis: basis for all operator fields to use
  • mesh::Mesh: mesh for operator

Keyword Arguments:

  • parameters::NamedTuple = ([wind = [1., 1.], τ = 1.0): named tuple of model parameters, defines wind speed and SUPG scaling

Returns:

  • SUPG mass matrix operator with basis on mesh

Example:

# supg mass matrix operator
mesh = Mesh2D(1.0, 1.0);
basis = TensorH1LagrangeBasis(3, 4, 1, mesh.dimension);
parameters = (wind = [1.0, 1.0], τ = 1.0);
supgmass = LFAToolkit.supgmassoperator(basis, mesh; parameters = parameters);

# verify
println(supgmass)

# output

finite element operator:
2d mesh:
    dx: 1.0
    dy: 1.0

2 inputs:
operator field:
  tensor product basis:
    numbernodes1d: 3
    numberquadraturepoints1d: 4
    numbercomponents: 1
    dimension: 2
  evaluation mode:
    interpolation
operator field:
  tensor product basis:
    numbernodes1d: 3
    numberquadraturepoints1d: 4
    numbercomponents: 1
    dimension: 2
  evaluation mode:
    quadratureweights

1 output:
operator field:
  tensor product basis:
    numbernodes1d: 3
    numberquadraturepoints1d: 4
    numbercomponents: 1
    dimension: 2
  evaluation modes:
    interpolation
    gradient
source