[P]arallel [Hi]gh-order [Li]brary for [P]DEs
Latest
Parallel High-Order Library for PDEs through hp-adaptive Discontinuous Galerkin methods
|
TargetFunctional base class. More...
#include <target_functional.h>
Public 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. | |
Public Member Functions | |
TargetFunctional (std::shared_ptr< DGBase< dim, real >> dg_input, const bool uses_solution_values=true, const bool uses_solution_gradient=true) | |
Constructor. More... | |
TargetFunctional (std::shared_ptr< DGBase< dim, real >> dg_input, const dealii::LinearAlgebra::distributed::Vector< real > &target_solution, const bool uses_solution_values=true, const bool uses_solution_gradient=true) | |
Constructor. More... | |
TargetFunctional (std::shared_ptr< DGBase< dim, real >> dg_input, std::shared_ptr< 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. More... | |
TargetFunctional (std::shared_ptr< DGBase< dim, real >> dg_input, const dealii::LinearAlgebra::distributed::Vector< real > &target_solution, std::shared_ptr< 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. More... | |
virtual real | evaluate_functional (const bool compute_dIdW=false, const bool compute_dIdX=false, const bool compute_d2I=false) override |
Evaluates the functional derivative with respect to the solution variable. More... | |
dealii::LinearAlgebra::distributed::Vector< real > | evaluate_dIdw_finiteDifferences (DGBase< dim, real > &dg, const PHiLiP::Physics::PhysicsBase< dim, nstate, real > &physics, const double stepsize) |
dealii::LinearAlgebra::distributed::Vector< real > | evaluate_dIdX_finiteDifferences (DGBase< dim, real > &dg, const PHiLiP::Physics::PhysicsBase< dim, nstate, real > &physics, const double stepsize) |
![]() | |
virtual | ~Functional ()=default |
Destructor. | |
Functional (std::shared_ptr< PHiLiP::DGBase< dim, real, dealii::parallel::distributed::Triangulation< dim > >> _dg, const bool _uses_solution_values=true, const bool _uses_solution_gradient=true) | |
Functional (std::shared_ptr< PHiLiP::DGBase< dim, real, dealii::parallel::distributed::Triangulation< dim > >> _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) |
dealii::LinearAlgebra::distributed::Vector< real > | evaluate_dIdw_finiteDifferences (DGBase< dim, real, dealii::parallel::distributed::Triangulation< dim > > &dg, const PHiLiP::Physics::PhysicsBase< dim, nstate, real > &physics, const double stepsize) |
dealii::LinearAlgebra::distributed::Vector< real > | evaluate_dIdX_finiteDifferences (DGBase< dim, real, dealii::parallel::distributed::Triangulation< dim > > &dg, const PHiLiP::Physics::PhysicsBase< dim, nstate, real > &physics, const double stepsize) |
Protected Member Functions | |
virtual real | evaluate_volume_cell_functional (const Physics::PhysicsBase< dim, nstate, real > &physics, const std::vector< real > &soln_coeff, const std::vector< real > &target_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 std::vector< real > &target_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. | |
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 std::vector< real > &target_soln_coeff, const dealii::FESystem< dim > &fe_solution, const std::vector< real > &coords_coeff, const dealii::FESystem< dim > &fe_metric, const dealii::Quadrature< dim-1 > &face_quadrature, const unsigned int face_number) const |
Corresponding real function to evaluate a cell's face 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 std::vector< real > &target_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-1 > &face_quadrature, const unsigned int face_number) const |
Corresponding FadFadType function to evaluate a cell's face functional. | |
virtual real | evaluate_volume_integrand (const PHiLiP::Physics::PhysicsBase< dim, nstate, real > &, const dealii::Point< dim, real > &, const std::array< real, nstate > &soln_at_q, const std::array< real, nstate > &target_soln_at_q, const std::array< dealii::Tensor< 1, dim, 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 > &soln_at_q, const std::array< real, nstate > &target_soln_at_q, const std::array< dealii::Tensor< 1, dim, 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 > &soln_at_q, const std::array< real, nstate > &target_soln_at_q, const std::array< dealii::Tensor< 1, dim, real >, nstate > &, const std::array< dealii::Tensor< 1, dim, real >, nstate > &) const |
Virtual function for computation of cell face 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 > &soln_at_q, const std::array< real, nstate > &target_soln_at_q, const std::array< dealii::Tensor< 1, dim, FadFadType >, nstate > &, const std::array< dealii::Tensor< 1, dim, FadFadType >, nstate > &) const |
Virtual function for Sacado computation of cell face functional term and derivatives. More... | |
![]() | |
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... | |
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. | |
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 | |
const dealii::LinearAlgebra::distributed::Vector< real > | target_solution |
Solution used to evaluate target functional. | |
![]() | |
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 |
Update flags needed at volume points. | |
const dealii::UpdateFlags | face_update_flags |
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 Member Functions | |
template<typename real2 > | |
real2 | evaluate_volume_cell_functional (const Physics::PhysicsBase< dim, nstate, real2 > &physics, const std::vector< real2 > &soln_coeff, const std::vector< real > &target_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. | |
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 std::vector< real > &target_soln_coeff, const dealii::FESystem< dim > &fe_solution, const std::vector< real2 > &coords_coeff, const dealii::FESystem< dim > &fe_metric, const dealii::Quadrature< dim-1 > &face_quadrature, const unsigned int face_number) const |
Templated function to evaluate a cell's face functional. | |
Additional Inherited Members | |
![]() | |
std::shared_ptr< DGBase< dim, real, dealii::parallel::distributed::Triangulation< dim > > > | 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. | |
TargetFunctional base class.
This base class is used to compute an inverse functional involving integration over the discretized volume and boundary of the domain. It differs from the Functional class by storing a target solution and having default volume and surface integrals corresponding to the L2-norm of the solution difference.
Often this is written in the form
\[ \mathcal{J}\left( \mathbf{u} \right) = \int_{\Omega} \left( p_{\Omega} \left( \mathbf{u} \right) - g_{\Omega} \left( \mathbf{u}_t \right) \right)^2 \mathrm{d} \Omega + \int_{\Gamma} \left( p_{\Gamma} \left( \mathbf{u} \right) - g_{\Gamma} \left( \mathbf{u}_t \right) \right)^2 \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 44 of file target_functional.h.
PHiLiP::TargetFunctional< dim, nstate, real >::TargetFunctional | ( | std::shared_ptr< DGBase< dim, real >> | dg_input, |
const bool | uses_solution_values = true , |
||
const bool | uses_solution_gradient = true |
||
) |
Constructor.
The target solution is initialized with the one currently within DGBase. 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.
Definition at line 51 of file target_functional.cpp.
PHiLiP::TargetFunctional< dim, nstate, real >::TargetFunctional | ( | std::shared_ptr< DGBase< dim, real >> | dg_input, |
const dealii::LinearAlgebra::distributed::Vector< real > & | target_solution, | ||
const bool | uses_solution_values = true , |
||
const bool | uses_solution_gradient = true |
||
) |
Constructor.
The target solution is provided instead of using the current solution in the DG object. 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.
Definition at line 64 of file target_functional.cpp.
PHiLiP::TargetFunctional< dim, nstate, real >::TargetFunctional | ( | std::shared_ptr< DGBase< dim, real >> | dg_input, |
std::shared_ptr< 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 base on DGBase. The target solution is initialized with the one currently within DGBase.
Definition at line 79 of file target_functional.cpp.
PHiLiP::TargetFunctional< dim, nstate, real >::TargetFunctional | ( | std::shared_ptr< DGBase< dim, real >> | dg_input, |
const dealii::LinearAlgebra::distributed::Vector< real > & | target_solution, | ||
std::shared_ptr< 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 base on DGBase The target solution is provided instead of using the current solution in the DG object.
Definition at line 89 of file target_functional.cpp.
|
inlineprotectedvirtual |
Virtual function for computation of cell face functional term.
Used only in the computation of evaluate_function(). If not overriden returns 0.
Reimplemented in PHiLiP::TargetWallPressure< dim, nstate, real >.
Definition at line 290 of file target_functional.h.
|
inlineprotectedvirtual |
Virtual function for Sacado computation of cell face functional term and derivatives.
Used only in the computation of evaluate_dIdw(). If not overriden returns 0.
Reimplemented in PHiLiP::TargetWallPressure< dim, nstate, real >.
Definition at line 309 of file target_functional.h.
dealii::LinearAlgebra::distributed::Vector< real > PHiLiP::TargetFunctional< dim, nstate, real >::evaluate_dIdw_finiteDifferences | ( | PHiLiP::DGBase< dim, real > & | dg, |
const PHiLiP::Physics::PhysicsBase< dim, nstate, real > & | physics, | ||
const double | stepsize | ||
) |
Finite difference evaluation of dIdW.
Definition at line 463 of file target_functional.cpp.
dealii::LinearAlgebra::distributed::Vector< real > PHiLiP::TargetFunctional< dim, nstate, real >::evaluate_dIdX_finiteDifferences | ( | PHiLiP::DGBase< dim, real > & | dg, |
const PHiLiP::Physics::PhysicsBase< dim, nstate, real > & | physics, | ||
const double | stepsize | ||
) |
Finite difference evaluation of dIdX.
Definition at line 580 of file target_functional.cpp.
|
overridevirtual |
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 from PHiLiP::Functional< dim, nstate, real >.
Reimplemented in PHiLiP::TargetWallPressure< dim, nstate, real >.
Definition at line 313 of file target_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::Tests::BoundaryInverseTarget< dim, nstate, real >, PHiLiP::TargetWallPressure< dim, nstate, real >, and PHiLiP::TargetBoundaryFunctional< dim, nstate, real >.
Definition at line 255 of file target_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::Tests::BoundaryInverseTarget< dim, nstate, real >, PHiLiP::TargetWallPressure< dim, nstate, real >, and PHiLiP::TargetBoundaryFunctional< dim, nstate, real >.
Definition at line 273 of file target_functional.h.