libMesh
Public Types | Public Member Functions | Public Attributes | Protected Member Functions | Private Member Functions | List of all members
libMesh::TransientSystem< Base > Class Template Reference

Manages storage and variables for transient systems. More...

#include <transient_system.h>

Inheritance diagram for libMesh::TransientSystem< Base >:
[legend]

Public Types

typedef TransientSystem< Base > sys_type
 The type of system. More...
 

Public Member Functions

 TransientSystem (EquationSystems &es, const std::string &name, const unsigned int number)
 Constructor. More...
 
 TransientSystem (const TransientSystem &)=delete
 Special functions. More...
 
TransientSystemoperator= (const TransientSystem &)=delete
 
 TransientSystem (TransientSystem &&)=default
 
TransientSystemoperator= (TransientSystem &&)=delete
 
virtual ~TransientSystem ()
 
sys_typesystem ()
 
virtual void clear () override
 Clear all the data structures associated with the system. More...
 
virtual std::string system_type () const override
 
Number old_solution (const dof_id_type global_dof_number) const
 
Number older_solution (const dof_id_type global_dof_number) const
 

Public Attributes

NumericVector< Number > * old_local_solution
 All the values I need to compute my contribution to the simulation at hand. More...
 
NumericVector< Number > * older_local_solution
 All the values I need to compute my contribution to the simulation at hand. More...
 

Protected Member Functions

virtual void re_update () override
 Re-update the local values when the mesh has changed. More...
 

Private Member Functions

virtual void add_old_vectors ()
 Helper function for (re-)adding old and older solution vectors. More...
 

Detailed Description

template<class Base>
class libMesh::TransientSystem< Base >

Manages storage and variables for transient systems.

This class is a specialized system for solving transient systems, e.g., systems with a time dependency. It provides appropriate storage, manages variables and ensures consistency for these systems.

This template class adds the functionality to manage, in addition to the solution vector, an old solution and an older solution vector. These vectors are useful to simulate transient systems.

Note
Additional vectors/matrices can be added via parent class interfaces.
Author
Benjamin S. Kirk
Date
2004 Used for solving transient systems of equations.

Definition at line 57 of file transient_system.h.

Member Typedef Documentation

◆ sys_type

template<class Base>
typedef TransientSystem<Base> libMesh::TransientSystem< Base >::sys_type

The type of system.

Definition at line 84 of file transient_system.h.

Constructor & Destructor Documentation

◆ TransientSystem() [1/3]

template<class Base >
libMesh::TransientSystem< Base >::TransientSystem ( EquationSystems es,
const std::string &  name,
const unsigned int  number 
)

Constructor.

Definition at line 45 of file transient_system.C.

47  :
48 
49  Base(es, name_in, number_in),
50  old_local_solution(nullptr),
51  older_local_solution(nullptr)
52 {
53  this->add_old_vectors();
54 }
virtual void add_old_vectors()
Helper function for (re-)adding old and older solution vectors.
NumericVector< Number > * old_local_solution
All the values I need to compute my contribution to the simulation at hand.
NumericVector< Number > * older_local_solution
All the values I need to compute my contribution to the simulation at hand.

◆ TransientSystem() [2/3]

template<class Base>
libMesh::TransientSystem< Base >::TransientSystem ( const TransientSystem< Base > &  )
delete

Special functions.

◆ TransientSystem() [3/3]

template<class Base>
libMesh::TransientSystem< Base >::TransientSystem ( TransientSystem< Base > &&  )
default

◆ ~TransientSystem()

template<class Base >
libMesh::TransientSystem< Base >::~TransientSystem ( )
virtualdefault

Member Function Documentation

◆ add_old_vectors()

template<class Base >
void libMesh::TransientSystem< Base >::add_old_vectors ( )
privatevirtual

Helper function for (re-)adding old and older solution vectors.

Definition at line 135 of file transient_system.C.

136 {
137  ParallelType type =
138 #ifdef LIBMESH_ENABLE_GHOSTED
139  GHOSTED;
140 #else
141  SERIAL;
142 #endif
143 
144  old_local_solution = &(this->add_vector("_transient_old_local_solution", true, type));
145  older_local_solution = &(this->add_vector("_transient_older_local_solution", true, type));
146 }
NumericVector< Number > * old_local_solution
All the values I need to compute my contribution to the simulation at hand.
NumericVector< Number > * older_local_solution
All the values I need to compute my contribution to the simulation at hand.
ParallelType
Defines an enum for parallel data structure types.

◆ clear()

template<class Base >
void libMesh::TransientSystem< Base >::clear ( )
overridevirtual

Clear all the data structures associated with the system.

Reimplemented in libMesh::TransientRBConstruction.

Definition at line 64 of file transient_system.C.

65 {
66  // clear the parent data
67  Base::clear();
68 
69  // Restore us to a "basic" state
70  this->add_old_vectors();
71 }
virtual void add_old_vectors()
Helper function for (re-)adding old and older solution vectors.

◆ old_solution()

template<class Base >
Number libMesh::TransientSystem< Base >::old_solution ( const dof_id_type  global_dof_number) const
Returns
The old solution (at the previous timestep) for the specified global DOF.

Definition at line 111 of file transient_system.C.

Referenced by assemble_cd(), and assemble_stokes().

112 {
113  // Check the sizes
114  libmesh_assert_less (global_dof_number, this->get_dof_map().n_dofs());
115  libmesh_assert_less (global_dof_number, old_local_solution->size());
116 
117  return (*old_local_solution)(global_dof_number);
118 }
virtual numeric_index_type size() const =0
NumericVector< Number > * old_local_solution
All the values I need to compute my contribution to the simulation at hand.

◆ older_solution()

template<class Base >
Number libMesh::TransientSystem< Base >::older_solution ( const dof_id_type  global_dof_number) const
Returns
The older solution (two timesteps ago) for the specified global DOF.

Definition at line 123 of file transient_system.C.

124 {
125  // Check the sizes
126  libmesh_assert_less (global_dof_number, this->get_dof_map().n_dofs());
127  libmesh_assert_less (global_dof_number, older_local_solution->size());
128 
129  return (*older_local_solution)(global_dof_number);
130 }
virtual numeric_index_type size() const =0
NumericVector< Number > * older_local_solution
All the values I need to compute my contribution to the simulation at hand.

◆ operator=() [1/2]

template<class Base>
TransientSystem& libMesh::TransientSystem< Base >::operator= ( const TransientSystem< Base > &  )
delete

◆ operator=() [2/2]

template<class Base>
TransientSystem& libMesh::TransientSystem< Base >::operator= ( TransientSystem< Base > &&  )
delete

◆ re_update()

template<class Base >
void libMesh::TransientSystem< Base >::re_update ( )
overrideprotectedvirtual

Re-update the local values when the mesh has changed.

This method takes the data updated by update() and makes it up-to-date on the current mesh.

Definition at line 76 of file transient_system.C.

77 {
78  // re_update the parent system
79  Base::re_update ();
80 
81  const std::vector<dof_id_type> & send_list = this->get_dof_map().get_send_list ();
82 
83  const dof_id_type first_local_dof = Base::get_dof_map().first_dof();
84  const dof_id_type end_local_dof = Base::get_dof_map().end_dof();
85 
86  // Check sizes
87  libmesh_assert_greater_equal (end_local_dof, first_local_dof);
88  libmesh_assert_greater_equal (older_local_solution->size(), send_list.size());
89  libmesh_assert_greater_equal (old_local_solution->size(), send_list.size());
90 
91  // Even if we don't have to do anything ourselves, localize() may
92  // use parallel_only tools
93  // if (first_local_dof == end_local_dof)
94  // return;
95 
96  // Update the old & older solutions with the send_list,
97  // which may have changed since their last update.
98  older_local_solution->localize (first_local_dof,
99  end_local_dof-1,
100  send_list);
101 
102  old_local_solution->localize (first_local_dof,
103  end_local_dof-1,
104  send_list);
105 }
virtual numeric_index_type size() const =0
NumericVector< Number > * old_local_solution
All the values I need to compute my contribution to the simulation at hand.
NumericVector< Number > * older_local_solution
All the values I need to compute my contribution to the simulation at hand.
uint8_t dof_id_type
Definition: id_types.h:67
virtual void localize(std::vector< T > &v_local) const =0
Creates a copy of the global vector in the local vector v_local.

◆ system()

template<class Base>
sys_type& libMesh::TransientSystem< Base >::system ( )
inline
Returns
A reference to *this.

Definition at line 89 of file transient_system.h.

89 { return *this; }

◆ system_type()

template<class Base >
std::string libMesh::TransientSystem< Base >::system_type ( ) const
inlineoverridevirtual
Returns
"Transient" prepended to T::system_type(). Helps in identifying the system type in an equation system file.

Definition at line 172 of file transient_system.h.

173 {
174  std::string type = "Transient";
175  type += Base::system_type ();
176 
177  return type;
178 }

Member Data Documentation

◆ old_local_solution

template<class Base>
NumericVector<Number>* libMesh::TransientSystem< Base >::old_local_solution

All the values I need to compute my contribution to the simulation at hand.

Think of this as the current solution with any ghost values needed from other processors.

Definition at line 126 of file transient_system.h.

Referenced by main(), run_timestepping(), SystemsTest::testProjectCube(), SystemsTest::testProjectLine(), SystemsTest::testProjectSquare(), and SystemsTest::tripleValueTest().

◆ older_local_solution

template<class Base>
NumericVector<Number>* libMesh::TransientSystem< Base >::older_local_solution

All the values I need to compute my contribution to the simulation at hand.

Think of this as the current solution with any ghost values needed from other processors.

Definition at line 134 of file transient_system.h.

Referenced by run_timestepping(), SystemsTest::testProjectCube(), SystemsTest::testProjectLine(), SystemsTest::testProjectSquare(), and SystemsTest::tripleValueTest().


The documentation for this class was generated from the following files: