libMesh
examples
adjoints
adjoints_ex4
L-shaped.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
LaplaceSystem
:
public
FEMSystem
12
{
13
public
:
14
// Constructor
15
LaplaceSystem
(
EquationSystems
& es,
16
const
std::string & name_in,
17
const
unsigned
int
number_in) :
18
FEMSystem
(es, name_in, number_in),
19
_fe_family(
"LAGRANGE"
),
20
_fe_order(1),
21
_analytic_jacobians(true)
22
{ this->init_qois(2); }
23
24
std::string &
fe_family
() {
return
_fe_family; }
25
unsigned
int
&
fe_order
() {
return
_fe_order; }
26
bool
&
analytic_jacobians
() {
return
_analytic_jacobians; }
27
28
// Postprocessing function which we are going to override for this application
29
30
virtual
void
postprocess();
31
32
Number
&
get_QoI_value
(std::string type,
unsigned
int
QoI_index)
33
{
34
if
(type ==
"exact"
)
35
return
exact_QoI[QoI_index];
36
else
37
return
computed_QoI[QoI_index];
38
}
39
40
protected
:
41
// System initialization
42
virtual
void
init_data ();
43
44
// Context initialization
45
virtual
void
init_context (
DiffContext
& context);
46
47
// Element residual and jacobian calculations
48
// Time dependent parts
49
virtual
bool
element_time_derivative (
bool
request_jacobian,
50
DiffContext
& context);
51
52
// Constraint parts
53
virtual
bool
side_constraint (
bool
request_jacobian,
54
DiffContext
& context);
55
56
// Overloading the postprocess function
57
58
virtual
void
element_postprocess(
DiffContext
& context);
59
60
virtual
void
side_postprocess(
DiffContext
& context);
61
62
// Overloading the qoi function on elements
63
64
virtual
void
element_qoi_derivative (
DiffContext
& context,
65
const
QoISet
& qois);
66
67
// Overloading the qoi function on sides
68
69
virtual
void
side_qoi_derivative (
DiffContext
& context,
70
const
QoISet
& qois);
71
72
Number
exact_solution
(
const
Point
&);
73
74
// Variables to hold the computed QoIs
75
76
Number
computed_QoI[2];
77
78
// Variables to read in the exact QoIs from l-shaped.in
79
80
Number
exact_QoI[2];
81
82
// The FE type to use
83
std::string _fe_family;
84
unsigned
int
_fe_order;
85
86
// Calculate Jacobians analytically or not?
87
bool
_analytic_jacobians;
88
};
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
exact_solution
Number(* exact_solution)(const Point &p, const Parameters &, const std::string &, const std::string &)
Definition:
adaptivity_ex4.C:137
libMesh::QoISet
Data structure for specifying which Quantities of Interest should be calculated in an adjoint or a pa...
Definition:
qoi_set.h:45
LaplaceSystem::get_QoI_value
Number & get_QoI_value(std::string type, unsigned int QoI_index)
Definition:
L-shaped.h:32
LaplaceSystem::fe_family
std::string & fe_family()
Definition:
L-shaped.h:24
LaplaceSystem
Definition:
L-shaped.h:11
LaplaceSystem::LaplaceSystem
LaplaceSystem(EquationSystems &es, const std::string &name_in, const unsigned int number_in)
Definition:
L-shaped.h:15
LaplaceSystem::analytic_jacobians
bool & analytic_jacobians()
Definition:
L-shaped.h:26
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
LaplaceSystem::fe_order
unsigned int & fe_order()
Definition:
L-shaped.h:25
libMesh::Number
Real Number
Definition:
libmesh_common.h:214
libMesh::Point
A Point defines a location in LIBMESH_DIM dimensional Real space.
Definition:
point.h:39
Generated on Thu Apr 18 2024 15:17:48 for libMesh by
1.8.14