[P]arallel [Hi]gh-order [Li]brary for [P]DEs
Latest
Parallel High-Order Library for PDEs through hp-adaptive Discontinuous Galerkin methods
|
#include <lift_drag.hpp>
Public Types | |
enum | Functional_types { lift, drag } |
Switch between lift and drag functional types. | |
Public Member Functions | |
LiftDragFunctional (std::shared_ptr< DGBase< dim, real, MeshType >> dg_input, const Functional_types functional_type) | |
Constructor. | |
real | evaluate_functional (const bool compute_dIdW=false, const bool compute_dIdX=false, const bool compute_d2I=false) override |
Destructor. | |
template<typename real2 > | |
real2 | evaluate_boundary_integrand (const PHiLiP::Physics::PhysicsBase< dim, nstate, real2 > &physics, const unsigned int boundary_id, const dealii::Point< dim, real2 > &, const dealii::Tensor< 1, dim, real2 > &normal, const std::array< real2, nstate > &soln_at_q, const std::array< dealii::Tensor< 1, dim, real2 >, nstate > &) const |
Virtual function for computation of cell boundary functional term. More... | |
virtual real | evaluate_boundary_integrand (const PHiLiP::Physics::PhysicsBase< dim, nstate, real > &physics, const unsigned int boundary_id, const dealii::Point< dim, real > &phys_coord, const dealii::Tensor< 1, dim, real > &normal, const std::array< real, nstate > &soln_at_q, const std::array< dealii::Tensor< 1, dim, real >, nstate > &soln_grad_at_q) const override |
Virtual function for computation of cell boundary functional term. More... | |
virtual FadFadType | evaluate_boundary_integrand (const PHiLiP::Physics::PhysicsBase< dim, nstate, FadFadType > &physics, const unsigned int boundary_id, const dealii::Point< dim, FadFadType > &phys_coord, const dealii::Tensor< 1, dim, FadFadType > &normal, const std::array< FadFadType, nstate > &soln_at_q, const std::array< dealii::Tensor< 1, dim, FadFadType >, nstate > &soln_grad_at_q) const override |
Virtual function for Sacado computation of cell boundary functional term and derivatives. More... | |
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 | ~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) |
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) |
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 | |
double | initialize_force_dimensionalization_factor () |
Compute force dimensionalization factor. | |
dealii::Tensor< 2, dim, double > | initialize_rotation_matrix (const double angle_of_attack) |
Initialize rotation matrix based on given angle of attack. | |
dealii::Tensor< 1, dim, double > | initialize_lift_vector (const dealii::Tensor< 2, dim, double > &rotation_matrix) |
Initialize lift vector with given rotation matrix based on angle of attack. More... | |
dealii::Tensor< 1, dim, double > | initialize_drag_vector (const dealii::Tensor< 2, dim, double > &rotation_matrix) |
Initialize drag vector with given rotation matrix based on angle of attack. More... | |
Private Attributes | |
const Functional_types | functional_type |
Switches between lift and drag. | |
const Physics::Euler< dim, dim+2, FadFadType > & | euler_fad_fad |
Casts DG's physics into an Euler physics reference. | |
const double | angle_of_attack |
Angle of attack retrieved from euler_fad_fad. | |
const dealii::Tensor< 2, dim, double > | rotation_matrix |
Rotation matrix based on angle of attack. | |
const dealii::Tensor< 1, dim, double > | lift_vector |
Lift force scaling based on the rotation matrix applied on a [0 1]^T vector. Assumes that the lift is the force in the positive y-direction. | |
const dealii::Tensor< 1, dim, double > | drag_vector |
Drag force scaling based on the rotation matrix applied on a [1 0]^T vector. Assumes that the drag is the force in the positive x-direction. | |
dealii::Tensor< 1, dim, double > | force_vector |
Used force scaling vector depending whether this functional represents lift or drag. | |
const double | force_dimensionalization_factor |
Pressure induced drag is given by. More... | |
Additional Inherited Members | |
![]() | |
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. | |
![]() | |
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. | |
![]() | |
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. | |
Target boundary values. Simply zero out the default volume contribution.
Definition at line 16 of file lift_drag.hpp.
|
inline |
Virtual function for computation of cell boundary functional term.
Used only in the computation of evaluate_function(). If not overriden returns 0.
Definition at line 92 of file lift_drag.hpp.
|
inlineoverridevirtual |
Virtual function for computation of cell boundary functional term.
Used only in the computation of evaluate_function(). If not overriden returns 0.
Reimplemented from PHiLiP::Functional< dim, nstate, real, MeshType >.
Definition at line 118 of file lift_drag.hpp.
|
inlineoverridevirtual |
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 from PHiLiP::Functional< dim, nstate, real, MeshType >.
Definition at line 137 of file lift_drag.hpp.
|
inlinevirtual |
Virtual function for computation of cell volume functional term.
Used only in the computation of evaluate_function(). If not overriden returns 0.
Reimplemented from PHiLiP::Functional< dim, nstate, real, MeshType >.
Definition at line 156 of file lift_drag.hpp.
|
inlinevirtual |
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 from PHiLiP::Functional< dim, nstate, real, MeshType >.
Definition at line 165 of file lift_drag.hpp.
|
private |
Initialize drag vector with given rotation matrix based on angle of attack.
Use convention that drag is in the positive x-direction.
Definition at line 85 of file lift_drag.cpp.
|
private |
Initialize lift vector with given rotation matrix based on angle of attack.
Use convention that lift is in the positive y-direction.
Definition at line 64 of file lift_drag.cpp.
|
private |
Pressure induced drag is given by.
\[ C_D = \frac{2}{L \rho_{\infty} V_{\infty}^2} \int_{S} p (\mathbf{n} \cdot \text{DragVector}) ds \]
\[ C_D = \text{force_dimensionalization_factor} \int_{S} p (\mathbf{n} \cdot \text{DragVector}) ds \]
Definition at line 61 of file lift_drag.hpp.