libMesh
Classes | Functions
libMesh::RBDataDeserialization Namespace Reference

Classes

class  RBEIMEvaluationDeserialization
 This class de-serializes a RBEIMEvaluation object using the Cap'n Proto library. More...
 
class  RBEvaluationDeserialization
 This class de-serializes an RBEvaluation object using the Cap'n Proto library. More...
 
class  RBSCMEvaluationDeserialization
 This class de-serializes a RBSCMEvaluation object using the Cap'n Proto library. More...
 
class  TransientRBEvaluationDeserialization
 This class de-serializes a TransientRBEvaluation object using the Cap'n Proto library. More...
 

Functions

void load_parameter_ranges (RBParametrized &rb_evaluation, RBData::ParameterRanges::Reader &parameter_ranges, RBData::DiscreteParameterList::Reader &discrete_parameters_list)
 Load parameter ranges and discrete parameter values into an RBEvaluation from the corresponding structure in the buffer. More...
 
template<typename RBEvaluationReaderNumber >
void load_rb_evaluation_data (RBEvaluation &rb_evaluation, RBEvaluationReaderNumber &rb_evaluation_reader, bool read_error_bound_data)
 Load an RB evaluation from a corresponding reader structure in the buffer. More...
 
template<typename RBEvaluationReaderNumber , typename TransRBEvaluationReaderNumber >
void load_transient_rb_evaluation_data (TransientRBEvaluation &trans_rb_eval, RBEvaluationReaderNumber &rb_evaluation_reader, TransRBEvaluationReaderNumber &trans_rb_eval_reader, bool read_error_bound_data)
 Load an RB evaluation from a corresponding reader structure in the buffer. More...
 
template<typename RBEvaluationReaderNumber , typename RBEIMEvaluationReaderNumber >
void load_rb_eim_evaluation_data (RBEIMEvaluation &rb_eim_eval, RBEvaluationReaderNumber &rb_evaluation_reader, RBEIMEvaluationReaderNumber &rb_eim_eval_reader)
 Load an EIM RB evaluation from a corresponding reader structure in the buffer. More...
 
void load_rb_scm_evaluation_data (RBSCMEvaluation &rb_scm_eval, RBData::RBSCMEvaluation::Reader &rb_scm_eval_reader)
 Load an SCM RB evaluation from a corresponding reader structure in the buffer. More...
 
void load_point (RBData::Point3D::Reader point_reader, Point &point)
 Helper function that loads point data. More...
 
void load_elem_into_mesh (RBData::MeshElem::Reader mesh_elem_reader, libMesh::Elem *elem, libMesh::ReplicatedMesh &mesh)
 Helper function that loads element data. More...
 

Function Documentation

void libMesh::RBDataDeserialization::load_elem_into_mesh ( RBData::MeshElem::Reader  mesh_elem_reader,
libMesh::Elem elem,
libMesh::ReplicatedMesh mesh 
)

Helper function that loads element data.

Definition at line 823 of file rb_data_deserialization.C.

References libMesh::ReplicatedMesh::add_elem(), libMesh::ReplicatedMesh::add_node(), libMesh::Elem::n_nodes(), libMesh::Elem::set_node(), and libMesh::Elem::subdomain_id().

Referenced by load_rb_eim_evaluation_data().

826 {
827  auto mesh_elem_point_list = mesh_elem_reader.getPoints();
828  unsigned int n_points = mesh_elem_point_list.size();
829 
830  if (n_points != elem->n_nodes())
831  libmesh_error_msg("Wrong number of nodes for element type");
832 
833  for (unsigned int i=0; i < n_points; ++i)
834  {
835  libMesh::Node * node = new libMesh::Node(mesh_elem_point_list[i].getX(),
836  mesh_elem_point_list[i].getY(),
837  mesh_elem_point_list[i].getZ());
838 
839  mesh.add_node(node);
840 
841  elem->set_node(i) = node;
842  }
843 
844  elem->subdomain_id() = mesh_elem_reader.getSubdomainId();
845 
846  mesh.add_elem(elem);
847 }
virtual Node *& set_node(const unsigned int i)
Definition: elem.h:1941
A Node is like a Point, but with more information.
Definition: node.h:52
virtual Elem * add_elem(Elem *e) libmesh_override
Add elem e to the end of the element array.
virtual Node * add_node(Node *n) libmesh_override
Add Node n to the end of the vertex array.
virtual unsigned int n_nodes() const =0
subdomain_id_type subdomain_id() const
Definition: elem.h:1951
void libMesh::RBDataDeserialization::load_parameter_ranges ( RBParametrized rb_evaluation,
RBData::ParameterRanges::Reader &  parameter_ranges,
RBData::DiscreteParameterList::Reader &  discrete_parameters_list 
)

Load parameter ranges and discrete parameter values into an RBEvaluation from the corresponding structure in the buffer.

Definition at line 239 of file rb_data_deserialization.C.

References libMesh::RBParametrized::initialize_parameters(), libMesh::Real, and libMesh::RBParameters::set_value().

Referenced by load_rb_evaluation_data(), and load_rb_scm_evaluation_data().

242 {
243  // Continuous parameters
244  RBParameters parameters_min;
245  RBParameters parameters_max;
246  {
247  unsigned int n_parameter_ranges = parameter_ranges.getNames().size();
248 
249  for (unsigned int i=0; i<n_parameter_ranges; ++i)
250  {
251  std::string parameter_name = parameter_ranges.getNames()[i];
252  Real min_value = parameter_ranges.getMinValues()[i];
253  Real max_value = parameter_ranges.getMaxValues()[i];
254 
255  parameters_min.set_value(parameter_name, min_value);
256  parameters_max.set_value(parameter_name, max_value);
257  }
258  }
259 
260  // Discrete parameters
261  std::map<std::string, std::vector<Real>> discrete_parameter_values;
262  {
263  unsigned int n_discrete_parameters = discrete_parameters_list.getNames().size();
264 
265  for (unsigned int i=0; i<n_discrete_parameters; ++i)
266  {
267  std::string parameter_name = discrete_parameters_list.getNames()[i];
268 
269  auto value_list = discrete_parameters_list.getValues()[i];
270  unsigned int n_values = value_list.size();
271  std::vector<Real> values(n_values);
272  for (unsigned int j=0; j<n_values; ++j)
273  values[j] = value_list[j];
274 
275  discrete_parameter_values[parameter_name] = values;
276  }
277  }
278 
279  rb_evaluation.initialize_parameters(parameters_min,
280  parameters_max,
281  discrete_parameter_values);
282 }
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real
void libMesh::RBDataDeserialization::load_point ( RBData::Point3D::Reader  point_reader,
Point point 
)

Helper function that loads point data.

Definition at line 811 of file rb_data_deserialization.C.

Referenced by load_rb_eim_evaluation_data().

812 {
813  point(0) = point_reader.getX();
814 
815  if (LIBMESH_DIM >= 2)
816  point(1) = point_reader.getY();
817 
818  if (LIBMESH_DIM >= 3)
819  point(2) = point_reader.getZ();
820 }
template<typename RBEvaluationReaderNumber , typename RBEIMEvaluationReaderNumber >
void libMesh::RBDataDeserialization::load_rb_eim_evaluation_data ( RBEIMEvaluation rb_eim_eval,
RBEvaluationReaderNumber &  rb_evaluation_reader,
RBEIMEvaluationReaderNumber &  rb_eim_eval_reader 
)

Load an EIM RB evaluation from a corresponding reader structure in the buffer.

Templated to deal with both Real and Complex numbers.

Definition at line 627 of file rb_data_deserialization.C.

References libMesh::Elem::build(), libMesh::ReplicatedMesh::clear(), libMesh::RBEIMEvaluation::get_interpolation_points_mesh(), libMesh::RBEvaluation::get_n_basis_functions(), libMesh::RBEIMEvaluation::interpolation_matrix, libMesh::RBEIMEvaluation::interpolation_points, libMesh::RBEIMEvaluation::interpolation_points_elem, libMesh::RBEIMEvaluation::interpolation_points_var, load_elem_into_mesh(), load_point(), and load_rb_evaluation_data().

Referenced by libMesh::RBDataDeserialization::RBEIMEvaluationDeserialization::read_from_file().

630 {
631  // We use read_error_bound_data=false here, since the RBEvaluation error bound data
632  // is not relevant to EIM.
633  load_rb_evaluation_data(rb_eim_evaluation,
634  rb_evaluation_reader,
635  /*read_error_bound_data*/ false);
636 
637  unsigned int n_bfs = rb_eim_evaluation.get_n_basis_functions();
638 
639  // EIM interpolation matrix
640  {
641  auto interpolation_matrix_list = rb_eim_evaluation_reader.getInterpolationMatrix();
642 
643  if (interpolation_matrix_list.size() != n_bfs*(n_bfs+1)/2)
644  libmesh_error_msg("Size error while reading the eim inner product matrix.");
645 
646  for (unsigned int i=0; i<n_bfs; ++i)
647  for (unsigned int j=0; j<=i; ++j)
648  {
649  unsigned int offset = i*(i+1)/2 + j;
650  rb_eim_evaluation.interpolation_matrix(i,j) =
651  load_scalar_value(interpolation_matrix_list[offset]);
652  }
653  }
654 
655  // Interpolation points
656  {
657  auto interpolation_points_list =
658  rb_eim_evaluation_reader.getInterpolationPoints();
659 
660  if (interpolation_points_list.size() != n_bfs)
661  libmesh_error_msg("Size error while reading the eim interpolation points.");
662 
663  rb_eim_evaluation.interpolation_points.resize(n_bfs);
664  for (unsigned int i=0; i<n_bfs; ++i)
665  {
666  load_point(interpolation_points_list[i],
667  rb_eim_evaluation.interpolation_points[i]);
668  }
669  }
670 
671  // Interpolation points variables
672  {
673  auto interpolation_points_var_list =
674  rb_eim_evaluation_reader.getInterpolationPointsVar();
675  rb_eim_evaluation.interpolation_points_var.resize(n_bfs);
676 
677  if (interpolation_points_var_list.size() != n_bfs)
678  libmesh_error_msg("Size error while reading the eim interpolation variables.");
679 
680  for (unsigned int i=0; i<n_bfs; ++i)
681  {
682  rb_eim_evaluation.interpolation_points_var[i] =
683  interpolation_points_var_list[i];
684  }
685  }
686 
687  // Interpolation elements
688  {
689  libMesh::dof_id_type elem_id = 0;
690  libMesh::ReplicatedMesh & interpolation_points_mesh =
691  rb_eim_evaluation.get_interpolation_points_mesh();
692  interpolation_points_mesh.clear();
693 
694  auto interpolation_points_elem_list =
695  rb_eim_evaluation_reader.getInterpolationPointsElems();
696  unsigned int n_interpolation_elems = interpolation_points_elem_list.size();
697 
698  if (n_interpolation_elems != n_bfs)
699  libmesh_error_msg("The number of elements should match the number of basis functions");
700 
701  rb_eim_evaluation.interpolation_points_elem.resize(n_interpolation_elems);
702 
703  for (unsigned int i=0; i<n_interpolation_elems; ++i)
704  {
705  auto mesh_elem_reader = interpolation_points_elem_list[i];
706  std::string elem_type_name = mesh_elem_reader.getType().cStr();
707  libMesh::ElemType elem_type =
708  libMesh::Utility::string_to_enum<libMesh::ElemType>(elem_type_name);
709 
710  libMesh::Elem * elem = libMesh::Elem::build(elem_type).release();
711  elem->set_id(elem_id++);
712  load_elem_into_mesh(mesh_elem_reader, elem, interpolation_points_mesh);
713 
714  rb_eim_evaluation.interpolation_points_elem[i] = elem;
715  }
716  }
717 }
The ReplicatedMesh class is derived from the MeshBase class, and is used to store identical copies of...
static UniquePtr< Elem > build(const ElemType type, Elem *p=libmesh_nullptr)
Definition: elem.C:238
void load_rb_evaluation_data(RBEvaluation &rb_evaluation, RBEvaluationReaderNumber &rb_evaluation_reader, bool read_error_bound_data)
Load an RB evaluation from a corresponding reader structure in the buffer.
ElemType
Defines an enum for geometric element types.
virtual void clear() libmesh_override
Clear all internal data.
This is the base class from which all geometric element types are derived.
Definition: elem.h:89
void load_point(RBData::Point3D::Reader point_reader, Point &point)
Helper function that loads point data.
void load_elem_into_mesh(RBData::MeshElem::Reader mesh_elem_reader, libMesh::Elem *elem, libMesh::ReplicatedMesh &mesh)
Helper function that loads element data.
uint8_t dof_id_type
Definition: id_types.h:64
template<typename RBEvaluationReaderNumber >
void libMesh::RBDataDeserialization::load_rb_evaluation_data ( RBEvaluation rb_evaluation,
RBEvaluationReaderNumber &  rb_evaluation_reader,
bool  read_error_bound_data 
)

Load an RB evaluation from a corresponding reader structure in the buffer.

Definition at line 285 of file rb_data_deserialization.C.

References libMesh::RBEvaluation::Aq_Aq_representor_innerprods, libMesh::RBEvaluation::compute_RB_inner_product, libMesh::RBEvaluation::Fq_Aq_representor_innerprods, libMesh::RBEvaluation::Fq_representor_innerprods, libMesh::RBThetaExpansion::get_n_A_terms(), libMesh::RBThetaExpansion::get_n_F_terms(), libMesh::RBThetaExpansion::get_n_output_terms(), libMesh::RBThetaExpansion::get_n_outputs(), libMesh::RBEvaluation::get_rb_theta_expansion(), load_parameter_ranges(), libMesh::RBEvaluation::output_dual_innerprods, libMesh::RBEvaluation::RB_Aq_vector, libMesh::RBEvaluation::RB_Fq_vector, libMesh::RBEvaluation::RB_inner_product_matrix, libMesh::RBEvaluation::RB_output_vectors, libMesh::RBEvaluation::resize_data_structures(), and libMesh::RBEvaluation::set_n_basis_functions().

Referenced by load_rb_eim_evaluation_data(), load_transient_rb_evaluation_data(), and libMesh::RBDataDeserialization::RBEvaluationDeserialization::read_from_file().

288 {
289  // Set number of basis functions
290  unsigned int n_bfs = rb_evaluation_reader.getNBfs();
291  rb_evaluation.set_n_basis_functions(n_bfs);
292 
293  rb_evaluation.resize_data_structures(n_bfs, read_error_bound_data);
294 
295  auto parameter_ranges =
296  rb_evaluation_reader.getParameterRanges();
297  auto discrete_parameters_list =
298  rb_evaluation_reader.getDiscreteParameters();
299 
300  load_parameter_ranges(rb_evaluation,
301  parameter_ranges,
302  discrete_parameters_list);
303 
304  const RBThetaExpansion & rb_theta_expansion = rb_evaluation.get_rb_theta_expansion();
305 
306  unsigned int n_F_terms = rb_theta_expansion.get_n_F_terms();
307  unsigned int n_A_terms = rb_theta_expansion.get_n_A_terms();
308 
309  if (read_error_bound_data)
310  {
311 
312  // Fq representor inner-product data
313  {
314  unsigned int Q_f_hat = n_F_terms*(n_F_terms+1)/2;
315 
316  auto fq_innerprods_list = rb_evaluation_reader.getFqInnerprods();
317  if (fq_innerprods_list.size() != Q_f_hat)
318  libmesh_error_msg("Size error while reading Fq representor norm data from buffer.");
319 
320  for (unsigned int i=0; i < Q_f_hat; ++i)
321  rb_evaluation.Fq_representor_innerprods[i] = load_scalar_value(fq_innerprods_list[i]);
322  }
323 
324  // Fq_Aq representor inner-product data
325  {
326  auto fq_aq_innerprods_list = rb_evaluation_reader.getFqAqInnerprods();
327  if (fq_aq_innerprods_list.size() != n_F_terms*n_A_terms*n_bfs)
328  libmesh_error_msg("Size error while reading Fq-Aq representor norm data from buffer.");
329 
330  for (unsigned int q_f=0; q_f<n_F_terms; ++q_f)
331  for (unsigned int q_a=0; q_a<n_A_terms; ++q_a)
332  for (unsigned int i=0; i<n_bfs; ++i)
333  {
334  unsigned int offset = q_f*n_A_terms*n_bfs + q_a*n_bfs + i;
335  rb_evaluation.Fq_Aq_representor_innerprods[q_f][q_a][i] =
336  load_scalar_value(fq_aq_innerprods_list[offset]);
337  }
338  }
339 
340  // Aq_Aq representor inner-product data
341  {
342  unsigned int Q_a_hat = n_A_terms*(n_A_terms+1)/2;
343  auto aq_aq_innerprods_list = rb_evaluation_reader.getAqAqInnerprods();
344  if (aq_aq_innerprods_list.size() != Q_a_hat*n_bfs*n_bfs)
345  libmesh_error_msg("Size error while reading Aq-Aq representor norm data from buffer.");
346 
347  for (unsigned int i=0; i<Q_a_hat; ++i)
348  for (unsigned int j=0; j<n_bfs; ++j)
349  for (unsigned int l=0; l<n_bfs; ++l)
350  {
351  unsigned int offset = i*n_bfs*n_bfs + j*n_bfs + l;
352  rb_evaluation.Aq_Aq_representor_innerprods[i][j][l] =
353  load_scalar_value(aq_aq_innerprods_list[offset]);
354  }
355  }
356 
357  // Output dual inner-product data
358  {
359  unsigned int n_outputs = rb_theta_expansion.get_n_outputs();
360  auto output_innerprod_outer = rb_evaluation_reader.getOutputDualInnerprods();
361 
362  if (output_innerprod_outer.size() != n_outputs)
363  libmesh_error_msg("Incorrect number of outputs detected in the buffer");
364 
365  for (unsigned int output_id=0; output_id<n_outputs; ++output_id)
366  {
367  unsigned int n_output_terms = rb_theta_expansion.get_n_output_terms(output_id);
368 
369  unsigned int Q_l_hat = n_output_terms*(n_output_terms+1)/2;
370  auto output_innerprod_inner = output_innerprod_outer[output_id];
371 
372  if (output_innerprod_inner.size() != Q_l_hat)
373  libmesh_error_msg("Incorrect number of output terms detected in the buffer");
374 
375  for (unsigned int q=0; q<Q_l_hat; ++q)
376  {
377  rb_evaluation.output_dual_innerprods[output_id][q] =
378  load_scalar_value(output_innerprod_inner[q]);
379  }
380  }
381  }
382  }
383 
384  // Output vectors
385  {
386  unsigned int n_outputs = rb_theta_expansion.get_n_outputs();
387  auto output_vector_outer = rb_evaluation_reader.getOutputVectors();
388 
389  if (output_vector_outer.size() != n_outputs)
390  libmesh_error_msg("Incorrect number of outputs detected in the buffer");
391 
392  for (unsigned int output_id=0; output_id<n_outputs; ++output_id)
393  {
394  unsigned int n_output_terms = rb_theta_expansion.get_n_output_terms(output_id);
395 
396  auto output_vector_middle = output_vector_outer[output_id];
397  if (output_vector_middle.size() != n_output_terms)
398  libmesh_error_msg("Incorrect number of output terms detected in the buffer");
399 
400  for (unsigned int q_l=0; q_l<n_output_terms; ++q_l)
401  {
402  auto output_vectors_inner_list = output_vector_middle[q_l];
403 
404  if (output_vectors_inner_list.size() != n_bfs)
405  libmesh_error_msg("Incorrect number of output terms detected in the buffer");
406 
407  for (unsigned int j=0; j<n_bfs; ++j)
408  {
409  rb_evaluation.RB_output_vectors[output_id][q_l](j) =
410  load_scalar_value(output_vectors_inner_list[j]);
411  }
412  }
413  }
414  }
415 
416  // Fq vectors and Aq matrices
417  {
418  auto rb_fq_vectors_outer_list = rb_evaluation_reader.getRbFqVectors();
419  if (rb_fq_vectors_outer_list.size() != n_F_terms)
420  libmesh_error_msg("Incorrect number of Fq vectors detected in the buffer");
421 
422  for (unsigned int q_f=0; q_f<n_F_terms; ++q_f)
423  {
424  auto rb_fq_vectors_inner_list = rb_fq_vectors_outer_list[q_f];
425  if (rb_fq_vectors_inner_list.size() != n_bfs)
426  libmesh_error_msg("Incorrect Fq vector size detected in the buffer");
427 
428  for (unsigned int i=0; i < n_bfs; ++i)
429  {
430  rb_evaluation.RB_Fq_vector[q_f](i) =
431  load_scalar_value(rb_fq_vectors_inner_list[i]);
432  }
433  }
434 
435  auto rb_Aq_matrices_outer_list = rb_evaluation_reader.getRbAqMatrices();
436  if (rb_Aq_matrices_outer_list.size() != n_A_terms)
437  libmesh_error_msg("Incorrect number of Aq matrices detected in the buffer");
438 
439  for (unsigned int q_a=0; q_a<n_A_terms; ++q_a)
440  {
441  auto rb_Aq_matrices_inner_list = rb_Aq_matrices_outer_list[q_a];
442  if (rb_Aq_matrices_inner_list.size() != n_bfs*n_bfs)
443  libmesh_error_msg("Incorrect Aq matrix size detected in the buffer");
444 
445  for (unsigned int i=0; i<n_bfs; ++i)
446  for (unsigned int j=0; j<n_bfs; ++j)
447  {
448  unsigned int offset = i*n_bfs+j;
449  rb_evaluation.RB_Aq_vector[q_a](i,j) =
450  load_scalar_value(rb_Aq_matrices_inner_list[offset]);
451  }
452  }
453  }
454 
455  // Inner-product matrix
456  if (rb_evaluation.compute_RB_inner_product)
457  {
458  auto rb_inner_product_matrix_list =
459  rb_evaluation_reader.getRbInnerProductMatrix();
460 
461  if (rb_inner_product_matrix_list.size() != n_bfs*n_bfs)
462  libmesh_error_msg("Size error while reading the inner product matrix.");
463 
464  for (unsigned int i=0; i<n_bfs; ++i)
465  for (unsigned int j=0; j<n_bfs; ++j)
466  {
467  unsigned int offset = i*n_bfs + j;
468  rb_evaluation.RB_inner_product_matrix(i,j) =
469  load_scalar_value(rb_inner_product_matrix_list[offset]);
470  }
471  }
472 }
void load_parameter_ranges(RBParametrized &rb_evaluation, RBData::ParameterRanges::Reader &parameter_ranges, RBData::DiscreteParameterList::Reader &discrete_parameters_list)
Load parameter ranges and discrete parameter values into an RBEvaluation from the corresponding struc...
void libMesh::RBDataDeserialization::load_rb_scm_evaluation_data ( RBSCMEvaluation rb_scm_eval,
RBData::RBSCMEvaluation::Reader &  rb_scm_eval_reader 
)

Load an SCM RB evaluation from a corresponding reader structure in the buffer.

Unlike the other functions above, this does not need to be templated because an RBSCMEvaluation only stores Real values, and hence doesn't depend on whether we're using complex numbers or not.

Definition at line 720 of file rb_data_deserialization.C.

References libMesh::RBSCMEvaluation::B_max, libMesh::RBSCMEvaluation::B_min, libMesh::RBSCMEvaluation::C_J, libMesh::RBSCMEvaluation::C_J_stability_vector, libMesh::RBThetaExpansion::get_n_A_terms(), libMesh::RBSCMEvaluation::get_rb_theta_expansion(), load_parameter_ranges(), libMesh::Real, and libMesh::RBSCMEvaluation::SCM_UB_vectors.

Referenced by libMesh::RBDataDeserialization::RBSCMEvaluationDeserialization::read_from_file().

722 {
723  auto parameter_ranges =
724  rb_scm_evaluation_reader.getParameterRanges();
725  auto discrete_parameters_list =
726  rb_scm_evaluation_reader.getDiscreteParameters();
727  load_parameter_ranges(rb_scm_eval,
728  parameter_ranges,
729  discrete_parameters_list);
730 
731  unsigned int n_A_terms = rb_scm_eval.get_rb_theta_expansion().get_n_A_terms();
732 
733  {
734  auto b_min_list = rb_scm_evaluation_reader.getBMin();
735 
736  if (b_min_list.size() != n_A_terms)
737  libmesh_error_msg("Size error while reading B_min");
738 
739  rb_scm_eval.B_min.clear();
740  for (unsigned int i=0; i<n_A_terms; i++)
741  rb_scm_eval.B_min.push_back(b_min_list[i]);
742  }
743 
744  {
745  auto b_max_list = rb_scm_evaluation_reader.getBMax();
746 
747  if (b_max_list.size() != n_A_terms)
748  libmesh_error_msg("Size error while reading B_max");
749 
750  rb_scm_eval.B_max.clear();
751  for (unsigned int i=0; i<n_A_terms; i++)
752  rb_scm_eval.B_max.push_back(b_max_list[i]);
753  }
754 
755  {
756  auto cJ_stability_vector =
757  rb_scm_evaluation_reader.getCJStabilityVector();
758 
759  rb_scm_eval.C_J_stability_vector.clear();
760  for (unsigned int i=0; i<cJ_stability_vector.size(); i++)
761  rb_scm_eval.C_J_stability_vector.push_back( cJ_stability_vector[i] );
762  }
763 
764  {
765  auto cJ_parameters_outer =
766  rb_scm_evaluation_reader.getCJ();
767 
768  rb_scm_eval.C_J.resize(cJ_parameters_outer.size());
769  for (unsigned int i=0; i<cJ_parameters_outer.size(); i++)
770  {
771  auto cJ_parameters_inner =
772  cJ_parameters_outer[i];
773 
774  for (unsigned int j=0; j<cJ_parameters_inner.size(); j++)
775  {
776  std::string param_name = cJ_parameters_inner[j].getName();
777  Real param_value = cJ_parameters_inner[j].getValue();
778  rb_scm_eval.C_J[i].set_value(param_name, param_value);
779  }
780  }
781  }
782 
783  {
784  auto scm_ub_vectors =
785  rb_scm_evaluation_reader.getScmUbVectors();
786 
787  // The number of UB vectors is the same as the number of C_J values
788  unsigned int n_C_J_values = rb_scm_evaluation_reader.getCJ().size();
789 
790  if (scm_ub_vectors.size() != n_C_J_values*n_A_terms)
791  libmesh_error_msg("Size mismatch in SCB UB vectors");
792 
793  rb_scm_eval.SCM_UB_vectors.resize( n_C_J_values );
794  for (unsigned int i=0; i<n_C_J_values; i++)
795  {
796  rb_scm_eval.SCM_UB_vectors[i].resize(n_A_terms);
797  for (unsigned int j=0; j<n_A_terms; j++)
798  {
799  unsigned int offset = i*n_A_terms + j;
800  rb_scm_eval.SCM_UB_vectors[i][j] = scm_ub_vectors[offset];
801 
802  }
803  }
804  }
805 
806 }
void load_parameter_ranges(RBParametrized &rb_evaluation, RBData::ParameterRanges::Reader &parameter_ranges, RBData::DiscreteParameterList::Reader &discrete_parameters_list)
Load parameter ranges and discrete parameter values into an RBEvaluation from the corresponding struc...
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real
template<typename RBEvaluationReaderNumber , typename TransRBEvaluationReaderNumber >
void libMesh::RBDataDeserialization::load_transient_rb_evaluation_data ( TransientRBEvaluation trans_rb_eval,
RBEvaluationReaderNumber &  rb_evaluation_reader,
TransRBEvaluationReaderNumber &  trans_rb_eval_reader,
bool  read_error_bound_data 
)

Load an RB evaluation from a corresponding reader structure in the buffer.

Templated to deal with both Real and Complex numbers.

Definition at line 475 of file rb_data_deserialization.C.

References libMesh::TransientRBEvaluation::Aq_Mq_representor_innerprods, libMesh::TransientRBEvaluation::Fq_Mq_representor_innerprods, libMesh::RBThetaExpansion::get_n_A_terms(), libMesh::RBThetaExpansion::get_n_F_terms(), libMesh::TransientRBThetaExpansion::get_n_M_terms(), libMesh::RBEvaluation::get_rb_theta_expansion(), libMesh::TransientRBEvaluation::initial_L2_error_all_N, load_rb_evaluation_data(), libMesh::TransientRBEvaluation::Mq_Mq_representor_innerprods, libMesh::TransientRBEvaluation::RB_initial_condition_all_N, libMesh::TransientRBEvaluation::RB_L2_matrix, libMesh::TransientRBEvaluation::RB_M_q_vector, libMesh::RBTemporalDiscretization::set_delta_t(), libMesh::RBTemporalDiscretization::set_euler_theta(), libMesh::RBTemporalDiscretization::set_n_time_steps(), and libMesh::RBTemporalDiscretization::set_time_step().

Referenced by libMesh::RBDataDeserialization::TransientRBEvaluationDeserialization::read_from_file().

479 {
480  load_rb_evaluation_data(trans_rb_eval,
481  rb_eval_reader,
482  read_error_bound_data);
483 
484  trans_rb_eval.set_delta_t( trans_rb_eval_reader.getDeltaT() );
485  trans_rb_eval.set_euler_theta( trans_rb_eval_reader.getEulerTheta() );
486  trans_rb_eval.set_n_time_steps( trans_rb_eval_reader.getNTimeSteps() );
487  trans_rb_eval.set_time_step( trans_rb_eval_reader.getTimeStep() );
488 
489  unsigned int n_bfs = rb_eval_reader.getNBfs();
490  unsigned int n_F_terms = trans_rb_eval.get_rb_theta_expansion().get_n_F_terms();
491  unsigned int n_A_terms = trans_rb_eval.get_rb_theta_expansion().get_n_A_terms();
492 
493  TransientRBThetaExpansion & trans_theta_expansion =
494  cast_ref<TransientRBThetaExpansion &>(trans_rb_eval.get_rb_theta_expansion());
495  unsigned int n_M_terms = trans_theta_expansion.get_n_M_terms();
496 
497  // L2 matrix
498  {
499  auto rb_L2_matrix_list =
500  trans_rb_eval_reader.getRbL2Matrix();
501 
502  if (rb_L2_matrix_list.size() != n_bfs*n_bfs)
503  libmesh_error_msg("Size error while reading the L2 matrix.");
504 
505  for (unsigned int i=0; i<n_bfs; ++i)
506  for (unsigned int j=0; j<n_bfs; ++j)
507  {
508  unsigned int offset = i*n_bfs + j;
509  trans_rb_eval.RB_L2_matrix(i,j) =
510  load_scalar_value(rb_L2_matrix_list[offset]);
511  }
512  }
513 
514  // Mq matrices
515  {
516  auto rb_Mq_matrices_outer_list = trans_rb_eval_reader.getRbMqMatrices();
517 
518  if (rb_Mq_matrices_outer_list.size() != n_M_terms)
519  libmesh_error_msg("Incorrect number of Mq matrices detected in the buffer");
520 
521  for (unsigned int q_m=0; q_m < n_M_terms; ++q_m)
522  {
523  auto rb_Mq_matrices_inner_list = rb_Mq_matrices_outer_list[q_m];
524  if (rb_Mq_matrices_inner_list.size() != n_bfs*n_bfs)
525  libmesh_error_msg("Incorrect Mq matrix size detected in the buffer");
526 
527  for (unsigned int i=0; i<n_bfs; ++i)
528  for (unsigned int j=0; j<n_bfs; ++j)
529  {
530  unsigned int offset = i*n_bfs+j;
531  trans_rb_eval.RB_M_q_vector[q_m](i,j) =
532  load_scalar_value(rb_Mq_matrices_inner_list[offset]);
533  }
534  }
535  }
536 
537  // The initial condition and L2 error at t=0.
538  {
539  auto initial_l2_errors_reader =
540  trans_rb_eval_reader.getInitialL2Errors();
541  if (initial_l2_errors_reader.size() != n_bfs)
542  libmesh_error_msg("Incorrect number of initial L2 error terms detected in the buffer");
543 
544  auto initial_conditions_outer_list =
545  trans_rb_eval_reader.getInitialConditions();
546  if (initial_conditions_outer_list.size() != n_bfs)
547  libmesh_error_msg("Incorrect number of outer initial conditions detected in the buffer");
548 
549  for (unsigned int i=0; i<n_bfs; i++)
550  {
551  trans_rb_eval.initial_L2_error_all_N[i] =
552  initial_l2_errors_reader[i];
553 
554  auto initial_conditions_inner_list = initial_conditions_outer_list[i];
555  if (initial_conditions_inner_list.size() != (i+1))
556  libmesh_error_msg("Incorrect number of inner initial conditions detected in the buffer");
557 
558  for (unsigned int j=0; j<=i; j++)
559  {
560  trans_rb_eval.RB_initial_condition_all_N[i](j) =
561  load_scalar_value(initial_conditions_inner_list[j]);
562  }
563  }
564  }
565 
566 
567  if (read_error_bound_data)
568  {
569  // Fq_Mq data
570  {
571  auto fq_mq_innerprods_list = trans_rb_eval_reader.getFqMqInnerprods();
572  if (fq_mq_innerprods_list.size() != n_F_terms*n_M_terms*n_bfs)
573  libmesh_error_msg("Size error while reading Fq-Mq representor data from buffer.");
574 
575  for (unsigned int q_f=0; q_f<n_F_terms; ++q_f)
576  for (unsigned int q_m=0; q_m<n_M_terms; ++q_m)
577  for (unsigned int i=0; i<n_bfs; ++i)
578  {
579  unsigned int offset = q_f*n_M_terms*n_bfs + q_m*n_bfs + i;
580  trans_rb_eval.Fq_Mq_representor_innerprods[q_f][q_m][i] =
581  load_scalar_value(fq_mq_innerprods_list[offset]);
582  }
583  }
584 
585 
586  // Mq_Mq representor inner-product data
587  {
588  unsigned int Q_m_hat = n_M_terms*(n_M_terms+1)/2;
589  auto mq_mq_innerprods_list = trans_rb_eval_reader.getMqMqInnerprods();
590  if (mq_mq_innerprods_list.size() != Q_m_hat*n_bfs*n_bfs)
591  libmesh_error_msg("Size error while reading Mq-Mq representor data from buffer.");
592 
593  for (unsigned int i=0; i<Q_m_hat; ++i)
594  for (unsigned int j=0; j<n_bfs; ++j)
595  for (unsigned int l=0; l<n_bfs; ++l)
596  {
597  unsigned int offset = i*n_bfs*n_bfs + j*n_bfs + l;
598  trans_rb_eval.Mq_Mq_representor_innerprods[i][j][l] =
599  load_scalar_value(mq_mq_innerprods_list[offset]);
600  }
601  }
602 
603  // Aq_Mq representor inner-product data
604  {
605  auto aq_mq_innerprods_list =
606  trans_rb_eval_reader.getAqMqInnerprods();
607  if (aq_mq_innerprods_list.size() != n_A_terms*n_M_terms*n_bfs*n_bfs)
608  libmesh_error_msg("Size error while reading Aq-Mq representor data from buffer.");
609 
610  for (unsigned int q_a=0; q_a<n_A_terms; q_a++)
611  for (unsigned int q_m=0; q_m<n_M_terms; q_m++)
612  for (unsigned int i=0; i<n_bfs; i++)
613  for (unsigned int j=0; j<n_bfs; j++)
614  {
615  unsigned int offset =
616  q_a*(n_M_terms*n_bfs*n_bfs) + q_m*(n_bfs*n_bfs) + i*n_bfs + j;
617 
618  trans_rb_eval.Aq_Mq_representor_innerprods[q_a][q_m][i][j] =
619  load_scalar_value(aq_mq_innerprods_list[offset]);
620  }
621  }
622 
623  }
624 }
void load_rb_evaluation_data(RBEvaluation &rb_evaluation, RBEvaluationReaderNumber &rb_evaluation_reader, bool read_error_bound_data)
Load an RB evaluation from a corresponding reader structure in the buffer.