libMesh
Classes | Functions
libMesh::RBDataSerialization Namespace Reference

Classes

class  RBEIMEvaluationSerialization
 This class serializes an RBEIMEvaluation object using the Cap'n Proto library. More...
 
class  RBEvaluationSerialization
 This class serializes an RBEvaluation object using the Cap'n Proto library. More...
 
class  RBSCMEvaluationSerialization
 This class serializes an RBSCMEvaluation object using the Cap'n Proto library. More...
 
class  TransientRBEvaluationSerialization
 This class serializes a TransientRBEvaluation object using the Cap'n Proto library. More...
 

Functions

void add_parameter_ranges_to_builder (const RBParametrized &rb_evaluation, RBData::ParameterRanges::Builder &parameter_ranges, RBData::DiscreteParameterList::Builder &discrete_parameters_list)
 Add parameter ranges for continuous and discrete parameters. More...
 
template<typename RBEvaluationBuilderNumber >
void add_rb_evaluation_data_to_builder (RBEvaluation &rb_eval, RBEvaluationBuilderNumber &rb_eval_builder)
 Add data for an RBEvaluation to the builder. More...
 
template<typename RBEvaluationBuilderNumber , typename TransRBEvaluationBuilderNumber >
void add_transient_rb_evaluation_data_to_builder (TransientRBEvaluation &trans_rb_eval, RBEvaluationBuilderNumber &rb_eval_builder, TransRBEvaluationBuilderNumber &trans_rb_eval_builder)
 Add data for a TransientRBEvaluation to the builder. More...
 
template<typename RBEvaluationBuilderNumber , typename RBEIMEvaluationBuilderNumber >
void add_rb_eim_evaluation_data_to_builder (RBEIMEvaluation &rb_eim_eval, RBEvaluationBuilderNumber &rb_eval_builder, RBEIMEvaluationBuilderNumber &rb_eim_eval_builder)
 Add data for an RBEIMEvaluation to the builder. More...
 
void add_rb_scm_evaluation_data_to_builder (RBSCMEvaluation &rb_scm_eval, RBData::RBSCMEvaluation::Builder &rb_scm_eval_builder)
 Add data for an RBSCMEvaluation to the builder. More...
 
void add_point_to_builder (const Point &point, RBData::Point3D::Builder point_builder)
 Helper function that adds point data. More...
 
void add_elem_to_builder (const Elem &elem, RBData::MeshElem::Builder mesh_elem_builder)
 Helper function that adds element data. More...
 

Function Documentation

void libMesh::RBDataSerialization::add_elem_to_builder ( const Elem elem,
RBData::MeshElem::Builder  mesh_elem_builder 
)

Helper function that adds element data.

Definition at line 751 of file rb_data_serialization.C.

References add_point_to_builder(), libMesh::Utility::enum_to_string(), libMesh::Elem::n_nodes(), libMesh::Elem::node_ref(), libMesh::Elem::subdomain_id(), and libMesh::Elem::type().

Referenced by add_rb_eim_evaluation_data_to_builder().

752 {
753  std::string elem_type_string = libMesh::Utility::enum_to_string(elem.type());
754 
755  mesh_elem_builder.setType(elem_type_string.c_str());
756  mesh_elem_builder.setSubdomainId(elem.subdomain_id());
757 
758  unsigned int n_points = elem.n_nodes();
759  auto mesh_elem_point_list = mesh_elem_builder.initPoints(n_points);
760 
761  for (unsigned int j=0; j < n_points; ++j)
762  {
763  add_point_to_builder(elem.node_ref(j), mesh_elem_point_list[j]);
764  }
765 }
void add_point_to_builder(const Point &point, RBData::Point3D::Builder point_builder)
Helper function that adds point data.
std::string enum_to_string(const T e)
void libMesh::RBDataSerialization::add_parameter_ranges_to_builder ( const RBParametrized rb_evaluation,
RBData::ParameterRanges::Builder &  parameter_ranges,
RBData::DiscreteParameterList::Builder &  discrete_parameters_list 
)

Add parameter ranges for continuous and discrete parameters.

Definition at line 262 of file rb_data_serialization.C.

References libMesh::RBParametrized::get_discrete_parameter_values(), libMesh::RBParametrized::get_n_continuous_params(), libMesh::RBParametrized::get_n_discrete_params(), libMesh::RBParametrized::get_parameter_names(), libMesh::RBParametrized::get_parameters_max(), libMesh::RBParametrized::get_parameters_min(), libMesh::RBParameters::get_value(), and libMesh::RBParametrized::is_discrete_parameter().

Referenced by add_rb_evaluation_data_to_builder(), and add_rb_scm_evaluation_data_to_builder().

265 {
266  // Continuous parameters
267  {
268  unsigned int n_continuous_parameters = rb_evaluation.get_n_continuous_params();
269  auto names = parameter_ranges_list.initNames(n_continuous_parameters);
270  auto mins = parameter_ranges_list.initMinValues(n_continuous_parameters);
271  auto maxs = parameter_ranges_list.initMaxValues(n_continuous_parameters);
272 
273  std::set<std::string> parameter_names = rb_evaluation.get_parameter_names();
274  const RBParameters & parameters_min = rb_evaluation.get_parameters_min();
275  const RBParameters & parameters_max = rb_evaluation.get_parameters_max();
276 
277  unsigned int count = 0;
278  for (const auto & parameter_name : parameter_names)
279  {
280  if (!rb_evaluation.is_discrete_parameter(parameter_name))
281  {
282  names.set(count, parameter_name);
283  mins.set(count, parameters_min.get_value(parameter_name));
284  maxs.set(count, parameters_max.get_value(parameter_name));
285 
286  ++count;
287  }
288  }
289 
290  if (count != n_continuous_parameters)
291  libmesh_error_msg("Mismatch in number of continuous parameters");
292  }
293 
294  // Discrete parameters
295  {
296  unsigned int n_discrete_parameters = rb_evaluation.get_n_discrete_params();
297  auto names = discrete_parameters_list.initNames(n_discrete_parameters);
298  auto values_outer = discrete_parameters_list.initValues(n_discrete_parameters);
299 
300  const std::map<std::string, std::vector<Real>> & discrete_parameters =
301  rb_evaluation.get_discrete_parameter_values();
302 
303  unsigned int count = 0;
304  for (const auto & discrete_parameter : discrete_parameters)
305  {
306  names.set(count, discrete_parameter.first);
307 
308  const std::vector<Real> & values = discrete_parameter.second;
309  unsigned int n_values = values.size();
310 
311  values_outer.init(count, n_values);
312  auto values_inner = values_outer[count];
313  for (unsigned int i=0; i<n_values; ++i)
314  {
315  values_inner.set(i, values[i]);
316  }
317 
318  ++count;
319  }
320 
321  if (count != n_discrete_parameters)
322  libmesh_error_msg("Mismatch in number of discrete parameters");
323  }
324 }
void libMesh::RBDataSerialization::add_point_to_builder ( const Point point,
RBData::Point3D::Builder  point_builder 
)

Helper function that adds point data.

Definition at line 740 of file rb_data_serialization.C.

Referenced by add_elem_to_builder(), and add_rb_eim_evaluation_data_to_builder().

741 {
742  point_builder.setX(point(0));
743 
744  if (LIBMESH_DIM >= 2)
745  point_builder.setY(point(1));
746 
747  if (LIBMESH_DIM >= 3)
748  point_builder.setZ(point(2));
749 }
template<typename RBEvaluationBuilderNumber , typename RBEIMEvaluationBuilderNumber >
void libMesh::RBDataSerialization::add_rb_eim_evaluation_data_to_builder ( RBEIMEvaluation rb_eim_eval,
RBEvaluationBuilderNumber &  rb_eval_builder,
RBEIMEvaluationBuilderNumber &  rb_eim_eval_builder 
)

Add data for an RBEIMEvaluation to the builder.

Templated to deal with both Real and Complex numbers.

Definition at line 601 of file rb_data_serialization.C.

References add_elem_to_builder(), add_point_to_builder(), add_rb_evaluation_data_to_builder(), libMesh::RBEvaluation::get_n_basis_functions(), libMesh::RBEIMEvaluation::interpolation_matrix, libMesh::RBEIMEvaluation::interpolation_points, libMesh::RBEIMEvaluation::interpolation_points_elem, and libMesh::RBEIMEvaluation::interpolation_points_var.

Referenced by libMesh::RBDataSerialization::RBEIMEvaluationSerialization::write_to_file().

604 {
605  add_rb_evaluation_data_to_builder(rb_eim_evaluation, rb_evaluation_builder);
606 
607  unsigned int n_bfs = rb_eim_evaluation.get_n_basis_functions();
608 
609  // EIM interpolation matrix
610  {
611  // We store the lower triangular part of an NxN matrix, the size of which is given by
612  // (N(N + 1))/2
613  unsigned int half_matrix_size = n_bfs*(n_bfs+1)/2;
614 
615  auto interpolation_matrix_list =
616  rb_eim_evaluation_builder.initInterpolationMatrix(half_matrix_size);
617  for (unsigned int i=0; i < n_bfs; ++i)
618  for (unsigned int j=0; j <= i; ++j)
619  {
620  unsigned int offset = i*(i+1)/2 + j;
621  set_scalar_in_list(interpolation_matrix_list,
622  offset,
623  rb_eim_evaluation.interpolation_matrix(i,j));
624  }
625  }
626 
627  // Interpolation points
628  {
629  auto interpolation_points_list =
630  rb_eim_evaluation_builder.initInterpolationPoints(n_bfs);
631  for (unsigned int i=0; i < n_bfs; ++i)
632  add_point_to_builder(rb_eim_evaluation.interpolation_points[i],
633  interpolation_points_list[i]);
634  }
635 
636  // Interpolation points variables
637  {
638  auto interpolation_points_var_list =
639  rb_eim_evaluation_builder.initInterpolationPointsVar(n_bfs);
640  for (unsigned int i=0; i<n_bfs; ++i)
641  interpolation_points_var_list.set(i,
642  rb_eim_evaluation.interpolation_points_var[i]);
643  }
644 
645  // Interpolation elements
646  {
647  unsigned int n_interpolation_elems =
648  rb_eim_evaluation.interpolation_points_elem.size();
649  auto interpolation_points_elem_list =
650  rb_eim_evaluation_builder.initInterpolationPointsElems(n_interpolation_elems);
651 
652  if (n_interpolation_elems != n_bfs)
653  libmesh_error_msg("The number of elements should match the number of basis functions");
654 
655  for (unsigned int i=0; i<n_interpolation_elems; ++i)
656  {
657  const libMesh::Elem & elem = *rb_eim_evaluation.interpolation_points_elem[i];
658  auto mesh_elem_builder = interpolation_points_elem_list[i];
659  add_elem_to_builder(elem, mesh_elem_builder);
660  }
661  }
662 }
This is the base class from which all geometric element types are derived.
Definition: elem.h:89
void add_rb_evaluation_data_to_builder(RBEvaluation &rb_eval, RBEvaluationBuilderNumber &rb_eval_builder)
Add data for an RBEvaluation to the builder.
void add_point_to_builder(const Point &point, RBData::Point3D::Builder point_builder)
Helper function that adds point data.
void add_elem_to_builder(const Elem &elem, RBData::MeshElem::Builder mesh_elem_builder)
Helper function that adds element data.
template<typename RBEvaluationBuilderNumber >
void libMesh::RBDataSerialization::add_rb_evaluation_data_to_builder ( RBEvaluation rb_eval,
RBEvaluationBuilderNumber &  rb_eval_builder 
)

Add data for an RBEvaluation to the builder.

Definition at line 327 of file rb_data_serialization.C.

References add_parameter_ranges_to_builder(), 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::RBEvaluation::get_n_basis_functions(), libMesh::RBThetaExpansion::get_n_F_terms(), libMesh::RBThetaExpansion::get_n_output_terms(), libMesh::RBThetaExpansion::get_n_outputs(), libMesh::RBEvaluation::get_rb_theta_expansion(), libMesh::RBEvaluation::output_dual_innerprods, libMesh::RBEvaluation::RB_Aq_vector, libMesh::RBEvaluation::RB_Fq_vector, libMesh::RBEvaluation::RB_inner_product_matrix, and libMesh::RBEvaluation::RB_output_vectors.

Referenced by add_rb_eim_evaluation_data_to_builder(), add_transient_rb_evaluation_data_to_builder(), and libMesh::RBDataSerialization::RBEvaluationSerialization::write_to_file().

329 {
330  const RBThetaExpansion & rb_theta_expansion = rb_eval.get_rb_theta_expansion();
331 
332  unsigned int n_F_terms = rb_theta_expansion.get_n_F_terms();
333  unsigned int n_A_terms = rb_theta_expansion.get_n_A_terms();
334 
335  // Number of basis functions
336  unsigned int n_bfs = rb_eval.get_n_basis_functions();
337  rb_evaluation_builder.setNBfs(n_bfs);
338 
339  // Fq representor inner-product data
340  {
341  unsigned int Q_f_hat = n_F_terms*(n_F_terms+1)/2;
342 
343  auto fq_innerprods_list = rb_evaluation_builder.initFqInnerprods(Q_f_hat);
344 
345  for (unsigned int i=0; i<Q_f_hat; i++)
346  set_scalar_in_list(fq_innerprods_list,
347  i,
348  rb_eval.Fq_representor_innerprods[i]);
349  }
350 
351  // FqAq representor inner-product data
352  {
353  auto fq_aq_innerprods_list =
354  rb_evaluation_builder.initFqAqInnerprods(n_F_terms*n_A_terms*n_bfs);
355 
356  for (unsigned int q_f=0; q_f < n_F_terms; ++q_f)
357  for (unsigned int q_a=0; q_a < n_A_terms; ++q_a)
358  for (unsigned int i=0; i < n_bfs; ++i)
359  {
360  unsigned int offset = q_f*n_A_terms*n_bfs + q_a*n_bfs + i;
361  set_scalar_in_list(
362  fq_aq_innerprods_list, offset,
363  rb_eval.Fq_Aq_representor_innerprods[q_f][q_a][i]);
364  }
365  }
366 
367  // AqAq representor inner-product data
368  {
369  unsigned int Q_a_hat = n_A_terms*(n_A_terms+1)/2;
370  auto aq_aq_innerprods_list =
371  rb_evaluation_builder.initAqAqInnerprods(Q_a_hat*n_bfs*n_bfs);
372 
373  for (unsigned int i=0; i < Q_a_hat; ++i)
374  for (unsigned int j=0; j < n_bfs; ++j)
375  for (unsigned int l=0; l < n_bfs; ++l)
376  {
377  unsigned int offset = i*n_bfs*n_bfs + j*n_bfs + l;
378  set_scalar_in_list(
379  aq_aq_innerprods_list,
380  offset,
381  rb_eval.Aq_Aq_representor_innerprods[i][j][l]);
382  }
383  }
384 
385  // Output dual inner-product data, and output vectors
386  {
387  unsigned int n_outputs = rb_theta_expansion.get_n_outputs();
388  auto output_innerprod_outer = rb_evaluation_builder.initOutputDualInnerprods(n_outputs);
389  auto output_vector_outer = rb_evaluation_builder.initOutputVectors(n_outputs);
390 
391  for (unsigned int output_id=0; output_id < n_outputs; ++output_id)
392  {
393  unsigned int n_output_terms = rb_theta_expansion.get_n_output_terms(output_id);
394 
395  {
396  unsigned int Q_l_hat = n_output_terms*(n_output_terms+1)/2;
397  auto output_innerprod_inner = output_innerprod_outer.init(output_id, Q_l_hat);
398  for (unsigned int q=0; q < Q_l_hat; ++q)
399  {
400  set_scalar_in_list(
401  output_innerprod_inner, q, rb_eval.output_dual_innerprods[output_id][q]);
402  }
403  }
404 
405  {
406  auto output_vector_middle = output_vector_outer.init(output_id, n_output_terms);
407  for (unsigned int q_l=0; q_l<n_output_terms; ++q_l)
408  {
409  auto output_vector_inner = output_vector_middle.init(q_l, n_bfs);
410  for (unsigned int j=0; j<n_bfs; ++j)
411  {
412  set_scalar_in_list(
413  output_vector_inner, j, rb_eval.RB_output_vectors[output_id][q_l](j));
414  }
415  }
416  }
417  }
418  }
419 
420  // Fq vectors and Aq matrices
421  {
422  unsigned int n_F_terms = rb_theta_expansion.get_n_F_terms();
423  unsigned int n_A_terms = rb_theta_expansion.get_n_A_terms();
424 
425  auto rb_fq_vectors_outer_list = rb_evaluation_builder.initRbFqVectors(n_F_terms);
426  for (unsigned int q_f=0; q_f < n_F_terms; ++q_f)
427  {
428  auto rb_fq_vectors_inner_list = rb_fq_vectors_outer_list.init(q_f, n_bfs);
429  for (unsigned int i=0; i<n_bfs; i++)
430  set_scalar_in_list(rb_fq_vectors_inner_list, i, rb_eval.RB_Fq_vector[q_f](i));
431  }
432 
433  auto rb_Aq_matrices_outer_list = rb_evaluation_builder.initRbAqMatrices(n_A_terms);
434  for (unsigned int q_a=0; q_a < n_A_terms; ++q_a)
435  {
436  auto rb_Aq_matrices_inner_list = rb_Aq_matrices_outer_list.init(q_a, n_bfs*n_bfs);
437  for (unsigned int i=0; i < n_bfs; ++i)
438  for (unsigned int j=0; j < n_bfs; ++j)
439  {
440  unsigned int offset = i*n_bfs+j;
441  set_scalar_in_list(rb_Aq_matrices_inner_list, offset, rb_eval.RB_Aq_vector[q_a](i,j));
442  }
443  }
444  }
445 
446  // Inner-product matrix
447  if (rb_eval.compute_RB_inner_product)
448  {
449  auto rb_inner_product_matrix_list =
450  rb_evaluation_builder.initRbInnerProductMatrix(n_bfs*n_bfs);
451 
452  for (unsigned int i=0; i < n_bfs; ++i)
453  for (unsigned int j=0; j < n_bfs; ++j)
454  {
455  unsigned int offset = i*n_bfs + j;
456  set_scalar_in_list(
457  rb_inner_product_matrix_list,
458  offset,
459  rb_eval.RB_inner_product_matrix(i,j) );
460  }
461  }
462 
463  auto parameter_ranges_list =
464  rb_evaluation_builder.initParameterRanges();
465  auto discrete_parameters_list =
466  rb_evaluation_builder.initDiscreteParameters();
468  parameter_ranges_list,
469  discrete_parameters_list);
470 }
void add_parameter_ranges_to_builder(const RBParametrized &rb_evaluation, RBData::ParameterRanges::Builder &parameter_ranges, RBData::DiscreteParameterList::Builder &discrete_parameters_list)
Add parameter ranges for continuous and discrete parameters.
void libMesh::RBDataSerialization::add_rb_scm_evaluation_data_to_builder ( RBSCMEvaluation rb_scm_eval,
RBData::RBSCMEvaluation::Builder &  rb_scm_eval_builder 
)

Add data for an RBSCMEvaluation to the builder.

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 665 of file rb_data_serialization.C.

References add_parameter_ranges_to_builder(), libMesh::RBSCMEvaluation::B_max, libMesh::RBSCMEvaluation::B_min, libMesh::RBSCMEvaluation::C_J, libMesh::RBSCMEvaluation::C_J_stability_vector, libMesh::RBSCMEvaluation::get_B_max(), libMesh::RBSCMEvaluation::get_B_min(), libMesh::RBSCMEvaluation::get_C_J_stability_constraint(), libMesh::RBThetaExpansion::get_n_A_terms(), libMesh::RBSCMEvaluation::get_rb_theta_expansion(), and libMesh::RBSCMEvaluation::get_SCM_UB_vector().

Referenced by libMesh::RBDataSerialization::RBSCMEvaluationSerialization::write_to_file().

667 {
668  auto parameter_ranges_list =
669  rb_scm_eval_builder.initParameterRanges();
670  auto discrete_parameters_list =
671  rb_scm_eval_builder.initDiscreteParameters();
673  parameter_ranges_list,
674  discrete_parameters_list);
675 
676  {
677  if (rb_scm_eval.B_min.size() != rb_scm_eval.get_rb_theta_expansion().get_n_A_terms())
678  libmesh_error_msg("Size error while writing B_min");
679  auto b_min_list = rb_scm_eval_builder.initBMin( rb_scm_eval.B_min.size() );
680  for (std::size_t i=0; i<rb_scm_eval.B_min.size(); i++)
681  b_min_list.set(i, rb_scm_eval.get_B_min(i));
682  }
683 
684  {
685  if (rb_scm_eval.B_max.size() != rb_scm_eval.get_rb_theta_expansion().get_n_A_terms())
686  libmesh_error_msg("Size error while writing B_max");
687 
688  auto b_max_list = rb_scm_eval_builder.initBMax( rb_scm_eval.B_max.size() );
689  for (std::size_t i=0; i<rb_scm_eval.B_max.size(); i++)
690  b_max_list.set(i, rb_scm_eval.get_B_max(i));
691  }
692 
693  {
694  auto cj_stability_vector =
695  rb_scm_eval_builder.initCJStabilityVector( rb_scm_eval.C_J_stability_vector.size() );
696  for (std::size_t i=0; i<rb_scm_eval.C_J_stability_vector.size(); i++)
697  cj_stability_vector.set(i, rb_scm_eval.get_C_J_stability_constraint(i));
698  }
699 
700  {
701  auto cj_parameters_outer =
702  rb_scm_eval_builder.initCJ( rb_scm_eval.C_J.size() );
703 
704  for (std::size_t i=0; i<rb_scm_eval.C_J.size(); i++)
705  {
706  auto cj_parameters_inner =
707  cj_parameters_outer.init(i, rb_scm_eval.C_J[i].n_parameters());
708 
709  RBParameters::const_iterator it = rb_scm_eval.C_J[i].begin();
710  RBParameters::const_iterator it_end = rb_scm_eval.C_J[i].end();
711 
712  unsigned int count = 0;
713  for ( ; it != it_end; ++it)
714  {
715  cj_parameters_inner[count].setName( it->first );
716  cj_parameters_inner[count].setValue( it->second );
717  count++;
718  }
719 
720  }
721  }
722 
723  {
724  unsigned int n_C_J_values = rb_scm_eval.C_J.size();
725  unsigned int n_A_terms = rb_scm_eval.get_rb_theta_expansion().get_n_A_terms();
726  unsigned int n_values = n_C_J_values*n_A_terms;
727  auto scm_ub_vectors =
728  rb_scm_eval_builder.initScmUbVectors( n_values );
729 
730  for (unsigned int i=0; i<n_C_J_values; i++)
731  for (unsigned int j=0; j<n_A_terms; j++)
732  {
733  unsigned int offset = i*n_A_terms + j;
734  scm_ub_vectors.set(offset, rb_scm_eval.get_SCM_UB_vector(i,j));
735  }
736  }
737 }
void add_parameter_ranges_to_builder(const RBParametrized &rb_evaluation, RBData::ParameterRanges::Builder &parameter_ranges, RBData::DiscreteParameterList::Builder &discrete_parameters_list)
Add parameter ranges for continuous and discrete parameters.
template<typename RBEvaluationBuilderNumber , typename TransRBEvaluationBuilderNumber >
void libMesh::RBDataSerialization::add_transient_rb_evaluation_data_to_builder ( TransientRBEvaluation trans_rb_eval,
RBEvaluationBuilderNumber &  rb_eval_builder,
TransRBEvaluationBuilderNumber &  trans_rb_eval_builder 
)

Add data for a TransientRBEvaluation to the builder.

Templated to deal with both Real and Complex numbers.

Definition at line 473 of file rb_data_serialization.C.

References add_rb_evaluation_data_to_builder(), libMesh::TransientRBEvaluation::Aq_Mq_representor_innerprods, libMesh::TransientRBEvaluation::Fq_Mq_representor_innerprods, libMesh::RBTemporalDiscretization::get_delta_t(), libMesh::RBTemporalDiscretization::get_euler_theta(), libMesh::RBEvaluation::get_n_basis_functions(), libMesh::TransientRBThetaExpansion::get_n_M_terms(), libMesh::RBTemporalDiscretization::get_n_time_steps(), libMesh::RBEvaluation::get_rb_theta_expansion(), libMesh::RBTemporalDiscretization::get_time_step(), libMesh::TransientRBEvaluation::initial_L2_error_all_N, libMesh::TransientRBEvaluation::Mq_Mq_representor_innerprods, libMesh::TransientRBEvaluation::RB_initial_condition_all_N, libMesh::TransientRBEvaluation::RB_L2_matrix, and libMesh::TransientRBEvaluation::RB_M_q_vector.

Referenced by libMesh::RBDataSerialization::TransientRBEvaluationSerialization::write_to_file().

476 {
477  add_rb_evaluation_data_to_builder(trans_rb_eval, rb_eval_builder);
478 
479  trans_rb_eval_builder.setDeltaT(trans_rb_eval.get_delta_t());
480  trans_rb_eval_builder.setEulerTheta(trans_rb_eval.get_euler_theta());
481  trans_rb_eval_builder.setNTimeSteps(trans_rb_eval.get_n_time_steps());
482  trans_rb_eval_builder.setTimeStep(trans_rb_eval.get_time_step());
483 
484  unsigned int n_bfs = trans_rb_eval.get_n_basis_functions();
485 
486  // L2-inner-product matrix
487  {
488  auto rb_L2_matrix_list =
489  trans_rb_eval_builder.initRbL2Matrix(n_bfs*n_bfs);
490 
491  for (unsigned int i=0; i<n_bfs; ++i)
492  for (unsigned int j=0; j<n_bfs; ++j)
493  {
494  unsigned int offset = i*n_bfs + j;
495  set_scalar_in_list(rb_L2_matrix_list,
496  offset,
497  trans_rb_eval.RB_L2_matrix(i,j));
498  }
499  }
500 
501  TransientRBThetaExpansion & trans_theta_expansion =
502  cast_ref<TransientRBThetaExpansion &>(trans_rb_eval.get_rb_theta_expansion());
503  unsigned int n_M_terms = trans_theta_expansion.get_n_M_terms();
504  // Mq matrices
505  {
506  auto rb_Mq_matrices_outer_list = trans_rb_eval_builder.initRbMqMatrices(n_M_terms);
507  for (unsigned int q_m=0; q_m < n_M_terms; ++q_m)
508  {
509  auto rb_Mq_matrices_inner_list = rb_Mq_matrices_outer_list.init(q_m, n_bfs*n_bfs);
510  for (unsigned int i=0; i < n_bfs; ++i)
511  for (unsigned int j=0; j < n_bfs; ++j)
512  {
513  unsigned int offset = i*n_bfs+j;
514  set_scalar_in_list(rb_Mq_matrices_inner_list,
515  offset,
516  trans_rb_eval.RB_M_q_vector[q_m](i,j));
517  }
518  }
519  }
520 
521  // The initial condition and L2 error at t=0.
522  // We store the values for each RB space of dimension (0,...,n_basis_functions).
523  {
524  auto initial_l2_errors_builder =
525  trans_rb_eval_builder.initInitialL2Errors(n_bfs);
526  auto initial_conditions_outer_list =
527  trans_rb_eval_builder.initInitialConditions(n_bfs);
528 
529  for (unsigned int i=0; i<n_bfs; i++)
530  {
531  initial_l2_errors_builder.set(i, trans_rb_eval.initial_L2_error_all_N[i]);
532 
533  auto initial_conditions_inner_list =
534  initial_conditions_outer_list.init(i, i+1);
535  for (unsigned int j=0; j<=i; j++)
536  {
537  set_scalar_in_list(initial_conditions_inner_list,
538  j,
539  trans_rb_eval.RB_initial_condition_all_N[i](j));
540  }
541  }
542  }
543 
544  // FqMq representor inner-product data
545  {
546  unsigned int n_F_terms = trans_theta_expansion.get_n_F_terms();
547  auto fq_mq_innerprods_list =
548  trans_rb_eval_builder.initFqMqInnerprods(n_F_terms*n_M_terms*n_bfs);
549 
550  for (unsigned int q_f=0; q_f<n_F_terms; ++q_f)
551  for (unsigned int q_m=0; q_m<n_M_terms; ++q_m)
552  for (unsigned int i=0; i<n_bfs; ++i)
553  {
554  unsigned int offset = q_f*n_M_terms*n_bfs + q_m*n_bfs + i;
555  set_scalar_in_list(fq_mq_innerprods_list,
556  offset,
557  trans_rb_eval.Fq_Mq_representor_innerprods[q_f][q_m][i]);
558  }
559  }
560 
561  // MqMq representor inner-product data
562  {
563  unsigned int Q_m_hat = n_M_terms*(n_M_terms+1)/2;
564  auto mq_mq_innerprods_list =
565  trans_rb_eval_builder.initMqMqInnerprods(Q_m_hat*n_bfs*n_bfs);
566 
567  for (unsigned int i=0; i < Q_m_hat; ++i)
568  for (unsigned int j=0; j < n_bfs; ++j)
569  for (unsigned int l=0; l < n_bfs; ++l)
570  {
571  unsigned int offset = i*n_bfs*n_bfs + j*n_bfs + l;
572  set_scalar_in_list(mq_mq_innerprods_list,
573  offset,
574  trans_rb_eval.Mq_Mq_representor_innerprods[i][j][l]);
575  }
576  }
577 
578  // AqMq representor inner-product data
579  {
580  unsigned int n_A_terms = trans_theta_expansion.get_n_A_terms();
581 
582  auto aq_mq_innerprods_list =
583  trans_rb_eval_builder.initAqMqInnerprods(n_A_terms*n_M_terms*n_bfs*n_bfs);
584 
585  for (unsigned int q_a=0; q_a<n_A_terms; q_a++)
586  for (unsigned int q_m=0; q_m<n_M_terms; q_m++)
587  for (unsigned int i=0; i<n_bfs; i++)
588  for (unsigned int j=0; j<n_bfs; j++)
589  {
590  unsigned int offset =
591  q_a*(n_M_terms*n_bfs*n_bfs) + q_m*(n_bfs*n_bfs) + i*n_bfs + j;
592  set_scalar_in_list(aq_mq_innerprods_list,
593  offset,
594  trans_rb_eval.Aq_Mq_representor_innerprods[q_a][q_m][i][j]);
595  }
596  }
597 
598 }
void add_rb_evaluation_data_to_builder(RBEvaluation &rb_eval, RBEvaluationBuilderNumber &rb_eval_builder)
Add data for an RBEvaluation to the builder.