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 RBEIMEvaluationBuilderNumber >
void add_rb_eim_evaluation_data_to_builder (RBEIMEvaluation &rb_eim_eval, 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...
 

Function Documentation

◆ add_parameter_ranges_to_builder()

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 246 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_parameters_max(), libMesh::RBParametrized::get_parameters_min(), libMesh::RBParameters::get_value(), and libMesh::RBParametrized::is_discrete_parameter().

Referenced by add_rb_eim_evaluation_data_to_builder(), add_rb_evaluation_data_to_builder(), and add_rb_scm_evaluation_data_to_builder().

249 {
250  // Continuous parameters
251  {
252  unsigned int n_continuous_parameters = rb_evaluation.get_n_continuous_params();
253  auto names = parameter_ranges_list.initNames(n_continuous_parameters);
254  auto mins = parameter_ranges_list.initMinValues(n_continuous_parameters);
255  auto maxs = parameter_ranges_list.initMaxValues(n_continuous_parameters);
256 
257  const RBParameters & parameters_min = rb_evaluation.get_parameters_min();
258  const RBParameters & parameters_max = rb_evaluation.get_parameters_max();
259 
260  // We could loop over either parameters_min or parameters_max, they should have the same keys.
261  unsigned int count = 0;
262  for (const auto & [key, val] : parameters_min)
263  if (!rb_evaluation.is_discrete_parameter(key))
264  {
265  names.set(count, key);
266  mins.set(count, parameters_min.get_value(key));
267  maxs.set(count, parameters_max.get_value(key));
268  ++count;
269  }
270 
271  libmesh_error_msg_if(count != n_continuous_parameters, "Mismatch in number of continuous parameters");
272  }
273 
274  // Discrete parameters
275  {
276  unsigned int n_discrete_parameters = rb_evaluation.get_n_discrete_params();
277  auto names = discrete_parameters_list.initNames(n_discrete_parameters);
278  auto values_outer = discrete_parameters_list.initValues(n_discrete_parameters);
279 
280  const std::map<std::string, std::vector<Real>> & discrete_parameters =
281  rb_evaluation.get_discrete_parameter_values();
282 
283  unsigned int count = 0;
284  for (const auto & discrete_parameter : discrete_parameters)
285  {
286  names.set(count, discrete_parameter.first);
287 
288  const std::vector<Real> & values = discrete_parameter.second;
289  unsigned int n_values = values.size();
290 
291  values_outer.init(count, n_values);
292  auto values_inner = values_outer[count];
293  for (unsigned int i=0; i<n_values; ++i)
294  {
295  values_inner.set(i, values[i]);
296  }
297 
298  ++count;
299  }
300 
301  libmesh_error_msg_if(count != n_discrete_parameters, "Mismatch in number of discrete parameters");
302  }
303 }

◆ add_point_to_builder()

void libMesh::RBDataSerialization::add_point_to_builder ( const Point point,
RBData::Point3D::Builder  point_builder 
)

Helper function that adds point data.

Definition at line 930 of file rb_data_serialization.C.

Referenced by add_rb_eim_evaluation_data_to_builder().

931 {
932  point_builder.setX(point(0));
933 
934  if (LIBMESH_DIM >= 2)
935  point_builder.setY(point(1));
936 
937  if (LIBMESH_DIM >= 3)
938  point_builder.setZ(point(2));
939 }

◆ add_rb_eim_evaluation_data_to_builder()

template<typename RBEIMEvaluationBuilderNumber >
void libMesh::RBDataSerialization::add_rb_eim_evaluation_data_to_builder ( RBEIMEvaluation rb_eim_eval,
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 580 of file rb_data_serialization.C.

References add_parameter_ranges_to_builder(), add_point_to_builder(), libMesh::RBEIMEvaluation::get_eim_solutions_for_training_set(), libMesh::RBEIMEvaluation::get_error_indicator_interpolation_row(), libMesh::RBEIMEvaluation::get_interpolation_matrix(), libMesh::RBEIMEvaluation::get_interpolation_points_boundary_id(), libMesh::RBEIMEvaluation::get_interpolation_points_comp(), libMesh::RBEIMEvaluation::get_interpolation_points_elem_id(), libMesh::RBEIMEvaluation::get_interpolation_points_elem_type(), libMesh::RBEIMEvaluation::get_interpolation_points_JxW_all_qp(), libMesh::RBEIMEvaluation::get_interpolation_points_node_id(), libMesh::RBEIMEvaluation::get_interpolation_points_phi_i_all_qp(), libMesh::RBEIMEvaluation::get_interpolation_points_phi_i_qp(), libMesh::RBEIMEvaluation::get_interpolation_points_qp(), libMesh::RBEIMEvaluation::get_interpolation_points_side_index(), libMesh::RBEIMEvaluation::get_interpolation_points_spatial_indices(), libMesh::RBEIMEvaluation::get_interpolation_points_subdomain_id(), libMesh::RBEIMEvaluation::get_interpolation_points_xyz(), libMesh::RBEIMEvaluation::get_interpolation_points_xyz_perturbations(), libMesh::RBEIMEvaluation::get_n_basis_functions(), libMesh::RBEIMEvaluation::get_n_interpolation_points(), libMesh::RBEIMEvaluation::get_n_interpolation_points_spatial_indices(), libMesh::RBEIMEvaluation::get_parametrized_function(), libMesh::index_range(), libMesh::RBEIMEvaluation::initialize_interpolation_points_spatial_indices(), libMesh::RBParametrizedFunction::is_lookup_table, libMesh::make_range(), libMesh::RBParametrizedFunction::on_mesh_nodes(), libMesh::RBParametrizedFunction::on_mesh_sides(), libMesh::DenseVector< T >::size(), and libMesh::RBEIMEvaluation::use_eim_error_indicator().

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

582 {
583  // Number of basis functions
584  unsigned int n_bfs = rb_eim_evaluation.get_n_basis_functions();
585  rb_eim_evaluation_builder.setNBfs(n_bfs);
586 
587  // EIM interpolation matrix
588  {
589  // We store the lower triangular part of an NxN matrix, the size of which is given by
590  // (N(N + 1))/2
591  unsigned int half_matrix_size = n_bfs*(n_bfs+1)/2;
592 
593  auto interpolation_matrix_list =
594  rb_eim_evaluation_builder.initInterpolationMatrix(half_matrix_size);
595  for (unsigned int i=0; i < n_bfs; ++i)
596  for (unsigned int j=0; j <= i; ++j)
597  {
598  unsigned int offset = i*(i+1)/2 + j;
599  set_scalar_in_list(interpolation_matrix_list,
600  offset,
601  rb_eim_evaluation.get_interpolation_matrix()(i,j));
602  }
603  }
604 
605  // We check use_eim_error_indicator() and the number of interpolation
606  // points in order to set use_error_indicator. This is because even if
607  // use_error_indicator() is true, we may not be using an error indicator,
608  // e.g. if there are no parameters in the model then we will not use an
609  // error indicator, and we can detect this by comparing the number of
610  // interpolation points with n_bfs.
611  bool use_error_indicator =
612  (rb_eim_evaluation.use_eim_error_indicator() &&
613  (rb_eim_evaluation.get_n_interpolation_points() > n_bfs));
614 
615  if (use_error_indicator)
616  {
617  auto error_indicator_data_list =
618  rb_eim_evaluation_builder.initEimErrorIndicatorInterpData(n_bfs);
619  const auto & error_indicator_row =
620  rb_eim_evaluation.get_error_indicator_interpolation_row();
621  for (unsigned int i=0; i < n_bfs; ++i)
622  {
623  set_scalar_in_list(error_indicator_data_list,
624  i,
625  error_indicator_row(i));
626  }
627  }
628 
629  // If we're using the EIM error indicator then we store one extra
630  // interpolation point and associated data, hence we increment n_bfs
631  // here so that we write out the extra data below.
632  //
633  // However the interpolation matrix and _error_indicator_interpolation_row
634  // does not include this extra point, which is why we write it out above
635  // before n_bfs is incremented.
636  if (use_error_indicator)
637  n_bfs++;
638 
639  libmesh_error_msg_if(n_bfs != rb_eim_evaluation.get_n_interpolation_points(),
640  "Number of basis functions should match number of interpolation points");
641 
642  auto parameter_ranges_list =
643  rb_eim_evaluation_builder.initParameterRanges();
644  auto discrete_parameters_list =
645  rb_eim_evaluation_builder.initDiscreteParameters();
646  add_parameter_ranges_to_builder(rb_eim_evaluation,
647  parameter_ranges_list,
648  discrete_parameters_list);
649 
650  // Interpolation points
651  {
652  auto interpolation_points_list =
653  rb_eim_evaluation_builder.initInterpolationXyz(n_bfs);
654  for (unsigned int i=0; i < n_bfs; ++i)
655  add_point_to_builder(rb_eim_evaluation.get_interpolation_points_xyz(i),
656  interpolation_points_list[i]);
657  }
658 
659  // Interpolation points comps
660  {
661  auto interpolation_points_comp_list =
662  rb_eim_evaluation_builder.initInterpolationComp(n_bfs);
663  for (unsigned int i=0; i<n_bfs; ++i)
664  interpolation_points_comp_list.set(i,
665  rb_eim_evaluation.get_interpolation_points_comp(i));
666  }
667 
668  // Interpolation points subdomain IDs
669  {
670  auto interpolation_points_subdomain_id_list =
671  rb_eim_evaluation_builder.initInterpolationSubdomainId(n_bfs);
672  for (unsigned int i=0; i<n_bfs; ++i)
673  interpolation_points_subdomain_id_list.set(i,
674  rb_eim_evaluation.get_interpolation_points_subdomain_id(i));
675  }
676 
677  // Interpolation points boundary IDs, relevant if the parametrized function is defined
678  // on mesh sides or nodesets
679  if (rb_eim_evaluation.get_parametrized_function().on_mesh_sides() ||
680  rb_eim_evaluation.get_parametrized_function().on_mesh_nodes())
681  {
682  auto interpolation_points_boundary_id_list =
683  rb_eim_evaluation_builder.initInterpolationBoundaryId(n_bfs);
684  for (unsigned int i=0; i<n_bfs; ++i)
685  interpolation_points_boundary_id_list.set(i,
686  rb_eim_evaluation.get_interpolation_points_boundary_id(i));
687  }
688 
689  // Interpolation points element IDs
690  {
691  auto interpolation_points_elem_id_list =
692  rb_eim_evaluation_builder.initInterpolationElemId(n_bfs);
693  for (unsigned int i=0; i<n_bfs; ++i)
694  interpolation_points_elem_id_list.set(i,
695  rb_eim_evaluation.get_interpolation_points_elem_id(i));
696  }
697 
698  // Interpolation points node IDs, relevant if the parametrized function is defined on mesh sides
699  if (rb_eim_evaluation.get_parametrized_function().on_mesh_nodes())
700  {
701  auto interpolation_points_node_id_list =
702  rb_eim_evaluation_builder.initInterpolationNodeId(n_bfs);
703  for (unsigned int i=0; i<n_bfs; ++i)
704  interpolation_points_node_id_list.set(i,
705  rb_eim_evaluation.get_interpolation_points_node_id(i));
706  }
707 
708  // Interpolation points side indices, relevant if the parametrized function is defined on mesh sides
709  if (rb_eim_evaluation.get_parametrized_function().on_mesh_sides())
710  {
711  auto interpolation_points_side_index_list =
712  rb_eim_evaluation_builder.initInterpolationSideIndex(n_bfs);
713  for (unsigned int i=0; i<n_bfs; ++i)
714  interpolation_points_side_index_list.set(i,
715  rb_eim_evaluation.get_interpolation_points_side_index(i));
716  }
717 
718  // Interpolation points quadrature point indices
719  {
720  auto interpolation_points_qp_list =
721  rb_eim_evaluation_builder.initInterpolationQp(n_bfs);
722  for (unsigned int i=0; i<n_bfs; ++i)
723  interpolation_points_qp_list.set(i,
724  rb_eim_evaluation.get_interpolation_points_qp(i));
725  }
726 
727  // Interpolation points perturbations
728  {
729  auto interpolation_points_list_outer =
730  rb_eim_evaluation_builder.initInterpolationXyzPerturb(n_bfs);
731  for (unsigned int i=0; i < n_bfs; ++i)
732  {
733  const std::vector<Point> & perturbs = rb_eim_evaluation.get_interpolation_points_xyz_perturbations(i);
734  auto interpolation_points_list_inner = interpolation_points_list_outer.init(i, perturbs.size());
735 
736  for (unsigned int j : index_range(perturbs))
737  {
738  add_point_to_builder(perturbs[j], interpolation_points_list_inner[j]);
739  }
740  }
741  }
742 
743  // Interpolation points JxW values at each qp
744  {
745  auto interpolation_points_list_outer =
746  rb_eim_evaluation_builder.initInterpolationJxWAllQp(n_bfs);
747  for (unsigned int i=0; i < n_bfs; ++i)
748  {
749  const std::vector<Real> & JxW = rb_eim_evaluation.get_interpolation_points_JxW_all_qp(i);
750  auto interpolation_points_list_inner = interpolation_points_list_outer.init(i, JxW.size());
751 
752  for (unsigned int j : index_range(JxW))
753  {
754  // Here we can use set() instead of set_scalar_in_list() because
755  // phi stores real-valued data only.
756  interpolation_points_list_inner.set(j, JxW[j]);
757  }
758  }
759  }
760 
761  // Interpolation points phi values at each qp
762  {
763  auto interpolation_points_list_outer =
764  rb_eim_evaluation_builder.initInterpolationPhiValuesAllQp(n_bfs);
765 
766  for (unsigned int i=0; i < n_bfs; ++i)
767  {
768  const auto & phi_i_all_qp = rb_eim_evaluation.get_interpolation_points_phi_i_all_qp(i);
769  auto interpolation_points_list_middle = interpolation_points_list_outer.init(i, phi_i_all_qp.size());
770 
771  for (unsigned int j : index_range(phi_i_all_qp))
772  {
773  auto interpolation_points_list_inner = interpolation_points_list_middle.init(j, phi_i_all_qp[j].size());
774 
775  for (auto k : index_range(phi_i_all_qp[j]))
776  interpolation_points_list_inner.set(k, phi_i_all_qp[j][k]);
777  }
778  }
779  }
780 
781  // Element type for the element that contains each interpolation point
782  {
783  auto interpolation_points_elem_type_list =
784  rb_eim_evaluation_builder.initInterpolationElemType(n_bfs);
785  for (unsigned int i=0; i<n_bfs; ++i)
786  interpolation_points_elem_type_list.set(i,
787  rb_eim_evaluation.get_interpolation_points_elem_type(i));
788  }
789 
790  // Optionally store EIM solutions for the training set
791  if (rb_eim_evaluation.get_parametrized_function().is_lookup_table)
792  {
793  const std::vector<DenseVector<Number>> & eim_solutions = rb_eim_evaluation.get_eim_solutions_for_training_set();
794 
795  auto eim_rhs_list_outer =
796  rb_eim_evaluation_builder.initEimSolutionsForTrainingSet(eim_solutions.size());
797  for (auto i : make_range(eim_solutions.size()))
798  {
799  const DenseVector<Number> & values = eim_solutions[i];
800  auto eim_rhs_list_inner = eim_rhs_list_outer.init(i, values.size());
801 
802  for (auto j : index_range(values))
803  {
804  set_scalar_in_list(eim_rhs_list_inner,
805  j,
806  values(j));
807  }
808  }
809  }
810 
811  // The shape function values at the interpolation points. This can be used to evaluate nodal data
812  // at EIM interpolation points, which are at quadrature points.
813  {
814  auto interpolation_points_list_outer =
815  rb_eim_evaluation_builder.initInterpolationPhiValues(n_bfs);
816  for (unsigned int i=0; i < n_bfs; ++i)
817  {
818  const std::vector<Real> & phi_i_qp_vec = rb_eim_evaluation.get_interpolation_points_phi_i_qp(i);
819  auto interpolation_points_list_inner = interpolation_points_list_outer.init(i, phi_i_qp_vec.size());
820 
821  for (unsigned int j : index_range(phi_i_qp_vec))
822  {
823  // Here we can use set() instead of set_scalar_in_list() because
824  // phi stores real-valued data only.
825  interpolation_points_list_inner.set(j, phi_i_qp_vec[j]);
826  }
827  }
828  }
829 
830  // Initialize EIM spatial indices for the interpolation points, and if there are
831  // any spatial indices then we store them in the buffer.
832  rb_eim_evaluation.initialize_interpolation_points_spatial_indices();
833  if (rb_eim_evaluation.get_n_interpolation_points_spatial_indices() > 0)
834  {
835  libmesh_error_msg_if(n_bfs != rb_eim_evaluation.get_n_interpolation_points_spatial_indices(),
836  "Error: Number of spatial indices should match number of EIM basis functions");
837 
838  auto interpolation_points_spatial_indices_builder =
839  rb_eim_evaluation_builder.initInterpolationSpatialIndices(n_bfs);
840  for (unsigned int i=0; i<n_bfs; ++i)
841  {
842  const std::vector<unsigned int> & spatial_indices =
843  rb_eim_evaluation.get_interpolation_points_spatial_indices(i);
844  unsigned int n_spatial_indices = spatial_indices.size();
845 
846  auto spatial_indices_builder =
847  interpolation_points_spatial_indices_builder.init(i, n_spatial_indices);
848 
849  for (auto j : make_range(n_spatial_indices))
850  {
851  spatial_indices_builder.set(j, spatial_indices[j]);
852  }
853  }
854  }
855 }
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 add_point_to_builder(const Point &point, RBData::Point3D::Builder point_builder)
Helper function that adds point data.
IntRange< T > make_range(T beg, T end)
The 2-parameter make_range() helper function returns an IntRange<T> when both input parameters are of...
Definition: int_range.h:134
auto index_range(const T &sizable)
Helper function that returns an IntRange<std::size_t> representing all the indices of the passed-in v...
Definition: int_range.h:111

◆ add_rb_evaluation_data_to_builder()

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 306 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_transient_rb_evaluation_data_to_builder(), and libMesh::RBDataSerialization::RBEvaluationSerialization::write_to_file().

308 {
309  const RBThetaExpansion & rb_theta_expansion = rb_eval.get_rb_theta_expansion();
310 
311  unsigned int n_F_terms = rb_theta_expansion.get_n_F_terms();
312  unsigned int n_A_terms = rb_theta_expansion.get_n_A_terms();
313 
314  // Number of basis functions
315  unsigned int n_bfs = rb_eval.get_n_basis_functions();
316  rb_evaluation_builder.setNBfs(n_bfs);
317 
318  // Fq representor inner-product data
319  {
320  unsigned int Q_f_hat = n_F_terms*(n_F_terms+1)/2;
321 
322  auto fq_innerprods_list = rb_evaluation_builder.initFqInnerprods(Q_f_hat);
323 
324  for (unsigned int i=0; i<Q_f_hat; i++)
325  set_scalar_in_list(fq_innerprods_list,
326  i,
327  rb_eval.Fq_representor_innerprods[i]);
328  }
329 
330  // FqAq representor inner-product data
331  {
332  auto fq_aq_innerprods_list =
333  rb_evaluation_builder.initFqAqInnerprods(n_F_terms*n_A_terms*n_bfs);
334 
335  for (unsigned int q_f=0; q_f < n_F_terms; ++q_f)
336  for (unsigned int q_a=0; q_a < n_A_terms; ++q_a)
337  for (unsigned int i=0; i < n_bfs; ++i)
338  {
339  unsigned int offset = q_f*n_A_terms*n_bfs + q_a*n_bfs + i;
340  set_scalar_in_list(
341  fq_aq_innerprods_list, offset,
342  rb_eval.Fq_Aq_representor_innerprods[q_f][q_a][i]);
343  }
344  }
345 
346  // AqAq representor inner-product data
347  {
348  unsigned int Q_a_hat = n_A_terms*(n_A_terms+1)/2;
349  auto aq_aq_innerprods_list =
350  rb_evaluation_builder.initAqAqInnerprods(Q_a_hat*n_bfs*n_bfs);
351 
352  for (unsigned int i=0; i < Q_a_hat; ++i)
353  for (unsigned int j=0; j < n_bfs; ++j)
354  for (unsigned int l=0; l < n_bfs; ++l)
355  {
356  unsigned int offset = i*n_bfs*n_bfs + j*n_bfs + l;
357  set_scalar_in_list(
358  aq_aq_innerprods_list,
359  offset,
360  rb_eval.Aq_Aq_representor_innerprods[i][j][l]);
361  }
362  }
363 
364  // Output dual inner-product data, and output vectors
365  {
366  unsigned int n_outputs = rb_theta_expansion.get_n_outputs();
367  auto output_innerprod_outer = rb_evaluation_builder.initOutputDualInnerprods(n_outputs);
368  auto output_vector_outer = rb_evaluation_builder.initOutputVectors(n_outputs);
369 
370  for (unsigned int output_id=0; output_id < n_outputs; ++output_id)
371  {
372  unsigned int n_output_terms = rb_theta_expansion.get_n_output_terms(output_id);
373 
374  {
375  unsigned int Q_l_hat = n_output_terms*(n_output_terms+1)/2;
376  auto output_innerprod_inner = output_innerprod_outer.init(output_id, Q_l_hat);
377  for (unsigned int q=0; q < Q_l_hat; ++q)
378  {
379  set_scalar_in_list(
380  output_innerprod_inner, q, rb_eval.output_dual_innerprods[output_id][q]);
381  }
382  }
383 
384  {
385  auto output_vector_middle = output_vector_outer.init(output_id, n_output_terms);
386  for (unsigned int q_l=0; q_l<n_output_terms; ++q_l)
387  {
388  auto output_vector_inner = output_vector_middle.init(q_l, n_bfs);
389  for (unsigned int j=0; j<n_bfs; ++j)
390  {
391  set_scalar_in_list(
392  output_vector_inner, j, rb_eval.RB_output_vectors[output_id][q_l](j));
393  }
394  }
395  }
396  }
397  }
398 
399  // Fq vectors and Aq matrices
400  {
401  unsigned int n_F_terms = rb_theta_expansion.get_n_F_terms();
402  unsigned int n_A_terms = rb_theta_expansion.get_n_A_terms();
403 
404  auto rb_fq_vectors_outer_list = rb_evaluation_builder.initRbFqVectors(n_F_terms);
405  for (unsigned int q_f=0; q_f < n_F_terms; ++q_f)
406  {
407  auto rb_fq_vectors_inner_list = rb_fq_vectors_outer_list.init(q_f, n_bfs);
408  for (unsigned int i=0; i<n_bfs; i++)
409  set_scalar_in_list(rb_fq_vectors_inner_list, i, rb_eval.RB_Fq_vector[q_f](i));
410  }
411 
412  auto rb_Aq_matrices_outer_list = rb_evaluation_builder.initRbAqMatrices(n_A_terms);
413  for (unsigned int q_a=0; q_a < n_A_terms; ++q_a)
414  {
415  auto rb_Aq_matrices_inner_list = rb_Aq_matrices_outer_list.init(q_a, n_bfs*n_bfs);
416  for (unsigned int i=0; i < n_bfs; ++i)
417  for (unsigned int j=0; j < n_bfs; ++j)
418  {
419  unsigned int offset = i*n_bfs+j;
420  set_scalar_in_list(rb_Aq_matrices_inner_list, offset, rb_eval.RB_Aq_vector[q_a](i,j));
421  }
422  }
423  }
424 
425  // Inner-product matrix
426  if (rb_eval.compute_RB_inner_product)
427  {
428  auto rb_inner_product_matrix_list =
429  rb_evaluation_builder.initRbInnerProductMatrix(n_bfs*n_bfs);
430 
431  for (unsigned int i=0; i < n_bfs; ++i)
432  for (unsigned int j=0; j < n_bfs; ++j)
433  {
434  unsigned int offset = i*n_bfs + j;
435  set_scalar_in_list(
436  rb_inner_product_matrix_list,
437  offset,
438  rb_eval.RB_inner_product_matrix(i,j) );
439  }
440  }
441 
442  auto parameter_ranges_list =
443  rb_evaluation_builder.initParameterRanges();
444  auto discrete_parameters_list =
445  rb_evaluation_builder.initDiscreteParameters();
447  parameter_ranges_list,
448  discrete_parameters_list);
449 }
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.

◆ add_rb_scm_evaluation_data_to_builder()

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 858 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(), libMesh::RBSCMEvaluation::get_SCM_UB_vector(), and libMesh::index_range().

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

860 {
861  auto parameter_ranges_list =
862  rb_scm_eval_builder.initParameterRanges();
863  auto discrete_parameters_list =
864  rb_scm_eval_builder.initDiscreteParameters();
866  parameter_ranges_list,
867  discrete_parameters_list);
868 
869  {
870  libmesh_error_msg_if(rb_scm_eval.B_min.size() != rb_scm_eval.get_rb_theta_expansion().get_n_A_terms(),
871  "Size error while writing B_min");
872  auto b_min_list = rb_scm_eval_builder.initBMin( rb_scm_eval.B_min.size() );
873  for (auto i : index_range(rb_scm_eval.B_min))
874  b_min_list.set(i, rb_scm_eval.get_B_min(i));
875  }
876 
877  {
878  libmesh_error_msg_if(rb_scm_eval.B_max.size() != rb_scm_eval.get_rb_theta_expansion().get_n_A_terms(),
879  "Size error while writing B_max");
880 
881  auto b_max_list = rb_scm_eval_builder.initBMax( rb_scm_eval.B_max.size() );
882  for (auto i : index_range(rb_scm_eval.B_max))
883  b_max_list.set(i, rb_scm_eval.get_B_max(i));
884  }
885 
886  {
887  auto cj_stability_vector =
888  rb_scm_eval_builder.initCJStabilityVector( rb_scm_eval.C_J_stability_vector.size() );
889  for (auto i : index_range(rb_scm_eval.C_J_stability_vector))
890  cj_stability_vector.set(i, rb_scm_eval.get_C_J_stability_constraint(i));
891  }
892 
893  {
894  auto cj_parameters_outer =
895  rb_scm_eval_builder.initCJ( rb_scm_eval.C_J.size() );
896 
897  for (auto i : index_range(rb_scm_eval.C_J))
898  {
899  auto cj_parameters_inner =
900  cj_parameters_outer.init(i, rb_scm_eval.C_J[i].n_parameters());
901 
902  unsigned int count = 0;
903  for (const auto & [key, val] : rb_scm_eval.C_J[i])
904  {
905  cj_parameters_inner[count].setName(key);
906  cj_parameters_inner[count].setValue(val);
907  count++;
908  }
909 
910  }
911  }
912 
913  {
914  unsigned int n_C_J_values = rb_scm_eval.C_J.size();
915  unsigned int n_A_terms = rb_scm_eval.get_rb_theta_expansion().get_n_A_terms();
916  unsigned int n_values = n_C_J_values*n_A_terms;
917  auto scm_ub_vectors =
918  rb_scm_eval_builder.initScmUbVectors( n_values );
919 
920  for (unsigned int i=0; i<n_C_J_values; i++)
921  for (unsigned int j=0; j<n_A_terms; j++)
922  {
923  unsigned int offset = i*n_A_terms + j;
924  scm_ub_vectors.set(offset, rb_scm_eval.get_SCM_UB_vector(i,j));
925  }
926  }
927 }
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.
auto index_range(const T &sizable)
Helper function that returns an IntRange<std::size_t> representing all the indices of the passed-in v...
Definition: int_range.h:111

◆ add_transient_rb_evaluation_data_to_builder()

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 452 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::RBThetaExpansion::get_n_A_terms(), libMesh::RBEvaluation::get_n_basis_functions(), libMesh::RBThetaExpansion::get_n_F_terms(), 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().

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