Line data Source code
1 : //* This file is part of the MOOSE framework
2 : //* https://mooseframework.inl.gov
3 : //*
4 : //* All rights reserved, see COPYRIGHT for full restrictions
5 : //* https://github.com/idaholab/moose/blob/master/COPYRIGHT
6 : //*
7 : //* Licensed under LGPL 2.1, please see LICENSE for details
8 : //* https://www.gnu.org/licenses/lgpl-2.1.html
9 :
10 : #include "RestartableEquationSystems.h"
11 :
12 : #include "DataIO.h"
13 :
14 : #include "libmesh/dof_map.h"
15 : #include "libmesh/dof_object.h"
16 : #include "libmesh/elem.h"
17 : #include "libmesh/node.h"
18 :
19 : const std::string RestartableEquationSystems::SystemHeader::system_solution_name =
20 : "SYSTEM_SOLUTION";
21 :
22 62498 : RestartableEquationSystems::RestartableEquationSystems(libMesh::MeshBase & mesh)
23 62498 : : _es(mesh), _load_all_vectors(true)
24 : {
25 62498 : }
26 :
27 : RestartableEquationSystems::EquationSystemsHeader
28 45214 : RestartableEquationSystems::buildHeader(
29 : const std::vector<const libMesh::DofObject *> & ordered_objects) const
30 : {
31 45214 : EquationSystemsHeader es_header;
32 :
33 : // Systems
34 136666 : for (const auto sys_num : make_range(_es.n_systems()))
35 : {
36 91452 : const auto & sys = _es.get_system(sys_num);
37 :
38 91452 : SystemHeader sys_header;
39 91452 : sys_header.name = sys.name();
40 91452 : sys_header.type = sys.system_type();
41 :
42 : // Variables in the system
43 183043 : for (const auto var_num : make_range(sys.n_vars()))
44 : {
45 91591 : const auto & var = sys.variable(var_num);
46 :
47 91591 : VariableHeader var_header;
48 91591 : var_header.name = var.name();
49 91591 : var_header.type = var.type();
50 91591 : var_header.size = 0;
51 91591 : var_header.variable = &var;
52 :
53 : mooseAssert(_es.comm().verify("sys_" + sys.name() + "_var_" + var.name()),
54 : "Out of order in parallel");
55 : mooseAssert(!sys_header.variables.count(var.name()), "Already inserted");
56 :
57 : // Non-SCALAR variable
58 91591 : if (var.type().family != SCALAR)
59 : {
60 18530584 : for (const auto & obj : ordered_objects)
61 18440684 : var_header.size += sizeof(Real) * obj->n_comp(sys.number(), var.number());
62 : }
63 : // SCALAR variable on the last rank
64 1691 : else if (_es.processor_id() == _es.n_processors() - 1)
65 : {
66 1654 : std::vector<dof_id_type> scalar_dofs;
67 1654 : sys.get_dof_map().SCALAR_dof_indices(scalar_dofs, var.number());
68 1654 : var_header.size += sizeof(Real) * scalar_dofs.size();
69 1654 : }
70 :
71 91591 : sys_header.variables.emplace(var.name(), var_header);
72 91591 : }
73 :
74 : // System vector
75 91452 : auto & sys_vec_header = sys_header.vectors[SystemHeader::system_solution_name];
76 91452 : sys_vec_header.name = SystemHeader::system_solution_name;
77 91452 : sys_vec_header.vector = sys.solution.get();
78 91452 : sys_vec_header.type = sys_vec_header.vector->type();
79 430719 : for (const auto vec_num : make_range(sys.n_vectors()))
80 : {
81 : mooseAssert(_es.comm().verify("sys_" + sys.name() + "_vec_" + sys.vector_name(vec_num)),
82 : "Out of order in parallel");
83 339267 : const auto & name = sys.vector_name(vec_num);
84 : mooseAssert(!sys_header.vectors.count(name), "Already inserted");
85 339267 : auto & vec_header = sys_header.vectors[name];
86 339267 : vec_header.name = sys.vector_name(vec_num);
87 339267 : vec_header.projections = sys.vector_preservation(vec_header.name);
88 339267 : vec_header.vector = &sys.get_vector(vec_header.name);
89 339267 : vec_header.type = vec_header.vector->type();
90 : }
91 :
92 : // System in this EquationSystems
93 : mooseAssert(!es_header.systems.count(sys.name()), "Already inserted");
94 91452 : es_header.systems.emplace(sys.name(), sys_header);
95 91452 : }
96 :
97 : // Setup the positions in each vector for easy access later
98 45214 : std::size_t offset = 0;
99 136666 : for (auto & sys_name_header_pair : es_header.systems)
100 : {
101 91452 : auto & sys_header = sys_name_header_pair.second;
102 522171 : for (auto & vec_name_header_pair : sys_header.vectors)
103 : {
104 430719 : auto & vec_header = vec_name_header_pair.second;
105 861875 : for (const auto & var_name_header_pair : sys_header.variables)
106 : {
107 431156 : const auto & var_header = var_name_header_pair.second;
108 : mooseAssert(!vec_header.variable_offset.count(var_header.name), "Already inserted");
109 431156 : vec_header.variable_offset[var_header.name] = offset;
110 431156 : offset += var_header.size;
111 : }
112 : }
113 : }
114 :
115 45214 : es_header.data_size = offset;
116 :
117 45214 : return es_header;
118 0 : }
119 :
120 : std::vector<const libMesh::DofObject *>
121 59455 : RestartableEquationSystems::orderDofObjects() const
122 : {
123 59455 : std::vector<const libMesh::DofObject *> objects;
124 356730 : auto add = [&objects](const auto begin, const auto end)
125 : {
126 118910 : std::set<const libMesh::DofObject *, libMesh::CompareDofObjectsByID> ordered(begin, end);
127 118910 : objects.insert(objects.end(), ordered.begin(), ordered.end());
128 178365 : };
129 :
130 59455 : const auto & mesh = _es.get_mesh();
131 59455 : add(mesh.local_elements_begin(), mesh.local_elements_end());
132 59455 : add(mesh.local_nodes_begin(), mesh.local_nodes_end());
133 :
134 118910 : return objects;
135 0 : }
136 :
137 : void
138 45214 : RestartableEquationSystems::store(std::ostream & stream) const
139 : {
140 : // Order objects (elements and then nodes) by ID for storing
141 45214 : const auto ordered_objects = orderDofObjects();
142 :
143 : // Store the header (systems, variables, vectors)
144 45214 : EquationSystemsHeader es_header = buildHeader(ordered_objects);
145 45214 : dataStore(stream, es_header, nullptr);
146 :
147 : // Store the ordered objects so we can do a sanity check on if we're
148 : // loading the same thing
149 : {
150 45214 : std::vector<dof_id_type> ordered_objects_ids(ordered_objects.size());
151 9802013 : for (const auto i : index_range(ordered_objects))
152 9756799 : ordered_objects_ids[i] = ordered_objects[i]->id();
153 45214 : dataStore(stream, ordered_objects_ids, nullptr);
154 45214 : }
155 :
156 : #ifndef NDEBUG
157 : const std::size_t data_initial_position = static_cast<std::size_t>(stream.tellp());
158 : #endif
159 :
160 : // Store each system
161 136666 : for (const auto & sys_name_header_pair : es_header.systems)
162 : {
163 91452 : const auto & sys_header = sys_name_header_pair.second;
164 91452 : const auto & sys = _es.get_system(sys_header.name);
165 :
166 : // Store each vector and variable
167 522171 : for (const auto & vec_name_header_pair : sys_header.vectors)
168 : {
169 430719 : const auto & vec_header = vec_name_header_pair.second;
170 430719 : const auto & vec = *vec_header.vector;
171 861875 : for (const auto & var_name_header_pair : sys_header.variables)
172 : {
173 431156 : const auto & var_header = var_name_header_pair.second;
174 431156 : const auto & var = *var_header.variable;
175 :
176 : #ifndef NDEBUG
177 : const std::size_t var_initial_position = stream.tellp();
178 : #endif
179 :
180 : // Non-SCALAR variable
181 431156 : if (var.type().family != SCALAR)
182 : {
183 : // Store for each component of each element and node
184 81862601 : for (const auto & obj : ordered_objects)
185 123034066 : for (const auto comp : make_range(obj->n_comp(sys.number(), var.number())))
186 : {
187 41594489 : auto val = vec(obj->dof_number(sys.number(), var.number(), comp));
188 41594489 : dataStore(stream, val, nullptr);
189 : }
190 : }
191 : // SCALAR variable on the last rank
192 8132 : else if (_es.processor_id() == _es.n_processors() - 1)
193 : {
194 8021 : const auto & dof_map = sys.get_dof_map();
195 8021 : std::vector<dof_id_type> scalar_dofs;
196 8021 : dof_map.SCALAR_dof_indices(scalar_dofs, var.number());
197 17408 : for (const auto dof : scalar_dofs)
198 : {
199 9387 : auto val = vec(dof);
200 9387 : dataStore(stream, val, nullptr);
201 : }
202 8021 : }
203 :
204 : #ifndef NDEBUG
205 : const std::size_t data_offset = var_initial_position - data_initial_position;
206 : mooseAssert(vec_header.variable_offset.at(var_header.name) == data_offset,
207 : "Invalid offset");
208 :
209 : const std::size_t current_position = static_cast<std::size_t>(stream.tellp());
210 : const std::size_t var_size = current_position - var_initial_position;
211 : mooseAssert(var_size == sys_header.variables.at(var.name()).size, "Incorrect assumed size");
212 : #endif
213 : }
214 : }
215 : }
216 :
217 : mooseAssert((data_initial_position + es_header.data_size) ==
218 : static_cast<std::size_t>(stream.tellp()),
219 : "Incorrect assumed size");
220 45214 : }
221 :
222 : void
223 14241 : RestartableEquationSystems::load(std::istream & stream)
224 : {
225 : // Load the header (systems, variables, vectors)
226 : // We do this first so that the loader can make informed decisions
227 : // on what to put where based on everything that is available
228 14241 : _loaded_header.systems.clear();
229 14241 : dataLoad(stream, _loaded_header, nullptr);
230 :
231 : // Order objects (elements and then node) by ID for storing
232 14241 : _loaded_ordered_objects = orderDofObjects();
233 :
234 : // Sanity check on if we're loading the same thing
235 : {
236 14241 : std::vector<dof_id_type> from_ordered_objects_ids;
237 14241 : dataLoad(stream, from_ordered_objects_ids, nullptr);
238 14241 : if (_loaded_ordered_objects.size() != from_ordered_objects_ids.size())
239 0 : mooseError("RestartableEquationSystems::load(): Number of previously stored elements/nodes (",
240 0 : _loaded_ordered_objects.size(),
241 : ") does not "
242 : "match the current number of elements/nodes (",
243 0 : from_ordered_objects_ids.size(),
244 : ")");
245 3667448 : for (const auto i : index_range(_loaded_ordered_objects))
246 3653207 : if (_loaded_ordered_objects[i]->id() != from_ordered_objects_ids[i])
247 0 : mooseError("RestartableEquationSystems::load(): Id of previously stored element/node (",
248 0 : _loaded_ordered_objects[i]->id(),
249 : ") does not "
250 : "match the current element/node id (",
251 0 : from_ordered_objects_ids[i],
252 : ")");
253 14241 : }
254 :
255 14241 : _loaded_stream_data_begin = static_cast<std::size_t>(stream.tellg());
256 :
257 : // Load everything that we have available at the moment
258 42840 : for (const auto & sys_name_header_pair : _loaded_header.systems)
259 : {
260 28599 : const auto & sys_header = sys_name_header_pair.second;
261 28599 : if (!_es.has_system(sys_header.name))
262 62 : continue;
263 28537 : auto & sys = _es.get_system(sys_header.name);
264 :
265 28537 : bool modified_sys = false;
266 :
267 148949 : for (const auto & vec_name_header_pair : sys_header.vectors)
268 : {
269 120412 : bool modified_vec = false;
270 :
271 120412 : const auto & vec_header = vec_name_header_pair.second;
272 120412 : const bool is_solution = vec_header.name == SystemHeader::system_solution_name;
273 :
274 120412 : if (!is_solution && !sys.have_vector(vec_header.name))
275 : {
276 102 : if (_load_all_vectors)
277 18 : sys.add_vector(vec_header.name, vec_header.projections, vec_header.type);
278 : else
279 84 : continue;
280 : }
281 :
282 120328 : auto & vec = is_solution ? *sys.solution : sys.get_vector(vec_header.name);
283 :
284 238918 : for (const auto & var_name_header_pair : sys_header.variables)
285 : {
286 118590 : const auto & var_header = var_name_header_pair.second;
287 118590 : if (!sys.has_variable(var_header.name))
288 14 : continue;
289 118576 : const auto & var = sys.variable(sys.variable_number(var_header.name));
290 118576 : if (var.type() != var_header.type)
291 0 : continue;
292 :
293 118576 : restore(sys_header, vec_header, var_header, sys, vec, var, stream);
294 118576 : modified_vec = true;
295 : }
296 :
297 120328 : if (modified_vec)
298 : {
299 96893 : vec.close();
300 96893 : modified_sys = true;
301 : }
302 : }
303 :
304 28537 : if (modified_sys)
305 21831 : sys.update();
306 : }
307 :
308 : // Move the stream to the end of our data so that we make RestartableDataReader happy
309 14241 : stream.seekg(_loaded_stream_data_begin + _loaded_header.data_size);
310 14241 : }
311 :
312 : void
313 118576 : RestartableEquationSystems::restore(const SystemHeader & from_sys_header,
314 : const VectorHeader & from_vec_header,
315 : const VariableHeader & from_var_header,
316 : const libMesh::System & to_sys,
317 : libMesh::NumericVector<libMesh::Number> & to_vec,
318 : const libMesh::Variable & to_var,
319 : std::istream & stream)
320 : {
321 : #ifndef NDEBUG
322 : const auto sys_it = _loaded_header.systems.find(from_sys_header.name);
323 : mooseAssert(sys_it != _loaded_header.systems.end(), "System does not exist");
324 : const auto & sys_header = sys_it->second;
325 : mooseAssert(sys_header == from_sys_header, "Not my system");
326 : const auto vec_it = sys_header.vectors.find(from_vec_header.name);
327 : mooseAssert(vec_it != sys_header.vectors.end(), "Vector does not exist");
328 : const auto & vec_header = vec_it->second;
329 : mooseAssert(vec_header == from_vec_header, "Not my vector");
330 : const auto var_it = sys_header.variables.find(from_var_header.name);
331 : mooseAssert(var_it != sys_header.variables.end(), "Variable does not exist");
332 : const auto & var_header = var_it->second;
333 : mooseAssert(var_header == from_var_header, "Not my variable");
334 : #endif
335 :
336 : const auto error =
337 0 : [&from_sys_header, &from_vec_header, &from_var_header, &to_sys, &to_var](auto... args)
338 : {
339 0 : mooseError("An error occured while restoring a system:\n",
340 : args...,
341 : "\n\nFrom system: ",
342 0 : from_sys_header.name,
343 : "\nFrom vector: ",
344 0 : from_vec_header.name,
345 : "\nFrom variable: ",
346 0 : from_var_header.name,
347 : "\nTo system: ",
348 0 : to_sys.name(),
349 : "\nTo variable: ",
350 0 : to_var.name());
351 118576 : };
352 :
353 118576 : if (from_var_header.type != to_var.type())
354 0 : error("Cannot restore to a variable of a different type");
355 :
356 118576 : const auto offset = from_vec_header.variable_offset.at(from_var_header.name);
357 118576 : stream.seekg(_loaded_stream_data_begin + offset);
358 :
359 : // Non-SCALAR variable
360 118576 : if (to_var.type().family != SCALAR)
361 : {
362 27280799 : for (const auto & obj : _loaded_ordered_objects)
363 41227670 : for (const auto comp : make_range(obj->n_comp(to_sys.number(), to_var.number())))
364 : {
365 : Real val;
366 14064752 : dataLoad(stream, val, nullptr);
367 14064752 : to_vec.set(obj->dof_number(to_sys.number(), to_var.number(), comp), val);
368 : }
369 : }
370 : // SCALAR variable on the last rank
371 695 : else if (_es.processor_id() == _es.n_processors() - 1)
372 : {
373 692 : const auto & dof_map = to_sys.get_dof_map();
374 692 : std::vector<dof_id_type> scalar_dofs;
375 692 : dof_map.SCALAR_dof_indices(scalar_dofs, to_var.number());
376 1627 : for (const auto dof : scalar_dofs)
377 : {
378 : Real val;
379 935 : dataLoad(stream, val, nullptr);
380 935 : to_vec.set(dof, val);
381 : }
382 692 : }
383 118576 : }
384 :
385 : void
386 45214 : dataStore(std::ostream & stream, RestartableEquationSystems & res, void *)
387 : {
388 45214 : res.store(stream);
389 45214 : }
390 :
391 : void
392 14241 : dataLoad(std::istream & stream, RestartableEquationSystems & res, void *)
393 : {
394 14241 : res.load(stream);
395 14241 : }
396 :
397 : void
398 45214 : dataStore(std::ostream & stream, RestartableEquationSystems::EquationSystemsHeader & header, void *)
399 : {
400 45214 : dataStore(stream, header.systems, nullptr);
401 45214 : dataStore(stream, header.data_size, nullptr);
402 45214 : }
403 :
404 : void
405 14241 : dataLoad(std::istream & stream, RestartableEquationSystems::EquationSystemsHeader & header, void *)
406 : {
407 14241 : dataLoad(stream, header.systems, nullptr);
408 14241 : dataLoad(stream, header.data_size, nullptr);
409 14241 : }
410 :
411 : void
412 91452 : dataStore(std::ostream & stream, RestartableEquationSystems::SystemHeader & header, void *)
413 : {
414 91452 : dataStore(stream, header.name, nullptr);
415 91452 : dataStore(stream, header.type, nullptr);
416 91452 : dataStore(stream, header.variables, nullptr);
417 91452 : dataStore(stream, header.vectors, nullptr);
418 91452 : }
419 :
420 : void
421 28599 : dataLoad(std::istream & stream, RestartableEquationSystems::SystemHeader & header, void *)
422 : {
423 28599 : dataLoad(stream, header.name, nullptr);
424 28599 : dataLoad(stream, header.type, nullptr);
425 28599 : dataLoad(stream, header.variables, nullptr);
426 28599 : dataLoad(stream, header.vectors, nullptr);
427 28599 : }
428 :
429 : void
430 91591 : dataStore(std::ostream & stream, RestartableEquationSystems::VariableHeader & header, void *)
431 : {
432 91591 : dataStore(stream, header.name, nullptr);
433 91591 : dataStore(stream, header.type, nullptr);
434 91591 : }
435 : void
436 29028 : dataLoad(std::istream & stream, RestartableEquationSystems::VariableHeader & header, void *)
437 : {
438 29028 : dataLoad(stream, header.name, nullptr);
439 29028 : dataLoad(stream, header.type, nullptr);
440 29028 : }
441 :
442 : void
443 430719 : dataStore(std::ostream & stream, RestartableEquationSystems::VectorHeader & header, void *)
444 : {
445 430719 : dataStore(stream, header.name, nullptr);
446 430719 : dataStore(stream, header.projections, nullptr);
447 430719 : dataStore(stream, header.type, nullptr);
448 430719 : dataStore(stream, header.variable_offset, nullptr);
449 430719 : }
450 : void
451 120540 : dataLoad(std::istream & stream, RestartableEquationSystems::VectorHeader & header, void *)
452 : {
453 120540 : dataLoad(stream, header.name, nullptr);
454 120540 : dataLoad(stream, header.projections, nullptr);
455 120540 : dataLoad(stream, header.type, nullptr);
456 120540 : dataLoad(stream, header.variable_offset, nullptr);
457 120540 : }
|
