[P]arallel [Hi]gh-order [Li]brary for [P]DEs
Latest
Parallel High-Order Library for PDEs through hp-adaptive Discontinuous Galerkin methods
|
Functional base class. More...
#include <functional.h>
Public Member Functions | |
virtual | ~Functional ()=default |
Destructor. | |
Functional (std::shared_ptr< PHiLiP::DGBase< dim, real, MeshType >> _dg, const bool _uses_solution_values=true, const bool _uses_solution_gradient=true) | |
Functional (std::shared_ptr< PHiLiP::DGBase< dim, real, MeshType >> _dg, std::shared_ptr< PHiLiP::Physics::PhysicsBase< dim, nstate, Sacado::Fad::DFad< Sacado::Fad::DFad< real >> >> _physics_fad_fad, const bool _uses_solution_values=true, const bool _uses_solution_gradient=true) | |
void | set_state (const dealii::LinearAlgebra::distributed::Vector< real > &solution_set) |
void | set_geom (const dealii::LinearAlgebra::distributed::Vector< real > &volume_nodes_set) |
virtual real | evaluate_functional (const bool compute_dIdW=false, const bool compute_dIdX=false, const bool compute_d2I=false) |
Evaluates the functional derivative with respect to the solution variable. More... | |
dealii::LinearAlgebra::distributed::Vector< real > | evaluate_dIdw_finiteDifferences (DGBase< dim, real, MeshType > &dg, const PHiLiP::Physics::PhysicsBase< dim, nstate, real > &physics, const double stepsize) |
dealii::LinearAlgebra::distributed::Vector< real > | evaluate_dIdX_finiteDifferences (DGBase< dim, real, MeshType > &dg, const PHiLiP::Physics::PhysicsBase< dim, nstate, real > &physics, const double stepsize) |
Public Attributes | |
std::shared_ptr< DGBase< dim, real, MeshType > > | dg |
Smart pointer to DGBase. | |
real | current_functional_value |
Store the functional value from the last time evaluate_functional() was called. | |
dealii::LinearAlgebra::distributed::Vector< real > | dIdw |
Vector for storing the derivatives with respect to each solution DoF. | |
dealii::LinearAlgebra::distributed::Vector< real > | dIdX |
Vector for storing the derivatives with respect to each grid DoF. | |
std::shared_ptr< dealii::TrilinosWrappers::SparseMatrix > | d2IdWdW |
Sparse matrix for storing the functional partial second derivatives. | |
std::shared_ptr< dealii::TrilinosWrappers::SparseMatrix > | d2IdWdX |
Sparse matrix for storing the functional partial second derivatives. | |
std::shared_ptr< dealii::TrilinosWrappers::SparseMatrix > | d2IdXdX |
Sparse matrix for storing the functional partial second derivatives. | |
Protected Member Functions | |
void | allocate_derivatives (const bool compute_dIdW, const bool compute_dIdX, const bool compute_d2I) |
Allocate and setup the derivative vectors/matrices. More... | |
void | allocate_dIdX (dealii::LinearAlgebra::distributed::Vector< real > &dIdX) const |
Allocate and setup the derivative dIdX vector. More... | |
void | set_derivatives (const bool compute_dIdW, const bool compute_dIdX, const bool compute_d2I, const Sacado::Fad::DFad< Sacado::Fad::DFad< real >> volume_local_sum, std::vector< dealii::types::global_dof_index > cell_soln_dofs_indices, std::vector< dealii::types::global_dof_index > cell_metric_dofs_indices) |
Set the derivative vectors/matrices. More... | |
void | need_compute (bool &compute_value, bool &compute_dIdW, bool &compute_dIdX, bool &compute_d2I) |
Checks which derivatives actually need to be recomputed. More... | |
template<typename real2 > | |
real2 | evaluate_volume_cell_functional (const Physics::PhysicsBase< dim, nstate, real2 > &physics, const std::vector< real2 > &soln_coeff, const dealii::FESystem< dim > &fe_solution, const std::vector< real2 > &coords_coeff, const dealii::FESystem< dim > &fe_metric, const dealii::Quadrature< dim > &volume_quadrature) const |
Templated function to evaluate a cell's volume functional. | |
virtual real | evaluate_volume_cell_functional (const Physics::PhysicsBase< dim, nstate, real > &physics, const std::vector< real > &soln_coeff, const dealii::FESystem< dim > &fe_solution, const std::vector< real > &coords_coeff, const dealii::FESystem< dim > &fe_metric, const dealii::Quadrature< dim > &volume_quadrature) const |
Corresponding real function to evaluate a cell's volume functional. | |
virtual Sacado::Fad::DFad< Sacado::Fad::DFad< real > > | evaluate_volume_cell_functional (const Physics::PhysicsBase< dim, nstate, Sacado::Fad::DFad< Sacado::Fad::DFad< real >>> &physics_fad_fad, const std::vector< Sacado::Fad::DFad< Sacado::Fad::DFad< real >> > &soln_coeff, const dealii::FESystem< dim > &fe_solution, const std::vector< Sacado::Fad::DFad< Sacado::Fad::DFad< real >> > &coords_coeff, const dealii::FESystem< dim > &fe_metric, const dealii::Quadrature< dim > &volume_quadrature) const |
Corresponding FadFadType function to evaluate a cell's volume functional. | |
template<typename real2 > | |
real2 | evaluate_boundary_cell_functional (const Physics::PhysicsBase< dim, nstate, real2 > &physics, const unsigned int boundary_id, const std::vector< real2 > &soln_coeff, const dealii::FESystem< dim > &fe_solution, const std::vector< real2 > &coords_coeff, const dealii::FESystem< dim > &fe_metric, const unsigned int face_number, const dealii::Quadrature< dim-1 > &face_quadrature) const |
Templated function to evaluate a cell's boundary functional. More... | |
virtual real | evaluate_boundary_cell_functional (const Physics::PhysicsBase< dim, nstate, real > &physics, const unsigned int boundary_id, const std::vector< real > &soln_coeff, const dealii::FESystem< dim > &fe_solution, const std::vector< real > &coords_coeff, const dealii::FESystem< dim > &fe_metric, const unsigned int face_number, const dealii::Quadrature< dim-1 > &face_quadrature) const |
Corresponding real function to evaluate a cell's boundary functional. | |
virtual Sacado::Fad::DFad< Sacado::Fad::DFad< real > > | evaluate_boundary_cell_functional (const Physics::PhysicsBase< dim, nstate, Sacado::Fad::DFad< Sacado::Fad::DFad< real >>> &physics_fad_fad, const unsigned int boundary_id, const std::vector< Sacado::Fad::DFad< Sacado::Fad::DFad< real >> > &soln_coeff, const dealii::FESystem< dim > &fe_solution, const std::vector< Sacado::Fad::DFad< Sacado::Fad::DFad< real >> > &coords_coeff, const dealii::FESystem< dim > &fe_metric, const unsigned int face_number, const dealii::Quadrature< dim-1 > &face_quadrature) const |
Corresponding FadFadType function to evaluate a cell's boundary functional. | |
virtual real | evaluate_volume_integrand (const PHiLiP::Physics::PhysicsBase< dim, nstate, real > &, const dealii::Point< dim, real > &, const std::array< real, nstate > &, const std::array< dealii::Tensor< 1, dim, real >, nstate > &) const |
Virtual function for computation of cell volume functional term. More... | |
virtual FadFadType | evaluate_volume_integrand (const PHiLiP::Physics::PhysicsBase< dim, nstate, FadFadType > &, const dealii::Point< dim, FadFadType > &, const std::array< FadFadType, nstate > &, const std::array< dealii::Tensor< 1, dim, FadFadType >, nstate > &) const |
Virtual function for Sacado computation of cell volume functional term and derivatives. More... | |
virtual real | evaluate_boundary_integrand (const PHiLiP::Physics::PhysicsBase< dim, nstate, real > &, const unsigned int, const dealii::Point< dim, real > &, const dealii::Tensor< 1, dim, real > &, const std::array< real, nstate > &, const std::array< dealii::Tensor< 1, dim, real >, nstate > &) const |
Virtual function for computation of cell boundary functional term. More... | |
virtual FadFadType | evaluate_boundary_integrand (const PHiLiP::Physics::PhysicsBase< dim, nstate, FadFadType > &, const unsigned int, const dealii::Point< dim, FadFadType > &, const dealii::Tensor< 1, dim, FadFadType > &, const std::array< FadFadType, nstate > &, const std::array< dealii::Tensor< 1, dim, FadFadType >, nstate > &) const |
Virtual function for Sacado computation of cell boundary functional term and derivatives. More... | |
Protected Attributes | |
std::shared_ptr< Physics::PhysicsBase< dim, nstate, FadFadType > > | physics_fad_fad |
Physics that should correspond to the one in DGBase. | |
const dealii::UpdateFlags | volume_update_flags = dealii::update_values | dealii::update_gradients | dealii::update_quadrature_points | dealii::update_JxW_values |
Update flags needed at volume points. | |
const dealii::UpdateFlags | face_update_flags = dealii::update_values | dealii::update_gradients | dealii::update_quadrature_points | dealii::update_JxW_values | dealii::update_normal_vectors |
Update flags needed at face points. | |
const bool | uses_solution_values |
Will evaluate solution values at quadrature points. | |
const bool | uses_solution_gradient |
Will evaluate solution gradient at quadrature points. | |
dealii::ConditionalOStream | pcout |
Parallel std::cout that only outputs on mpi_rank==0. | |
Private Types | |
using | FadType = Sacado::Fad::DFad< real > |
Sacado AD type for first derivatives. | |
using | FadFadType = Sacado::Fad::DFad< FadType > |
Sacado AD type that allows 2nd derivatives. | |
Private Member Functions | |
void | init_vectors () |
Private Attributes | |
dealii::LinearAlgebra::distributed::Vector< double > | solution_value |
dealii::LinearAlgebra::distributed::Vector< double > | volume_nodes_value |
dealii::LinearAlgebra::distributed::Vector< double > | solution_dIdW |
dealii::LinearAlgebra::distributed::Vector< double > | volume_nodes_dIdW |
dealii::LinearAlgebra::distributed::Vector< double > | solution_dIdX |
dealii::LinearAlgebra::distributed::Vector< double > | volume_nodes_dIdX |
dealii::LinearAlgebra::distributed::Vector< double > | solution_d2I |
dealii::LinearAlgebra::distributed::Vector< double > | volume_nodes_d2I |
Functional base class.
This base class is used to compute a functional of interest (for example, lift or drag) involving integration over the discretized volume and boundary of the domain. Often this is written in the form
\[ \mathcal{J}\left( \mathbf{u} \right) = \int_{\Omega} g_{\Omega} \left( \mathbf{u} \right) \mathrm{d} \Omega + \int_{\Gamma} g_{\Gamma} \left( \mathbf{u} \right) \mathrm{d} \Gamma \]
where the cellwise or boundary edgewise functions \( g_{\Omega} \left( \mathbf{u} \right) \) and \( g_{\Gamma} \left( \mathbf{u} \right) \) are to be overridden in the derived class. Also computes the functional derivatives which are involved in the computation of the adjoint. If derivatives are needed, the Sacado versions of these functions must also be defined.
Definition at line 43 of file functional.h.
|
explicit |
Constructor. Since we don't have access to the Physics through DGBase, we recreate a Physics based on the parameter file of DGBase. However, this will not work if the physics have been overriden through DGWeak::set_physics() as seen in the diffusion_exact_adjoint test case.
PHiLiP::Functional< dim, nstate, real, MeshType >::Functional | ( | std::shared_ptr< PHiLiP::DGBase< dim, real, MeshType >> | _dg, |
std::shared_ptr< PHiLiP::Physics::PhysicsBase< dim, nstate, Sacado::Fad::DFad< Sacado::Fad::DFad< real >> >> | _physics_fad_fad, | ||
const bool | _uses_solution_values = true , |
||
const bool | _uses_solution_gradient = true |
||
) |
Constructor. Uses provided physics instead of creating a new one based on DGBase
Definition at line 365 of file functional.cpp.
|
protected |
Allocate and setup the derivative vectors/matrices.
Helper function to simplify the evaluate_functional
Definition at line 400 of file functional.cpp.
|
protected |
Allocate and setup the derivative dIdX vector.
Helper function to simplify the evaluate_functional
Definition at line 388 of file functional.cpp.
|
protected |
Templated function to evaluate a cell's boundary functional.
Used only in the computation of evaluate_function(). If not overriden returns 0. Used only in the computation of evaluate_dIdw(). If not overriden returns 0.
Definition at line 565 of file functional.cpp.
|
inlineprotectedvirtual |
Virtual function for computation of cell boundary functional term.
Used only in the computation of evaluate_function(). If not overriden returns 0.
Reimplemented in PHiLiP::OutletPressureIntegral< dim, nstate, real, MeshType >, PHiLiP::FunctionalErrorNormLpBoundary< dim, nstate, real, MeshType >, PHiLiP::FunctionalWeightedIntegralBoundary< dim, nstate, real, MeshType >, PHiLiP::FunctionalNormLpBoundary< dim, nstate, real, MeshType >, L2_Norm_Functional< dim, nstate, real >, L2_Norm_Functional< dim, nstate, real >, PHiLiP::LiftDragFunctional< dim, nstate, real, MeshType >, L2_Norm_Functional< dim, nstate, real >, and L2_Norm_Functional< dim, nstate, real >.
Definition at line 285 of file functional.h.
|
inlineprotectedvirtual |
Virtual function for Sacado computation of cell boundary functional term and derivatives.
Used only in the computation of evaluate_dIdw(). If not overriden returns 0.
Reimplemented in PHiLiP::OutletPressureIntegral< dim, nstate, real, MeshType >, PHiLiP::FunctionalErrorNormLpBoundary< dim, nstate, real, MeshType >, PHiLiP::FunctionalWeightedIntegralBoundary< dim, nstate, real, MeshType >, PHiLiP::FunctionalNormLpBoundary< dim, nstate, real, MeshType >, L2_Norm_Functional< dim, nstate, real >, L2_Norm_Functional< dim, nstate, real >, PHiLiP::LiftDragFunctional< dim, nstate, real, MeshType >, L2_Norm_Functional< dim, nstate, real >, and L2_Norm_Functional< dim, nstate, real >.
Definition at line 296 of file functional.h.
dealii::LinearAlgebra::distributed::Vector< real > PHiLiP::Functional< dim, nstate, real, MeshType >::evaluate_dIdw_finiteDifferences | ( | PHiLiP::DGBase< dim, real, MeshType > & | dg, |
const PHiLiP::Physics::PhysicsBase< dim, nstate, real > & | physics, | ||
const double | stepsize | ||
) |
Finite difference evaluation of dIdW to verify against analytical.
Definition at line 984 of file functional.cpp.
dealii::LinearAlgebra::distributed::Vector< real > PHiLiP::Functional< dim, nstate, real, MeshType >::evaluate_dIdX_finiteDifferences | ( | PHiLiP::DGBase< dim, real, MeshType > & | dg, |
const PHiLiP::Physics::PhysicsBase< dim, nstate, real > & | physics, | ||
const double | stepsize | ||
) |
Finite difference evaluation of dIdX to verify against analytical.
Definition at line 1107 of file functional.cpp.
|
virtual |
Evaluates the functional derivative with respect to the solution variable.
Loops over the discretized domain and determines the sensitivity of the functional value to each solution node. Computed from
\[ \left. \frac{\partial \mathcal{J}_h}{\partial \mathbf{u}} \right|_{\mathbf{u}_h} = \sum_{k=1}^{N_e} \left. \frac{\partial}{\partial \mathbf{u}} \int_{\Omega_h^k} g_{\Omega} \left( \mathbf{u}_h^k \right) \mathrm{d} \Omega \right|_{\mathbf{u}_h} + \sum_{k=1}^{N_b} \left. \frac{\partial}{\partial \mathbf{u}} \int_{\Gamma_h^k} g_{\Gamma} \left( \mathbf{u}_h^k \right) \mathrm{d} \Gamma \right|_{\mathbf{u}_h} \]
Calls the functions evaluate_volume_integrand() and evaluate_cell_boundary() to be overridden
Reimplemented in PHiLiP::TargetFunctional< dim, nstate, real >, PHiLiP::LiftDragFunctional< dim, nstate, real, MeshType >, and PHiLiP::TargetWallPressure< dim, nstate, real >.
Definition at line 795 of file functional.cpp.
|
inlineprotectedvirtual |
Virtual function for computation of cell volume functional term.
Used only in the computation of evaluate_function(). If not overriden returns 0.
Reimplemented in PHiLiP::SolutionIntegral< dim, nstate, real, MeshType >, PHiLiP::FunctionalErrorNormLpVolume< dim, nstate, real, MeshType >, PHiLiP::FunctionalWeightedIntegralVolume< dim, nstate, real, MeshType >, PHiLiP::FunctionalNormLpVolume< dim, nstate, real, MeshType >, L2_Norm_Functional< dim, nstate, real >, PHiLiP::Tests::DiffusionFunctional< dim, nstate, real >, L2_Norm_Functional< dim, nstate, real >, PHiLiP::LiftDragFunctional< dim, nstate, real, MeshType >, L2_Norm_Functional< dim, nstate, real >, L2_Norm_Functional< dim, nstate, real >, and PHiLiP::Tests::L2normError< dim, nstate, real >.
Definition at line 268 of file functional.h.
|
inlineprotectedvirtual |
Virtual function for Sacado computation of cell volume functional term and derivatives.
Used only in the computation of evaluate_dIdw(). If not overriden returns 0.
Reimplemented in PHiLiP::SolutionIntegral< dim, nstate, real, MeshType >, PHiLiP::FunctionalErrorNormLpVolume< dim, nstate, real, MeshType >, PHiLiP::FunctionalWeightedIntegralVolume< dim, nstate, real, MeshType >, PHiLiP::FunctionalNormLpVolume< dim, nstate, real, MeshType >, L2_Norm_Functional< dim, nstate, real >, PHiLiP::Tests::DiffusionFunctional< dim, nstate, real >, L2_Norm_Functional< dim, nstate, real >, PHiLiP::LiftDragFunctional< dim, nstate, real, MeshType >, L2_Norm_Functional< dim, nstate, real >, L2_Norm_Functional< dim, nstate, real >, and PHiLiP::Tests::L2normError< dim, nstate, real >.
Definition at line 277 of file functional.h.
|
private |
Initializes/allocates the vectors used to check if we recompute the functional or its derivatives.
Definition at line 341 of file functional.cpp.
|
protected |
Checks which derivatives actually need to be recomputed.
If the stored solution and mesh are the same as the one used to previously compute the derivative, then we do not need to recompute them.
Definition at line 690 of file functional.cpp.
|
protected |
Set the derivative vectors/matrices.
Helper function to simplify the evaluate_functional
Definition at line 436 of file functional.cpp.
void PHiLiP::Functional< dim, nstate, real, MeshType >::set_geom | ( | const dealii::LinearAlgebra::distributed::Vector< real > & | volume_nodes_set | ) |
Set the associated DGBase's HighOrderGrid's volume_nodes to volume_nodes_set
.
Definition at line 382 of file functional.cpp.
void PHiLiP::Functional< dim, nstate, real, MeshType >::set_state | ( | const dealii::LinearAlgebra::distributed::Vector< real > & | solution_set | ) |
Set the associated DGBase's solution to solution_set
.
Definition at line 376 of file functional.cpp.
|
private |
Modal coefficients of the solution used to compute d2I last Will be used to avoid recomputing d2I.
Definition at line 157 of file functional.h.
|
private |
Modal coefficients of the solution used to compute dIdW last Will be used to avoid recomputing dIdW.
Definition at line 143 of file functional.h.
|
private |
Modal coefficients of the solution used to compute dIdX last Will be used to avoid recomputing dIdX.
Definition at line 150 of file functional.h.
|
private |
Modal coefficients of the solution used to compute dIdW last Will be used to avoid recomputing dIdW.
Definition at line 136 of file functional.h.
|
private |
Modal coefficients of the grid nodes used to compute d2I last Will be used to avoid recomputing d2I.
Definition at line 160 of file functional.h.
|
private |
Modal coefficients of the grid nodes used to compute dIdW last Will be used to avoid recomputing dIdW.
Definition at line 146 of file functional.h.
|
private |
Modal coefficients of the grid nodes used to compute dIdX last Will be used to avoid recomputing dIdX.
Definition at line 153 of file functional.h.
|
private |
Modal coefficients of the grid nodes used to compute dIdW last Will be used to avoid recomputing dIdW.
Definition at line 139 of file functional.h.