[P]arallel [Hi]gh-order [Li]brary for [P]DEs
Latest
Parallel High-Order Library for PDEs through hp-adaptive Discontinuous Galerkin methods
|
Sum Factorization derived class. More...
#include <operators.h>
Public Member Functions | |
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... | |
![]() | |
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 | |
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. | |
![]() | |
const unsigned int | max_degree |
Max polynomial degree. | |
const unsigned int | max_grid_degree |
Max grid degree. | |
const int | nstate |
Number of states. | |
Additional Inherited Members | |
![]() | |
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. | |
Sum Factorization derived class.
Definition at line 131 of file operators.h.
void PHiLiP::OPERATOR::SumFactorizedOperators< dim, n_faces, real >::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.
Often, we compute a dot product in dim, where each matrix multiplictaion uses sum factorization. Example, consider taking the reference divergence of the reference flux:
\[ \nabla^r \cdot f^r = \frac{\partial f^r}{\partial \xi} + \frac{\partial f^r}{\partial \eta} + \frac{\partial f^r}{\partial \zeta} = \left( \frac{d \mathbf{\chi}(\mathbf{\xi})}{d\xi} \otimes \mathbf{\chi}(\mathbf{\eta}) \otimes \mathbf{\chi}(\mathbf{\zeta}) \right) \left(\hat{\mathbf{f}^r}\right)^T + \left( \mathbf{\chi}(\mathbf{\xi}) \otimes \frac{d\mathbf{\chi}(\mathbf{\eta})}{d\eta} \otimes \mathbf{\chi}(\mathbf{\zeta}) \right) \left(\hat{\mathbf{f}^r}\right)^T + \left( \mathbf{\chi}(\mathbf{\xi}) \otimes \mathbf{\chi}(\mathbf{\eta} \otimes \frac{d\mathbf{\chi}(\mathbf{\zeta})}{d\zeta})\right) \left(\hat{\mathbf{f}^r}\right)^T, \]
where we use sum factorization to evaluate each matrix-vector multiplication in each dim direction.
Definition at line 381 of file operators.cpp.
void PHiLiP::OPERATOR::SumFactorizedOperators< dim, n_faces, real >::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.
Note that we have the factor of 2.0 definied in this function. We also make use of the structure of the flux basis to get the matrix vector product after the Hadamard product to be \( \mathcal{O}(n^{d+1})\).
Definition at line 536 of file operators.cpp.
void PHiLiP::OPERATOR::SumFactorizedOperators< dim, n_faces, real >::Hadamard_product | ( | const dealii::FullMatrix< real > & | input_mat1, |
const dealii::FullMatrix< real > & | input_mat2, | ||
dealii::FullMatrix< real > & | output_mat | ||
) |
Computes a single Hadamard product.
For input mat1 \( A \) and input mat2 \( B \), this computes \( A \circ B = C \implies \left( C \right)_{ij} = \left( A \right)_{ij}\left( B \right)_{ij}\).
Definition at line 795 of file operators.cpp.
|
overridevirtual |
Computes the inner product between a matrix and a vector multiplied by some weight function.
That is, we compute \( \int Awu d\mathbf{\Omega}_r = \mathbf{A}^T \text{diag}(w) \mathbf{u}^T \). When using this function, pass \( \mathbf{A} \) and NOT it's transpose–the function transposes it in the first few lines.
Implements PHiLiP::OPERATOR::OperatorsBase< dim, n_faces, real >.
Definition at line 460 of file operators.cpp.
|
overridevirtual |
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.
Uses sum-factorization with BLAS techniques to solve the the matrix-vector multiplication, where the matrix is the tensor product of three one-dimensional matrices. We use the standard notation that x runs the fastest, then y, and z runs the slowest. For an operator \(\mathbf{A}\) of size \(n^d\) with \(n\) the one dimensional dofs and \( d\) the dim, for ouput row vector \( v \) and input row vector \( u \), we compute \(v^T=\mathbf{A}u^T\). Lastly, the adding allows the result to add onto the previous output_vect scaled by "factor".
Implements PHiLiP::OPERATOR::OperatorsBase< dim, n_faces, real >.
Definition at line 196 of file operators.cpp.
void PHiLiP::OPERATOR::SumFactorizedOperators< dim, n_faces, real >::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.
This is for the case where the operator of size dim is the dyadic product of the same 1D operator in each direction
Definition at line 308 of file operators.cpp.
void PHiLiP::OPERATOR::SumFactorizedOperators< dim, n_faces, real >::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.
Often times we have to interpolate to a surface, where in multiple dimensions, that's the tensor product of a surface operator with volume operators. This simplifies the function call. Explicitly, this passes basis_surf in the direction by face_number, and basis_vol in all other directions.
Definition at line 319 of file operators.cpp.
void PHiLiP::OPERATOR::SumFactorizedOperators< dim, n_faces, real >::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.
The Hadamard product comes up naturally in calculcating 2-point fluxes, so we needed an efficient way to compute it. Using the commutative property of Hadamard products: \( (A \otimes B) \circ ( C \otimes D) = (A\circ C) \otimes (B\circ D) \), we can find a "sum-factorization" type expression for \( A \circ U \), where here \( A = A_x \otimes A_y \otimes A_z \) and \( U \) is an \( n \times n \) matrix.
We make use of the flux basis being collocated on flux nodes, so the directions that aren't the derivative are identity. (Note that weights can be a diagonal matrix not necessarily identity). This results in the Hadamard product only needing \( \mathcal{O}(n^{d+1})\) flops to compute non-zero entries.
This is NOT for GENERAL Hadamard products since those are \( \mathcal{O}(n^{2d})\) .
Definition at line 682 of file operators.cpp.
std::array<dealii::FullMatrix<double>,2> PHiLiP::OPERATOR::SumFactorizedOperators< dim, n_faces, real >::oneD_surf_operator |
Stores the one dimensional surface operator.
Note that in 1D there are only 2 faces
Definition at line 360 of file operators.h.