diff_physics.h
Go to the documentation of this file.
1 
2 // The libMesh Finite Element Library.
3 // Copyright (C) 2002-2014 Benjamin S. Kirk, John W. Peterson, Roy H. Stogner
4 
5 // This library is free software; you can redistribute it and/or
6 // modify it under the terms of the GNU Lesser General Public
7 // License as published by the Free Software Foundation; either
8 // version 2.1 of the License, or (at your option) any later version.
9 
10 // This library is distributed in the hope that it will be useful,
11 // but WITHOUT ANY WARRANTY; without even the implied warranty of
12 // MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
13 // Lesser General Public License for more details.
14 
15 // You should have received a copy of the GNU Lesser General Public
16 // License along with this library; if not, write to the Free Software
17 // Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA
18 
19 
20 
21 #ifndef LIBMESH_DIFF_PHYSICS_H
22 #define LIBMESH_DIFF_PHYSICS_H
23 
24 // Local Includes
25 #include "libmesh/libmesh.h"
26 #include "libmesh/auto_ptr.h"
27 #include "libmesh/diff_context.h"
28 
29 // C++ includes
30 #include <vector>
31 
32 namespace libMesh
33 {
34 
35 // Forward Declarations
36 class System;
37 class TimeSolver;
38 
39 template <typename T> class NumericVector;
40 
53 // ------------------------------------------------------------
54 // DifferentiableSystem class definition
55 
57 {
58 public:
59 
65  compute_internal_sides (false),
66  _mesh_sys (NULL),
67  _mesh_x_var (libMesh::invalid_uint),
68  _mesh_y_var (libMesh::invalid_uint),
69  _mesh_z_var (libMesh::invalid_uint)
70  {}
71 
75  virtual ~DifferentiablePhysics ();
76 
81 
85  virtual void clear_physics ();
86 
90  virtual void init_physics (const System& sys);
91 
105  virtual bool element_time_derivative (bool request_jacobian,
106  DiffContext &) {
107  return request_jacobian;
108  }
109 
123  virtual bool element_constraint (bool request_jacobian,
124  DiffContext &) {
125  return request_jacobian;
126  }
127 
135 
152  virtual bool side_time_derivative (bool request_jacobian,
153  DiffContext &) {
154  return request_jacobian;
155  }
156 
172  virtual bool side_constraint (bool request_jacobian,
173  DiffContext &) {
174  return request_jacobian;
175  }
176 
188  virtual void time_evolving (unsigned int var) {
189  if (_time_evolving.size() <= var)
190  _time_evolving.resize(var+1, false);
191  _time_evolving[var] = true;
192  }
193 
201  bool is_time_evolving (unsigned int var) const {
202  return _time_evolving[var];
203  }
204 
213  virtual bool eulerian_residual (bool request_jacobian,
214  DiffContext &) {
215  return request_jacobian;
216  }
217 
229  virtual bool mass_residual (bool request_jacobian,
230  DiffContext &) {
231  return request_jacobian;
232  }
233 
246  virtual bool side_mass_residual (bool request_jacobian,
247  DiffContext &) {
248  return request_jacobian;
249  }
250 
251  /*
252  * Prepares the result of a build_context() call for use.
253  *
254  * Most FEMSystem-based problems will need to reimplement this in order to
255  * call FE::get_*() as their particular physics requires.
256  */
257  virtual void init_context(DiffContext &) {}
258 
282  virtual void set_mesh_system(System* sys);
283 
289  const System* get_mesh_system() const;
290 
296 
311  virtual void set_mesh_x_var(unsigned int var);
312 
317  unsigned int get_mesh_x_var() const;
318 
323  virtual void set_mesh_y_var(unsigned int var);
324 
329  unsigned int get_mesh_y_var() const;
330 
335  virtual void set_mesh_z_var(unsigned int var);
336 
341  unsigned int get_mesh_z_var() const;
342 
343 
344 protected:
345 
350 
355 
360  std::vector<bool> _time_evolving;
361 };
362 
363 
364 
365 // ------------------------------------------------------------
366 // DifferentiablePhysics inline methods
367 
368 
369 inline
371 {
372  // For now we assume that we're doing fully coupled mesh motion
373 // if (sys && sys != this)
374 // libmesh_not_implemented();
375 
376  // For the foreseeable future we'll assume that we keep these
377  // Systems in the same EquationSystems
378  // libmesh_assert_equal_to (&this->get_equation_systems(),
379  // &sys->get_equation_systems());
380 
381  // And for the immediate future this code may not even work
382  libmesh_experimental();
383 
384  _mesh_sys = sys;
385 }
386 
387 
388 
389 inline
391 {
392  _mesh_x_var = var;
393 }
394 
395 
396 
397 inline
399 {
400  _mesh_y_var = var;
401 }
402 
403 
404 
405 inline
407 {
408  _mesh_z_var = var;
409 }
410 
411 
412 
413 inline
415 {
416  return _mesh_sys;
417 }
418 
419 inline
421 {
422  return _mesh_sys;
423 }
424 
425 inline
427 {
428  return _mesh_x_var;
429 }
430 
431 inline
433 {
434  return _mesh_y_var;
435 }
436 
437 inline
439 {
440  return _mesh_z_var;
441 }
442 
443 
444 
445 } // namespace libMesh
446 
447 
448 #endif // LIBMESH_DIFF_PHYSICS_H

Site Created By: libMesh Developers
Last modified: February 07 2014 16:57:04 UTC

Hosted By:
SourceForge.net Logo