1 #ifndef __WEAK_DISCONTINUOUSGALERKIN_H__ 2 #define __WEAK_DISCONTINUOUSGALERKIN_H__ 4 #include "dg_base_state.hpp" 5 #include "solution/local_solution.hpp" 12 #if PHILIP_DIM==1 // dealii::parallel::distributed::Triangulation<dim> does not work for 1D 13 template <
int dim,
int nstate,
typename real,
typename MeshType = dealii::Triangulation<dim>>
15 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);
35 typename dealii::DoFHandler<dim>::active_cell_iterator cell,
36 const dealii::types::global_dof_index current_cell_index,
37 const std::vector<dealii::types::global_dof_index> &cell_dofs_indices,
38 const std::vector<dealii::types::global_dof_index> &metric_dof_indices,
39 const unsigned int poly_degree,
40 const unsigned int grid_degree,
48 std::array<std::vector<real>,dim> &,
49 dealii::hp::FEValues<dim,dim> &fe_values_collection_volume,
50 dealii::hp::FEValues<dim,dim> &fe_values_collection_volume_lagrange,
51 const dealii::FESystem<dim,dim> ¤t_fe_ref,
52 dealii::Vector<real> &local_rhs_int_cell,
53 std::vector<dealii::Tensor<1,dim,real>> &,
55 const bool compute_dRdW,
const bool compute_dRdX,
const bool compute_d2R);
59 typename dealii::DoFHandler<dim>::active_cell_iterator cell,
60 const dealii::types::global_dof_index current_cell_index,
61 const unsigned int iface,
62 const unsigned int boundary_id,
64 const std::vector<dealii::types::global_dof_index> &cell_dofs_indices,
65 const std::vector<dealii::types::global_dof_index> &metric_dof_indices,
66 const unsigned int poly_degree,
67 const unsigned int grid_degree,
75 std::array<std::vector<real>,dim> &,
76 dealii::hp::FEFaceValues<dim,dim> &fe_values_collection_face_int,
77 const dealii::FESystem<dim,dim> ¤t_fe_ref,
78 dealii::Vector<real> &local_rhs_int_cell,
79 std::vector<dealii::Tensor<1,dim,real>> &,
81 const bool compute_dRdW,
const bool compute_dRdX,
const bool compute_d2R);
85 typename dealii::DoFHandler<dim>::active_cell_iterator cell,
86 typename dealii::DoFHandler<dim>::active_cell_iterator neighbor_cell,
87 const dealii::types::global_dof_index current_cell_index,
88 const dealii::types::global_dof_index neighbor_cell_index,
89 const unsigned int iface,
90 const unsigned int neighbor_iface,
92 const std::vector<dealii::types::global_dof_index> ¤t_dofs_indices,
93 const std::vector<dealii::types::global_dof_index> &neighbor_dofs_indices,
94 const std::vector<dealii::types::global_dof_index> ¤t_metric_dofs_indices,
95 const std::vector<dealii::types::global_dof_index> &neighbor_metric_dofs_indices,
110 std::array<std::vector<real>,dim> &,
111 dealii::hp::FEFaceValues<dim,dim> &fe_values_collection_face_int,
112 dealii::hp::FEFaceValues<dim,dim> &fe_values_collection_face_ext,
113 dealii::Vector<real> ¤t_cell_rhs,
114 dealii::Vector<real> &neighbor_cell_rhs,
115 std::vector<dealii::Tensor<1,dim,real>> &,
116 dealii::LinearAlgebra::distributed::Vector<double> &rhs,
117 std::array<dealii::LinearAlgebra::distributed::Vector<double>,dim> &,
119 const bool compute_dRdW,
const bool compute_dRdX,
const bool compute_d2R);
123 typename dealii::DoFHandler<dim>::active_cell_iterator cell,
124 typename dealii::DoFHandler<dim>::active_cell_iterator neighbor_cell,
125 const dealii::types::global_dof_index current_cell_index,
126 const dealii::types::global_dof_index neighbor_cell_index,
127 const unsigned int iface,
128 const unsigned int neighbor_iface,
129 const unsigned int neighbor_i_subface,
131 const std::vector<dealii::types::global_dof_index> ¤t_dofs_indices,
132 const std::vector<dealii::types::global_dof_index> &neighbor_dofs_indices,
133 const std::vector<dealii::types::global_dof_index> ¤t_metric_dofs_indices,
134 const std::vector<dealii::types::global_dof_index> &neighbor_metric_dofs_indices,
149 std::array<std::vector<real>,dim> &,
150 dealii::hp::FEFaceValues<dim,dim> &fe_values_collection_face_int,
151 dealii::hp::FESubfaceValues<dim,dim> &fe_values_collection_subface,
152 dealii::Vector<real> ¤t_cell_rhs,
153 dealii::Vector<real> &neighbor_cell_rhs,
154 std::vector<dealii::Tensor<1,dim,real>> &,
155 dealii::LinearAlgebra::distributed::Vector<double> &rhs,
156 std::array<dealii::LinearAlgebra::distributed::Vector<double>,dim> &,
158 const bool compute_dRdW,
const bool compute_dRdX,
const bool compute_d2R);
169 template <
typename real2>
171 typename dealii::DoFHandler<dim>::active_cell_iterator cell,
172 const dealii::types::global_dof_index current_cell_index,
175 const std::vector<real> &local_dual,
176 const dealii::Quadrature<dim> &quadrature,
178 std::vector<real2> &rhs, real2 &dual_dot_residual,
179 const bool compute_metric_derivatives,
180 const dealii::FEValues<dim,dim> &fe_values_vol);
185 template <
typename real2>
187 typename dealii::DoFHandler<dim>::active_cell_iterator cell,
188 const dealii::types::global_dof_index current_cell_index,
191 const std::vector< real > &local_dual,
192 const unsigned int face_number,
193 const unsigned int boundary_id,
197 const dealii::FEFaceValuesBase<dim,dim> &fe_values_boundary,
199 const dealii::Quadrature<dim-1> &quadrature,
200 std::vector<real2> &rhs,
201 real2 &dual_dot_residual,
202 const bool compute_metric_derivatives);
207 template <
typename real2>
209 typename dealii::DoFHandler<dim>::active_cell_iterator cell,
210 const dealii::types::global_dof_index current_cell_index,
211 const dealii::types::global_dof_index neighbor_cell_index,
216 const std::vector< double > &dual_int,
217 const std::vector< double > &dual_ext,
218 const std::pair<unsigned int, int> face_subface_int,
219 const std::pair<unsigned int, int> face_subface_ext,
220 const typename dealii::QProjector<dim>::DataSetDescriptor face_data_set_int,
221 const typename dealii::QProjector<dim>::DataSetDescriptor face_data_set_ext,
225 const dealii::FEFaceValuesBase<dim,dim> &fe_values_int,
226 const dealii::FEFaceValuesBase<dim,dim> &fe_values_ext,
228 const dealii::Quadrature<dim-1> &face_quadrature,
229 std::vector<real2> &rhs_int,
230 std::vector<real2> &rhs_ext,
231 real2 &dual_dot_residual,
232 const bool compute_dRdW,
const bool compute_dRdX,
const bool compute_d2R);
240 template <
typename real2>
242 typename dealii::DoFHandler<dim>::active_cell_iterator cell,
243 const dealii::types::global_dof_index current_cell_index,
244 const dealii::FEValues<dim,dim> &fe_values_vol,
245 const dealii::FESystem<dim,dim> &fe_soln,
246 const dealii::Quadrature<dim> &quadrature,
247 const std::vector<dealii::types::global_dof_index> &metric_dof_indices,
248 const std::vector<dealii::types::global_dof_index> &soln_dof_indices,
249 dealii::Vector<real> &local_rhs_cell,
250 const dealii::FEValues<dim,dim> &fe_values_lagrange,
252 const bool compute_dRdW,
const bool compute_dRdX,
const bool compute_d2R);
258 template <
typename adtype>
260 typename dealii::DoFHandler<dim>::active_cell_iterator cell,
261 const dealii::types::global_dof_index current_cell_index,
262 const unsigned int face_number,
263 const unsigned int boundary_id,
264 const dealii::FEFaceValuesBase<dim,dim> &fe_values_boundary,
266 const dealii::FESystem<dim,dim> &fe_soln,
267 const dealii::Quadrature<dim-1> &quadrature,
268 const std::vector<dealii::types::global_dof_index> &metric_dof_indices,
269 const std::vector<dealii::types::global_dof_index> &soln_dof_indices,
273 dealii::Vector<real> &local_rhs_cell,
274 const bool compute_dRdW,
const bool compute_dRdX,
const bool compute_d2R);
282 template <
typename adtype>
284 typename dealii::DoFHandler<dim>::active_cell_iterator cell,
285 const dealii::types::global_dof_index current_cell_index,
286 const dealii::types::global_dof_index neighbor_cell_index,
287 const std::pair<unsigned int, int> face_subface_int,
288 const std::pair<unsigned int, int> face_subface_ext,
289 const typename dealii::QProjector<dim>::DataSetDescriptor face_data_set_int,
290 const typename dealii::QProjector<dim>::DataSetDescriptor face_data_set_ext,
291 const dealii::FEFaceValuesBase<dim,dim> &,
292 const dealii::FEFaceValuesBase<dim,dim> &,
294 const dealii::FESystem<dim,dim> &fe_int,
295 const dealii::FESystem<dim,dim> &fe_ext,
296 const dealii::Quadrature<dim-1> &face_quadrature,
297 const std::vector<dealii::types::global_dof_index> &metric_dof_indices_int,
298 const std::vector<dealii::types::global_dof_index> &metric_dof_indices_ext,
299 const std::vector<dealii::types::global_dof_index> &soln_dof_indices_int,
300 const std::vector<dealii::types::global_dof_index> &soln_dof_indices_ext,
304 dealii::Vector<real> &local_rhs_int_cell,
305 dealii::Vector<real> &local_rhs_ext_cell,
306 const bool compute_dRdW,
const bool compute_dRdX,
const bool compute_d2R);
314 typename dealii::DoFHandler<dim>::active_cell_iterator cell,
315 const dealii::types::global_dof_index current_cell_index,
316 const dealii::FEValues<dim,dim> &,
317 const dealii::FESystem<dim,dim> &fe,
318 const dealii::Quadrature<dim> &quadrature,
319 const std::vector<dealii::types::global_dof_index> &metric_dof_indices,
320 const std::vector<dealii::types::global_dof_index> &soln_dof_indices,
321 dealii::Vector<real> &local_rhs_cell,
322 const dealii::FEValues<dim,dim> &,
323 const bool compute_dRdW,
const bool compute_dRdX,
const bool compute_d2R);
329 typename dealii::DoFHandler<dim>::active_cell_iterator cell,
330 const dealii::types::global_dof_index current_cell_index,
331 const unsigned int face_number,
332 const unsigned int boundary_id,
333 const dealii::FEFaceValuesBase<dim,dim> &fe_values_boundary,
335 const dealii::FESystem<dim,dim> &fe,
336 const dealii::Quadrature<dim-1> &quadrature,
337 const std::vector<dealii::types::global_dof_index> &metric_dof_indices,
338 const std::vector<dealii::types::global_dof_index> &soln_dof_indices,
339 dealii::Vector<real> &local_rhs_cell,
340 const bool compute_dRdW,
const bool compute_dRdX,
const bool compute_d2R);
348 typename dealii::DoFHandler<dim>::active_cell_iterator cell,
349 const dealii::types::global_dof_index current_cell_index,
350 const dealii::types::global_dof_index neighbor_cell_index,
351 const std::pair<unsigned int, int> face_subface_int,
352 const std::pair<unsigned int, int> face_subface_ext,
353 const typename dealii::QProjector<dim>::DataSetDescriptor face_data_set_int,
354 const typename dealii::QProjector<dim>::DataSetDescriptor face_data_set_ext,
355 const dealii::FEFaceValuesBase<dim,dim> &,
356 const dealii::FEFaceValuesBase<dim,dim> &,
358 const dealii::FESystem<dim,dim> &fe_int,
359 const dealii::FESystem<dim,dim> &fe_ext,
360 const dealii::Quadrature<dim-1> &face_quadrature,
361 const std::vector<dealii::types::global_dof_index> &metric_dof_indices_int,
362 const std::vector<dealii::types::global_dof_index> &metric_dof_indices_ext,
363 const std::vector<dealii::types::global_dof_index> &soln_dof_indices_int,
364 const std::vector<dealii::types::global_dof_index> &soln_dof_indices_ext,
365 dealii::Vector<real> &local_rhs_int_cell,
366 dealii::Vector<real> &local_rhs_ext_cell,
367 const bool compute_dRdW,
const bool compute_dRdX,
const bool compute_d2R);
373 typename dealii::DoFHandler<dim>::active_cell_iterator cell,
374 const dealii::types::global_dof_index current_cell_index,
375 const dealii::FEValues<dim,dim> &fe_values_vol,
376 const dealii::FESystem<dim,dim> &fe_soln,
377 const dealii::Quadrature<dim> &quadrature,
378 const std::vector<dealii::types::global_dof_index> &metric_dof_indices,
379 const std::vector<dealii::types::global_dof_index> &soln_dof_indices,
380 dealii::Vector<real> &local_rhs_cell,
381 const dealii::FEValues<dim,dim> &fe_values_lagrange,
383 const bool compute_dRdW,
const bool compute_dRdX,
const bool compute_d2R);
388 typename dealii::DoFHandler<dim>::active_cell_iterator cell,
389 const dealii::types::global_dof_index current_cell_index,
390 const unsigned int face_number,
391 const unsigned int boundary_id,
392 const dealii::FEFaceValuesBase<dim,dim> &fe_values_boundary,
394 const dealii::FESystem<dim,dim> &fe_soln,
395 const dealii::Quadrature<dim-1> &quadrature,
396 const std::vector<dealii::types::global_dof_index> &metric_dof_indices,
397 const std::vector<dealii::types::global_dof_index> &soln_dof_indices,
401 dealii::Vector<real> &local_rhs_cell,
402 const bool compute_dRdW,
const bool compute_dRdX,
const bool compute_d2R);
407 typename dealii::DoFHandler<dim>::active_cell_iterator cell,
408 const dealii::types::global_dof_index current_cell_index,
409 const dealii::types::global_dof_index neighbor_cell_index,
410 const std::pair<unsigned int, int> face_subface_int,
411 const std::pair<unsigned int, int> face_subface_ext,
412 const typename dealii::QProjector<dim>::DataSetDescriptor face_data_set_int,
413 const typename dealii::QProjector<dim>::DataSetDescriptor face_data_set_ext,
414 const dealii::FEFaceValuesBase<dim,dim> &,
415 const dealii::FEFaceValuesBase<dim,dim> &,
417 const dealii::FESystem<dim,dim> &fe_int,
418 const dealii::FESystem<dim,dim> &fe_ext,
419 const dealii::Quadrature<dim-1> &face_quadrature,
420 const std::vector<dealii::types::global_dof_index> &metric_dof_indices_int,
421 const std::vector<dealii::types::global_dof_index> &metric_dof_indices_ext,
422 const std::vector<dealii::types::global_dof_index> &soln_dof_indices_int,
423 const std::vector<dealii::types::global_dof_index> &soln_dof_indices_ext,
427 dealii::Vector<real> &local_rhs_int_cell,
428 dealii::Vector<real> &local_rhs_ext_cell,
429 const bool compute_dRdW,
const bool compute_dRdX,
const bool compute_d2R);
433 typename dealii::DoFHandler<dim>::active_cell_iterator cell,
434 const dealii::types::global_dof_index current_cell_index,
435 const dealii::FEValues<dim,dim> &fe_values_volume,
436 const std::vector<dealii::types::global_dof_index> ¤t_dofs_indices,
437 const std::vector<dealii::types::global_dof_index> &metric_dof_indices,
438 const unsigned int poly_degree,
439 const unsigned int grid_degree,
440 dealii::Vector<real> ¤t_cell_rhs,
441 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. ...
Base class from which Advection, Diffusion, ConvectionDiffusion, and Euler is derived.
Class to store local solution coefficients and provide evaluation functions.
Base class of numerical flux associated with dissipation.
Files for the baseline physics.
void assemble_boundary_term(typename dealii::DoFHandler< dim >::active_cell_iterator cell, const dealii::types::global_dof_index current_cell_index, const LocalSolution< real2, dim, nstate > &local_solution, const LocalSolution< real2, dim, dim > &local_metric, const std::vector< real > &local_dual, const unsigned int face_number, const unsigned int boundary_id, const Physics::PhysicsBase< dim, nstate, real2 > &physics, const NumericalFlux::NumericalFluxConvective< dim, nstate, real2 > &conv_num_flux, const NumericalFlux::NumericalFluxDissipative< dim, nstate, real2 > &diss_num_flux, const dealii::FEFaceValuesBase< dim, dim > &fe_values_boundary, const real penalty, const dealii::Quadrature< dim-1 > &quadrature, std::vector< real2 > &rhs, real2 &dual_dot_residual, const bool compute_metric_derivatives)
Main function responsible for evaluating the boundary integral and the specified derivatives.
DGWeak class templated on the number of state variables.
Main parameter class that contains the various other sub-parameter classes.
void assemble_auxiliary_residual()
Assembles the auxiliary equations' residuals and solves for the auxiliary variables.
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)
Evaluate the integral over the cell volume and the specified derivatives.
void allocate_dual_vector()
Allocate the dual vector for optimization.
void assemble_face_residual(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, const Physics::PhysicsBase< dim, nstate, real > &physics, const NumericalFlux::NumericalFluxConvective< dim, nstate, real > &conv_num_flux, const NumericalFlux::NumericalFluxDissipative< dim, nstate, real > &diss_num_flux, 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 face.
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 > &, OPERATOR::basis_functions< dim, 2 *dim, real > &, OPERATOR::local_basis_stiffness< dim, 2 *dim, real > &, OPERATOR::vol_projection_operator< dim, 2 *dim, real > &, OPERATOR::vol_projection_operator< dim, 2 *dim, real > &, OPERATOR::metric_operators< real, dim, 2 *dim > &, OPERATOR::mapping_shape_functions< dim, 2 *dim, real > &, std::array< std::vector< real >, dim > &, 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 >> &, const bool, const bool compute_dRdW, const bool compute_dRdX, const bool compute_d2R)
Builds the necessary fe values and assembles volume residual.
Base metric operators class that stores functions used in both the volume and on surface.
void assemble_volume_term(typename dealii::DoFHandler< dim >::active_cell_iterator cell, const dealii::types::global_dof_index current_cell_index, const LocalSolution< real2, dim, nstate > &local_solution, const LocalSolution< real2, dim, dim > &local_metric, const std::vector< real > &local_dual, const dealii::Quadrature< dim > &quadrature, const Physics::PhysicsBase< dim, nstate, real2 > &physics, std::vector< real2 > &rhs, real2 &dual_dot_residual, const bool compute_metric_derivatives, const dealii::FEValues< dim, dim > &fe_values_vol)
Main function responsible for evaluating the integral over the cell volume and the specified derivati...
void assemble_face_codi_taped_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, const Physics::PhysicsBase< dim, nstate, adtype > &physics, const NumericalFlux::NumericalFluxConvective< dim, nstate, adtype > &conv_num_flux, const NumericalFlux::NumericalFluxDissipative< dim, nstate, adtype > &diss_num_flux, 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)
Preparation of CoDiPack taping for internal cell faces integrals, and derivative evaluation.
The mapping shape functions evaluated at the desired nodes (facet set included in volume grid nodes f...
Base class of numerical flux associated with convection.
void assemble_boundary_codi_taped_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_soln, 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, const Physics::PhysicsBase< dim, nstate, adtype > &physics, const NumericalFlux::NumericalFluxConvective< dim, nstate, adtype > &conv_num_flux, const NumericalFlux::NumericalFluxDissipative< dim, nstate, adtype > &diss_num_flux, dealii::Vector< real > &local_rhs_cell, const bool compute_dRdW, const bool compute_dRdX, const bool compute_d2R)
Preparation of CoDiPack taping for boundary integral, and derivative evaluation.
void assemble_volume_residual(typename dealii::DoFHandler< dim >::active_cell_iterator cell, const dealii::types::global_dof_index current_cell_index, const dealii::FEValues< dim, dim > &fe_values_vol, const dealii::FESystem< dim, dim > &fe_soln, 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 > &fe_values_lagrange, const Physics::PhysicsBase< dim, nstate, real > &physics, const bool compute_dRdW, const bool compute_dRdX, const bool compute_d2R)
Evaluate the integral over the cell volume.
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, const unsigned int, const unsigned int, const unsigned int, OPERATOR::basis_functions< dim, 2 *dim, real > &, OPERATOR::basis_functions< dim, 2 *dim, real > &, OPERATOR::basis_functions< dim, 2 *dim, real > &, OPERATOR::basis_functions< dim, 2 *dim, real > &, OPERATOR::local_basis_stiffness< dim, 2 *dim, real > &, OPERATOR::vol_projection_operator< dim, 2 *dim, real > &, OPERATOR::vol_projection_operator< dim, 2 *dim, real > &, OPERATOR::metric_operators< real, dim, 2 *dim > &, OPERATOR::metric_operators< real, dim, 2 *dim > &, OPERATOR::mapping_shape_functions< dim, 2 *dim, real > &, std::array< std::vector< real >, dim > &, 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 >> &, dealii::LinearAlgebra::distributed::Vector< double > &rhs, std::array< dealii::LinearAlgebra::distributed::Vector< double >, dim > &, const bool, const bool compute_dRdW, const bool compute_dRdX, const bool compute_d2R)
Builds the necessary fe values and assembles subface residual.
void assemble_face_term(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 LocalSolution< real2, dim, nstate > &soln_int, const LocalSolution< real2, dim, nstate > &soln_ext, const LocalSolution< real2, dim, dim > &metric_int, const LocalSolution< real2, dim, dim > &metric_ext, const std::vector< double > &dual_int, const std::vector< double > &dual_ext, 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 Physics::PhysicsBase< dim, nstate, real2 > &physics, const NumericalFlux::NumericalFluxConvective< dim, nstate, real2 > &conv_num_flux, const NumericalFlux::NumericalFluxDissipative< dim, nstate, real2 > &diss_num_flux, const dealii::FEFaceValuesBase< dim, dim > &fe_values_int, const dealii::FEFaceValuesBase< dim, dim > &fe_values_ext, const real penalty, const dealii::Quadrature< dim-1 > &face_quadrature, std::vector< real2 > &rhs_int, std::vector< real2 > &rhs_ext, real2 &dual_dot_residual, const bool compute_dRdW, const bool compute_dRdX, const bool compute_d2R)
Main function responsible for evaluating the internal face integral and the specified derivatives...
DGWeak(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.
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, const unsigned int, const unsigned int, const unsigned int, OPERATOR::basis_functions< dim, 2 *dim, real > &, OPERATOR::basis_functions< dim, 2 *dim, real > &, OPERATOR::basis_functions< dim, 2 *dim, real > &, OPERATOR::basis_functions< dim, 2 *dim, real > &, OPERATOR::local_basis_stiffness< dim, 2 *dim, real > &, OPERATOR::vol_projection_operator< dim, 2 *dim, real > &, OPERATOR::vol_projection_operator< dim, 2 *dim, real > &, OPERATOR::metric_operators< real, dim, 2 *dim > &, OPERATOR::metric_operators< real, dim, 2 *dim > &, OPERATOR::mapping_shape_functions< dim, 2 *dim, real > &, std::array< std::vector< real >, dim > &, 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 >> &, dealii::LinearAlgebra::distributed::Vector< double > &rhs, std::array< dealii::LinearAlgebra::distributed::Vector< double >, dim > &, const bool, const bool compute_dRdW, const bool compute_dRdX, const bool compute_d2R)
Builds the necessary fe values and assembles face residual.
typename DGBase< dim, real, MeshType >::Triangulation Triangulation
Alias to base class Triangulation.
void assemble_boundary_residual(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_soln, 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, const Physics::PhysicsBase< dim, nstate, real > &physics, const NumericalFlux::NumericalFluxConvective< dim, nstate, real > &conv_num_flux, const NumericalFlux::NumericalFluxDissipative< dim, nstate, real > &diss_num_flux, dealii::Vector< real > &local_rhs_cell, const bool compute_dRdW, const bool compute_dRdX, const bool compute_d2R)
Evaluate the integral over the boundary.
Abstract class templated on the number of state variables.
void assemble_volume_codi_taped_derivatives(typename dealii::DoFHandler< dim >::active_cell_iterator cell, const dealii::types::global_dof_index current_cell_index, const dealii::FEValues< dim, dim > &fe_values_vol, const dealii::FESystem< dim, dim > &fe_soln, 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 > &fe_values_lagrange, const Physics::PhysicsBase< dim, nstate, real2 > &physics, const bool compute_dRdW, const bool compute_dRdX, const bool compute_d2R)
Preparation of CoDiPack taping for volume integral, and derivative evaluation.
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.
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 > &, OPERATOR::basis_functions< dim, 2 *dim, real > &, OPERATOR::local_basis_stiffness< dim, 2 *dim, real > &, OPERATOR::vol_projection_operator< dim, 2 *dim, real > &, OPERATOR::vol_projection_operator< dim, 2 *dim, real > &, OPERATOR::metric_operators< real, dim, 2 *dim > &, OPERATOR::mapping_shape_functions< dim, 2 *dim, real > &, std::array< std::vector< real >, dim > &, 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 >> &, const bool, const bool compute_dRdW, const bool compute_dRdX, const bool compute_d2R)
Builds the necessary fe values and assembles boundary residual.
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)
Evaluate the integral over the cell edges that are on domain boundaries and the specified derivatives...
DGBase is independent of the number of state variables.
Projection operator corresponding to basis functions onto M-norm (L2).
Local stiffness matrix without jacobian dependence.