[P]arallel [Hi]gh-order [Li]brary for [P]DEs  Latest
Parallel High-Order Library for PDEs through hp-adaptive Discontinuous Galerkin methods
PHiLiP::OPERATOR::metric_operators< real, dim, n_faces > Class Template Reference

Base metric operators class that stores functions used in both the volume and on surface. More...

#include <operators.h>

Inheritance diagram for PHiLiP::OPERATOR::metric_operators< real, dim, n_faces >:
Collaboration diagram for PHiLiP::OPERATOR::metric_operators< real, dim, n_faces >:

Public Member Functions

 metric_operators (const int nstate_input, const unsigned int max_degree_input, const unsigned int grid_degree_input, const bool store_vol_flux_nodes_input=false, const bool store_surf_flux_nodes_input=false, const bool store_Jacobian_input=false)
 Constructor.
 
void transform_physical_to_reference (const dealii::Tensor< 1, dim, real > &phys, const dealii::Tensor< 2, dim, real > &metric_cofactor, dealii::Tensor< 1, dim, real > &ref)
 Given a physical tensor, return the reference tensor.
 
void transform_reference_to_physical (const dealii::Tensor< 1, dim, real > &ref, const dealii::Tensor< 2, dim, real > &metric_cofactor, dealii::Tensor< 1, dim, real > &phys)
 Given a reference tensor, return the physical tensor.
 
void transform_physical_to_reference_vector (const dealii::Tensor< 1, dim, std::vector< real >> &phys, const dealii::Tensor< 2, dim, std::vector< real >> &metric_cofactor, dealii::Tensor< 1, dim, std::vector< real >> &ref)
 Given a physical tensor of vector of points, return the reference tensor of vector.
 
void transform_reference_unit_normal_to_physical_unit_normal (const unsigned int n_quad_pts, const dealii::Tensor< 1, dim, real > &ref, const dealii::Tensor< 2, dim, std::vector< real >> &metric_cofactor, std::vector< dealii::Tensor< 1, dim, real >> &phys)
 Given a reference tensor, return the physical tensor.
 
void build_determinant_volume_metric_Jacobian (const unsigned int n_quad_pts, const unsigned int n_metric_dofs, const std::array< std::vector< real >, dim > &mapping_support_points, mapping_shape_functions< dim, n_faces, real > &mapping_basis)
 Builds just the determinant of the volume metric determinant.
 
void build_volume_metric_operators (const unsigned int n_quad_pts, const unsigned int n_metric_dofs, const std::array< std::vector< real >, dim > &mapping_support_points, mapping_shape_functions< dim, n_faces, real > &mapping_basis, const bool use_invariant_curl_form=false)
 Builds the volume metric operators. More...
 
void build_facet_metric_operators (const unsigned int iface, const unsigned int n_quad_pts, const unsigned int n_metric_dofs, const std::array< std::vector< real >, dim > &mapping_support_points, mapping_shape_functions< dim, n_faces, real > &mapping_basis, const bool use_invariant_curl_form=false)
 Builds the facet metric operators. More...
 
- Public Member Functions inherited from PHiLiP::OPERATOR::SumFactorizedOperators< dim, n_faces, real >
 SumFactorizedOperators (const int nstate_input, const unsigned int max_degree_input, const unsigned int grid_degree_input)
 Precompute 1D operator in constructor.
 
void matrix_vector_mult (const std::vector< real > &input_vect, std::vector< real > &output_vect, const dealii::FullMatrix< double > &basis_x, const dealii::FullMatrix< double > &basis_y, const dealii::FullMatrix< double > &basis_z, const bool adding=false, const double factor=1.0) override
 Computes a matrix-vector product using sum-factorization. Pass the one-dimensional basis, where x runs the fastest, then y, and z runs the slowest. Also, assume each one-dimensional basis is the same size. More...
 
void divergence_matrix_vector_mult (const dealii::Tensor< 1, dim, std::vector< real >> &input_vect, std::vector< real > &output_vect, const dealii::FullMatrix< double > &basis_x, const dealii::FullMatrix< double > &basis_y, const dealii::FullMatrix< double > &basis_z, const dealii::FullMatrix< double > &gradient_basis_x, const dealii::FullMatrix< double > &gradient_basis_y, const dealii::FullMatrix< double > &gradient_basis_z)
 Computes the divergence using the sum factorization matrix-vector multiplication. More...
 
void divergence_matrix_vector_mult_1D (const dealii::Tensor< 1, dim, std::vector< real >> &input_vect, std::vector< real > &output_vect, const dealii::FullMatrix< double > &basis, const dealii::FullMatrix< double > &gradient_basis)
 Computes the divergence using sum-factorization where the basis are the same in each direction.
 
void gradient_matrix_vector_mult (const std::vector< real > &input_vect, dealii::Tensor< 1, dim, std::vector< real >> &output_vect, const dealii::FullMatrix< double > &basis_x, const dealii::FullMatrix< double > &basis_y, const dealii::FullMatrix< double > &basis_z, const dealii::FullMatrix< double > &gradient_basis_x, const dealii::FullMatrix< double > &gradient_basis_y, const dealii::FullMatrix< double > &gradient_basis_z)
 Computes the gradient of a scalar using sum-factorization.
 
void gradient_matrix_vector_mult_1D (const std::vector< real > &input_vect, dealii::Tensor< 1, dim, std::vector< real >> &output_vect, const dealii::FullMatrix< double > &basis, const dealii::FullMatrix< double > &gradient_basis)
 Computes the gradient of a scalar using sum-factorization where the basis are the same in each direction.
 
void inner_product (const std::vector< real > &input_vect, const std::vector< real > &weight_vect, std::vector< real > &output_vect, const dealii::FullMatrix< double > &basis_x, const dealii::FullMatrix< double > &basis_y, const dealii::FullMatrix< double > &basis_z, const bool adding=false, const double factor=1.0) override
 Computes the inner product between a matrix and a vector multiplied by some weight function. More...
 
void divergence_two_pt_flux_Hadamard_product (const dealii::Tensor< 1, dim, dealii::FullMatrix< real >> &input_mat, std::vector< real > &output_vect, const std::vector< real > &weights, const dealii::FullMatrix< double > &basis, const double scaling=2.0)
 Computes the divergence of the 2pt flux Hadamard products, then sums the rows. More...
 
void surface_two_pt_flux_Hadamard_product (const dealii::FullMatrix< real > &input_mat, std::vector< real > &output_vect_vol, std::vector< real > &output_vect_surf, const std::vector< real > &weights, const std::array< dealii::FullMatrix< double >, 2 > &surf_basis, const unsigned int iface, const unsigned int dim_not_zero, const double scaling=2.0)
 Computes the surface cross Hadamard products for skew-symmetric form from Eq. (15) in Chan, Jesse. "Skew-symmetric entropy stable modal discontinuous Galerkin formulations." Journal of Scientific Computing 81.1 (2019): 459-485.
 
void two_pt_flux_Hadamard_product (const dealii::FullMatrix< real > &input_mat, dealii::FullMatrix< real > &output_mat, const dealii::FullMatrix< double > &basis, const std::vector< real > &weights, const int direction)
 Computes the Hadamard product ONLY for 2pt flux calculations. More...
 
void sum_factorized_Hadamard_sparsity_pattern (const unsigned int rows_size, const unsigned int columns_size, std::vector< std::array< unsigned int, dim >> &rows, std::vector< std::array< unsigned int, dim >> &columns)
 Computes the rows and columns vectors with non-zero indices for sum-factorized Hadamard products.
 
void sum_factorized_Hadamard_basis_assembly (const unsigned int rows_size_1D, const unsigned int columns_size_1D, const std::vector< std::array< unsigned int, dim >> &rows, const std::vector< std::array< unsigned int, dim >> &columns, const dealii::FullMatrix< double > &basis, const std::vector< double > &weights, std::array< dealii::FullMatrix< double >, dim > &basis_sparse)
 Constructs the \( n^d \times n\) basis operator storing all non-zero entries for a "sum-factorized" Hadamard product.
 
void sum_factorized_Hadamard_surface_sparsity_pattern (const unsigned int rows_size, const unsigned int columns_size, std::vector< unsigned int > &rows, std::vector< unsigned int > &columns, const int dim_not_zero)
 Computes the rows and columns vectors with non-zero indices for surface sum-factorized Hadamard products.
 
void sum_factorized_Hadamard_surface_basis_assembly (const unsigned int rows_size, const unsigned int columns_size_1D, const std::vector< unsigned int > &rows, const std::vector< unsigned int > &columns, const dealii::FullMatrix< double > &basis, const std::vector< double > &weights, dealii::FullMatrix< double > &basis_sparse, const int dim_not_zero)
 Constructs the \( n^{d-1} \times n\) basis operator storing all non-zero entries for a "sum-factorized" surface Hadamard product.
 
void matrix_vector_mult_1D (const std::vector< real > &input_vect, std::vector< real > &output_vect, const dealii::FullMatrix< double > &basis_x, const bool adding=false, const double factor=1.0)
 Apply the matrix vector operation using the 1D operator in each direction. More...
 
void inner_product_1D (const std::vector< real > &input_vect, const std::vector< real > &weight_vect, std::vector< real > &output_vect, const dealii::FullMatrix< double > &basis_x, const bool adding=false, const double factor=1.0)
 Apply the inner product operation using the 1D operator in each direction.
 
void matrix_vector_mult_surface_1D (const unsigned int face_number, const std::vector< real > &input_vect, std::vector< real > &output_vect, const std::array< dealii::FullMatrix< double >, 2 > &basis_surf, const dealii::FullMatrix< double > &basis_vol, const bool adding=false, const double factor=1.0)
 Apply sum-factorization matrix vector multiplication on a surface. More...
 
void inner_product_surface_1D (const unsigned int face_number, const std::vector< real > &input_vect, const std::vector< real > &weight_vect, std::vector< real > &output_vect, const std::array< dealii::FullMatrix< double >, 2 > &basis_surf, const dealii::FullMatrix< double > &basis_vol, const bool adding=false, const double factor=1.0)
 Apply sum-factorization inner product on a surface.
 
void Hadamard_product (const dealii::FullMatrix< real > &input_mat1, const dealii::FullMatrix< real > &input_mat2, dealii::FullMatrix< real > &output_mat)
 Computes a single Hadamard product. More...
 
- Public Member Functions inherited from PHiLiP::OPERATOR::OperatorsBase< dim, n_faces, real >
virtual ~OperatorsBase ()=default
 Destructor.
 
 OperatorsBase (const int nstate_input, const unsigned int max_degree_input, const unsigned int grid_degree_input)
 Constructor.
 
dealii::FullMatrix< double > tensor_product (const dealii::FullMatrix< double > &basis_x, const dealii::FullMatrix< double > &basis_y, const dealii::FullMatrix< double > &basis_z)
 Returns the tensor product of matrices passed.
 
dealii::FullMatrix< double > tensor_product_state (const int nstate, const dealii::FullMatrix< double > &basis_x, const dealii::FullMatrix< double > &basis_y, const dealii::FullMatrix< double > &basis_z)
 Returns the tensor product of matrices passed, but makes it sparse diagonal by state. More...
 
double compute_factorial (double n)
 Standard function to compute factorial of a number.
 

Public Attributes

const bool store_Jacobian
 Flag if store metric Jacobian at flux nodes.
 
const bool store_vol_flux_nodes
 Flag if store metric Jacobian at flux nodes.
 
const bool store_surf_flux_nodes
 Flag if store metric Jacobian at flux nodes.
 
dealii::Tensor< 2, dim, std::vector< real > > metric_cofactor_vol
 The volume metric cofactor matrix.
 
dealii::Tensor< 2, dim, std::vector< real > > metric_cofactor_surf
 The facet metric cofactor matrix, for ONE face.
 
std::vector< real > det_Jac_vol
 The determinant of the metric Jacobian at volume cubature nodes.
 
std::vector< real > det_Jac_surf
 The determinant of the metric Jacobian at facet cubature nodes.
 
dealii::Tensor< 2, dim, std::vector< real > > metric_Jacobian_vol_cubature
 Stores the metric Jacobian at flux nodes.
 
dealii::Tensor< 1, dim, std::vector< real > > flux_nodes_vol
 Stores the physical volume flux nodes.
 
std::array< dealii::Tensor< 1, dim, std::vector< real > >, n_faces > flux_nodes_surf
 Stores the physical facet flux nodes.
 
- Public Attributes inherited from PHiLiP::OPERATOR::SumFactorizedOperators< dim, n_faces, real >
dealii::FullMatrix< double > oneD_vol_operator
 Stores the one dimensional volume operator.
 
std::array< dealii::FullMatrix< double >, 2 > oneD_surf_operator
 Stores the one dimensional surface operator. More...
 
dealii::FullMatrix< double > oneD_grad_operator
 Stores the one dimensional gradient operator.
 
std::array< dealii::FullMatrix< double >, 2 > oneD_surf_grad_operator
 Stores the one dimensional surface gradient operator.
 
- Public Attributes inherited from PHiLiP::OPERATOR::OperatorsBase< dim, n_faces, real >
const unsigned int max_degree
 Max polynomial degree.
 
const unsigned int max_grid_degree
 Max grid degree.
 
const int nstate
 Number of states.
 

Protected Member Functions

void build_metric_Jacobian (const unsigned int n_quad_pts, const std::array< std::vector< real >, dim > &mapping_support_points, const dealii::FullMatrix< double > &basis_x_flux_nodes, const dealii::FullMatrix< double > &basis_y_flux_nodes, const dealii::FullMatrix< double > &basis_z_flux_nodes, const dealii::FullMatrix< double > &grad_basis_x_flux_nodes, const dealii::FullMatrix< double > &grad_basis_y_flux_nodes, const dealii::FullMatrix< double > &grad_basis_z_flux_nodes, std::vector< dealii::Tensor< 2, dim, real >> &local_Jac)
 Builds the metric Jacobian evaluated at a vector of points. More...
 
void build_determinant_metric_Jacobian (const unsigned int n_quad_pts, const std::array< std::vector< real >, dim > &mapping_support_points, const dealii::FullMatrix< double > &basis_x_flux_nodes, const dealii::FullMatrix< double > &basis_y_flux_nodes, const dealii::FullMatrix< double > &basis_z_flux_nodes, const dealii::FullMatrix< double > &grad_basis_x_flux_nodes, const dealii::FullMatrix< double > &grad_basis_y_flux_nodes, const dealii::FullMatrix< double > &grad_basis_z_flux_nodes, std::vector< real > &det_metric_Jac)
 Assembles the determinant of metric Jacobian. More...
 
void build_local_metric_cofactor_matrix (const unsigned int n_quad_pts, const unsigned int n_metric_dofs, const std::array< std::vector< real >, dim > &mapping_support_points, const dealii::FullMatrix< double > &basis_x_grid_nodes, const dealii::FullMatrix< double > &basis_y_grid_nodes, const dealii::FullMatrix< double > &basis_z_grid_nodes, const dealii::FullMatrix< double > &basis_x_flux_nodes, const dealii::FullMatrix< double > &basis_y_flux_nodes, const dealii::FullMatrix< double > &basis_z_flux_nodes, const dealii::FullMatrix< double > &grad_basis_x_grid_nodes, const dealii::FullMatrix< double > &grad_basis_y_grid_nodes, const dealii::FullMatrix< double > &grad_basis_z_grid_nodes, const dealii::FullMatrix< double > &grad_basis_x_flux_nodes, const dealii::FullMatrix< double > &grad_basis_y_flux_nodes, const dealii::FullMatrix< double > &grad_basis_z_flux_nodes, dealii::Tensor< 2, dim, std::vector< real >> &metric_cofactor, const bool use_invariant_curl_form=false)
 Called on the fly and returns the metric cofactor at cubature nodes.
 
void compute_local_3D_cofactor (const unsigned int n_metric_dofs, const unsigned int n_quad_pts, const std::array< std::vector< real >, dim > &mapping_support_points, const dealii::FullMatrix< double > &basis_x_grid_nodes, const dealii::FullMatrix< double > &basis_y_grid_nodes, const dealii::FullMatrix< double > &basis_z_grid_nodes, const dealii::FullMatrix< double > &basis_x_flux_nodes, const dealii::FullMatrix< double > &basis_y_flux_nodes, const dealii::FullMatrix< double > &basis_z_flux_nodes, const dealii::FullMatrix< double > &grad_basis_x_grid_nodes, const dealii::FullMatrix< double > &grad_basis_y_grid_nodes, const dealii::FullMatrix< double > &grad_basis_z_grid_nodes, const dealii::FullMatrix< double > &grad_basis_x_flux_nodes, const dealii::FullMatrix< double > &grad_basis_y_flux_nodes, const dealii::FullMatrix< double > &grad_basis_z_flux_nodes, dealii::Tensor< 2, dim, std::vector< real >> &metric_cofactor, const bool use_invariant_curl_form=false)
 Computes local 3D cofactor matrix. More...
 

Additional Inherited Members

- Protected Attributes inherited from PHiLiP::OPERATOR::OperatorsBase< dim, n_faces, real >
unsigned int max_grid_degree_check
 Check to see if the metrics used are a higher order then the initialized grid.
 
const MPI_Comm mpi_communicator
 MPI communicator.
 
dealii::ConditionalOStream pcout
 Parallel std::cout that only outputs on mpi_rank==0.
 

Detailed Description

template<typename real, int dim, int n_faces>
class PHiLiP::OPERATOR::metric_operators< real, dim, n_faces >

Base metric operators class that stores functions used in both the volume and on surface.

Definition at line 1086 of file operators.h.

Member Function Documentation

◆ build_determinant_metric_Jacobian()

template<typename real , int dim, int n_faces>
void PHiLiP::OPERATOR::metric_operators< real, dim, n_faces >::build_determinant_metric_Jacobian ( const unsigned int  n_quad_pts,
const std::array< std::vector< real >, dim > &  mapping_support_points,
const dealii::FullMatrix< double > &  basis_x_flux_nodes,
const dealii::FullMatrix< double > &  basis_y_flux_nodes,
const dealii::FullMatrix< double > &  basis_z_flux_nodes,
const dealii::FullMatrix< double > &  grad_basis_x_flux_nodes,
const dealii::FullMatrix< double > &  grad_basis_y_flux_nodes,
const dealii::FullMatrix< double > &  grad_basis_z_flux_nodes,
std::vector< real > &  det_metric_Jac 
)
protected

Assembles the determinant of metric Jacobian.

\( \|J^\Omega \| = \mathbf{a}_1 \cdot (\mathbf{a}_2 \otimes \mathbf{a}_3)\), where \(\mathbf{a}_i = \mathbf{\nabla}^r x_i \) are the physical vector bases. Pass the 1D mapping shape functions evaluated at flux nodes, and the 1D gradient of mapping shape functions evaluated at flux nodes.

Definition at line 2479 of file operators.cpp.

◆ build_facet_metric_operators()

template<typename real , int dim, int n_faces>
void PHiLiP::OPERATOR::metric_operators< real, dim, n_faces >::build_facet_metric_operators ( const unsigned int  iface,
const unsigned int  n_quad_pts,
const unsigned int  n_metric_dofs,
const std::array< std::vector< real >, dim > &  mapping_support_points,
mapping_shape_functions< dim, n_faces, real > &  mapping_basis,
const bool  use_invariant_curl_form = false 
)

Builds the facet metric operators.

Builds and stores facet metric cofactor and determinant of metric Jacobian at the facet cubature nodes (one face). If passed flag to store Jacobian when metric operators is constructed, will also store the JAcobian at flux nodes.

Definition at line 2351 of file operators.cpp.

◆ build_metric_Jacobian()

template<typename real , int dim, int n_faces>
void PHiLiP::OPERATOR::metric_operators< real, dim, n_faces >::build_metric_Jacobian ( const unsigned int  n_quad_pts,
const std::array< std::vector< real >, dim > &  mapping_support_points,
const dealii::FullMatrix< double > &  basis_x_flux_nodes,
const dealii::FullMatrix< double > &  basis_y_flux_nodes,
const dealii::FullMatrix< double > &  basis_z_flux_nodes,
const dealii::FullMatrix< double > &  grad_basis_x_flux_nodes,
const dealii::FullMatrix< double > &  grad_basis_y_flux_nodes,
const dealii::FullMatrix< double > &  grad_basis_z_flux_nodes,
std::vector< dealii::Tensor< 2, dim, real >> &  local_Jac 
)
protected

Builds the metric Jacobian evaluated at a vector of points.

\( \mathbf{J} = [\mathbf{a}_1, \mathbf{a}_2, \mathbf{a}_3]^T \) where \(\mathbf{a}_i = \mathbf{\nabla}^r x_i \) are the physical vector bases.

Definition at line 2441 of file operators.cpp.

◆ build_volume_metric_operators()

template<typename real , int dim, int n_faces>
void PHiLiP::OPERATOR::metric_operators< real, dim, n_faces >::build_volume_metric_operators ( const unsigned int  n_quad_pts,
const unsigned int  n_metric_dofs,
const std::array< std::vector< real >, dim > &  mapping_support_points,
mapping_shape_functions< dim, n_faces, real > &  mapping_basis,
const bool  use_invariant_curl_form = false 
)

Builds the volume metric operators.

Builds and stores volume metric cofactor and determinant of metric Jacobian at the volume cubature nodes. If passed flag to store Jacobian when metric operators is constructed, will also store the JAcobian at flux nodes.

Definition at line 2294 of file operators.cpp.

◆ compute_local_3D_cofactor()

template<typename real , int dim, int n_faces>
void PHiLiP::OPERATOR::metric_operators< real, dim, n_faces >::compute_local_3D_cofactor ( const unsigned int  n_metric_dofs,
const unsigned int  n_quad_pts,
const std::array< std::vector< real >, dim > &  mapping_support_points,
const dealii::FullMatrix< double > &  basis_x_grid_nodes,
const dealii::FullMatrix< double > &  basis_y_grid_nodes,
const dealii::FullMatrix< double > &  basis_z_grid_nodes,
const dealii::FullMatrix< double > &  basis_x_flux_nodes,
const dealii::FullMatrix< double > &  basis_y_flux_nodes,
const dealii::FullMatrix< double > &  basis_z_flux_nodes,
const dealii::FullMatrix< double > &  grad_basis_x_grid_nodes,
const dealii::FullMatrix< double > &  grad_basis_y_grid_nodes,
const dealii::FullMatrix< double > &  grad_basis_z_grid_nodes,
const dealii::FullMatrix< double > &  grad_basis_x_flux_nodes,
const dealii::FullMatrix< double > &  grad_basis_y_flux_nodes,
const dealii::FullMatrix< double > &  grad_basis_z_flux_nodes,
dealii::Tensor< 2, dim, std::vector< real >> &  metric_cofactor,
const bool  use_invariant_curl_form = false 
)
protected

Computes local 3D cofactor matrix.

We compute the metric cofactor matrix \(\mathbf{C}_m\) via the conservative curl form of Abe 2014 and Kopriva 2006 by default. Can use invariant curl form by passing flag. To ensure consistent normals, we consider the two cubature sets, grid nodes (mapping-support-points), and flux nodes (quadrature nodes). The metric cofactor matrix is thus:

\[ (\mathbf{C})_{ni} = J(\mathbf{a}^i)_n= -\hat{\mathbf{e}}_i \cdot \nabla^r\times\mathbf{\Theta}(\mathbf{\xi}_{\text{flux nodes}}^r)\Big[ \mathbf{\Theta}(\mathbf{\xi}_{\text{grid nodes}}^r)\hat{\mathbf{x}}_l^{c^T} \nabla^r \mathbf{\Theta}(\mathbf{\xi}_{\text{grid nodes}}^r)\hat{\mathbf{x}}_m^{c^T} \Big] \text{, }\\i=1,2,3\text{, }n=1,2,3\text{ }(n,m,l)\text{ cyclic,} \]

for the conservative curl form, and

\[ (\mathbf{C})_{ni} = J(\mathbf{a}^i)_n= -\frac{1}{2}\hat{\mathbf{e}}_i \cdot \nabla^r\times\mathbf{\Theta}(\mathbf{\xi}_{\text{flux nodes}}^r)\Big[ \mathbf{\Theta}(\mathbf{\xi}_{\text{grid nodes}}^r)\hat{\mathbf{x}}_l^{c^T} \nabla^r \mathbf{\Theta}(\mathbf{\xi}_{\text{grid nodes}}^r)\hat{\mathbf{x}}_m^{c^T} - \mathbf{\Theta}(\mathbf{\xi}_{\text{grid nodes}}^r)\hat{\mathbf{x}}_m^{c^T} \nabla^r \mathbf{\Theta}(\mathbf{\xi}_{\text{grid nodes}}^r)\hat{\mathbf{x}}_l^{c^T}\Big] \text{, }\\i=1,2,3\text{, }n=1,2,3\text{ }(n,m,l)\text{ cyclic,} \]

for the invariant curl form.
We let \(\mathbf{\Theta}(\mathbf{\xi}^r)\) represent the mapping shape functions.

Definition at line 2588 of file operators.cpp.


The documentation for this class was generated from the following files: