[P]arallel [Hi]gh-order [Li]brary for [P]DEs  Latest
Parallel High-Order Library for PDEs through hp-adaptive Discontinuous Galerkin methods
strong_dg.hpp
1 #ifndef __STRONG_DISCONTINUOUSGALERKIN_H__
2 #define __STRONG_DISCONTINUOUSGALERKIN_H__
3 
4 #include "dg_base_state.hpp"
5 
6 namespace PHiLiP {
7 
9 /* Contains the functions that need to be templated on the number of state variables.
10  */
11 #if PHILIP_DIM==1 // dealii::parallel::distributed::Triangulation<dim> does not work for 1D
12 template <int dim, int nstate, typename real, typename MeshType = dealii::Triangulation<dim>>
13 #else
14 template <int dim, int nstate, typename real, typename MeshType = dealii::parallel::distributed::Triangulation<dim>>
15 #endif
16 class DGStrong: public DGBaseState<dim, nstate, real, MeshType>
17 {
18 protected:
21 
22 public:
24  DGStrong(
25  const Parameters::AllParameters *const parameters_input,
26  const unsigned int degree,
27  const unsigned int max_degree_input,
28  const unsigned int grid_degree_input,
29  const std::shared_ptr<Triangulation> triangulation_input);
30 
32 
37 
39  void allocate_dual_vector ();
40 
41 private:
43  template<typename DoFCellAccessorType1, typename DoFCellAccessorType2>
45  const DoFCellAccessorType1 &current_cell,
46  const DoFCellAccessorType2 &current_metric_cell,
47  std::vector<dealii::LinearAlgebra::distributed::Vector<double>> &rhs);
48 
49 protected:
52  typename dealii::DoFHandler<dim>::active_cell_iterator cell,
53  const dealii::types::global_dof_index current_cell_index,
54  const std::vector<dealii::types::global_dof_index> &cell_dofs_indices,
55  const std::vector<dealii::types::global_dof_index> &metric_dof_indices,
56  const unsigned int poly_degree,
57  const unsigned int grid_degree,
61  OPERATOR::vol_projection_operator<dim,2*dim,real> &soln_basis_projection_oper_int,
62  OPERATOR::vol_projection_operator<dim,2*dim,real> &soln_basis_projection_oper_ext,
65  std::array<std::vector<real>,dim> &mapping_support_points,
66  dealii::hp::FEValues<dim,dim> &/*fe_values_collection_volume*/,
67  dealii::hp::FEValues<dim,dim> &/*fe_values_collection_volume_lagrange*/,
68  const dealii::FESystem<dim,dim> &/*current_fe_ref*/,
69  dealii::Vector<real> &local_rhs_int_cell,
70  std::vector<dealii::Tensor<1,dim,real>> &local_auxiliary_RHS,
71  const bool compute_auxiliary_right_hand_side,
72  const bool /*compute_dRdW*/, const bool /*compute_dRdX*/, const bool /*compute_d2R*/);
73 
76  typename dealii::DoFHandler<dim>::active_cell_iterator /*cell*/,
77  const dealii::types::global_dof_index current_cell_index,
78  const unsigned int iface,
79  const unsigned int boundary_id,
80  const real penalty,
81  const std::vector<dealii::types::global_dof_index> &cell_dofs_indices,
82  const std::vector<dealii::types::global_dof_index> &metric_dof_indices,
83  const unsigned int poly_degree,
84  const unsigned int grid_degree,
88  OPERATOR::vol_projection_operator<dim,2*dim,real> &soln_basis_projection_oper_int,
89  OPERATOR::vol_projection_operator<dim,2*dim,real> &soln_basis_projection_oper_ext,
92  std::array<std::vector<real>,dim> &mapping_support_points,
93  dealii::hp::FEFaceValues<dim,dim> &/*fe_values_collection_face_int*/,
94  const dealii::FESystem<dim,dim> &/*current_fe_ref*/,
95  dealii::Vector<real> &local_rhs_int_cell,
96  std::vector<dealii::Tensor<1,dim,real>> &local_auxiliary_RHS,
97  const bool compute_auxiliary_right_hand_side,
98  const bool /*compute_dRdW*/, const bool /*compute_dRdX*/, const bool /*compute_d2R*/);
99 
102  typename dealii::DoFHandler<dim>::active_cell_iterator /*cell*/,
103  typename dealii::DoFHandler<dim>::active_cell_iterator /*neighbor_cell*/,
104  const dealii::types::global_dof_index current_cell_index,
105  const dealii::types::global_dof_index neighbor_cell_index,
106  const unsigned int iface,
107  const unsigned int neighbor_iface,
108  const real penalty,
109  const std::vector<dealii::types::global_dof_index> &current_dofs_indices,
110  const std::vector<dealii::types::global_dof_index> &neighbor_dofs_indices,
111  const std::vector<dealii::types::global_dof_index> &current_metric_dofs_indices,
112  const std::vector<dealii::types::global_dof_index> &neighbor_metric_dofs_indices,
113  const unsigned int poly_degree_int,
114  const unsigned int poly_degree_ext,
115  const unsigned int grid_degree_int,
116  const unsigned int grid_degree_ext,
122  OPERATOR::vol_projection_operator<dim,2*dim,real> &soln_basis_projection_oper_int,
123  OPERATOR::vol_projection_operator<dim,2*dim,real> &soln_basis_projection_oper_ext,
127  std::array<std::vector<real>,dim> &mapping_support_points,
128  dealii::hp::FEFaceValues<dim,dim> &/*fe_values_collection_face_int*/,
129  dealii::hp::FEFaceValues<dim,dim> &/*fe_values_collection_face_ext*/,
130  dealii::Vector<real> &current_cell_rhs,
131  dealii::Vector<real> &neighbor_cell_rhs,
132  std::vector<dealii::Tensor<1,dim,real>> &current_cell_rhs_aux,
133  dealii::LinearAlgebra::distributed::Vector<double> &rhs,
134  std::array<dealii::LinearAlgebra::distributed::Vector<double>,dim> &rhs_aux,
135  const bool compute_auxiliary_right_hand_side,
136  const bool compute_dRdW, const bool compute_dRdX, const bool compute_d2R);
137 
139 
142  typename dealii::DoFHandler<dim>::active_cell_iterator cell,
143  typename dealii::DoFHandler<dim>::active_cell_iterator neighbor_cell,
144  const dealii::types::global_dof_index current_cell_index,
145  const dealii::types::global_dof_index neighbor_cell_index,
146  const unsigned int iface,
147  const unsigned int neighbor_iface,
148  const unsigned int /*neighbor_i_subface*/,
149  const real penalty,
150  const std::vector<dealii::types::global_dof_index> &current_dofs_indices,
151  const std::vector<dealii::types::global_dof_index> &neighbor_dofs_indices,
152  const std::vector<dealii::types::global_dof_index> &current_metric_dofs_indices,
153  const std::vector<dealii::types::global_dof_index> &neighbor_metric_dofs_indices,
154  const unsigned int poly_degree_int,
155  const unsigned int poly_degree_ext,
156  const unsigned int grid_degree_int,
157  const unsigned int grid_degree_ext,
163  OPERATOR::vol_projection_operator<dim,2*dim,real> &soln_basis_projection_oper_int,
164  OPERATOR::vol_projection_operator<dim,2*dim,real> &soln_basis_projection_oper_ext,
168  std::array<std::vector<real>,dim> &mapping_support_points,
169  dealii::hp::FEFaceValues<dim,dim> &fe_values_collection_face_int,
170  dealii::hp::FESubfaceValues<dim,dim> &/*fe_values_collection_subface*/,
171  dealii::Vector<real> &current_cell_rhs,
172  dealii::Vector<real> &neighbor_cell_rhs,
173  std::vector<dealii::Tensor<1,dim,real>> &current_cell_rhs_aux,
174  dealii::LinearAlgebra::distributed::Vector<double> &rhs,
175  std::array<dealii::LinearAlgebra::distributed::Vector<double>,dim> &rhs_aux,
176  const bool compute_auxiliary_right_hand_side,
177  const bool compute_dRdW, const bool compute_dRdX, const bool compute_d2R);
178 
179 public:
181 
186  const std::vector<dealii::types::global_dof_index> &current_dofs_indices,
187  const unsigned int poly_degree,
191  std::vector<dealii::Tensor<1,dim,real>> &local_auxiliary_RHS);
192 
193 protected:
196  const unsigned int iface,
197  const dealii::types::global_dof_index current_cell_index,
198  const unsigned int poly_degree,
199  const unsigned int boundary_id,
200  const std::vector<dealii::types::global_dof_index> &dofs_indices,
203  std::vector<dealii::Tensor<1,dim,real>> &local_auxiliary_RHS);
204 
205 public:
207 
213  const unsigned int iface,
214  const unsigned int neighbor_iface,
215  const dealii::types::global_dof_index current_cell_index,
216  const dealii::types::global_dof_index neighbor_cell_index,
217  const unsigned int poly_degree_int,
218  const unsigned int poly_degree_ext,
219  const std::vector<dealii::types::global_dof_index> &dof_indices_int,
220  const std::vector<dealii::types::global_dof_index> &dof_indices_ext,
224  std::vector<dealii::Tensor<1,dim,real>> &local_auxiliary_RHS_int,
225  std::vector<dealii::Tensor<1,dim,real>> &local_auxiliary_RHS_ext);
226 
227 protected:
229 
247  typename dealii::DoFHandler<dim>::active_cell_iterator cell,
248  const dealii::types::global_dof_index current_cell_index,
249  const std::vector<dealii::types::global_dof_index> &cell_dofs_indices,
250  const unsigned int poly_degree,
254  OPERATOR::vol_projection_operator<dim,2*dim,real> &soln_basis_projection_oper,
256  dealii::Vector<real> &local_rhs_int_cell);
257 
260  const unsigned int iface,
261  const dealii::types::global_dof_index current_cell_index,
262  const unsigned int boundary_id,
263  const unsigned int poly_degree,
264  const real penalty,
265  const std::vector<dealii::types::global_dof_index> &dof_indices,
268  OPERATOR::vol_projection_operator<dim,2*dim,real> &soln_basis_projection_oper,
270  dealii::Vector<real> &local_rhs_cell);
271 
273 
282  const unsigned int iface,
283  const unsigned int neighbor_iface,
284  const dealii::types::global_dof_index current_cell_index,
285  const dealii::types::global_dof_index neighbor_cell_index,
286  const unsigned int poly_degree_int,
287  const unsigned int poly_degree_ext,
288  const real penalty,
289  const std::vector<dealii::types::global_dof_index> &dof_indices_int,
290  const std::vector<dealii::types::global_dof_index> &dof_indices_ext,
295  OPERATOR::vol_projection_operator<dim,2*dim,real> &soln_basis_projection_oper_int,
296  OPERATOR::vol_projection_operator<dim,2*dim,real> &soln_basis_projection_oper_ext,
299  dealii::Vector<real> &local_rhs_int_cell,
300  dealii::Vector<real> &local_rhs_ext_cell);
301 
302 protected:
304 
306  typename dealii::DoFHandler<dim>::active_cell_iterator cell,
307  const dealii::types::global_dof_index current_cell_index,
308  const dealii::FEValues<dim,dim> &,//fe_values_vol,
309  const dealii::FESystem<dim,dim> &fe,
310  const dealii::Quadrature<dim> &quadrature,
311  const std::vector<dealii::types::global_dof_index> &metric_dof_indices,
312  const std::vector<dealii::types::global_dof_index> &soln_dof_indices,
313  dealii::Vector<real> &local_rhs_cell,
314  const dealii::FEValues<dim,dim> &/*fe_values_lagrange*/,
315  const bool compute_dRdW, const bool compute_dRdX, const bool compute_d2R);
316 
319  typename dealii::DoFHandler<dim>::active_cell_iterator cell,
320  const dealii::types::global_dof_index current_cell_index,
321  const unsigned int face_number,
322  const unsigned int boundary_id,
323  const dealii::FEFaceValuesBase<dim,dim> &fe_values_boundary,
324  const real penalty,
325  const dealii::FESystem<dim,dim> &fe,
326  const dealii::Quadrature<dim-1> &quadrature,
327  const std::vector<dealii::types::global_dof_index> &metric_dof_indices,
328  const std::vector<dealii::types::global_dof_index> &soln_dof_indices,
329  dealii::Vector<real> &local_rhs_cell,
330  const bool compute_dRdW, const bool compute_dRdX, const bool compute_d2R);
331 
333 
337  typename dealii::DoFHandler<dim>::active_cell_iterator cell,
338  const dealii::types::global_dof_index current_cell_index,
339  const dealii::types::global_dof_index neighbor_cell_index,
340  const std::pair<unsigned int, int> face_subface_int,
341  const std::pair<unsigned int, int> face_subface_ext,
342  const typename dealii::QProjector<dim>::DataSetDescriptor face_data_set_int,
343  const typename dealii::QProjector<dim>::DataSetDescriptor face_data_set_ext,
344  const dealii::FEFaceValuesBase<dim,dim> &,//fe_values_int,
345  const dealii::FEFaceValuesBase<dim,dim> &,//fe_values_ext,
346  const real penalty,
347  const dealii::FESystem<dim,dim> &fe_int,
348  const dealii::FESystem<dim,dim> &fe_ext,
349  const dealii::Quadrature<dim-1> &face_quadrature,
350  const std::vector<dealii::types::global_dof_index> &metric_dof_indices_int,
351  const std::vector<dealii::types::global_dof_index> &metric_dof_indices_ext,
352  const std::vector<dealii::types::global_dof_index> &soln_dof_indices_int,
353  const std::vector<dealii::types::global_dof_index> &soln_dof_indices_ext,
354  dealii::Vector<real> &local_rhs_int_cell,
355  dealii::Vector<real> &local_rhs_ext_cell,
356  const bool compute_dRdW, const bool compute_dRdX, const bool compute_d2R);
357 
360  typename dealii::DoFHandler<dim>::active_cell_iterator cell,
361  const dealii::types::global_dof_index current_cell_index,
362  const dealii::FEValues<dim,dim> &fe_values_volume,
363  const std::vector<dealii::types::global_dof_index> &current_dofs_indices,
364  const std::vector<dealii::types::global_dof_index> &metric_dof_indices,
365  const unsigned int poly_degree,
366  const unsigned int grid_degree,
367  dealii::Vector<real> &current_cell_rhs,
368  const dealii::FEValues<dim,dim> &fe_values_lagrange);
369 
370 
372 
373 }; // end of DGStrong class
374 
375 } // PHiLiP namespace
376 
377 #endif
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)
Evaluate the integral over the internal cell edges and its specified derivatives. ...
Definition: strong_dg.cpp:2916
void assemble_boundary_term_strong(const unsigned int iface, const dealii::types::global_dof_index current_cell_index, const unsigned int boundary_id, const unsigned int poly_degree, const real penalty, const std::vector< dealii::types::global_dof_index > &dof_indices, OPERATOR::basis_functions< dim, 2 *dim, real > &soln_basis, OPERATOR::basis_functions< dim, 2 *dim, real > &flux_basis, OPERATOR::vol_projection_operator< dim, 2 *dim, real > &soln_basis_projection_oper, OPERATOR::metric_operators< real, dim, 2 *dim > &metric_oper, dealii::Vector< real > &local_rhs_cell)
Strong form primary equation&#39;s boundary right-hand-side.
Definition: strong_dg.cpp:1312
void assemble_auxiliary_residual()
Assembles the auxiliary equations&#39; residuals and solves for the auxiliary variables.
Definition: strong_dg.cpp:398
void assemble_cell_auxiliary_residual(const DoFCellAccessorType1 &current_cell, const DoFCellAccessorType2 &current_metric_cell, std::vector< dealii::LinearAlgebra::distributed::Vector< double >> &rhs)
Assembles the auxiliary equations&#39; cell residuals.
DGStrong(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)
Constructor.
Definition: strong_dg.cpp:23
Files for the baseline physics.
Definition: ADTypes.hpp:10
void assemble_volume_term_strong(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 unsigned int poly_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, OPERATOR::metric_operators< real, dim, 2 *dim > &metric_oper, dealii::Vector< real > &local_rhs_int_cell)
Strong form primary equation&#39;s volume right-hand-side.
Definition: strong_dg.cpp:879
void assemble_face_term_strong(const unsigned int iface, const unsigned int neighbor_iface, const dealii::types::global_dof_index current_cell_index, const dealii::types::global_dof_index neighbor_cell_index, const unsigned int poly_degree_int, const unsigned int poly_degree_ext, const real penalty, const std::vector< dealii::types::global_dof_index > &dof_indices_int, const std::vector< dealii::types::global_dof_index > &dof_indices_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::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, dealii::Vector< real > &local_rhs_int_cell, dealii::Vector< real > &local_rhs_ext_cell)
Strong form primary equation&#39;s facet right-hand-side.
Definition: strong_dg.cpp:1808
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)
Evaluate the integral over the cell volume and the specified derivatives.
Definition: strong_dg.cpp:2789
Main parameter class that contains the various other sub-parameter classes.
void assemble_face_term_auxiliary_equation(const unsigned int iface, const unsigned int neighbor_iface, const dealii::types::global_dof_index current_cell_index, const dealii::types::global_dof_index neighbor_cell_index, const unsigned int poly_degree_int, const unsigned int poly_degree_ext, const std::vector< dealii::types::global_dof_index > &dof_indices_int, const std::vector< dealii::types::global_dof_index > &dof_indices_ext, OPERATOR::basis_functions< dim, 2 *dim, real > &soln_basis_int, OPERATOR::basis_functions< dim, 2 *dim, real > &soln_basis_ext, OPERATOR::metric_operators< real, dim, 2 *dim > &metric_oper_int, std::vector< dealii::Tensor< 1, dim, real >> &local_auxiliary_RHS_int, std::vector< dealii::Tensor< 1, dim, real >> &local_auxiliary_RHS_ext)
Evaluate the facet RHS for the auxiliary equation.
Definition: strong_dg.cpp:721
DGStrong class templated on the number of state variables.
Definition: strong_dg.hpp:16
Base metric operators class that stores functions used in both the volume and on surface.
Definition: operators.h:1086
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)
Assemble boundary term derivatives.
Definition: strong_dg.cpp:2637
The mapping shape functions evaluated at the desired nodes (facet set included in volume grid nodes f...
Definition: operators.h:1026
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, const real penalty, const std::vector< dealii::types::global_dof_index > &current_dofs_indices, const std::vector< dealii::types::global_dof_index > &neighbor_dofs_indices, const std::vector< dealii::types::global_dof_index > &current_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 > &, dealii::Vector< real > &current_cell_rhs, dealii::Vector< real > &neighbor_cell_rhs, std::vector< dealii::Tensor< 1, dim, real >> &current_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)
Builds the necessary operators and assembles subface residual.
Definition: strong_dg.cpp:313
void assemble_volume_term_auxiliary_equation(const std::vector< dealii::types::global_dof_index > &current_dofs_indices, const unsigned int poly_degree, OPERATOR::basis_functions< dim, 2 *dim, real > &soln_basis, OPERATOR::basis_functions< dim, 2 *dim, real > &flux_basis, OPERATOR::metric_operators< real, dim, 2 *dim > &metric_oper, std::vector< dealii::Tensor< 1, dim, real >> &local_auxiliary_RHS)
Evaluate the volume RHS for the auxiliary equation.
Definition: strong_dg.cpp:501
typename DGBase< dim, real, MeshType >::Triangulation Triangulation
Alias to base class Triangulation.
void assemble_boundary_term_auxiliary_equation(const unsigned int iface, const dealii::types::global_dof_index current_cell_index, const unsigned int poly_degree, const unsigned int boundary_id, const std::vector< dealii::types::global_dof_index > &dofs_indices, OPERATOR::basis_functions< dim, 2 *dim, real > &soln_basis, OPERATOR::metric_operators< real, dim, 2 *dim > &metric_oper, std::vector< dealii::Tensor< 1, dim, real >> &local_auxiliary_RHS)
Evaluate the boundary RHS for the auxiliary equation.
Definition: strong_dg.cpp:575
Abstract class templated on the number of state variables.
void assemble_face_term_and_build_operators(typename dealii::DoFHandler< dim >::active_cell_iterator, typename dealii::DoFHandler< dim >::active_cell_iterator, 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 > &current_dofs_indices, const std::vector< dealii::types::global_dof_index > &neighbor_dofs_indices, const std::vector< dealii::types::global_dof_index > &current_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 > &, dealii::hp::FEFaceValues< dim, dim > &, dealii::Vector< real > &current_cell_rhs, dealii::Vector< real > &neighbor_cell_rhs, std::vector< dealii::Tensor< 1, dim, real >> &current_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)
Builds the necessary operators and assembles face residual.
Definition: strong_dg.cpp:186
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 > &, dealii::hp::FEValues< dim, dim > &, const dealii::FESystem< dim, dim > &, 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, const bool, const bool)
Builds the necessary operators and assembles volume residual for either primary or auxiliary...
Definition: strong_dg.cpp:38
MeshType Triangulation
Definition: dg_base.hpp:89
void assemble_boundary_term_and_build_operators(typename dealii::DoFHandler< dim >::active_cell_iterator, 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 > &, const dealii::FESystem< dim, dim > &, 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, const bool, const bool)
Builds the necessary operators and assembles boundary residual for either primary or auxiliary...
Definition: strong_dg.cpp:125
DGBase is independent of the number of state variables.
Definition: dg_base.hpp:82
void allocate_dual_vector()
Allocate the dual vector for optimization.
Definition: strong_dg.cpp:3136
Projection operator corresponding to basis functions onto M-norm (L2).
Definition: operators.h:694
Local stiffness matrix without jacobian dependence.
Definition: operators.h:472
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 > &current_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 > &current_cell_rhs, const dealii::FEValues< dim, dim > &fe_values_lagrange)
Evaluate the integral over the cell volume.
Definition: strong_dg.cpp:3120