[P]arallel [Hi]gh-order [Li]brary for [P]DEs
Latest
Parallel High-Order Library for PDEs through hp-adaptive Discontinuous Galerkin methods
|
Abstract class templated on the number of state variables. More...
#include <dg_base_state.hpp>
Public Member Functions | |
DGBaseState (const Parameters::AllParameters *const parameters_input, const unsigned int degree, const unsigned int max_degree_input, const unsigned int grid_degree_input, const std::shared_ptr< Triangulation > triangulation_input) | |
< Input parameters. More... | |
void | set_physics (std::shared_ptr< Physics::PhysicsBase< dim, nstate, real > > pde_physics_double_input, std::shared_ptr< Physics::PhysicsBase< dim, nstate, FadType > > pde_physics_fad_input, std::shared_ptr< Physics::PhysicsBase< dim, nstate, RadType > > pde_physics_rad_input, std::shared_ptr< Physics::PhysicsBase< dim, nstate, FadFadType > > pde_physics_fad_fad_input, std::shared_ptr< Physics::PhysicsBase< dim, nstate, RadFadType > > pde_physics_rad_fad_input) |
void | allocate_model_variables () |
Allocate the necessary variables declared in src/physics/model.h. | |
void | update_model_variables () |
Update the necessary variables declared in src/physics/model.h. | |
void | set_use_auxiliary_eq () |
Set use_auxiliary_eq flag. | |
![]() | |
virtual | ~DGBase ()=default |
Destructor. | |
DGBase (const int nstate_input, const Parameters::AllParameters *const parameters_input, const unsigned int degree, const unsigned int max_degree_input, const unsigned int grid_degree_input, const std::shared_ptr< Triangulation > triangulation_input) | |
Principal constructor that will call delegated constructor. More... | |
void | reinit () |
Reinitializes the DG object after a change of triangulation. More... | |
DGBase (const int nstate_input, const Parameters::AllParameters *const parameters_input, const unsigned int degree, const unsigned int max_degree_input, const unsigned int grid_degree_input, const std::shared_ptr< Triangulation > triangulation_input, const MassiveCollectionTuple collection_tuple) | |
Delegated constructor that initializes collections. More... | |
void | set_high_order_grid (std::shared_ptr< HighOrderGrid< dim, real, MeshType >> new_high_order_grid) |
Sets the associated high order grid with the provided one. | |
void | set_all_cells_fe_degree (const unsigned int degree) |
Refers to a collection Mappings, which represents the high-order grid. More... | |
unsigned int | get_max_fe_degree () |
Gets the maximum value of currently active FE degree. | |
unsigned int | get_min_fe_degree () |
Gets the minimum value of currently active FE degree. | |
dealii::Point< dim > | coordinates_of_highest_refined_cell (bool check_for_p_refined_cell=false) |
Returns the coordinates of the most refined cell. | |
virtual void | allocate_system (const bool compute_dRdW=true, const bool compute_dRdX=true, const bool compute_d2R=true) |
Allocates the system. More... | |
void | time_scale_solution_update (dealii::LinearAlgebra::distributed::Vector< double > &solution_update, const real CFL) const |
Scales a solution update with the appropriate maximum time step. More... | |
void | time_scaled_mass_matrices (const real scale) |
void | reinit_operators_for_cell_residual_loop (const unsigned int poly_degree_int, const unsigned int poly_degree_ext, const unsigned int grid_degree, OPERATOR::basis_functions< dim, 2 *dim, real > &soln_basis_int, OPERATOR::basis_functions< dim, 2 *dim, real > &soln_basis_ext, OPERATOR::basis_functions< dim, 2 *dim, real > &flux_basis_int, OPERATOR::basis_functions< dim, 2 *dim, real > &flux_basis_ext, OPERATOR::local_basis_stiffness< dim, 2 *dim, real > &flux_basis_stiffness, OPERATOR::vol_projection_operator< dim, 2 *dim, real > &soln_basis_projection_oper_int, OPERATOR::vol_projection_operator< dim, 2 *dim, real > &soln_basis_projection_oper_ext, OPERATOR::mapping_shape_functions< dim, 2 *dim, real > &mapping_basis) |
Builds needed operators for cell residual loop. | |
void | reinit_operators_for_mass_matrix (const bool Cartesian_element, const unsigned int poly_degree, const unsigned int grid_degree, OPERATOR::mapping_shape_functions< dim, 2 *dim, real > &mapping_basis, OPERATOR::basis_functions< dim, 2 *dim, real > &basis, OPERATOR::local_mass< dim, 2 *dim, real > &reference_mass_matrix, OPERATOR::local_Flux_Reconstruction_operator< dim, 2 *dim, real > &reference_FR, OPERATOR::local_Flux_Reconstruction_operator_aux< dim, 2 *dim, real > &reference_FR_aux, OPERATOR::derivative_p< dim, 2 *dim, real > &deriv_p) |
Builds needed operators to compute mass matrices/inverses efficiently. | |
void | evaluate_mass_matrices (bool do_inverse_mass_matrix=false) |
Allocates and evaluates the mass matrices for the entire grid. | |
void | evaluate_local_metric_dependent_mass_matrix_and_set_in_global_mass_matrix (const bool Cartesian_element, const bool do_inverse_mass_matrix, const unsigned int poly_degree, const unsigned int curr_grid_degree, const unsigned int n_quad_pts, const unsigned int n_dofs_cell, const std::vector< dealii::types::global_dof_index > dofs_indices, OPERATOR::metric_operators< real, dim, 2 *dim > &metric_oper, OPERATOR::basis_functions< dim, 2 *dim, real > &basis, OPERATOR::local_mass< dim, 2 *dim, real > &reference_mass_matrix, OPERATOR::local_Flux_Reconstruction_operator< dim, 2 *dim, real > &reference_FR, OPERATOR::local_Flux_Reconstruction_operator_aux< dim, 2 *dim, real > &reference_FR_aux, OPERATOR::derivative_p< dim, 2 *dim, real > &deriv_p) |
Evaluates the metric dependent local mass matrices and inverses, then sets them in the global matrices. | |
void | apply_inverse_global_mass_matrix (const dealii::LinearAlgebra::distributed::Vector< double > &input_vector, dealii::LinearAlgebra::distributed::Vector< double > &output_vector, const bool use_auxiliary_eq=false) |
Applies the inverse of the local metric dependent mass matrices when the global is not stored. More... | |
void | apply_global_mass_matrix (const dealii::LinearAlgebra::distributed::Vector< double > &input_vector, dealii::LinearAlgebra::distributed::Vector< double > &output_vector, const bool use_auxiliary_eq=false, const bool use_unmodified_mass_matrix=false) |
Applies the local metric dependent mass matrices when the global is not stored. More... | |
std::vector< real > | evaluate_time_steps (const bool exact_time_stepping) |
Evaluates the maximum stable time step. More... | |
void | add_mass_matrices (const real scale) |
Add mass matrices to the system scaled by a factor (likely time-step) More... | |
void | add_time_scaled_mass_matrices () |
Add time scaled mass matrices to the system. More... | |
double | get_residual_l2norm () const |
Returns the L2-norm of the right_hand_side vector. | |
double | get_residual_linfnorm () const |
Returns the Linf-norm of the right_hand_side vector. | |
unsigned int | n_dofs () const |
Number of degrees of freedom. | |
void | set_anisotropic_flags () |
Set anisotropic flags based on jump indicator. More... | |
template<typename real2 > | |
real2 | discontinuity_sensor (const dealii::Quadrature< dim > &volume_quadrature, const std::vector< real2 > &soln_coeff_high, const dealii::FiniteElement< dim, dim > &fe_high, const std::vector< real2 > &jac_det) |
void | set_dual (const dealii::LinearAlgebra::distributed::Vector< real > &dual_input) |
Sets the stored dual variables used to compute the dual dotted with the residual Hessians. | |
dealii::SparsityPattern | get_dRdX_sparsity_pattern () |
Evaluate SparsityPattern of dRdX. | |
dealii::SparsityPattern | get_dRdW_sparsity_pattern () |
Evaluate SparsityPattern of dRdW. | |
dealii::SparsityPattern | get_d2RdWdW_sparsity_pattern () |
Evaluate SparsityPattern of the residual Hessian dual.d2RdWdW. | |
dealii::SparsityPattern | get_d2RdXdX_sparsity_pattern () |
Evaluate SparsityPattern of the residual Hessian dual.d2RdXdX. | |
dealii::SparsityPattern | get_d2RdWdX_sparsity_pattern () |
Evaluate SparsityPattern of the residual Hessian dual.d2RdXdW. | |
dealii::SparsityPattern | get_dRdXs_sparsity_pattern () |
Evaluate SparsityPattern of dRdXs. | |
dealii::SparsityPattern | get_d2RdXsdXs_sparsity_pattern () |
Evaluate SparsityPattern of the residual Hessian dual.d2RdXsdXs. | |
dealii::SparsityPattern | get_d2RdWdXs_sparsity_pattern () |
Evaluate SparsityPattern of the residual Hessian dual.d2RdXsdW. | |
dealii::TrilinosWrappers::SparseMatrix | get_dRdX_finite_differences (dealii::SparsityPattern dRdX_sparsity_pattern) |
Evaluate dRdX using finite-differences. | |
void | initialize_manufactured_solution () |
Virtual function defined in DG. | |
void | output_results_vtk (const unsigned int cycle, const double current_time=0.0) |
Output solution. | |
void | output_face_results_vtk (const unsigned int cycle, const double current_time=0.0) |
Output Euler face solution. | |
void | assemble_residual (const bool compute_dRdW=false, const bool compute_dRdX=false, const bool compute_d2R=false, const double CFL_mass=0.0) |
Main loop of the DG class. More... | |
template<typename DoFCellAccessorType1 , typename DoFCellAccessorType2 > | |
void | assemble_cell_residual (const DoFCellAccessorType1 ¤t_cell, const DoFCellAccessorType2 ¤t_metric_cell, const bool compute_dRdW, const bool compute_dRdX, const bool compute_d2R, dealii::hp::FEValues< dim, dim > &fe_values_collection_volume, dealii::hp::FEFaceValues< dim, dim > &fe_values_collection_face_int, dealii::hp::FEFaceValues< dim, dim > &fe_values_collection_face_ext, dealii::hp::FESubfaceValues< dim, dim > &fe_values_collection_subface, dealii::hp::FEValues< dim, dim > &fe_values_collection_volume_lagrange, OPERATOR::basis_functions< dim, 2 *dim, real > &soln_basis_int, OPERATOR::basis_functions< dim, 2 *dim, real > &soln_basis_ext, OPERATOR::basis_functions< dim, 2 *dim, real > &flux_basis_int, OPERATOR::basis_functions< dim, 2 *dim, real > &flux_basis_ext, OPERATOR::local_basis_stiffness< dim, 2 *dim, real > &flux_basis_stiffness, OPERATOR::vol_projection_operator< dim, 2 *dim, real > &soln_basis_projection_oper_int, OPERATOR::vol_projection_operator< dim, 2 *dim, real > &soln_basis_projection_oper_ext, OPERATOR::mapping_shape_functions< dim, 2 *dim, real > &mapping_basis, const bool compute_auxiliary_right_hand_side, dealii::LinearAlgebra::distributed::Vector< double > &rhs, std::array< dealii::LinearAlgebra::distributed::Vector< double >, dim > &rhs_aux) |
Used in assemble_residual(). More... | |
void | set_current_time (const real current_time_input) |
Sets the current time within DG to be used for unsteady source terms. | |
void | allocate_auxiliary_equation () |
Allocates the auxiliary equations' variables and right hand side (primarily for Strong form diffusive) | |
virtual void | assemble_auxiliary_residual ()=0 |
Asembles the auxiliary equations' residuals and solves. | |
virtual void | allocate_dual_vector ()=0 |
Allocate the dual vector for optimization. More... | |
void | update_artificial_dissipation_discontinuity_sensor () |
Update discontinuity sensor. More... | |
Public Attributes | |
std::shared_ptr< Physics::PhysicsBase< dim, nstate, real > > | pde_physics_double |
Contains the physics of the PDE with real type. | |
std::shared_ptr< Physics::ModelBase< dim, nstate, real > > | pde_model_double |
Contains the model terms of the PDEType == PhysicsModel with real type. | |
std::unique_ptr< NumericalFlux::NumericalFluxConvective< dim, nstate, real > > | conv_num_flux_double |
Convective numerical flux with real type. | |
std::unique_ptr< NumericalFlux::NumericalFluxDissipative< dim, nstate, real > > | diss_num_flux_double |
Dissipative numerical flux with real type. | |
std::shared_ptr< ArtificialDissipationBase< dim, nstate > > | artificial_dissip |
Link to Artificial dissipation class (with three dissipation types, depending on the input). | |
std::shared_ptr< Physics::PhysicsBase< dim, nstate, FadType > > | pde_physics_fad |
Contains the physics of the PDE with FadType. | |
std::shared_ptr< Physics::ModelBase< dim, nstate, FadType > > | pde_model_fad |
Contains the model terms of the PDEType == PhysicsModel with FadType. | |
std::unique_ptr< NumericalFlux::NumericalFluxConvective< dim, nstate, FadType > > | conv_num_flux_fad |
Convective numerical flux with FadType. | |
std::unique_ptr< NumericalFlux::NumericalFluxDissipative< dim, nstate, FadType > > | diss_num_flux_fad |
Dissipative numerical flux with FadType. | |
std::shared_ptr< Physics::PhysicsBase< dim, nstate, RadType > > | pde_physics_rad |
Contains the physics of the PDE with RadType. | |
std::shared_ptr< Physics::ModelBase< dim, nstate, RadType > > | pde_model_rad |
Contains the model terms of the PDEType == PhysicsModel with RadType. | |
std::unique_ptr< NumericalFlux::NumericalFluxConvective< dim, nstate, RadType > > | conv_num_flux_rad |
Convective numerical flux with RadType. | |
std::unique_ptr< NumericalFlux::NumericalFluxDissipative< dim, nstate, RadType > > | diss_num_flux_rad |
Dissipative numerical flux with RadType. | |
std::shared_ptr< Physics::PhysicsBase< dim, nstate, FadFadType > > | pde_physics_fad_fad |
Contains the physics of the PDE with FadFadType. | |
std::shared_ptr< Physics::ModelBase< dim, nstate, FadFadType > > | pde_model_fad_fad |
Contains the model terms of the PDEType == PhysicsModel with FadFadType. | |
std::unique_ptr< NumericalFlux::NumericalFluxConvective< dim, nstate, FadFadType > > | conv_num_flux_fad_fad |
Convective numerical flux with FadFadType. | |
std::unique_ptr< NumericalFlux::NumericalFluxDissipative< dim, nstate, FadFadType > > | diss_num_flux_fad_fad |
Dissipative numerical flux with FadFadType. | |
std::shared_ptr< Physics::PhysicsBase< dim, nstate, RadFadType > > | pde_physics_rad_fad |
Contains the physics of the PDE with RadFadDtype. | |
std::shared_ptr< Physics::ModelBase< dim, nstate, RadFadType > > | pde_model_rad_fad |
Contains the model terms of the PDEType == PhysicsModel with RadFadType. | |
std::unique_ptr< NumericalFlux::NumericalFluxConvective< dim, nstate, RadFadType > > | conv_num_flux_rad_fad |
Convective numerical flux with RadFadDtype. | |
std::unique_ptr< NumericalFlux::NumericalFluxDissipative< dim, nstate, RadFadType > > | diss_num_flux_rad_fad |
Dissipative numerical flux with RadFadDtype. | |
![]() | |
const Parameters::AllParameters *const | all_parameters |
Pointer to all parameters. | |
const int | nstate |
Number of state variables. More... | |
const unsigned int | initial_degree |
Initial polynomial degree assigned during constructor. | |
const unsigned int | max_degree |
Maximum degree used for p-refi1nement. More... | |
const unsigned int | max_grid_degree |
Maximum grid degree used for hp-refi1nement. More... | |
std::shared_ptr< Triangulation > | triangulation |
Mesh. | |
dealii::SparsityPattern | sparsity_pattern |
Sparsity pattern used on the system_matrix. More... | |
dealii::SparsityPattern | mass_sparsity_pattern |
Sparsity pattern used on the system_matrix. More... | |
dealii::TrilinosWrappers::SparseMatrix | time_scaled_global_mass_matrix |
Global mass matrix divided by the time scales. More... | |
dealii::TrilinosWrappers::SparseMatrix | global_mass_matrix |
Global mass matrix. More... | |
dealii::TrilinosWrappers::SparseMatrix | global_inverse_mass_matrix |
Global inverser mass matrix. More... | |
dealii::TrilinosWrappers::SparseMatrix | global_mass_matrix_auxiliary |
Global auxiliary mass matrix. More... | |
dealii::TrilinosWrappers::SparseMatrix | global_inverse_mass_matrix_auxiliary |
Global inverse of the auxiliary mass matrix. | |
dealii::TrilinosWrappers::SparseMatrix | system_matrix |
dealii::TrilinosWrappers::SparseMatrix | system_matrix_transpose |
std::unique_ptr< Epetra_RowMatrixTransposer > | epetra_rowmatrixtransposer_dRdW |
Epetra_RowMatrixTransposer used to transpose the system_matrix. | |
dealii::TrilinosWrappers::SparseMatrix | dRdXv |
dealii::TrilinosWrappers::SparseMatrix | d2RdWdW |
dealii::TrilinosWrappers::SparseMatrix | d2RdXdX |
dealii::TrilinosWrappers::SparseMatrix | d2RdWdX |
dealii::LinearAlgebra::distributed::Vector< double > | right_hand_side |
Residual of the current solution. More... | |
dealii::IndexSet | locally_owned_dofs |
Locally own degrees of freedom. | |
dealii::IndexSet | ghost_dofs |
Locally relevant ghost degrees of freedom. | |
dealii::IndexSet | locally_relevant_dofs |
Union of locally owned degrees of freedom and relevant ghost degrees of freedom. | |
dealii::IndexSet | locally_owned_dofs_grid |
Locally own degrees of freedom for the grid. | |
dealii::IndexSet | ghost_dofs_grid |
Locally relevant ghost degrees of freedom for the grid. | |
dealii::IndexSet | locally_relevant_dofs_grid |
dealii::LinearAlgebra::distributed::Vector< double > | solution |
Current modal coefficients of the solution. More... | |
std::array< dealii::LinearAlgebra::distributed::Vector< double >, dim > | auxiliary_right_hand_side |
The auxiliary equations' right hand sides. | |
std::array< dealii::LinearAlgebra::distributed::Vector< double >, dim > | auxiliary_solution |
The auxiliary equations' solution. | |
dealii::Vector< double > | cell_volume |
Time it takes for the maximum wavespeed to cross the cell domain. More... | |
dealii::Vector< double > | max_dt_cell |
Time it takes for the maximum wavespeed to cross the cell domain. More... | |
dealii::Vector< double > | reduced_mesh_weights |
dealii::Vector< double > | artificial_dissipation_coeffs |
Artificial dissipation in each cell. | |
dealii::Vector< double > | artificial_dissipation_se |
Artificial dissipation error ratio sensor in each cell. | |
dealii::LinearAlgebra::distributed::Vector< real > | dual |
Current optimization dual variables corresponding to the residual constraints also known as the adjoint. More... | |
bool | update_artificial_diss |
const dealii::hp::FECollection< dim > | fe_collection |
Finite Element Collection for p-finite-element to represent the solution. More... | |
dealii::hp::QCollection< dim > | volume_quadrature_collection |
Finite Element Collection to represent the high-order grid. More... | |
dealii::hp::QCollection< dim-1 > | face_quadrature_collection |
Quadrature used to evaluate face integrals. | |
const dealii::hp::FECollection< dim > | fe_collection_lagrange |
Lagrange basis used in strong form. More... | |
const dealii::hp::FECollection< 1 > | oneD_fe_collection |
1D Finite Element Collection for p-finite-element to represent the solution More... | |
const dealii::hp::FECollection< 1 > | oneD_fe_collection_1state |
1D Finite Element Collection for p-finite-element to represent the solution for a single state. More... | |
const dealii::hp::FECollection< 1 > | oneD_fe_collection_flux |
1D collocated flux basis used in strong form More... | |
dealii::hp::QCollection< 1 > | oneD_quadrature_collection |
1D quadrature to generate Lagrange polynomials for the sake of flux interpolation. | |
dealii::QGauss< 0 > | oneD_face_quadrature |
1D surface quadrature is always one single point for all poly degrees. | |
dealii::DoFHandler< dim > | dof_handler |
Finite Element Collection to represent the high-order grid. More... | |
std::shared_ptr< HighOrderGrid< dim, real, MeshType > > | high_order_grid |
High order grid that will provide the MappingFEField. | |
double | assemble_residual_time |
Computational time for assembling residual. | |
bool | freeze_artificial_dissipation |
Flag to freeze artificial dissipation. | |
double | max_artificial_dissipation_coeff |
Stores maximum artificial dissipation while assembling the residual. | |
bool | use_auxiliary_eq |
Flag for using the auxiliary equation. | |
Protected Types | |
using | Triangulation = typename DGBase< dim, real, MeshType >::Triangulation |
Alias to base class Triangulation. | |
Protected Member Functions | |
real | evaluate_CFL (std::vector< std::array< real, nstate > > soln_at_q, const real artificial_dissipation, const real cell_diameter, const unsigned int cell_degree) |
Evaluate the time it takes for the maximum wavespeed to cross the cell domain. More... | |
void | reset_numerical_fluxes () |
Reinitializes the numerical fluxes based on the current physics. More... | |
![]() | |
virtual void | assemble_volume_term_and_build_operators (typename dealii::DoFHandler< dim >::active_cell_iterator cell, const dealii::types::global_dof_index current_cell_index, const std::vector< dealii::types::global_dof_index > &cell_dofs_indices, const std::vector< dealii::types::global_dof_index > &metric_dof_indices, const unsigned int poly_degree, const unsigned int grid_degree, OPERATOR::basis_functions< dim, 2 *dim, real > &soln_basis, OPERATOR::basis_functions< dim, 2 *dim, real > &flux_basis, OPERATOR::local_basis_stiffness< dim, 2 *dim, real > &flux_basis_stiffness, OPERATOR::vol_projection_operator< dim, 2 *dim, real > &soln_basis_projection_oper_int, OPERATOR::vol_projection_operator< dim, 2 *dim, real > &soln_basis_projection_oper_ext, OPERATOR::metric_operators< real, dim, 2 *dim > &metric_oper, OPERATOR::mapping_shape_functions< dim, 2 *dim, real > &mapping_basis, std::array< std::vector< real >, dim > &mapping_support_points, dealii::hp::FEValues< dim, dim > &fe_values_collection_volume, dealii::hp::FEValues< dim, dim > &fe_values_collection_volume_lagrange, const dealii::FESystem< dim, dim > ¤t_fe_ref, dealii::Vector< real > &local_rhs_int_cell, std::vector< dealii::Tensor< 1, dim, real >> &local_auxiliary_RHS, const bool compute_auxiliary_right_hand_side, const bool compute_dRdW, const bool compute_dRdX, const bool compute_d2R)=0 |
Builds the necessary operators/fe values and assembles volume residual. | |
virtual void | assemble_boundary_term_and_build_operators (typename dealii::DoFHandler< dim >::active_cell_iterator cell, const dealii::types::global_dof_index current_cell_index, const unsigned int iface, const unsigned int boundary_id, const real penalty, const std::vector< dealii::types::global_dof_index > &cell_dofs_indices, const std::vector< dealii::types::global_dof_index > &metric_dof_indices, const unsigned int poly_degree, const unsigned int grid_degree, OPERATOR::basis_functions< dim, 2 *dim, real > &soln_basis, OPERATOR::basis_functions< dim, 2 *dim, real > &flux_basis, OPERATOR::local_basis_stiffness< dim, 2 *dim, real > &flux_basis_stiffness, OPERATOR::vol_projection_operator< dim, 2 *dim, real > &soln_basis_projection_oper_int, OPERATOR::vol_projection_operator< dim, 2 *dim, real > &soln_basis_projection_oper_ext, OPERATOR::metric_operators< real, dim, 2 *dim > &metric_oper, OPERATOR::mapping_shape_functions< dim, 2 *dim, real > &mapping_basis, std::array< std::vector< real >, dim > &mapping_support_points, dealii::hp::FEFaceValues< dim, dim > &fe_values_collection_face_int, const dealii::FESystem< dim, dim > ¤t_fe_ref, dealii::Vector< real > &local_rhs_int_cell, std::vector< dealii::Tensor< 1, dim, real >> &local_auxiliary_RHS, const bool compute_auxiliary_right_hand_side, const bool compute_dRdW, const bool compute_dRdX, const bool compute_d2R)=0 |
Builds the necessary operators/fe values and assembles boundary residual. | |
virtual void | assemble_face_term_and_build_operators (typename dealii::DoFHandler< dim >::active_cell_iterator cell, typename dealii::DoFHandler< dim >::active_cell_iterator neighbor_cell, const dealii::types::global_dof_index current_cell_index, const dealii::types::global_dof_index neighbor_cell_index, const unsigned int iface, const unsigned int neighbor_iface, const real penalty, const std::vector< dealii::types::global_dof_index > ¤t_dofs_indices, const std::vector< dealii::types::global_dof_index > &neighbor_dofs_indices, const std::vector< dealii::types::global_dof_index > ¤t_metric_dofs_indices, const std::vector< dealii::types::global_dof_index > &neighbor_metric_dofs_indices, const unsigned int poly_degree_int, const unsigned int poly_degree_ext, const unsigned int grid_degree_int, const unsigned int grid_degree_ext, OPERATOR::basis_functions< dim, 2 *dim, real > &soln_basis_int, OPERATOR::basis_functions< dim, 2 *dim, real > &soln_basis_ext, OPERATOR::basis_functions< dim, 2 *dim, real > &flux_basis_int, OPERATOR::basis_functions< dim, 2 *dim, real > &flux_basis_ext, OPERATOR::local_basis_stiffness< dim, 2 *dim, real > &flux_basis_stiffness, OPERATOR::vol_projection_operator< dim, 2 *dim, real > &soln_basis_projection_oper_int, OPERATOR::vol_projection_operator< dim, 2 *dim, real > &soln_basis_projection_oper_ext, OPERATOR::metric_operators< real, dim, 2 *dim > &metric_oper_int, OPERATOR::metric_operators< real, dim, 2 *dim > &metric_oper_ext, OPERATOR::mapping_shape_functions< dim, 2 *dim, real > &mapping_basis, std::array< std::vector< real >, dim > &mapping_support_points, dealii::hp::FEFaceValues< dim, dim > &fe_values_collection_face_int, dealii::hp::FEFaceValues< dim, dim > &fe_values_collection_face_ext, dealii::Vector< real > ¤t_cell_rhs, dealii::Vector< real > &neighbor_cell_rhs, std::vector< dealii::Tensor< 1, dim, real >> ¤t_cell_rhs_aux, dealii::LinearAlgebra::distributed::Vector< double > &rhs, std::array< dealii::LinearAlgebra::distributed::Vector< double >, dim > &rhs_aux, const bool compute_auxiliary_right_hand_side, const bool compute_dRdW, const bool compute_dRdX, const bool compute_d2R)=0 |
Builds the necessary operators/fe values and assembles face residual. | |
virtual void | assemble_subface_term_and_build_operators (typename dealii::DoFHandler< dim >::active_cell_iterator cell, typename dealii::DoFHandler< dim >::active_cell_iterator neighbor_cell, const dealii::types::global_dof_index current_cell_index, const dealii::types::global_dof_index neighbor_cell_index, const unsigned int iface, const unsigned int neighbor_iface, const unsigned int neighbor_i_subface, const real penalty, const std::vector< dealii::types::global_dof_index > ¤t_dofs_indices, const std::vector< dealii::types::global_dof_index > &neighbor_dofs_indices, const std::vector< dealii::types::global_dof_index > ¤t_metric_dofs_indices, const std::vector< dealii::types::global_dof_index > &neighbor_metric_dofs_indices, const unsigned int poly_degree_int, const unsigned int poly_degree_ext, const unsigned int grid_degree_int, const unsigned int grid_degree_ext, OPERATOR::basis_functions< dim, 2 *dim, real > &soln_basis_int, OPERATOR::basis_functions< dim, 2 *dim, real > &soln_basis_ext, OPERATOR::basis_functions< dim, 2 *dim, real > &flux_basis_int, OPERATOR::basis_functions< dim, 2 *dim, real > &flux_basis_ext, OPERATOR::local_basis_stiffness< dim, 2 *dim, real > &flux_basis_stiffness, OPERATOR::vol_projection_operator< dim, 2 *dim, real > &soln_basis_projection_oper_int, OPERATOR::vol_projection_operator< dim, 2 *dim, real > &soln_basis_projection_oper_ext, OPERATOR::metric_operators< real, dim, 2 *dim > &metric_oper_int, OPERATOR::metric_operators< real, dim, 2 *dim > &metric_oper_ext, OPERATOR::mapping_shape_functions< dim, 2 *dim, real > &mapping_basis, std::array< std::vector< real >, dim > &mapping_support_points, dealii::hp::FEFaceValues< dim, dim > &fe_values_collection_face_int, dealii::hp::FESubfaceValues< dim, dim > &fe_values_collection_subface, dealii::Vector< real > ¤t_cell_rhs, dealii::Vector< real > &neighbor_cell_rhs, std::vector< dealii::Tensor< 1, dim, real >> ¤t_cell_rhs_aux, dealii::LinearAlgebra::distributed::Vector< double > &rhs, std::array< dealii::LinearAlgebra::distributed::Vector< double >, dim > &rhs_aux, const bool compute_auxiliary_right_hand_side, const bool compute_dRdW, const bool compute_dRdX, const bool compute_d2R)=0 |
Builds the necessary operators/fe values and assembles subface residual. | |
virtual void | assemble_volume_term_derivatives (typename dealii::DoFHandler< dim >::active_cell_iterator cell, const dealii::types::global_dof_index current_cell_index, const dealii::FEValues< dim, dim > &, const dealii::FESystem< dim, dim > &fe, const dealii::Quadrature< dim > &quadrature, const std::vector< dealii::types::global_dof_index > &metric_dof_indices, const std::vector< dealii::types::global_dof_index > &soln_dof_indices, dealii::Vector< real > &local_rhs_cell, const dealii::FEValues< dim, dim > &, const bool compute_dRdW, const bool compute_dRdX, const bool compute_d2R)=0 |
Evaluate the integral over the cell volume and the specified derivatives. More... | |
virtual void | assemble_boundary_term_derivatives (typename dealii::DoFHandler< dim >::active_cell_iterator cell, const dealii::types::global_dof_index current_cell_index, const unsigned int face_number, const unsigned int boundary_id, const dealii::FEFaceValuesBase< dim, dim > &fe_values_boundary, const real penalty, const dealii::FESystem< dim, dim > &fe, const dealii::Quadrature< dim-1 > &quadrature, const std::vector< dealii::types::global_dof_index > &metric_dof_indices, const std::vector< dealii::types::global_dof_index > &soln_dof_indices, dealii::Vector< real > &local_rhs_cell, const bool compute_dRdW, const bool compute_dRdX, const bool compute_d2R)=0 |
Evaluate the integral over the cell edges that are on domain boundaries and the specified derivatives. More... | |
virtual void | assemble_face_term_derivatives (typename dealii::DoFHandler< dim >::active_cell_iterator cell, const dealii::types::global_dof_index current_cell_index, const dealii::types::global_dof_index neighbor_cell_index, const std::pair< unsigned int, int > face_subface_int, const std::pair< unsigned int, int > face_subface_ext, const typename dealii::QProjector< dim >::DataSetDescriptor face_data_set_int, const typename dealii::QProjector< dim >::DataSetDescriptor face_data_set_ext, const dealii::FEFaceValuesBase< dim, dim > &, const dealii::FEFaceValuesBase< dim, dim > &, const real penalty, const dealii::FESystem< dim, dim > &fe_int, const dealii::FESystem< dim, dim > &fe_ext, const dealii::Quadrature< dim-1 > &face_quadrature, const std::vector< dealii::types::global_dof_index > &metric_dof_indices_int, const std::vector< dealii::types::global_dof_index > &metric_dof_indices_ext, const std::vector< dealii::types::global_dof_index > &soln_dof_indices_int, const std::vector< dealii::types::global_dof_index > &soln_dof_indices_ext, dealii::Vector< real > &local_rhs_int_cell, dealii::Vector< real > &local_rhs_ext_cell, const bool compute_dRdW, const bool compute_dRdX, const bool compute_d2R)=0 |
Evaluate the integral over the internal cell edges and its specified derivatives. More... | |
virtual void | assemble_volume_term_explicit (typename dealii::DoFHandler< dim >::active_cell_iterator cell, const dealii::types::global_dof_index current_cell_index, const dealii::FEValues< dim, dim > &fe_values_volume, const std::vector< dealii::types::global_dof_index > ¤t_dofs_indices, const std::vector< dealii::types::global_dof_index > &metric_dof_indices, const unsigned int poly_degree, const unsigned int grid_degree, dealii::Vector< real > ¤t_cell_rhs, const dealii::FEValues< dim, dim > &fe_values_lagrange)=0 |
Evaluate the integral over the cell volume. | |
Additional Inherited Members | |
![]() | |
using | Triangulation = MeshType |
using | MassiveCollectionTuple = std::tuple< dealii::hp::FECollection< dim >, dealii::hp::QCollection< dim >, dealii::hp::QCollection< dim-1 >, dealii::hp::FECollection< dim >, dealii::hp::FECollection< 1 >, dealii::hp::FECollection< 1 >, dealii::hp::FECollection< 1 >, dealii::hp::QCollection< 1 > > |
Makes for cleaner doxygen documentation. | |
![]() | |
real | current_time |
The current time set in set_current_time() | |
const dealii::FE_Q< dim > | fe_q_artificial_dissipation |
Continuous distribution of artificial dissipation. | |
dealii::DoFHandler< dim > | dof_handler_artificial_dissipation |
Degrees of freedom handler for C0 artificial dissipation. | |
dealii::LinearAlgebra::distributed::Vector< double > | artificial_dissipation_c0 |
Artificial dissipation coefficients. | |
const dealii::UpdateFlags | volume_update_flags |
Update flags needed at volume points. More... | |
const dealii::UpdateFlags | face_update_flags |
Update flags needed at face points. More... | |
const dealii::UpdateFlags | neighbor_face_update_flags = dealii::update_values | dealii::update_gradients | dealii::update_quadrature_points | dealii::update_JxW_values |
Update flags needed at neighbor' face points. More... | |
MPI_Comm | mpi_communicator |
MPI communicator. | |
dealii::ConditionalOStream | pcout |
Parallel std::cout that only outputs on mpi_rank==0. | |
Abstract class templated on the number of state variables.
Definition at line 18 of file dg_base_state.hpp.
PHiLiP::DGBaseState< dim, nstate, real, MeshType >::DGBaseState | ( | const Parameters::AllParameters *const | parameters_input, |
const unsigned int | degree, | ||
const unsigned int | max_degree_input, | ||
const unsigned int | grid_degree_input, | ||
const std::shared_ptr< Triangulation > | triangulation_input | ||
) |
|
protected |
Evaluate the time it takes for the maximum wavespeed to cross the cell domain.
Currently only uses the convective eigenvalues. Future changes would take in account the maximum diffusivity and take the minimum time between dx/conv_eig and dx*dx/max_visc to determine the minimum travel time of information.
Furthermore, a more robust implementation would convert the values to a Bezier basis where the maximum and minimum values would be bounded by the Bernstein modal coefficients.
Definition at line 187 of file dg_base_state.cpp.
|
protected |
Reinitializes the numerical fluxes based on the current physics.
Usually called after setting physics.
Definition at line 41 of file dg_base_state.cpp.
void PHiLiP::DGBaseState< dim, nstate, real, MeshType >::set_physics | ( | std::shared_ptr< Physics::PhysicsBase< dim, nstate, real > > | pde_physics_double_input, |
std::shared_ptr< Physics::PhysicsBase< dim, nstate, FadType > > | pde_physics_fad_input, | ||
std::shared_ptr< Physics::PhysicsBase< dim, nstate, RadType > > | pde_physics_rad_input, | ||
std::shared_ptr< Physics::PhysicsBase< dim, nstate, FadFadType > > | pde_physics_fad_fad_input, | ||
std::shared_ptr< Physics::PhysicsBase< dim, nstate, RadFadType > > | pde_physics_rad_fad_input | ||
) |
Change the physics object. Must provide all the AD types to ensure that the derivatives are consistent.
Definition at line 75 of file dg_base_state.cpp.