1 #ifndef __STRONG_DISCONTINUOUSGALERKIN_H__ 2 #define __STRONG_DISCONTINUOUSGALERKIN_H__ 4 #include "dg_base_state.hpp" 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>>
14 template <
int dim,
int nstate,
typename real,
typename MeshType = dealii::parallel::distributed::Triangulation<dim>>
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);
43 template<
typename DoFCellAccessorType1,
typename DoFCellAccessorType2>
45 const DoFCellAccessorType1 ¤t_cell,
46 const DoFCellAccessorType2 ¤t_metric_cell,
47 std::vector<dealii::LinearAlgebra::distributed::Vector<double>> &rhs);
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,
65 std::array<std::vector<real>,dim> &mapping_support_points,
66 dealii::hp::FEValues<dim,dim> &,
67 dealii::hp::FEValues<dim,dim> &,
68 const dealii::FESystem<dim,dim> &,
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 ,
const bool ,
const bool );
76 typename dealii::DoFHandler<dim>::active_cell_iterator ,
77 const dealii::types::global_dof_index current_cell_index,
78 const unsigned int iface,
79 const unsigned int boundary_id,
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,
92 std::array<std::vector<real>,dim> &mapping_support_points,
93 dealii::hp::FEFaceValues<dim,dim> &,
94 const dealii::FESystem<dim,dim> &,
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 ,
const bool ,
const bool );
102 typename dealii::DoFHandler<dim>::active_cell_iterator ,
103 typename dealii::DoFHandler<dim>::active_cell_iterator ,
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,
109 const std::vector<dealii::types::global_dof_index> ¤t_dofs_indices,
110 const std::vector<dealii::types::global_dof_index> &neighbor_dofs_indices,
111 const std::vector<dealii::types::global_dof_index> ¤t_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,
127 std::array<std::vector<real>,dim> &mapping_support_points,
128 dealii::hp::FEFaceValues<dim,dim> &,
129 dealii::hp::FEFaceValues<dim,dim> &,
130 dealii::Vector<real> ¤t_cell_rhs,
131 dealii::Vector<real> &neighbor_cell_rhs,
132 std::vector<dealii::Tensor<1,dim,real>> ¤t_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);
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,
150 const std::vector<dealii::types::global_dof_index> ¤t_dofs_indices,
151 const std::vector<dealii::types::global_dof_index> &neighbor_dofs_indices,
152 const std::vector<dealii::types::global_dof_index> ¤t_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,
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> &,
171 dealii::Vector<real> ¤t_cell_rhs,
172 dealii::Vector<real> &neighbor_cell_rhs,
173 std::vector<dealii::Tensor<1,dim,real>> ¤t_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);
186 const std::vector<dealii::types::global_dof_index> ¤t_dofs_indices,
187 const unsigned int poly_degree,
191 std::vector<dealii::Tensor<1,dim,real>> &local_auxiliary_RHS);
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);
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);
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,
256 dealii::Vector<real> &local_rhs_int_cell);
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,
265 const std::vector<dealii::types::global_dof_index> &dof_indices,
270 dealii::Vector<real> &local_rhs_cell);
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,
289 const std::vector<dealii::types::global_dof_index> &dof_indices_int,
290 const std::vector<dealii::types::global_dof_index> &dof_indices_ext,
299 dealii::Vector<real> &local_rhs_int_cell,
300 dealii::Vector<real> &local_rhs_ext_cell);
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> &,
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> &,
315 const bool compute_dRdW,
const bool compute_dRdX,
const bool compute_d2R);
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,
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);
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> &,
345 const dealii::FEFaceValuesBase<dim,dim> &,
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);
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> ¤t_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> ¤t_cell_rhs,
368 const dealii::FEValues<dim,dim> &fe_values_lagrange);
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. ...
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's boundary right-hand-side.
void assemble_auxiliary_residual()
Assembles the auxiliary equations' residuals and solves for the auxiliary variables.
void assemble_cell_auxiliary_residual(const DoFCellAccessorType1 ¤t_cell, const DoFCellAccessorType2 ¤t_metric_cell, std::vector< dealii::LinearAlgebra::distributed::Vector< double >> &rhs)
Assembles the auxiliary equations' 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.
Files for the baseline physics.
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's volume right-hand-side.
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's facet right-hand-side.
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.
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.
DGStrong class templated on the number of state variables.
Base metric operators class that stores functions used in both the volume and on surface.
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.
The mapping shape functions evaluated at the desired nodes (facet set included in volume grid nodes f...
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 > ¤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 > &, 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)
Builds the necessary operators and assembles subface residual.
void assemble_volume_term_auxiliary_equation(const std::vector< dealii::types::global_dof_index > ¤t_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.
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.
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 > ¤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 > &, dealii::hp::FEFaceValues< dim, dim > &, 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)
Builds the necessary operators and assembles face residual.
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...
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...
DGBase is independent of the number of state variables.
void allocate_dual_vector()
Allocate the dual vector for optimization.
Projection operator corresponding to basis functions onto M-norm (L2).
Local stiffness matrix without jacobian dependence.
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)
Evaluate the integral over the cell volume.