libMesh
diff_physics.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_PHYSICS_H
21 #define LIBMESH_DIFF_PHYSICS_H
22 
23 // Local Includes
24 #include "libmesh/libmesh.h"
25 #include "libmesh/auto_ptr.h"
26 
27 // C++ includes
28 #include <vector>
29 
30 namespace libMesh
31 {
32 
33 // Forward Declarations
34 class System;
35 class DiffContext;
36 
75 {
76 public:
77 
83  compute_internal_sides (false),
88  {}
89 
93  virtual ~DifferentiablePhysics ();
94 
99 
103  virtual void clear_physics ();
104 
108  virtual void init_physics (const System & sys);
109 
123  virtual bool element_time_derivative (bool request_jacobian,
124  DiffContext &) {
125  return request_jacobian;
126  }
127 
141  virtual bool element_constraint (bool request_jacobian,
142  DiffContext &) {
143  return request_jacobian;
144  }
145 
153 
170  virtual bool side_time_derivative (bool request_jacobian,
171  DiffContext &) {
172  return request_jacobian;
173  }
174 
190  virtual bool side_constraint (bool request_jacobian,
191  DiffContext &) {
192  return request_jacobian;
193  }
194 
208  virtual bool nonlocal_time_derivative (bool request_jacobian,
209  DiffContext &) {
210  return request_jacobian;
211  }
212 
226  virtual bool nonlocal_constraint (bool request_jacobian,
227  DiffContext &) {
228  return request_jacobian;
229  }
230 
231 
247 #ifdef LIBMESH_ENABLE_DEPRECATED
248  virtual void time_evolving (unsigned int var)
249  {
250  libmesh_deprecated();
251  this->time_evolving(var,1);
252  }
253 #endif
254 
267  virtual void time_evolving (unsigned int var, unsigned int order);
268 
276  bool is_time_evolving (unsigned int var) const
277  {
278  libmesh_assert_less(var,_time_evolving.size());
279  libmesh_assert( _time_evolving[var] == 0 ||
280  _time_evolving[var] == 1 ||
281  _time_evolving[var] == 2 );
282  return _time_evolving[var];
283  }
284 
293  virtual bool eulerian_residual (bool request_jacobian,
294  DiffContext &) {
295  return request_jacobian;
296  }
297 
317  virtual bool mass_residual (bool request_jacobian,
318  DiffContext &) {
319  return request_jacobian;
320  }
321 
334  virtual bool side_mass_residual (bool request_jacobian,
335  DiffContext &) {
336  return request_jacobian;
337  }
338 
355  virtual bool nonlocal_mass_residual (bool request_jacobian,
356  DiffContext & c);
357 
373  virtual bool damping_residual (bool request_jacobian,
374  DiffContext &) {
375  return request_jacobian;
376  }
377 
390  virtual bool side_damping_residual (bool request_jacobian,
391  DiffContext &) {
392  return request_jacobian;
393  }
394 
405  virtual bool nonlocal_damping_residual (bool request_jacobian,
406  DiffContext &) {
407  return request_jacobian;
408  }
409 
410  /*
411  * Prepares the result of a build_context() call for use.
412  *
413  * Most FEMSystem-based problems will need to reimplement this in order to
414  * call FE::get_*() as their particular physics requires.
415  */
416  virtual void init_context(DiffContext &) {}
417 
441  virtual void set_mesh_system(System * sys);
442 
448  const System * get_mesh_system() const;
449 
455 
470  virtual void set_mesh_x_var(unsigned int var);
471 
476  unsigned int get_mesh_x_var() const;
477 
482  virtual void set_mesh_y_var(unsigned int var);
483 
488  unsigned int get_mesh_y_var() const;
489 
494  virtual void set_mesh_z_var(unsigned int var);
495 
500  unsigned int get_mesh_z_var() const;
501 
507  bool _eulerian_time_deriv (bool request_jacobian,
508  DiffContext &);
509 
511  { return !_first_order_vars.empty(); }
512 
516  const std::set<unsigned int> & get_first_order_vars() const
517  { return _first_order_vars; }
518 
519  bool is_first_order_var( unsigned int var ) const
520  { return _first_order_vars.find(var) != _first_order_vars.end(); }
521 
522 
524  { return !_second_order_vars.empty(); }
525 
529  const std::set<unsigned int> & get_second_order_vars() const
530  { return _second_order_vars; }
531 
532  bool is_second_order_var( unsigned int var ) const
533  { return _second_order_vars.find(var) != _second_order_vars.end(); }
534 
535 
536 protected:
537 
542 
547 
553  std::vector<unsigned int> _time_evolving;
554 
558  std::set<unsigned int> _first_order_vars;
559 
563  std::set<unsigned int> _second_order_vars;
564 
570  std::map<unsigned int,unsigned int> _second_order_dot_vars;
571 
572 };
573 
574 // ------------------------------------------------------------
575 // DifferentiablePhysics inline methods
576 
577 
578 inline
580 {
581  // For now we assume that we're doing fully coupled mesh motion
582  // if (sys && sys != this)
583  // libmesh_not_implemented();
584 
585  // For the foreseeable future we'll assume that we keep these
586  // Systems in the same EquationSystems
587  // libmesh_assert_equal_to (&this->get_equation_systems(),
588  // &sys->get_equation_systems());
589 
590  // And for the immediate future this code may not even work
591  libmesh_experimental();
592 
593  _mesh_sys = sys;
594 }
595 
596 
597 
598 inline
600 {
601  _mesh_x_var = var;
602 }
603 
604 
605 
606 inline
608 {
609  _mesh_y_var = var;
610 }
611 
612 
613 
614 inline
616 {
617  _mesh_z_var = var;
618 }
619 
620 
621 
622 inline
624 {
625  return _mesh_sys;
626 }
627 
628 inline
630 {
631  return _mesh_sys;
632 }
633 
634 inline
636 {
637  return _mesh_x_var;
638 }
639 
640 inline
642 {
643  return _mesh_y_var;
644 }
645 
646 inline
648 {
649  return _mesh_z_var;
650 }
651 
652 
653 
654 } // namespace libMesh
655 
656 
657 #endif // LIBMESH_DIFF_PHYSICS_H
virtual bool side_constraint(bool request_jacobian, DiffContext &)
Adds the constraint contribution on side of elem to elem_residual.
Definition: diff_physics.h:190
const std::set< unsigned int > & get_second_order_vars() const
Definition: diff_physics.h:529
virtual void clear_physics()
Clear any data structures associated with the physics.
Definition: diff_physics.C:33
unsigned int _mesh_x_var
Variables from which to acquire moving mesh information.
Definition: diff_physics.h:546
This class provides all data required for a physics package (e.g.
Definition: diff_context.h:54
const unsigned int invalid_uint
A number which is used quite often to represent an invalid or uninitialized value.
Definition: libmesh.h:184
virtual ~DifferentiablePhysics()
Destructor.
Definition: diff_physics.C:26
virtual bool element_time_derivative(bool request_jacobian, DiffContext &)
Adds the time derivative contribution on elem to elem_residual.
Definition: diff_physics.h:123
virtual bool eulerian_residual(bool request_jacobian, DiffContext &)
Adds a pseudo-convection contribution on elem to elem_residual, if the nodes of elem are being transl...
Definition: diff_physics.h:293
virtual bool side_damping_residual(bool request_jacobian, DiffContext &)
Subtracts a damping vector contribution on side of elem from elem_residual.
Definition: diff_physics.h:390
ImplicitSystem & sys
std::set< unsigned int > _second_order_vars
Variable indices for those variables that are second order in time.
Definition: diff_physics.h:563
const class libmesh_nullptr_t libmesh_nullptr
virtual bool nonlocal_time_derivative(bool request_jacobian, DiffContext &)
Adds any nonlocal time derivative contributions (e.g.
Definition: diff_physics.h:208
virtual bool damping_residual(bool request_jacobian, DiffContext &)
Subtracts a damping vector contribution on elem from elem_residual.
Definition: diff_physics.h:373
bool is_time_evolving(unsigned int var) const
Definition: diff_physics.h:276
The libMesh namespace provides an interface to certain functionality in the library.
std::map< unsigned int, unsigned int > _second_order_dot_vars
If the user adds any second order variables, then we need to also cache the map to their correspondin...
Definition: diff_physics.h:570
virtual bool nonlocal_constraint(bool request_jacobian, DiffContext &)
Adds any nonlocal constraint contributions (e.g.
Definition: diff_physics.h:226
bool have_second_order_vars() const
Definition: diff_physics.h:523
bool _eulerian_time_deriv(bool request_jacobian, DiffContext &)
This method simply combines element_time_derivative() and eulerian_residual(), which makes its addres...
Definition: diff_physics.C:102
virtual void set_mesh_x_var(unsigned int var)
Tells the DifferentiablePhysics that variable var from the mesh system should be used to update the x...
Definition: diff_physics.h:599
virtual bool mass_residual(bool request_jacobian, DiffContext &)
Subtracts a mass vector contribution on elem from elem_residual.
Definition: diff_physics.h:317
libmesh_assert(j)
std::unique_ptr< T > UniquePtr
Definition: auto_ptr.h:46
const System * get_mesh_system() const
Definition: diff_physics.h:623
unsigned int get_mesh_x_var() const
Definition: diff_physics.h:635
virtual void init_physics(const System &sys)
Initialize any data structures associated with the physics.
Definition: diff_physics.C:40
bool is_first_order_var(unsigned int var) const
Definition: diff_physics.h:519
virtual void set_mesh_y_var(unsigned int var)
Tells the DifferentiablePhysics that variable var from the mesh system should be used to update the y...
Definition: diff_physics.h:607
This is the base class for classes which contain information related to any physical process that mig...
Definition: system.h:76
System * _mesh_sys
System from which to acquire moving mesh information.
Definition: diff_physics.h:541
DifferentiablePhysics()
Constructor.
Definition: diff_physics.h:82
virtual void time_evolving(unsigned int var)
Tells the DiffSystem that variable var is evolving with respect to time.
Definition: diff_physics.h:248
This class provides a specific system class.
Definition: diff_physics.h:74
virtual void init_context(DiffContext &)
Definition: diff_physics.h:416
bool compute_internal_sides
compute_internal_sides is false by default, indicating that side_* computations will only be done on ...
Definition: diff_physics.h:152
bool is_second_order_var(unsigned int var) const
Definition: diff_physics.h:532
virtual bool nonlocal_damping_residual(bool request_jacobian, DiffContext &)
Subtracts any nonlocal damping vector contributions (e.g.
Definition: diff_physics.h:405
virtual UniquePtr< DifferentiablePhysics > clone_physics()=0
Copy of this object.
unsigned int get_mesh_z_var() const
Definition: diff_physics.h:647
virtual void set_mesh_z_var(unsigned int var)
Tells the DifferentiablePhysics that variable var from the mesh system should be used to update the z...
Definition: diff_physics.h:615
unsigned int get_mesh_y_var() const
Definition: diff_physics.h:641
const std::set< unsigned int > & get_first_order_vars() const
Definition: diff_physics.h:516
virtual bool side_mass_residual(bool request_jacobian, DiffContext &)
Subtracts a mass vector contribution on side of elem from elem_residual.
Definition: diff_physics.h:334
std::set< unsigned int > _first_order_vars
Variable indices for those variables that are first order in time.
Definition: diff_physics.h:558
virtual bool side_time_derivative(bool request_jacobian, DiffContext &)
Adds the time derivative contribution on side of elem to elem_residual.
Definition: diff_physics.h:170
virtual void set_mesh_system(System *sys)
Tells the DifferentiablePhysics that system sys contains the isoparametric Lagrangian variables which...
Definition: diff_physics.h:579
virtual bool element_constraint(bool request_jacobian, DiffContext &)
Adds the constraint contribution on elem to elem_residual.
Definition: diff_physics.h:141
std::vector< unsigned int > _time_evolving
Stores unsigned int to tell us which variables are evolving as first order in time (1)...
Definition: diff_physics.h:553
virtual bool nonlocal_mass_residual(bool request_jacobian, DiffContext &c)
Subtracts any nonlocal mass vector contributions (e.g.
Definition: diff_physics.C:63