libMesh
examples
adjoints
adjoints_ex6
poisson.h
Go to the documentation of this file.
1
#include "libmesh/enum_fe_family.h"
2
#include "libmesh/fem_system.h"
3
#include "libmesh/qoi_set.h"
4
#include "libmesh/system.h"
5
6
// Bring in everything from the libMesh namespace
7
using namespace
libMesh
;
8
9
// FEMSystem, TimeSolver and NewtonSolver will handle most tasks,
10
// but we must specify element residuals
11
class
PoissonSystem
:
public
FEMSystem
12
{
13
public
:
14
PoissonSystem
(
EquationSystems
& es,
15
const
std::string & name_in,
16
const
unsigned
int
number_in)
17
:
FEMSystem
(es, name_in, number_in),
18
_fe_family(
"LAGRANGE"
), _fe_order(2),
19
_analytic_jacobians(true) { this->init_qois(1); computed_QoI[0] = 0.0; }
20
21
std::string &
fe_family
() {
return
_fe_family; }
22
unsigned
int
&
fe_order
() {
return
_fe_order; }
23
bool
&
analytic_jacobians
() {
return
_analytic_jacobians; }
24
25
// Postprocessing function which we are going to override for this application
26
27
virtual
void
postprocess(
void
);
28
29
Number
&
get_QoI_value
(std::string type,
unsigned
int
QoI_index)
30
{
31
if
(type ==
"exact"
)
32
{
33
return
exact_QoI[QoI_index];
34
}
35
else
36
{
37
return
computed_QoI[QoI_index];
38
}
39
}
40
41
protected
:
42
// System initialization
43
virtual
void
init_data ();
44
45
// Context initialization
46
virtual
void
init_context (
DiffContext
& context);
47
48
// Element residual and jacobian calculations
49
// Time dependent parts
50
virtual
bool
element_time_derivative (
bool
request_jacobian,
51
DiffContext
& context);
52
53
// Overloading the postprocess function
54
55
virtual
void
element_postprocess(
DiffContext
& context);
56
57
// Forcing function coefficient
58
Real
alpha
;
59
60
// Indices for each variable
61
unsigned
int
T_var
;
62
63
// Variables to hold the computed QoIs
64
65
Number
computed_QoI[1];
66
67
// Variables to read in the exact QoIs from poisson.in
68
69
Number
exact_QoI[1];
70
71
// The FE type to use
72
std::string
_fe_family
;
73
unsigned
int
_fe_order
;
74
75
// Calculate Jacobians analytically or not?
76
bool
_analytic_jacobians
;
77
};
PoissonSystem::T_var
unsigned int T_var
Definition:
poisson.h:61
libMesh::EquationSystems
This is the EquationSystems class.
Definition:
equation_systems.h:68
libMesh::DiffContext
This class provides all data required for a physics package (e.g.
Definition:
diff_context.h:55
libMesh
The libMesh namespace provides an interface to certain functionality in the library.
libMesh::FEMSystem
This class provides a specific system class.
Definition:
fem_system.h:53
PoissonSystem::PoissonSystem
PoissonSystem(EquationSystems &es, const std::string &name_in, const unsigned int number_in)
Definition:
poisson.h:14
PoissonSystem::_fe_order
unsigned int _fe_order
Definition:
poisson.h:73
PoissonSystem
Definition:
poisson.h:11
PoissonSystem::_fe_family
std::string _fe_family
Definition:
poisson.h:72
PoissonSystem::fe_family
std::string & fe_family()
Definition:
poisson.h:21
PoissonSystem::fe_order
unsigned int & fe_order()
Definition:
poisson.h:22
libMesh::Real
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real
Definition:
libmesh_common.h:143
PoissonSystem::_analytic_jacobians
bool _analytic_jacobians
Definition:
poisson.h:76
PoissonSystem::analytic_jacobians
bool & analytic_jacobians()
Definition:
poisson.h:23
libMesh::Number
Real Number
Definition:
libmesh_common.h:214
PoissonSystem::alpha
Real alpha
Definition:
poisson.h:58
PoissonSystem::get_QoI_value
Number & get_QoI_value(std::string type, unsigned int QoI_index)
Definition:
poisson.h:29
Generated on Wed Mar 27 2024 17:39:10 for libMesh by
1.8.14