[P]arallel [Hi]gh-order [Li]brary for [P]DEs  Latest
Parallel High-Order Library for PDEs through hp-adaptive Discontinuous Galerkin methods
adjoint.h
1 #ifndef __ADJOINT_H__
2 #define __ADJOINT_H__
3 
4 /* includes */
5 #include <deal.II/distributed/solution_transfer.h>
6 #include <deal.II/fe/fe_q.h>
7 #include <deal.II/fe/fe_values.h>
8 #include <deal.II/lac/la_parallel_vector.h>
9 
10 #include <iostream>
11 #include <vector>
12 
13 #include "dg/dg_base.hpp"
14 #include "functional.h"
15 #include "parameters/all_parameters.h"
16 #include "physics/physics.h"
17 
18 namespace PHiLiP {
19 
21 
34  #if PHILIP_DIM==1 // dealii::parallel::distributed::Triangulation<dim> does not work for 1D
35 template <int dim, int nstate, typename real, typename MeshType = dealii::Triangulation<dim>>
36 #else
37 template <int dim, int nstate, typename real, typename MeshType = dealii::parallel::distributed::Triangulation<dim>>
38 #endif
39 class Adjoint
40 {
41 public:
42 
46  fine,
47  };
48 
50 
54  Adjoint(
55  std::shared_ptr< DGBase<dim,real,MeshType> > _dg,
56  std::shared_ptr< Functional<dim, nstate, real, MeshType> > _functional,
57  std::shared_ptr< Physics::PhysicsBase<dim,nstate,Sacado::Fad::DFad<real>> > _physics);
58 
60 
63  void reinit();
64  // to reinitialize with other pointers, just create a new class
65 
67 
71 
73 
76  void coarse_to_fine();
77 
79 
82  void fine_to_coarse();
83 
85 
92  dealii::LinearAlgebra::distributed::Vector<real> fine_grid_adjoint();
93 
95 
101  dealii::LinearAlgebra::distributed::Vector<real> coarse_grid_adjoint();
102 
104 
111  dealii::Vector<real> dual_weighted_residual();
112 
114 
118  void output_results_vtk(const unsigned int cycle);
119 
121  std::shared_ptr< DGBase<dim,real,MeshType> > dg;
123  std::shared_ptr< Functional<dim, nstate, real, MeshType> > functional;
125  std::shared_ptr< Physics::PhysicsBase<dim,nstate,Sacado::Fad::DFad<real>> > physics;
126 
128  std::shared_ptr<MeshType> triangulation;
130  dealii::LinearAlgebra::distributed::Vector<real> solution_coarse;
132  dealii::LinearAlgebra::distributed::Vector<real> dIdw_fine;
134  dealii::LinearAlgebra::distributed::Vector<real> dIdw_coarse;
136  dealii::LinearAlgebra::distributed::Vector<real> adjoint_fine;
138  dealii::LinearAlgebra::distributed::Vector<real> adjoint_coarse;
140 
142  dealii::Vector<real> dual_weighted_residual_fine;
143 
145  dealii::Vector<real> coarse_fe_index;
146 
149 
150 protected:
151  MPI_Comm mpi_communicator;
152  dealii::ConditionalOStream pcout;
153 
154 }; // Adjoint class
155 
156 
157 } // PHiLiP namespace
158 
159 #endif // __ADJOINT_H__
dealii::LinearAlgebra::distributed::Vector< real > dIdw_coarse
functional derivative (on the coarse grid)
Definition: adjoint.h:134
std::shared_ptr< Functional< dim, nstate, real, MeshType > > functional
Functional class pointer.
Definition: adjoint.h:123
void convert_to_state(AdjointStateEnum state)
Converts the adjoint to specified state.
Definition: adjoint.cpp:78
AdjointStateEnum adjoint_state
Current adjoint state.
Definition: adjoint.h:148
dealii::LinearAlgebra::distributed::Vector< real > fine_grid_adjoint()
Computes the fine grid adjoint.
Definition: adjoint.cpp:157
void fine_to_coarse()
Return the problem to the original solution and polynomial distribution.
Definition: adjoint.cpp:136
dealii::LinearAlgebra::distributed::Vector< real > adjoint_coarse
coarse grid adjoint ( )
Definition: adjoint.h:138
Base class from which Advection, Diffusion, ConvectionDiffusion, and Euler is derived.
Definition: physics.h:34
Adjoint class.
Definition: adjoint.h:39
std::shared_ptr< MeshType > triangulation
Grid.
Definition: adjoint.h:128
Adjoint(std::shared_ptr< DGBase< dim, real, MeshType > > _dg, std::shared_ptr< Functional< dim, nstate, real, MeshType > > _functional, std::shared_ptr< Physics::PhysicsBase< dim, nstate, Sacado::Fad::DFad< real >> > _physics)
Constructor.
Definition: adjoint.cpp:30
std::shared_ptr< Physics::PhysicsBase< dim, nstate, Sacado::Fad::DFad< real > > > physics
Problem physics (for calling the functional class)
Definition: adjoint.h:125
Files for the baseline physics.
Definition: ADTypes.hpp:10
AdjointStateEnum
For storing the current state in the adjoint.
Definition: adjoint.h:44
dealii::Vector< real > dual_weighted_residual()
compute the Dual Weighted Residual (DWR)
Definition: adjoint.cpp:216
dealii::LinearAlgebra::distributed::Vector< real > coarse_grid_adjoint()
Computes the coarse grid adjoint.
Definition: adjoint.cpp:187
Refined state.
Definition: adjoint.h:46
dealii::LinearAlgebra::distributed::Vector< real > dIdw_fine
functional derivative (on the fine grid)
Definition: adjoint.h:132
void reinit()
Reinitialize Adjoint with the same pointers.
Definition: adjoint.cpp:53
dealii::Vector< real > coarse_fe_index
Original FE_index distribution.
Definition: adjoint.h:145
MPI_Comm mpi_communicator
MPI communicator.
Definition: adjoint.h:151
dealii::LinearAlgebra::distributed::Vector< real > solution_coarse
original solution
Definition: adjoint.h:130
std::shared_ptr< DGBase< dim, real, MeshType > > dg
DG class pointer.
Definition: adjoint.h:121
void output_results_vtk(const unsigned int cycle)
Outputs the current solution and adjoint values.
Definition: adjoint.cpp:249
dealii::Vector< real > dual_weighted_residual_fine
dual weighted residual
Definition: adjoint.h:142
Initial state.
Definition: adjoint.h:45
void coarse_to_fine()
Projects the problem to a p-enriched space.
Definition: adjoint.cpp:93
Functional base class.
Definition: functional.h:43
DGBase is independent of the number of state variables.
Definition: dg_base.hpp:82
dealii::LinearAlgebra::distributed::Vector< real > adjoint_fine
fine grid adjoint ( )
Definition: adjoint.h:136
dealii::ConditionalOStream pcout
Parallel std::cout that only outputs on mpi_rank==0.
Definition: adjoint.h:152