transient_system.C
Go to the documentation of this file.
1 // The libMesh Finite Element Library.
2 // Copyright (C) 2002-2018 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 // C++ includes
21 
22 // Local includes
27 #include "libmesh/dof_map.h"
28 #include "libmesh/numeric_vector.h"
30 #include "libmesh/eigen_system.h"
31 
32 namespace libMesh
33 {
34 
35 
36 // ------------------------------------------------------------
37 // TransientSystem implementation
38 template <class Base>
40  const std::string & name_in,
41  const unsigned int number_in) :
42 
43  Base (es, name_in, number_in)
44 {
45  this->add_old_vectors();
46 }
47 
48 
49 
50 template <class Base>
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 }
61 
62 
63 
64 template <class Base>
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 }
83 
84 
85 
86 template <class Base>
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 }
117 
118 
119 
120 
121 template <class Base>
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 }
130 
131 
132 
133 template <class Base>
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 }
142 
143 
144 
145 template <class Base>
147 {
148 #ifdef LIBMESH_ENABLE_GHOSTED
149  old_local_solution =
150  std::unique_ptr<NumericVector<Number>>
151  (&(this->add_vector("_transient_old_local_solution", true, GHOSTED)));
152  older_local_solution =
153  std::unique_ptr<NumericVector<Number>>
154  (&(this->add_vector("_transient_older_local_solution", true, GHOSTED)));
155 #else
156  old_local_solution =
157  std::unique_ptr<NumericVector<Number>>
158  (&(this->add_vector("_transient_old_local_solution", true, SERIAL)));
159  older_local_solution =
160  std::unique_ptr<NumericVector<Number>>
161  (&(this->add_vector("_transient_older_local_solution", true, SERIAL)));
162 #endif
163 }
164 
165 // ------------------------------------------------------------
166 // TransientSystem instantiations
169 template class TransientSystem<ExplicitSystem>;
170 template class TransientSystem<System>;
171 template class TransientSystem<RBConstruction>;
172 #ifdef LIBMESH_HAVE_SLEPC
173 template class TransientSystem<EigenSystem>;
174 #endif
175 
176 } // namespace libMesh
Manages multiples systems of equations.
Manages storage and variables for transient systems.
TransientSystem(EquationSystems &es, const std::string &name, const unsigned int number)
uint8_t dof_id_type
Definition: id_types.h:64