https://mooseframework.inl.gov
DataIO.C
Go to the documentation of this file.
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 "DenseMatrix.h"
11 #include "DataIO.h"
12 #include "MooseMesh.h"
13 #include "FEProblemBase.h"
14 #include "NonlinearSystemBase.h"
15 
16 #include "libmesh/vector_value.h"
17 #include "libmesh/tensor_value.h"
18 
19 #include "libmesh/elem.h"
20 #include "libmesh/petsc_vector.h"
21 #include "libmesh/enum_solver_package.h"
22 #include "libmesh/petsc_solver_exception.h"
23 
24 using namespace libMesh;
25 
26 template <>
27 void
28 dataStore(std::ostream & stream, Real & v, void * /*context*/)
29 {
30  stream.write((char *)&v, sizeof(v));
31 }
32 
33 template <>
34 void
35 dataStore(std::ostream & stream, std::string & v, void * /*context*/)
36 {
37  // Write the size of the string
38  unsigned int size = v.size();
39  stream.write((char *)&size, sizeof(size));
40 
41  // Write the string (Do not store the null byte)
42  stream.write(v.c_str(), sizeof(char) * size);
43 }
44 
45 template <>
46 void
47 dataStore(std::ostream & stream, VariableName & v, void * context)
48 {
49  auto & name = static_cast<std::string &>(v);
50  dataStore(stream, name, context);
51 }
52 
53 template <>
54 void
55 dataStore(std::ostream & stream, UserObjectName & v, void * context)
56 {
57  auto & name = static_cast<std::string &>(v);
58  dataStore(stream, name, context);
59 }
60 
61 template <>
62 void
63 dataStore(std::ostream & stream, bool & v, void * /*context*/)
64 {
65  stream.write((char *)&v, sizeof(v));
66 }
67 
68 template <>
69 void
70 dataStore(std::ostream & stream, std::vector<bool> & v, void * context)
71 {
72  for (bool b : v)
73  dataStore(stream, b, context);
74 }
75 
76 template <>
77 void
78 dataStore(std::ostream & stream, RankTwoTensor & rtt, void * context)
79 {
80  dataStore(stream, rtt._coords, context);
81 }
82 
83 template <>
84 void
85 dataStore(std::ostream & stream, RankThreeTensor & rtht, void * context)
86 {
87  dataStore(stream, rtht._vals, context);
88 }
89 
90 template <>
91 void
92 dataStore(std::ostream & stream, RankFourTensor & rft, void * context)
93 {
94  dataStore(stream, rft._vals, context);
95 }
96 
97 template <>
98 void
99 dataStore(std::ostream & stream, ADReal & dn, void * context)
100 {
101  dataStore(stream, dn.value(), context);
102 
103  if (ADReal::do_derivatives)
104  {
105  auto & derivatives = dn.derivatives();
106  std::size_t size = derivatives.size();
107  dataStore(stream, size, context);
108  for (MooseIndex(size) i = 0; i < size; ++i)
109  {
110  dataStore(stream, derivatives.raw_index(i), context);
111  dataStore(stream, derivatives.raw_at(i), context);
112  }
113  }
114 }
115 
116 template <>
117 void
118 dataStore(std::ostream & stream, const Elem *& e, void * context)
119 {
120  // TODO: Write out the unique ID of this elem
122 
123  if (e)
124  {
125  id = e->id();
127  mooseError("Can't output Elems with invalid ids!");
128  }
129 
130  storeHelper(stream, id, context);
131 }
132 
133 template <>
134 void
135 dataStore(std::ostream & stream, const Node *& n, void * context)
136 {
137  // TODO: Write out the unique ID of this node
139 
140  if (n)
141  {
142  id = n->id();
144  mooseError("Can't output Nodes with invalid ids!");
145  }
146 
147  storeHelper(stream, id, context);
148 }
149 
150 template <>
151 void
152 dataStore(std::ostream & stream, Elem *& e, void * context)
153 {
154  // TODO: Write out the unique ID of this elem
156 
157  if (e)
158  {
159  id = e->id();
161  mooseError("Can't output Elems with invalid ids!");
162  }
163 
164  storeHelper(stream, id, context);
165 }
166 
167 template <>
168 void
169 dataStore(std::ostream & stream, Node *& n, void * context)
170 {
171  // TODO: Write out the unique ID of this node
173 
174  if (n)
175  {
176  id = n->id();
178  mooseError("Can't output Nodes with invalid ids!");
179  }
180 
181  storeHelper(stream, id, context);
182 }
183 
184 template <>
185 void
186 dataStore(std::ostream & stream, std::stringstream & s, void * /* context */)
187 {
188  const std::string & s_str = s.str();
189 
190  size_t s_size = s_str.size();
191  stream.write((char *)&s_size, sizeof(s_size));
192 
193  stream.write(s_str.c_str(), sizeof(char) * (s_str.size()));
194 }
195 
196 template <>
197 void
198 dataStore(std::ostream & stream, RealEigenVector & v, void * context)
199 {
200  unsigned int m = v.size();
201  stream.write((char *)&m, sizeof(m));
202  for (unsigned int i = 0; i < v.size(); i++)
203  {
204  Real r = v(i);
205  dataStore(stream, r, context);
206  }
207 }
208 
209 template <>
210 void
211 dataStore(std::ostream & stream, RealEigenMatrix & v, void * context)
212 {
213  unsigned int m = v.rows();
214  stream.write((char *)&m, sizeof(m));
215  unsigned int n = v.cols();
216  stream.write((char *)&n, sizeof(n));
217  for (unsigned int i = 0; i < m; i++)
218  for (unsigned int j = 0; j < n; j++)
219  {
220  Real r = v(i, j);
221  dataStore(stream, r, context);
222  }
223 }
224 
225 template <typename T>
226 void
227 dataStore(std::ostream & stream, TensorValue<T> & v, void * context)
228 {
229  for (const auto i : make_range(Moose::dim))
230  for (const auto j : make_range(Moose::dim))
231  {
232  T r = v(i, j);
233  dataStore(stream, r, context);
234  }
235 }
236 
237 template void dataStore(std::ostream & stream, TensorValue<Real> & v, void * context);
238 template void dataStore(std::ostream & stream, TensorValue<ADReal> & v, void * context);
239 
240 template <typename T>
241 void
242 dataStore(std::ostream & stream, DenseMatrix<T> & v, void * context)
243 {
244  unsigned int m = v.m();
245  unsigned int n = v.n();
246  stream.write((char *)&m, sizeof(m));
247  stream.write((char *)&n, sizeof(n));
248  for (unsigned int i = 0; i < m; i++)
249  for (unsigned int j = 0; j < n; j++)
250  {
251  T r = v(i, j);
252  dataStore(stream, r, context);
253  }
254 }
255 
256 template void dataStore(std::ostream & stream, DenseMatrix<Real> & v, void * context);
257 template void dataStore(std::ostream & stream, DenseMatrix<ADReal> & v, void * context);
258 
259 template <typename T>
260 void
261 dataStore(std::ostream & stream, VectorValue<T> & v, void * context)
262 {
263  // Obviously if someone loads data with different LIBMESH_DIM than was used for saving them, it
264  // won't work.
265  for (const auto i : make_range(Moose::dim))
266  {
267  T r = v(i);
268  dataStore(stream, r, context);
269  }
270 }
271 
272 template void dataStore(std::ostream & stream, VectorValue<Real> & v, void * context);
273 template void dataStore(std::ostream & stream, VectorValue<ADReal> & v, void * context);
274 
275 void
276 dataStore(std::ostream & stream, Point & p, void * context)
277 {
278  for (const auto i : make_range(Moose::dim))
279  {
280  Real r = p(i);
281  dataStore(stream, r, context);
282  }
283 }
284 
285 template <>
286 void
287 dataStore(std::ostream & stream, libMesh::Parameters & p, void * context)
288 {
289  // First store the size of the map
290  unsigned int size = p.n_parameters();
291  stream.write((char *)&size, sizeof(size));
292 
293  auto it = p.begin();
294  auto end = p.end();
295 
296  for (; it != end; ++it)
297  {
298  auto & key = const_cast<std::string &>(it->first);
299  auto type = it->second->type();
300 
301  storeHelper(stream, key, context);
302  storeHelper(stream, type, context);
303 
304 #define storescalar(ptype) \
305  else if (it->second->type() == demangle(typeid(ptype).name())) storeHelper( \
306  stream, \
307  (dynamic_cast<libMesh::Parameters::Parameter<ptype> *>(MooseUtils::get(it->second)))->get(), \
308  context)
309 
310  if (false)
311  ;
312  storescalar(Real);
313  storescalar(short);
314  storescalar(int);
315  storescalar(long);
316  storescalar(unsigned short);
317  storescalar(unsigned int);
318  storescalar(unsigned long);
319 
320 #undef storescalar
321  }
322 }
323 
324 template <>
325 void
326 dataStore(std::ostream & stream,
327  std::unique_ptr<libMesh::NumericVector<Number>> & v,
328  void * context)
329 {
330  // Classes may declare unique pointers to vectors as restartable data and never actually create
331  // vector instances. This happens for example in the `TimeIntegrator` class where subvector
332  // instances are only created if multiple time integrators are present
333  bool have_vector = v.get();
334  dataStore(stream, have_vector, context);
335  if (!have_vector)
336  return;
337 
338  mooseAssert(context, "Needs a context of the communicator");
339  const auto & comm = *static_cast<const libMesh::Parallel::Communicator *>(context);
340  mooseAssert(&comm == &v->comm(), "Inconsistent communicator");
341 
342  if (v->type() == GHOSTED)
343  mooseError("Cannot store ghosted numeric vectors");
344 
345  // Store the communicator size for sanity checking later
346  unsigned int comm_size = comm.size();
347  dataStore(stream, comm_size, nullptr);
348 
349  // Store the solver package so that we know what vector type to construct
350  libMesh::SolverPackage solver_package;
351  if (dynamic_cast<libMesh::PetscVector<Number> *>(v.get()))
352  solver_package = PETSC_SOLVERS;
353  else
354  mooseError("Can only store unique_ptrs of PetscVectors");
355  int solver_package_int = solver_package;
356  dataStore(stream, solver_package_int, nullptr);
357 
358  // Store the sizes
359  dof_id_type size = v->size();
360  dataStore(stream, size, nullptr);
361  dof_id_type local_size = v->local_size();
362  dataStore(stream, local_size, nullptr);
363 
364  // Store the vector itself
365  dataStore(stream, *v, nullptr);
366 }
367 
368 // global load functions
369 
370 template <>
371 void
372 dataLoad(std::istream & stream, Real & v, void * /*context*/)
373 {
374  stream.read((char *)&v, sizeof(v));
375 }
376 
377 template <>
378 void
379 dataLoad(std::istream & stream, std::string & v, void * /*context*/)
380 {
381  // Read the size of the string
382  unsigned int size = 0;
383  stream.read((char *)&size, sizeof(size));
384 
385  // Resize the string data
386  v.resize(size);
387 
388  // Read the string
389  stream.read(&v[0], sizeof(char) * size);
390 }
391 
392 template <>
393 void
394 dataLoad(std::istream & stream, VariableName & v, void * context)
395 {
396  auto & name = static_cast<std::string &>(v);
397  dataLoad(stream, name, context);
398 }
399 
400 template <>
401 void
402 dataLoad(std::istream & stream, UserObjectName & v, void * context)
403 {
404  auto & name = static_cast<std::string &>(v);
405  dataLoad(stream, name, context);
406 }
407 
408 template <>
409 void
410 dataLoad(std::istream & stream, bool & v, void * /*context*/)
411 {
412  stream.read((char *)&v, sizeof(v));
413 }
414 
415 template <>
416 void
417 dataLoad(std::istream & stream, std::vector<bool> & v, void * context)
418 {
419  for (bool b : v)
420  dataLoad(stream, b, context);
421 }
422 
423 template <>
424 void
425 dataLoad(std::istream & stream, ADReal & dn, void * context)
426 {
427  dataLoad(stream, dn.value(), context);
428 
429  if (ADReal::do_derivatives)
430  {
431  auto & derivatives = dn.derivatives();
432  std::size_t size = 0;
433  stream.read((char *)&size, sizeof(size));
434  derivatives.resize(size);
435 
436  for (MooseIndex(derivatives) i = 0; i < derivatives.size(); ++i)
437  {
438  dataLoad(stream, derivatives.raw_index(i), context);
439  dataLoad(stream, derivatives.raw_at(i), context);
440  }
441  }
442 }
443 
444 template <>
445 void
446 dataLoad(std::istream & stream, const Elem *& e, void * context)
447 {
448  if (!context)
449  mooseError("Can only load Elem objects using a MooseMesh context!");
450 
451  MooseMesh * mesh = static_cast<MooseMesh *>(context);
452 
453  // TODO: Write out the unique ID of this element
455 
456  loadHelper(stream, id, context);
457 
459  e = mesh->elemPtr(id);
460  else
461  e = NULL;
462 }
463 
464 template <>
465 void
466 dataLoad(std::istream & stream, const Node *& n, void * context)
467 {
468  if (!context)
469  mooseError("Can only load Node objects using a MooseMesh context!");
470 
471  MooseMesh * mesh = static_cast<MooseMesh *>(context);
472 
473  // TODO: Write out the unique ID of this nodeent
475 
476  loadHelper(stream, id, context);
477 
479  n = mesh->nodePtr(id);
480  else
481  n = NULL;
482 }
483 
484 template <>
485 void
486 dataLoad(std::istream & stream, Elem *& e, void * context)
487 {
488  if (!context)
489  mooseError("Can only load Elem objects using a MooseMesh context!");
490 
491  MooseMesh * mesh = static_cast<MooseMesh *>(context);
492 
493  // TODO: Write out the unique ID of this element
495 
496  loadHelper(stream, id, context);
497 
499  e = mesh->elemPtr(id);
500  else
501  e = NULL;
502 }
503 
504 template <>
505 void
506 dataLoad(std::istream & stream, Node *& n, void * context)
507 {
508  if (!context)
509  mooseError("Can only load Node objects using a MooseMesh context!");
510 
511  MooseMesh * mesh = static_cast<MooseMesh *>(context);
512 
513  // TODO: Write out the unique ID of this nodeent
515 
516  loadHelper(stream, id, context);
517 
519  n = mesh->nodePtr(id);
520  else
521  n = NULL;
522 }
523 
524 template <>
525 void
526 dataLoad(std::istream & stream, std::stringstream & s, void * /* context */)
527 {
528  size_t s_size = 0;
529  stream.read((char *)&s_size, sizeof(s_size));
530 
531  std::unique_ptr<char[]> s_s = std::make_unique<char[]>(s_size);
532  stream.read(s_s.get(), s_size);
533 
534  // Clear the stringstream before loading new data into it.
535  s.str(std::string());
536  s.write(s_s.get(), s_size);
537 }
538 
539 template <>
540 void
541 dataLoad(std::istream & stream, RealEigenVector & v, void * context)
542 {
543  unsigned int n = 0;
544  stream.read((char *)&n, sizeof(n));
545  v.resize(n);
546  for (unsigned int i = 0; i < n; i++)
547  {
548  Real r = 0;
549  dataLoad(stream, r, context);
550  v(i) = r;
551  }
552 }
553 
554 template <>
555 void
556 dataLoad(std::istream & stream, RealEigenMatrix & v, void * context)
557 {
558  unsigned int m = 0;
559  stream.read((char *)&m, sizeof(m));
560  unsigned int n = 0;
561  stream.read((char *)&n, sizeof(n));
562  v.resize(m, n);
563  for (unsigned int i = 0; i < m; i++)
564  for (unsigned int j = 0; j < n; j++)
565  {
566  Real r = 0;
567  dataLoad(stream, r, context);
568  v(i, j) = r;
569  }
570 }
571 
572 template <typename T>
573 void
574 dataLoad(std::istream & stream, TensorValue<T> & v, void * context)
575 {
576  // Obviously if someone loads data with different LIBMESH_DIM than was used for saving them, it
577  // won't work.
578  for (const auto i : make_range(Moose::dim))
579  for (const auto j : make_range(Moose::dim))
580  {
581  T r = 0;
582  dataLoad(stream, r, context);
583  v(i, j) = r;
584  }
585 }
586 
587 template void dataLoad(std::istream & stream, TensorValue<Real> & v, void * context);
588 template void dataLoad(std::istream & stream, TensorValue<ADReal> & v, void * context);
589 
590 template <typename T>
591 void
592 dataLoad(std::istream & stream, DenseMatrix<T> & v, void * context)
593 {
594  unsigned int m = 0, n = 0;
595  stream.read((char *)&m, sizeof(m));
596  stream.read((char *)&n, sizeof(n));
597  v.resize(m, n);
598  for (unsigned int i = 0; i < m; i++)
599  for (unsigned int j = 0; j < n; j++)
600  {
601  T r = 0;
602  dataLoad(stream, r, context);
603  v(i, j) = r;
604  }
605 }
606 
607 template void dataLoad(std::istream & stream, DenseMatrix<Real> & v, void * context);
608 template void dataLoad(std::istream & stream, DenseMatrix<ADReal> & v, void * context);
609 
610 template <typename T>
611 void
612 dataLoad(std::istream & stream, VectorValue<T> & v, void * context)
613 {
614  // Obviously if someone loads data with different LIBMESH_DIM than was used for saving them, it
615  // won't work.
616  for (const auto i : make_range(Moose::dim))
617  {
618  T r = 0;
619  dataLoad(stream, r, context);
620  v(i) = r;
621  }
622 }
623 
624 template void dataLoad(std::istream & stream, VectorValue<Real> & v, void * context);
625 template void dataLoad(std::istream & stream, VectorValue<ADReal> & v, void * context);
626 
627 void
628 dataLoad(std::istream & stream, Point & p, void * context)
629 {
630  for (const auto i : make_range(Moose::dim))
631  {
632  Real r = 0;
633  dataLoad(stream, r, context);
634  p(i) = r;
635  }
636 }
637 
638 template <>
639 void
640 dataLoad(std::istream & stream, libMesh::Parameters & p, void * context)
641 {
642  p.clear();
643 
644  // First read the size of the map
645  unsigned int size = 0;
646  stream.read((char *)&size, sizeof(size));
647 
648  for (unsigned int i = 0; i < size; i++)
649  {
650  std::string key, type;
651  loadHelper(stream, key, context);
652  loadHelper(stream, type, context);
653 
654 #define loadscalar(ptype) \
655  else if (type == demangle(typeid(ptype).name())) do \
656  { \
657  ptype & value = p.set<ptype>(key); \
658  loadHelper(stream, value, context); \
659  } \
660  while (0)
661 
662  if (false)
663  ;
664  loadscalar(Real);
665  loadscalar(short);
666  loadscalar(int);
667  loadscalar(long);
668  loadscalar(unsigned short);
669  loadscalar(unsigned int);
670  loadscalar(unsigned long);
671 
672 #undef loadscalar
673  }
674 }
675 
676 template <>
677 void
678 dataLoad(std::istream & stream, std::unique_ptr<libMesh::NumericVector<Number>> & v, void * context)
679 {
680  bool have_vector;
681  dataLoad(stream, have_vector, context);
682 
683  if (!have_vector)
684  return;
685 
686  mooseAssert(context, "Needs a context of the communicator");
687  const auto & comm = *static_cast<const libMesh::Parallel::Communicator *>(context);
688  if (v)
689  mooseAssert(&comm == &v->comm(), "Inconsistent communicator");
690 
691  // Load the communicator size for consistency checks
692  unsigned int comm_size;
693  dataLoad(stream, comm_size, nullptr);
694  mooseAssert(comm.size() == comm_size, "Inconsistent communicator size");
695 
696  // Load the solver package to build the vector
697  int solver_package_int;
698  dataLoad(stream, solver_package_int, nullptr);
699  libMesh::SolverPackage solver_package = static_cast<libMesh::SolverPackage>(solver_package_int);
700 
701  // Load the sizes
702  dof_id_type size, local_size;
703  dataLoad(stream, size, nullptr);
704  dataLoad(stream, local_size, nullptr);
705 
706  // Construct the vector given the type, only if we need to. v could be non-null here
707  // if we're advancing back and loading a backup
708  if (!v)
709  {
710  v = NumericVector<Number>::build(comm, solver_package);
711  v->init(size, local_size);
712  }
713  else
714  mooseAssert(v->type() != GHOSTED, "Cannot be ghosted");
715 
716  // Make sure that the sizes are consistent; this will happen if we're calling this
717  // on a vector that has already been loaded previously
718  mooseAssert(v->size() == size, "Inconsistent size");
719  mooseAssert(v->local_size() == local_size, "Inconsistent local size");
720 
721  // Now that we have an initialized vector, fill the entries
722  dataLoad(stream, *v, nullptr);
723 }
724 
725 template <>
726 void
727 dataLoad(std::istream & stream, Vec & v, void * context)
728 {
729  PetscInt local_size;
730  LibmeshPetscCallA(PETSC_COMM_WORLD, VecGetLocalSize(v, &local_size));
731  PetscScalar * array;
732  LibmeshPetscCallA(PETSC_COMM_WORLD, VecGetArray(v, &array));
733  for (PetscInt i = 0; i < local_size; i++)
734  dataLoad(stream, array[i], context);
735 
736  LibmeshPetscCallA(PETSC_COMM_WORLD, VecRestoreArray(v, &array));
737 }
738 
739 template <>
740 void
741 dataStore(std::ostream & stream, Vec & v, void * context)
742 {
743  PetscInt local_size;
744  LibmeshPetscCallA(PETSC_COMM_WORLD, VecGetLocalSize(v, &local_size));
745  PetscScalar * array;
746  LibmeshPetscCallA(PETSC_COMM_WORLD, VecGetArray(v, &array));
747  for (PetscInt i = 0; i < local_size; i++)
748  dataStore(stream, array[i], context);
749 
750  LibmeshPetscCallA(PETSC_COMM_WORLD, VecRestoreArray(v, &array));
751 }
std::string name(const ElemQuality q)
RankFourTensorTempl is designed to handle any N-dimensional fourth order tensor, C.
void mooseError(Args &&... args)
Emit an error message with the given stringified, concatenated args and terminate the application...
Definition: MooseError.h:302
MeshBase & mesh
static constexpr std::size_t dim
This is the dimension of all vector and tensor datastructures used in MOOSE.
Definition: Moose.h:154
unsigned int m() const
The following methods are specializations for using the libMesh::Parallel::packed_range_* routines fo...
T _coords[LIBMESH_DIM *LIBMESH_DIM]
DualNumber< Real, DNDerivativeType, true > ADReal
Definition: ADRealForward.h:46
RankThreeTensor is designed to handle any N-dimensional third order tensor, r.
void storeHelper(std::ostream &stream, P &data, void *context)
Scalar helper routine.
Definition: DataIO.h:893
dof_id_type id() const
void dataLoad(std::istream &stream, Real &v, void *)
Definition: DataIO.C:372
MooseMesh wraps a libMesh::Mesh object and enhances its capabilities by caching additional data and s...
Definition: MooseMesh.h:88
static const dof_id_type invalid_id
T _vals[N3]
The values of the rank-three tensor stored by index=((i * LIBMESH_DIM + j) * LIBMESH_DIM + k) ...
T _vals[N4]
The values of the rank-four tensor stored by index=(((i * LIBMESH_DIM + j) * LIBMESH_DIM + k) * LIBME...
Eigen::Matrix< Real, Eigen::Dynamic, Eigen::Dynamic > RealEigenMatrix
Definition: MooseTypes.h:149
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real
virtual void clear()
void resize(const unsigned int new_m, const unsigned int new_n)
IntRange< T > make_range(T beg, T end)
static std::unique_ptr< NumericVector< T > > build(const Parallel::Communicator &comm, SolverPackage solver_package=libMesh::default_solver_package(), ParallelType parallel_type=AUTOMATIC)
Eigen::Matrix< Real, Eigen::Dynamic, 1 > RealEigenVector
Definition: MooseTypes.h:146
std::size_t n_parameters() const
unsigned int n() const
void loadHelper(std::istream &stream, P &data, void *context)
Scalar helper routine.
Definition: DataIO.h:985
void dataStore(std::ostream &stream, Real &v, void *)
Definition: DataIO.C:28
uint8_t dof_id_type