libMesh
heatsystem.h
Go to the documentation of this file.
1 // The libMesh Finite Element Library.
2 // Copyright (C) 2002-2017 Benjamin S. Kirk, John W. Peterson, Roy H. Stogner
3 
4 // This library is free software; you can redistribute it and/or
5 // modify it under the terms of the GNU Lesser General Public
6 // License as published by the Free Software Foundation; either
7 // version 2.1 of the License, or (at your option) any later version.
8 
9 // This library is distributed in the hope that it will be useful,
10 // but WITHOUT ANY WARRANTY; without even the implied warranty of
11 // MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
12 // Lesser General Public License for more details.
13 
14 // You should have received a copy of the GNU Lesser General Public
15 // License along with this library; if not, write to the Free Software
16 // Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA
17 
18 
19 
20 #include "libmesh/enum_fe_family.h"
21 #include "libmesh/fem_system.h"
22 #include "libmesh/parameter_pointer.h"
23 #include "libmesh/parameter_vector.h"
24 
25 using namespace libMesh;
26 
27 // FEMSystem, TimeSolver and NewtonSolver will handle most tasks,
28 // but we must specify element residuals
29 class HeatSystem : public FEMSystem
30 {
31 public:
32  // Constructor
34  const std::string & name_in,
35  const unsigned int number_in)
36  : FEMSystem(es, name_in, number_in),
37  _k(1.0),
38  _fe_family("LAGRANGE"),
39  _fe_order(1),
40  _analytic_jacobians(true),
41  R_plus_dp(0.0),
42  R_minus_dp(0.0),
43  dp(1.e-6)
44  { qoi.resize(1); }
45 
46  std::string & fe_family() { return _fe_family; }
47  unsigned int & fe_order() { return _fe_order; }
48  Real & k() { return _k; }
49  bool & analytic_jacobians() { return _analytic_jacobians; }
50 
51  // A function to compute and accumulate residuals
52  void perturb_accumulate_residuals(ParameterVector & parameters);
53 
54  // Sensitivity Calculation
56  {
57  final_sensitivity = -(R_plus_dp - R_minus_dp)/(2*dp);
58 
59  return final_sensitivity;
60  }
61 
62  void set_tf(Real val)
63  {
64  tf = val;
65  }
66 
68  {
69  if (!parameter_vector.size())
70  for (std::size_t i = 0; i != parameters.size(); ++i)
71  parameter_vector.push_back
73  (new ParameterPointer<Number>(&parameters[i])));
74 
75  return parameter_vector;
76  }
77 
78  Number & get_QoI_value(unsigned int QoI_index)
79  {
80  return computed_QoI[QoI_index];
81  }
82 
83 protected:
84  // System initialization
85  virtual void init_data ();
86 
87  // Context initialization
88  virtual void init_context (DiffContext & context);
89 
90  // Element residual and jacobian calculations
91  // Time dependent parts
92  virtual bool element_time_derivative (bool request_jacobian,
93  DiffContext & context);
94 
95  // Constraint parts
96  // virtual bool side_constraint (bool request_jacobian,
97  // DiffContext & context);
98 
99  // RHS for adjoint problem
100  virtual void element_qoi_derivative (DiffContext & context,
101  const QoISet & /* qois */);
102 
103  //virtual void element_qoi (DiffContext & context, const QoISet & qois);
104 
105  // Parameters associated with the system
106  std::vector<Number> parameters;
107 
108  // The ParameterVector object that will contain pointers to
109  // the system parameters
111 
112  // The parameters to solve for
114 
115  // The final time parameter
117 
118  // Variables to hold the computed QoIs
119  Number computed_QoI[1];
120 
121  // The FE type to use
122  std::string _fe_family;
123  unsigned int _fe_order;
124 
125  // Index for T variable
126  unsigned int T_var;
127 
128  // Calculate Jacobians analytically or not?
130 
131  // Variables to hold the perturbed residuals
134 
135  // Perturbation parameter
137 
138  // The final computed sensitivity
140 };
This is the EquationSystems class.
Number & compute_final_sensitivity()
Definition: heatsystem.h:55
Real & k()
Definition: heatsystem.h:48
void push_back(UniquePtr< ParameterAccessor< Number >> new_accessor)
Adds an additional parameter accessor to the end of the vector.
This class provides all data required for a physics package (e.g.
Definition: diff_context.h:54
Number R_plus_dp
Definition: heatsystem.h:132
Data structure for specifying which Parameters should be independent variables in a parameter sensiti...
Accessor object allowing reading and modification of the independent variables in a parameter sensiti...
Data structure for specifying which Quantities of Interest should be calculated in an adjoint or a pa...
Definition: qoi_set.h:45
Number R_minus_dp
Definition: heatsystem.h:133
std::vector< Number > parameters
Definition: heatsystem.h:106
The libMesh namespace provides an interface to certain functionality in the library.
void set_tf(Real val)
Definition: heatsystem.h:62
This class provides a specific system class.
Definition: fem_system.h:53
unsigned int T_var
Definition: heatsystem.h:126
std::unique_ptr< T > UniquePtr
Definition: auto_ptr.h:46
ParameterVector & get_parameter_vector()
Definition: heatsystem.h:67
bool _analytic_jacobians
Definition: heatsystem.h:129
std::size_t size() const
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real
std::string _fe_family
Definition: heatsystem.h:122
Number final_sensitivity
Definition: heatsystem.h:139
unsigned int _fe_order
Definition: heatsystem.h:123
std::string & fe_family()
Definition: heatsystem.h:46
bool & analytic_jacobians()
Definition: heatsystem.h:49
Number & get_QoI_value(unsigned int QoI_index)
Definition: heatsystem.h:78
HeatSystem(EquationSystems &es, const std::string &name_in, const unsigned int number_in)
Definition: heatsystem.h:33
Accessor object allowing reading and modification of the independent variables in a parameter sensiti...
unsigned int & fe_order()
Definition: heatsystem.h:47
ParameterVector parameter_vector
Definition: heatsystem.h:110