libMesh
diff_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_DIFF_SYSTEM_H
21 #define LIBMESH_DIFF_SYSTEM_H
22 
23 // Local Includes
24 #include "libmesh/auto_ptr.h"
25 #include "libmesh/diff_context.h"
26 #include "libmesh/diff_physics.h"
27 #include "libmesh/diff_qoi.h"
28 #include "libmesh/implicit_system.h"
29 #include "libmesh/time_solver.h"
30 
31 // C++ includes
32 
33 namespace libMesh
34 {
35 
36 // Forward Declarations
37 class TimeSolver;
38 
39 template <typename T> class NumericVector;
40 
54  public virtual DifferentiablePhysics,
55  public virtual DifferentiableQoI
56 {
57 public:
58 
64  const std::string & name,
65  const unsigned int number);
66 
70  virtual ~DifferentiableSystem ();
71 
76 
81 
86  virtual void clear () libmesh_override;
87 
92  virtual void reinit () libmesh_override;
93 
98  virtual void assemble () libmesh_override;
99 
104  virtual LinearSolver<Number> * get_linear_solver() const libmesh_override;
105 
111  virtual std::pair<unsigned int, Real>
112  get_linear_solve_parameters() const libmesh_override;
113 
118  virtual void release_linear_solver(LinearSolver<Number> *) const libmesh_override;
119 
124  virtual void assembly (bool get_residual,
125  bool get_jacobian,
126  bool apply_heterogeneous_constraints = false,
127  bool apply_no_constraints = false) libmesh_override = 0;
128 
134  virtual void solve () libmesh_override;
135 
140  virtual std::pair<unsigned int, Real>
141  adjoint_solve (const QoISet & qoi_indices = QoISet()) libmesh_override;
142 
146  virtual UniquePtr<DifferentiablePhysics> clone_physics() libmesh_override
147  {
148  libmesh_not_implemented();
149  // dummy
151  }
152 
156  virtual UniquePtr<DifferentiableQoI> clone() libmesh_override
157  {
158  libmesh_not_implemented();
159  // dummy
160  return UniquePtr<DifferentiableQoI>(this);
161  }
162 
170  { return this->_diff_physics; }
171 
179  { return this->_diff_physics; }
180 
185  { this->_diff_physics = (physics_in->clone_physics()).release();
186  this->_diff_physics->init_physics(*this);}
187 
192 
198  const DifferentiableQoI * get_qoi() const
199  { return this->diff_qoi; }
200 
207  { return this->diff_qoi; }
208 
212  void attach_qoi( DifferentiableQoI * qoi_in )
213  { this->diff_qoi = (qoi_in->clone()).release();
214  // User needs to resize qoi system qoi accordingly
215  this->diff_qoi->init_qoi( this->qoi );}
216 
222 
230  {
231  time_solver.reset(_time_solver.release());
232  }
233 
238 
242  const TimeSolver & get_time_solver() const;
243 
249 
259 
264  virtual void postprocess () {}
265 
269  virtual void element_postprocess (DiffContext &) {}
270 
275  virtual void side_postprocess (DiffContext &) {}
276 
286  unsigned int get_second_order_dot_var( unsigned int var ) const;
287 
291  bool have_first_order_scalar_vars() const;
292 
296  bool have_second_order_scalar_vars() const;
297 
303 
309 
315 
320 
325 
330 
335 
340 
345 
350 
351 protected:
352 
359 
366 
378  virtual void init_data () libmesh_override;
379 
386 
395  void add_dot_var_dirichlet_bcs( unsigned int var_idx, unsigned int dot_var_idx);
396 
397 };
398 
399 // --------------------------------------------------------------
400 // DifferentiableSystem inline methods
401 inline
402 TimeSolver & DifferentiableSystem::get_time_solver()
403 {
404  libmesh_assert(time_solver.get());
405  libmesh_assert_equal_to (&(time_solver->system()), this);
406  return *time_solver;
407 }
408 
409 inline
411 {
412  libmesh_assert(time_solver.get());
413  libmesh_assert_equal_to (&(time_solver->system()), this);
414  return *time_solver;
415 }
416 
417 } // namespace libMesh
418 
419 
420 #endif // LIBMESH_DIFF_SYSTEM_H
void attach_physics(DifferentiablePhysics *physics_in)
Attach external Physics object.
Definition: diff_system.h:184
UniquePtr< TimeSolver > time_solver
A pointer to the solver object we&#39;re going to use.
Definition: diff_system.h:221
This is the EquationSystems class.
This class provides all data required for a physics package (e.g.
Definition: diff_context.h:54
const DifferentiableQoI * get_qoi() const
Definition: diff_system.h:198
virtual void assemble() libmesh_override
Prepares matrix and rhs for matrix assembly.
Definition: diff_system.C:145
This is a generic class that defines a solver to handle time integration of DifferentiableSystems.
Definition: time_solver.h:58
Data structure for specifying which Quantities of Interest should be calculated in an adjoint or a pa...
Definition: qoi_set.h:45
bool print_jacobian_norms
Set print_jacobian_norms to true to print |J| whenever it is assembled.
Definition: diff_system.h:329
virtual std::pair< unsigned int, Real > adjoint_solve(const QoISet &qoi_indices=QoISet()) libmesh_override
This function sets the _is_adjoint boolean member of TimeSolver to true and then calls the adjoint_so...
Definition: diff_system.C:164
virtual void init_data() libmesh_override
Initializes the member data fields associated with the system, so that, e.g., assemble() may be used...
Definition: diff_system.C:108
virtual UniquePtr< DiffContext > build_context()
Builds a DiffContext object with enough information to do evaluations on each element.
Definition: diff_system.C:137
DifferentiableQoI * get_qoi()
Definition: diff_system.h:206
void add_second_order_dot_vars()
Helper function to add "velocity" variables that are cousins to second order-in-time variables in the...
Definition: diff_system.C:198
bool postprocess_sides
If postprocess_sides is true (it is false by default), the postprocessing loop will loop over all sid...
Definition: diff_system.h:302
DifferentiablePhysics * get_physics()
Definition: diff_system.h:178
void add_dot_var_dirichlet_bcs(unsigned int var_idx, unsigned int dot_var_idx)
Helper function to and Dirichlet boundary conditions to "dot" variable cousins of second order variab...
Definition: diff_system.C:228
The libMesh namespace provides an interface to certain functionality in the library.
virtual LinearSolver< Number > * get_linear_solver() const libmesh_override
Definition: diff_system.C:175
bool print_jacobians
Set print_jacobians to true to print J whenever it is assembled.
Definition: diff_system.h:334
virtual void clear() libmesh_override
Clear all the data structures associated with the system.
Definition: diff_system.C:69
bool have_second_order_scalar_vars() const
Check for any second order vars that are also belong to FEFamily::SCALAR.
Definition: diff_system.C:348
const std::string & name() const
Definition: system.h:1998
virtual void solve() libmesh_override
Invokes the solver associated with the system.
Definition: diff_system.C:152
virtual UniquePtr< DifferentiableQoI > clone() libmesh_override
We don&#39;t allow systems to be attached to each other.
Definition: diff_system.h:156
bool print_element_residuals
Set print_element_residuals to true to print each R_elem contribution.
Definition: diff_system.h:344
libmesh_assert(j)
virtual void assembly(bool get_residual, bool get_jacobian, bool apply_heterogeneous_constraints=false, bool apply_no_constraints=false) libmesh_override=0
Assembles a residual in rhs and/or a jacobian in matrix, as requested.
DifferentiableQoI * diff_qoi
Pointer to object to use for quantity of interest assembly evaluations.
Definition: diff_system.h:365
std::unique_ptr< T > UniquePtr
Definition: auto_ptr.h:46
virtual void init_physics(const System &sys)
Initialize any data structures associated with the physics.
Definition: diff_physics.C:40
bool print_element_solutions
Set print_element_solutions to true to print each U_elem input.
Definition: diff_system.h:339
This class provides a specific system class.
Definition: diff_system.h:53
unsigned int get_second_order_dot_var(unsigned int var) const
For a given second order (in time) variable var, this method will return the index to the correspondi...
Definition: diff_system.C:313
This base class can be inherited from to provide interfaces to linear solvers from different packages...
Definition: linear_solver.h:58
bool print_solution_norms
Set print_residual_norms to true to print |U| whenever it is used in an assembly() call...
Definition: diff_system.h:308
void swap_physics(DifferentiablePhysics *&swap_physics)
Swap current physics object with external object.
Definition: diff_system.C:368
DifferentiablePhysics * _diff_physics
Pointer to object to use for physics assembly evaluations.
Definition: diff_system.h:358
virtual UniquePtr< DifferentiableQoI > clone()=0
Copy of this object.
bool print_residual_norms
Set print_residual_norms to true to print |F| whenever it is assembled.
Definition: diff_system.h:319
virtual void side_postprocess(DiffContext &)
Does any work that needs to be done on side of elem in a postprocessing loop.
Definition: diff_system.h:275
std::vector< Number > qoi
Values of the quantities of interest.
Definition: system.h:1553
void set_time_solver(UniquePtr< TimeSolver > _time_solver)
Sets the time_solver FIXME: This code is a little dangerous as it transfers ownership from the TimeSo...
Definition: diff_system.h:229
Real deltat
For time-dependent problems, this is the amount delta t to advance the solution in time...
Definition: diff_system.h:248
virtual void release_linear_solver(LinearSolver< Number > *) const libmesh_override
Releases a pointer to a linear solver acquired by this->get_linear_solver()
Definition: diff_system.C:194
bool have_first_order_scalar_vars() const
Check for any first order vars that are also belong to FEFamily::SCALAR.
Definition: diff_system.C:330
bool print_residuals
Set print_residuals to true to print F whenever it is assembled.
Definition: diff_system.h:324
virtual std::pair< unsigned int, Real > get_linear_solve_parameters() const libmesh_override
Definition: diff_system.C:184
DifferentiableSystem sys_type
The type of system.
Definition: diff_system.h:75
virtual UniquePtr< DifferentiablePhysics > clone_physics() libmesh_override
We don&#39;t allow systems to be attached to each other.
Definition: diff_system.h:146
void attach_qoi(DifferentiableQoI *qoi_in)
Attach external QoI object.
Definition: diff_system.h:212
const DifferentiablePhysics * get_physics() const
Definition: diff_system.h:169
virtual void element_postprocess(DiffContext &)
Does any work that needs to be done on elem in a postprocessing loop.
Definition: diff_system.h:269
bool print_solutions
Set print_solutions to true to print U whenever it is used in an assembly() call. ...
Definition: diff_system.h:314
This class provides a specific system class.
Definition: diff_physics.h:74
DifferentiableSystem(EquationSystems &es, const std::string &name, const unsigned int number)
Constructor.
Definition: diff_system.C:34
This class provides a specific system class.
Definition: diff_qoi.h:49
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real
unsigned int number() const
Definition: system.h:2006
ImplicitSystem Parent
The type of the parent.
Definition: diff_system.h:80
virtual UniquePtr< DifferentiablePhysics > clone_physics()=0
Copy of this object.
bool print_element_jacobians
Set print_element_jacobians to true to print each J_elem contribution.
Definition: diff_system.h:349
virtual void postprocess()
Executes a postprocessing loop over all elements, and if postprocess_sides is true over all sides...
Definition: diff_system.h:264
virtual ~DifferentiableSystem()
Destructor.
Definition: diff_system.C:56
virtual void reinit() libmesh_override
Reinitializes the member data fields associated with the system, so that, e.g., assemble() may be use...
Definition: diff_system.C:96
TimeSolver & get_time_solver()
Definition: diff_system.h:402
virtual void init_qoi(std::vector< Number > &)
Initialize system qoi.
Definition: diff_qoi.h:68
This class provides a specific system class.