1 #ifndef __ANISOTROPIC_MESH_ADAPTATION__ 2 #define __ANISOTROPIC_MESH_ADAPTATION__ 4 #include "dg/dg_base.hpp" 5 #include "functional/functional.h" 16 #if PHILIP_DIM==1 // dealii::parallel::distributed::Triangulation<dim> does not work for 1D 17 template <
int dim,
int nstate,
typename real,
typename MeshType = dealii::Triangulation<dim>>
19 template <
int dim,
int nstate,
typename real,
typename MeshType = dealii::parallel::distributed::Triangulation<dim>>
24 using VectorType = dealii::LinearAlgebra::distributed::Vector<double>;
31 const real _complexity,
32 const bool _use_goal_oriented_approach =
false);
75 std::vector<std::array<dealii::Tensor<1,dim,real>,nstate>> &flux_at_support_pts,
76 const dealii::FEValues<dim,dim> &fe_values_input,
77 const std::vector<dealii::types::global_dof_index> &dof_indices,
78 typename dealii::DoFHandler<dim>::active_cell_iterator cell)
const;
81 std::shared_ptr<DGBase<dim,real,MeshType>>
dg;
99 dealii::ConditionalOStream
pcout;
114 std::shared_ptr< Functional<dim, nstate, real, MeshType> >
functional;
std::shared_ptr< Functional< dim, nstate, real, MeshType > > functional
Functional to evaluate the adjoint for goal oriented anisotropic meshing.
int mpi_rank
Processor# of current processor.
void compute_feature_based_hessian()
Computes feature based hessian (i.e. hessian of the solution).
std::shared_ptr< DGBase< dim, real, MeshType > > dg
Shared pointer to DGBase.
void initialize_cellwise_metric_and_hessians()
Initializes cellwise metric and hessian to zero tensors.
dealii::ConditionalOStream pcout
std::cout only on processor #0.
std::shared_ptr< Physics::PhysicsBase< dim, nstate, real > > pde_physics_double
Contains the physics of the PDE with real type.
Files for the baseline physics.
void adapt_mesh()
Function which adapts mesh and loads in new mesh.
std::vector< dealii::Tensor< 2, dim, real > > cellwise_hessian
Stores hessian in each cell.
void get_flux_at_support_pts(std::vector< std::array< dealii::Tensor< 1, dim, real >, nstate >> &flux_at_support_pts, const dealii::FEValues< dim, dim > &fe_values_input, const std::vector< dealii::types::global_dof_index > &dof_indices, typename dealii::DoFHandler< dim >::active_cell_iterator cell) const
Returns flux coeffs by evaluating flux at support points of fe.
const real normLp
Lp Norm w.r.t. which the anlytical optimization is done.
dealii::LinearAlgebra::distributed::Vector< double > VectorType
Alias for dealii's parallel distributed vector.
unsigned int get_iquad_near_cellcenter(const dealii::Quadrature< dim > &volume_quadrature) const
Returns quadrature number of a point which is closest to the reference cell's center.
void compute_goal_oriented_hessian()
Computes pseudo Hessian for the goal oriented approach.
std::vector< dealii::Tensor< 2, dim, real > > cellwise_optimal_metric
Stores optimal metric in each cell.
const unsigned int max_dofs_per_cell
Stores max dofs per cell to initialize dof_indices.
void change_p_degree_and_interpolate_solution(const unsigned int poly_degree)
Change the polynomial order and interpolate solution.
int n_mpi
Total no. of processors.
AnisotropicMeshAdaptation(std::shared_ptr< DGBase< dim, real, MeshType > > dg_input, const real _norm_Lp, const real _complexity, const bool _use_goal_oriented_approach=false)
Constructor.
dealii::Tensor< 2, dim, real > get_positive_definite_tensor(const dealii::Tensor< 2, dim, real > &input_tensor) const
Returns positive tensor from an input tensor by taking absolute of eigenvalues.
void reconstruct_p2_solution()
Reconstructs p2 solution after interpolation.
const unsigned int initial_poly_degree
Stores initial polynomial degree.
const bool use_goal_oriented_approach
Flag to use goal oriented approach. It is set to false by default.
void compute_cellwise_optimal_metric()
Computes optimal metric depending on goal oriented or feature based approach.
void compute_abs_hessian()
Computes hessian using the input coefficients, which can be a solution sensor or (for goal oriented a...
DGBase is independent of the number of state variables.
MPI_Comm mpi_communicator
Alias for MPI_COMM_WORLD.
const real complexity
Analogue of number of vertices/elements in continuous mesh framework.