libMesh::TransientSystem< Base > Class Template Reference

Used for solving transient systems of equations. More...

#include <transient_system.h>

Inheritance diagram for libMesh::TransientSystem< Base >:

Public Types

typedef TransientSystem< Base > sys_type
 

Public Member Functions

 TransientSystem (EquationSystems &es, const std::string &name, const unsigned int number)
 
virtual ~TransientSystem ()
 
sys_typesystem ()
 
virtual void clear () libmesh_override
 
virtual std::string system_type () const libmesh_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

std::unique_ptr< NumericVector< Number > > old_local_solution
 
std::unique_ptr< NumericVector< Number > > older_local_solution
 

Protected Member Functions

virtual void re_update () libmesh_override
 

Private Member Functions

virtual void add_old_vectors ()
 

Detailed Description

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

Used for solving transient systems of equations.

This class provides a specific system class. It aims at transient systems, offering nothing more than just the essentials needed to solve a system.

Note
Additional vectors/matrices can be added via parent class interfaces.
Author
Benjamin S. Kirk
Date
2004

Definition at line 51 of file transient_system.h.

Member Typedef Documentation

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

The type of system.

Definition at line 71 of file transient_system.h.

Constructor & Destructor Documentation

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

Constructor. Initializes required data structures.

Definition at line 39 of file transient_system.C.

References libMesh::TransientSystem< Base >::add_old_vectors().

41  :
42 
43  Base (es, name_in, number_in)
44 {
45  this->add_old_vectors();
46 }
virtual void add_old_vectors()
template<class Base >
libMesh::TransientSystem< Base >::~TransientSystem ( )
virtual

Destructor.

Definition at line 51 of file transient_system.C.

References libMesh::TransientSystem< Base >::clear(), libMesh::TransientSystem< Base >::old_local_solution, and libMesh::TransientSystem< Base >::older_local_solution.

52 {
53  this->clear();
54 
55  // We still have std::unique_ptrs for API compatibility, but
56  // now that we're System::add_vector()ing these, we can trust
57  // the base class to handle memory management
58  old_local_solution.release();
59  older_local_solution.release();
60 }
std::unique_ptr< NumericVector< Number > > old_local_solution
std::unique_ptr< NumericVector< Number > > older_local_solution
virtual void clear() libmesh_override

Member Function Documentation

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

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

Definition at line 146 of file transient_system.C.

References libMesh::GHOSTED, libMesh::TransientSystem< Base >::old_local_solution, libMesh::TransientSystem< Base >::older_local_solution, and libMesh::SERIAL.

Referenced by libMesh::TransientSystem< Base >::clear(), and libMesh::TransientSystem< Base >::TransientSystem().

147 {
148 #ifdef LIBMESH_ENABLE_GHOSTED
150  std::unique_ptr<NumericVector<Number>>
151  (&(this->add_vector("_transient_old_local_solution", true, GHOSTED)));
153  std::unique_ptr<NumericVector<Number>>
154  (&(this->add_vector("_transient_older_local_solution", true, GHOSTED)));
155 #else
157  std::unique_ptr<NumericVector<Number>>
158  (&(this->add_vector("_transient_old_local_solution", true, SERIAL)));
160  std::unique_ptr<NumericVector<Number>>
161  (&(this->add_vector("_transient_older_local_solution", true, SERIAL)));
162 #endif
163 }
std::unique_ptr< NumericVector< Number > > old_local_solution
std::unique_ptr< NumericVector< Number > > older_local_solution
template<class Base >
void libMesh::TransientSystem< Base >::clear ( )
virtual

Clear all the data structures associated with the system.

Reimplemented in libMesh::TransientRBConstruction.

Definition at line 65 of file transient_system.C.

References libMesh::TransientSystem< Base >::add_old_vectors(), libMesh::TransientSystem< Base >::old_local_solution, and libMesh::TransientSystem< Base >::older_local_solution.

Referenced by libMesh::TransientSystem< RBConstruction >::system(), and libMesh::TransientSystem< Base >::~TransientSystem().

66 {
67  // clear the parent data
68  Base::clear();
69 
70  // the old & older local solutions
71  // are now deleted by System!
72  // old_local_solution->clear();
73  // older_local_solution->clear();
74 
75  // FIXME: This preserves maximum backwards compatibility,
76  // but is probably grossly unnecessary:
77  old_local_solution.release();
78  older_local_solution.release();
79 
80  // Restore us to a "basic" state
81  this->add_old_vectors();
82 }
virtual void add_old_vectors()
std::unique_ptr< NumericVector< Number > > old_local_solution
std::unique_ptr< NumericVector< Number > > older_local_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 122 of file transient_system.C.

References libMesh::TransientSystem< Base >::old_local_solution.

Referenced by libMesh::TransientSystem< RBConstruction >::system().

123 {
124  // Check the sizes
125  libmesh_assert_less (global_dof_number, this->get_dof_map().n_dofs());
126  libmesh_assert_less (global_dof_number, old_local_solution->size());
127 
128  return (*old_local_solution)(global_dof_number);
129 }
std::unique_ptr< NumericVector< Number > > old_local_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 134 of file transient_system.C.

References libMesh::TransientSystem< Base >::older_local_solution.

Referenced by libMesh::TransientSystem< RBConstruction >::system().

135 {
136  // Check the sizes
137  libmesh_assert_less (global_dof_number, this->get_dof_map().n_dofs());
138  libmesh_assert_less (global_dof_number, older_local_solution->size());
139 
140  return (*older_local_solution)(global_dof_number);
141 }
std::unique_ptr< NumericVector< Number > > older_local_solution
template<class Base >
void libMesh::TransientSystem< Base >::re_update ( )
protectedvirtual

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 87 of file transient_system.C.

References libMesh::TransientSystem< Base >::old_local_solution, and libMesh::TransientSystem< Base >::older_local_solution.

88 {
89  // re_update the parent system
90  Base::re_update ();
91 
92  const std::vector<dof_id_type> & send_list = this->get_dof_map().get_send_list ();
93 
94  const dof_id_type first_local_dof = Base::get_dof_map().first_dof();
95  const dof_id_type end_local_dof = Base::get_dof_map().end_dof();
96 
97  // Check sizes
98  libmesh_assert_greater_equal (end_local_dof, first_local_dof);
99  libmesh_assert_greater_equal (older_local_solution->size(), send_list.size());
100  libmesh_assert_greater_equal (old_local_solution->size(), send_list.size());
101 
102  // Even if we don't have to do anything ourselves, localize() may
103  // use parallel_only tools
104  // if (first_local_dof == end_local_dof)
105  // return;
106 
107  // Update the old & older solutions with the send_list,
108  // which may have changed since their last update.
109  older_local_solution->localize (first_local_dof,
110  end_local_dof-1,
111  send_list);
112 
113  old_local_solution->localize (first_local_dof,
114  end_local_dof-1,
115  send_list);
116 }
std::unique_ptr< NumericVector< Number > > old_local_solution
std::unique_ptr< NumericVector< Number > > older_local_solution
uint8_t dof_id_type
Definition: id_types.h:64
template<class Base>
sys_type& libMesh::TransientSystem< Base >::system ( )
inline
Returns
A reference to *this.

Definition at line 76 of file transient_system.h.

76 { return *this; }
template<class Base >
std::string libMesh::TransientSystem< Base >::system_type ( ) const
inlinevirtual
Returns
"Transient" prepended to T::system_type(). Helps in identifying the system type in an equation system file.

Definition at line 160 of file transient_system.h.

Referenced by libMesh::TransientSystem< RBConstruction >::system().

161 {
162  std::string type = "Transient";
163  type += Base::system_type ();
164 
165  return type;
166 }

Member Data Documentation

template<class Base>
std::unique_ptr<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 113 of file transient_system.h.

Referenced by libMesh::TransientSystem< Base >::add_old_vectors(), libMesh::TransientSystem< Base >::clear(), libMesh::TransientSystem< Base >::old_solution(), libMesh::TransientSystem< Base >::re_update(), and libMesh::TransientSystem< Base >::~TransientSystem().

template<class Base>
std::unique_ptr<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 121 of file transient_system.h.

Referenced by libMesh::TransientSystem< Base >::add_old_vectors(), libMesh::TransientSystem< Base >::clear(), libMesh::TransientSystem< Base >::older_solution(), libMesh::TransientSystem< Base >::re_update(), and libMesh::TransientSystem< Base >::~TransientSystem().


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