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

Manufactured solution used for grid studies to check convergence orders. More...

#include <manufactured_solution.h>

Inheritance diagram for PHiLiP::ManufacturedSolutionFunction< dim, real >:
Collaboration diagram for PHiLiP::ManufacturedSolutionFunction< dim, real >:

Public Member Functions

 ManufacturedSolutionFunction (const unsigned int nstate=1)
 Constructor that initializes base_values, amplitudes, frequencies. More...
 
virtual real value (const dealii::Point< dim, real > &point, const unsigned int istate=0) const =0
 Manufactured solution exact value. More...
 
virtual dealii::Tensor< 1, dim, real > gradient (const dealii::Point< dim, real > &point, const unsigned int istate=0) const =0
 Gradient of the exact manufactured solution. More...
 
dealii::Tensor< 1, dim, real > gradient_fd (const dealii::Point< dim, real > &point, const unsigned int istate=0) const
 Uses finite-difference to evaluate the gradient.
 
virtual dealii::SymmetricTensor< 2, dim, real > hessian (const dealii::Point< dim, real > &point, const unsigned int istate=0) const =0
 Hessian of the exact manufactured solution. More...
 
dealii::SymmetricTensor< 2, dim, real > hessian_fd (const dealii::Point< dim, real > &point, const unsigned int istate=0) const
 Uses finite-difference to evaluate the hessian.
 
std::vector< real > stdvector_values (const dealii::Point< dim, real > &point) const
 Same as Function::values() except it returns it into a std::vector format.
 
void vector_gradient (const dealii::Point< dim, real > &p, std::vector< dealii::Tensor< 1, dim, real > > &gradients) const
 See dealii::Function<dim,real>::vector_gradient.
 

Public Attributes

const unsigned int nstate
 

Protected Attributes

std::vector< double > base_values
 
std::vector< double > amplitudes
 
std::vector< dealii::Tensor< 1, dim, real > > frequencies
 

Detailed Description

template<int dim, typename real>
class PHiLiP::ManufacturedSolutionFunction< dim, real >

Manufactured solution used for grid studies to check convergence orders.

This class also provides derivatives necessary to evaluate source terms.

Definition at line 22 of file manufactured_solution.h.

Constructor & Destructor Documentation

◆ ManufacturedSolutionFunction()

template<int dim, typename real >
PHiLiP::ManufacturedSolutionFunction< dim, real >::ManufacturedSolutionFunction ( const unsigned int  nstate = 1)
explicit

Constructor that initializes base_values, amplitudes, frequencies.

Calls the Function(const unsigned int n_components) constructor in deal.II This sets the public attribute n_components = nstate, which can then be accessed by all the other functions

Definition at line 1108 of file manufactured_solution.cpp.

Member Function Documentation

◆ gradient()

template<int dim, typename real >
virtual dealii::Tensor<1,dim,real> PHiLiP::ManufacturedSolutionFunction< dim, real >::gradient ( const dealii::Point< dim, real > &  point,
const unsigned int  istate = 0 
) const
pure virtual

◆ hessian()

template<int dim, typename real >
virtual dealii::SymmetricTensor<2,dim,real> PHiLiP::ManufacturedSolutionFunction< dim, real >::hessian ( const dealii::Point< dim, real > &  point,
const unsigned int  istate = 0 
) const
pure virtual

Hessian of the exact manufactured solution.

For example

hess_u[s][0][0] = -A[s]*freq[s][0]*freq[s][0]*sin(freq[s][0]*x)*sin(freq[s][1]*y)*sin(freq[s][2]*z);
hess_u[s][0][1] = A[s]*freq[s][0]*freq[s][1]*cos(freq[s][0]*x)*cos(freq[s][1]*y)*sin(freq[s][2]*z);
hess_u[s][0][2] = A[s]*freq[s][0]*freq[s][2]*cos(freq[s][0]*x)*sin(freq[s][1]*y)*cos(freq[s][2]*z);
hess_u[s][1][0] = A[s]*freq[s][1]*freq[s][0]*cos(freq[s][0]*x)*cos(freq[s][1]*y)*sin(freq[s][2]*z);
hess_u[s][1][1] = -A[s]*freq[s][1]*freq[s][1]*sin(freq[s][0]*x)*sin(freq[s][1]*y)*sin(freq[s][2]*z);
hess_u[s][1][2] = A[s]*freq[s][1]*freq[s][2]*sin(freq[s][0]*x)*cos(freq[s][1]*y)*cos(freq[s][2]*z);
hess_u[s][2][0] = A[s]*freq[s][2]*freq[s][0]*cos(freq[s][0]*x)*sin(freq[s][1]*y)*cos(freq[s][2]*z);
hess_u[s][2][1] = A[s]*freq[s][2]*freq[s][1]*sin(freq[s][0]*x)*cos(freq[s][1]*y)*cos(freq[s][2]*z);
hess_u[s][2][2] = -A[s]*freq[s][2]*freq[s][2]*sin(freq[s][0]*x)*sin(freq[s][1]*y)*sin(freq[s][2]*z);

Note that this term is symmetric since \(\frac{\partial u }{\partial x \partial y} = \frac{\partial u }{\partial y \partial x} \)

Implemented in PHiLiP::ManufacturedSolutionNavahBase< dim, real >, PHiLiP::ManufacturedSolutionExample< dim, real >, PHiLiP::ManufacturedSolutionQuadratic< dim, real >, PHiLiP::ManufacturedSolutionSShock< dim, real >, PHiLiP::ManufacturedSolutionBoundaryLayer< dim, real >, PHiLiP::ManufacturedSolutionAtan< dim, real >, PHiLiP::ManufacturedSolutionEvenPoly< dim, real >, PHiLiP::ManufacturedSolutionPoly< dim, real >, PHiLiP::ManufacturedSolutionExp< dim, real >, PHiLiP::ManufacturedSolutionAdd< dim, real >, PHiLiP::ManufacturedSolutionCosine< dim, real >, PHiLiP::ManufacturedSolutionSine< dim, real >, PHiLiP::ManufacturedSolutionZero< dim, real >, PHiLiP::Tests::Shocked1D1State< dim, real >, PHiLiP::Tests::ManufacturedSolutionV< dim, real >, and PHiLiP::Tests::ManufacturedSolutionU< dim, real >.

◆ value()

Member Data Documentation

◆ amplitudes

template<int dim, typename real >
std::vector<double> PHiLiP::ManufacturedSolutionFunction< dim, real >::amplitudes
protected

Constants used to manufactured solution.

Definition at line 132 of file manufactured_solution.h.

◆ base_values

template<int dim, typename real >
std::vector<double> PHiLiP::ManufacturedSolutionFunction< dim, real >::base_values
protected

Constants used to manufactured solution.

Definition at line 131 of file manufactured_solution.h.

◆ frequencies

template<int dim, typename real >
std::vector<dealii::Tensor<1,dim,real> > PHiLiP::ManufacturedSolutionFunction< dim, real >::frequencies
protected

Constants used to manufactured solution.

Definition at line 133 of file manufactured_solution.h.

◆ nstate

template<int dim, typename real >
const unsigned int PHiLiP::ManufacturedSolutionFunction< dim, real >::nstate

Corresponds to n_components in the dealii::Function

Definition at line 36 of file manufactured_solution.h.


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