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

TargetFunctional base class. More...

#include <target_functional.h>

Inheritance diagram for PHiLiP::TargetFunctional< dim, nstate, real >:
Collaboration diagram for PHiLiP::TargetFunctional< dim, nstate, real >:

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)
 
- Public Member Functions inherited from PHiLiP::Functional< dim, nstate, real >
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...
 
- Protected Member Functions inherited from PHiLiP::Functional< dim, nstate, real >
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.
 
- Protected Attributes inherited from PHiLiP::Functional< dim, nstate, real >
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

- Public Attributes inherited from PHiLiP::Functional< dim, nstate, real >
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.
 

Detailed Description

template<int dim, int nstate, typename real>
class PHiLiP::TargetFunctional< dim, nstate, real >

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.

Constructor & Destructor Documentation

◆ TargetFunctional() [1/4]

template<int dim, int nstate, typename real >
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.

◆ TargetFunctional() [2/4]

template<int dim, int nstate, typename real >
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.

◆ TargetFunctional() [3/4]

template<int dim, int nstate, typename real >
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.

◆ TargetFunctional() [4/4]

template<int dim, int nstate, typename real >
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.

Member Function Documentation

◆ evaluate_boundary_integrand() [1/2]

template<int dim, int nstate, typename real >
virtual real PHiLiP::TargetFunctional< dim, nstate, 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
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.

◆ evaluate_boundary_integrand() [2/2]

template<int dim, int nstate, typename real >
virtual FadFadType PHiLiP::TargetFunctional< dim, nstate, real >::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
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.

◆ evaluate_dIdw_finiteDifferences()

template<int dim, int nstate, typename real >
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.

◆ evaluate_dIdX_finiteDifferences()

template<int dim, int nstate, typename real >
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.

◆ evaluate_functional()

template<int dim, int nstate, typename real >
real PHiLiP::TargetFunctional< dim, nstate, real >::evaluate_functional ( const bool  compute_dIdW = false,
const bool  compute_dIdX = false,
const bool  compute_d2I = false 
)
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.

◆ evaluate_volume_integrand() [1/2]

template<int dim, int nstate, typename real >
virtual real PHiLiP::TargetFunctional< dim, nstate, 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
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.

◆ evaluate_volume_integrand() [2/2]

template<int dim, int nstate, typename real >
virtual FadFadType PHiLiP::TargetFunctional< dim, nstate, real >::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
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.


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