fem_system.h
Go to the documentation of this file.
1 // The libMesh Finite Element Library.
2 // Copyright (C) 2002-2018 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 
92  virtual void assembly (bool get_residual,
93  bool get_jacobian,
94  bool apply_heterogeneous_constraints = false,
95  bool apply_no_constraints = false) override;
96 
105  virtual void solve () override;
106 
111  void mesh_position_get();
112 
117  void mesh_position_set();
118 
127  virtual std::unique_ptr<DiffContext> build_context() override;
128 
129  /*
130  * Prepares the result of a build_context() call for use.
131  *
132  * Most FEMSystem-based problems will need to reimplement this in order to
133  * call FE::get_*() as their particular physics requires.
134  */
135  virtual void init_context(DiffContext &) override;
136 
141  virtual void postprocess () override;
142 
151  virtual void assemble_qoi (const QoISet & indices = QoISet()) override;
152 
160  virtual void assemble_qoi_derivative (const QoISet & qoi_indices = QoISet(),
161  bool include_liftfunc = true,
162  bool apply_constraints = true) override;
163 
171 
182 
192  Real numerical_jacobian_h_for_var(unsigned int var_num) const;
193 
194  void set_numerical_jacobian_h_for_var(unsigned int var_num, Real new_h);
195 
210 
214  typedef bool (TimeSolver::*TimeSolverResPtr)(bool, DiffContext &);
215 
221  FEMContext & context) const;
222 
228  void numerical_elem_jacobian (FEMContext & context) const;
229 
235  void numerical_side_jacobian (FEMContext & context) const;
236 
242  void numerical_nonlocal_jacobian (FEMContext & context) const;
243 
244 protected:
249  virtual void init_data () override;
250 
251 private:
252  std::vector<Real> _numerical_jacobian_h_for_var;
253 };
254 
255 // --------------------------------------------------------------
256 // FEMSystem inline methods
257 inline
258 Real
259 FEMSystem::numerical_jacobian_h_for_var(unsigned int var_num) const
260 {
261  if ((var_num >= _numerical_jacobian_h_for_var.size()) ||
262  _numerical_jacobian_h_for_var[var_num] == Real(0))
263  return numerical_jacobian_h;
264 
265  return _numerical_jacobian_h_for_var[var_num];
266 }
267 
268 inline
270  Real new_h)
271 {
272  if (_numerical_jacobian_h_for_var.size() <= var_num)
273  _numerical_jacobian_h_for_var.resize(var_num+1,Real(0));
274 
275  libmesh_assert_greater(new_h, 0);
276 
277  _numerical_jacobian_h_for_var[var_num] = new_h;
278 }
279 
280 } // namespace libMesh
281 
282 
283 #endif // LIBMESH_FEM_SYSTEM_H
void numerical_elem_jacobian(FEMContext &context) const
Definition: fem_system.C:1294
Manages multiples systems of equations.
Real verify_analytic_jacobians
Definition: fem_system.h:209
FEMSystem sys_type
Definition: fem_system.h:74
Used to specify quantities of interest in a simulation.
Definition: qoi_set.h:45
virtual void assemble_qoi(const QoISet &indices=QoISet()) override
Definition: fem_system.C:1121
void numerical_jacobian(TimeSolverResPtr res, FEMContext &context) const
Definition: fem_system.C:1182
virtual ~FEMSystem()
Definition: fem_system.C:841
void numerical_side_jacobian(FEMContext &context) const
Definition: fem_system.C:1302
bool(TimeSolver::* TimeSolverResPtr)(bool, DiffContext &)
Definition: fem_system.h:214
void mesh_position_set()
Definition: fem_system.C:1061
virtual void postprocess() override
Definition: fem_system.C:1101
unsigned int number() const
Definition: system.h:2025
virtual void init_data() override
Definition: fem_system.C:847
bool fe_reinit_during_postprocess
Definition: fem_system.h:170
DifferentiableSystem Parent
Definition: fem_system.h:79
void numerical_nonlocal_jacobian(FEMContext &context) const
Definition: fem_system.C:1310
virtual void solve() override
Definition: fem_system.C:1050
std::vector< Real > _numerical_jacobian_h_for_var
Definition: fem_system.h:252
void set_numerical_jacobian_h_for_var(unsigned int var_num, Real new_h)
Definition: fem_system.h:269
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real
virtual void init_context(DiffContext &) override
Definition: fem_system.C:1342
virtual void assemble_qoi_derivative(const QoISet &qoi_indices=QoISet(), bool include_liftfunc=true, bool apply_constraints=true) override
Definition: fem_system.C:1151
virtual std::unique_ptr< DiffContext > build_context() override
Definition: fem_system.C:1318
const std::string & name() const
Definition: system.h:2017
FEMSystem(EquationSystems &es, const std::string &name, const unsigned int number)
Definition: fem_system.C:830
void mesh_position_get()
Definition: fem_system.C:1389
virtual void assembly(bool get_residual, bool get_jacobian, bool apply_heterogeneous_constraints=false, bool apply_no_constraints=false) override
Definition: fem_system.C:854
Real numerical_jacobian_h_for_var(unsigned int var_num) const
Definition: fem_system.h:259