libMesh
fem_system.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 #ifndef LIBMESH_FEM_SYSTEM_H
21 #define LIBMESH_FEM_SYSTEM_H
22 
23 // Local Includes
24 #include "libmesh/diff_system.h"
25 #include "libmesh/fem_physics.h"
26 
27 // C++ includes
28 #include <cstddef>
29 
30 namespace libMesh
31 {
32 
33 // Forward Declarations
34 class DiffContext;
35 class FEMContext;
36 
37 
54  public FEMPhysics
55 {
56 public:
57 
63  const std::string & name,
64  const unsigned int number);
65 
69  virtual ~FEMSystem ();
70 
75 
80 
86  virtual void assembly (bool get_residual,
87  bool get_jacobian,
88  bool apply_heterogeneous_constraints = false,
89  bool apply_no_constraints = false) libmesh_override;
90 
99  virtual void solve () libmesh_override;
100 
105  void mesh_position_get();
106 
111  void mesh_position_set();
112 
121  virtual UniquePtr<DiffContext> build_context() libmesh_override;
122 
123  /*
124  * Prepares the result of a build_context() call for use.
125  *
126  * Most FEMSystem-based problems will need to reimplement this in order to
127  * call FE::get_*() as their particular physics requires.
128  */
129  virtual void init_context(DiffContext &) libmesh_override;
130 
135  virtual void postprocess () libmesh_override;
136 
145  virtual void assemble_qoi (const QoISet & indices = QoISet()) libmesh_override;
146 
154  virtual void assemble_qoi_derivative (const QoISet & qoi_indices = QoISet(),
155  bool include_liftfunc = true,
156  bool apply_constraints = true) libmesh_override;
157 
165 
176 
186  Real numerical_jacobian_h_for_var(unsigned int var_num) const;
187 
188  void set_numerical_jacobian_h_for_var(unsigned int var_num, Real new_h);
189 
204 
208  typedef bool (TimeSolver::*TimeSolverResPtr)(bool, DiffContext &);
209 
215  FEMContext & context) const;
216 
222  void numerical_elem_jacobian (FEMContext & context) const;
223 
229  void numerical_side_jacobian (FEMContext & context) const;
230 
236  void numerical_nonlocal_jacobian (FEMContext & context) const;
237 
238 protected:
243  virtual void init_data () libmesh_override;
244 
245 private:
247 };
248 
249 // --------------------------------------------------------------
250 // FEMSystem inline methods
251 inline
252 Real
253 FEMSystem::numerical_jacobian_h_for_var(unsigned int var_num) const
254 {
255  if ((var_num >= _numerical_jacobian_h_for_var.size()) ||
256  _numerical_jacobian_h_for_var[var_num] == Real(0))
257  return numerical_jacobian_h;
258 
259  return _numerical_jacobian_h_for_var[var_num];
260 }
261 
262 inline
264  Real new_h)
265 {
266  if (_numerical_jacobian_h_for_var.size() <= var_num)
267  _numerical_jacobian_h_for_var.resize(var_num+1,Real(0));
268 
269  libmesh_assert_greater(new_h, 0);
270 
271  _numerical_jacobian_h_for_var[var_num] = new_h;
272 }
273 
274 } // namespace libMesh
275 
276 
277 #endif // LIBMESH_FEM_SYSTEM_H
This is the EquationSystems class.
void numerical_side_jacobian(FEMContext &context) const
Uses the results of multiple side_residual() calls to numerically differentiate the corresponding jac...
Definition: fem_system.C:1305
Real verify_analytic_jacobians
If verify_analytic_jacobian is equal to zero (as it is by default), no numeric jacobians will be calc...
Definition: fem_system.h:203
This class provides all data required for a physics package (e.g.
Definition: diff_context.h:54
Real numerical_jacobian_h_for_var(unsigned int var_num) const
If numerical_jacobian_h_for_var(var_num) is changed from its default value (numerical_jacobian_h), the FEMSystem will perturb solution vector entries for variable var_num by that amount when calculating finite differences with respect to that variable.
Definition: fem_system.h:253
FEMSystem sys_type
The type of system.
Definition: fem_system.h:74
This is a generic class that defines a solver to handle time integration of DifferentiableSystems.
Definition: time_solver.h:58
This class provides a specific system class.
Definition: fem_physics.h:44
Data structure for specifying which Quantities of Interest should be calculated in an adjoint or a pa...
Definition: qoi_set.h:45
Real numerical_jacobian_h
If calculating numeric jacobians is required, the FEMSystem will perturb each solution vector entry b...
Definition: fem_system.h:175
bool(TimeSolver::* TimeSolverResPtr)(bool, DiffContext &)
Syntax sugar to make numerical_jacobian() declaration easier.
Definition: fem_system.h:208
The libMesh namespace provides an interface to certain functionality in the library.
This class provides a specific system class.
Definition: fem_system.h:53
void set_numerical_jacobian_h_for_var(unsigned int var_num, Real new_h)
Definition: fem_system.h:263
const std::string & name() const
Definition: system.h:1998
virtual UniquePtr< DiffContext > build_context() libmesh_override
Builds a FEMContext object with enough information to do evaluations on each element.
Definition: fem_system.C:1321
std::unique_ptr< T > UniquePtr
Definition: auto_ptr.h:46
void numerical_nonlocal_jacobian(FEMContext &context) const
Uses the results of multiple side_residual() calls to numerically differentiate the corresponding jac...
Definition: fem_system.C:1313
This class provides a specific system class.
Definition: diff_system.h:53
virtual void solve() libmesh_override
Invokes the solver associated with the system.
Definition: fem_system.C:1056
virtual void assemble_qoi_derivative(const QoISet &qoi_indices=QoISet(), bool include_liftfunc=true, bool apply_constraints=true) libmesh_override
Runs a qoi derivative assembly loop over all elements, and if assemble_qoi_sides is true over all sid...
Definition: fem_system.C:1157
bool fe_reinit_during_postprocess
If fe_reinit_during_postprocess is true (it is true by default), FE objects will be reinit()ed with t...
Definition: fem_system.h:164
DifferentiableSystem Parent
The type of the parent.
Definition: fem_system.h:79
This class provides all data required for a physics package (e.g.
Definition: fem_context.h:61
FEMSystem(EquationSystems &es, const std::string &name, const unsigned int number)
Constructor.
Definition: fem_system.C:828
std::vector< Real > _numerical_jacobian_h_for_var
Definition: fem_system.h:246
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real
unsigned int number() const
Definition: system.h:2006
void mesh_position_set()
Tells the FEMSystem to set the mesh nodal coordinates which should correspond to degree of freedom co...
Definition: fem_system.C:1067
void numerical_elem_jacobian(FEMContext &context) const
Uses the results of multiple element_residual() calls to numerically differentiate the corresponding ...
Definition: fem_system.C:1297
virtual void assembly(bool get_residual, bool get_jacobian, bool apply_heterogeneous_constraints=false, bool apply_no_constraints=false) libmesh_override
Prepares matrix or rhs for matrix assembly.
Definition: fem_system.C:852
void numerical_jacobian(TimeSolverResPtr res, FEMContext &context) const
Uses the results of multiple res calls to numerically differentiate the corresponding jacobian...
Definition: fem_system.C:1188
virtual ~FEMSystem()
Destructor.
Definition: fem_system.C:839
virtual void init_data() libmesh_override
Initializes the member data fields associated with the system, so that, e.g., assemble() may be used...
Definition: fem_system.C:845
void mesh_position_get()
Tells the FEMSystem to set the degree of freedom coefficients which should correspond to mesh nodal c...
Definition: fem_system.C:1392
virtual void assemble_qoi(const QoISet &indices=QoISet()) libmesh_override
Runs a qoi assembly loop over all elements, and if assemble_qoi_sides is true over all sides...
Definition: fem_system.C:1127
virtual void postprocess() libmesh_override
Runs a postprocessing loop over all elements, and if postprocess_sides is true over all sides...
Definition: fem_system.C:1107
virtual void init_context(DiffContext &) libmesh_override
Definition: fem_system.C:1345