20 #include "libmesh/libmesh_config.h" 21 #if defined(LIBMESH_HAVE_CAPNPROTO) 24 #include "libmesh/rb_data_serialization.h" 25 #include "libmesh/rb_eim_evaluation.h" 26 #include "libmesh/enum_to_string.h" 27 #include "libmesh/transient_rb_theta_expansion.h" 28 #include "libmesh/rb_evaluation.h" 29 #include "libmesh/transient_rb_evaluation.h" 30 #include "libmesh/rb_eim_evaluation.h" 31 #include "libmesh/rb_scm_evaluation.h" 32 #include "libmesh/elem.h" 33 #include "libmesh/int_range.h" 34 #include "libmesh/rb_parametrized_function.h" 35 #include "libmesh/libmesh_logging.h" 38 #include <capnp/serialize.h> 39 #include <capnp/serialize-packed.h> 44 #ifdef LIBMESH_HAVE_UNISTD_H 59 template <
typename T,
typename U>
60 void set_scalar_in_list(T list,
unsigned int i, U
value)
62 #ifdef LIBMESH_USE_COMPLEX_NUMBERS 72 namespace RBDataSerialization
87 LOG_SCOPE(
"write_to_file()",
"RBEvaluationSerialization");
91 capnp::MallocMessageBuilder message;
93 #ifndef LIBMESH_USE_COMPLEX_NUMBERS 94 RBData::RBEvaluationReal::Builder rb_eval_builder =
95 message.initRoot<RBData::RBEvaluationReal>();
97 RBData::RBEvaluationComplex::Builder rb_eval_builder =
98 message.initRoot<RBData::RBEvaluationComplex>();
103 int fd = open(path.c_str(), O_WRONLY | O_CREAT | O_TRUNC, 0664);
104 libmesh_error_msg_if(!fd,
"Error opening a write-only file descriptor to " + path);
107 capnp::writePackedMessageToFd(fd, message);
109 capnp::writeMessageToFd(fd, message);
111 int error = close(fd);
112 libmesh_error_msg_if(error,
"Error closing a write-only file descriptor to " + path);
123 _trans_rb_eval(trans_rb_eval)
131 LOG_SCOPE(
"write_to_file()",
"TransientRBEvaluationSerialization");
135 capnp::MallocMessageBuilder message;
137 #ifndef LIBMESH_USE_COMPLEX_NUMBERS 138 RBData::TransientRBEvaluationReal::Builder trans_rb_eval_builder =
139 message.initRoot<RBData::TransientRBEvaluationReal>();
140 RBData::RBEvaluationReal::Builder rb_eval_builder =
141 trans_rb_eval_builder.initRbEvaluation();
143 RBData::TransientRBEvaluationComplex::Builder trans_rb_eval_builder =
144 message.initRoot<RBData::TransientRBEvaluationComplex>();
145 RBData::RBEvaluationComplex::Builder rb_eval_builder =
146 trans_rb_eval_builder.initRbEvaluation();
151 trans_rb_eval_builder);
153 int fd = open(path.c_str(), O_WRONLY | O_CREAT | O_TRUNC, 0664);
154 libmesh_error_msg_if(!fd,
"Error opening a write-only file descriptor to " + path);
157 capnp::writePackedMessageToFd(fd, message);
159 capnp::writeMessageToFd(fd, message);
161 int error = close(fd);
162 libmesh_error_msg_if(error,
"Error closing a write-only file descriptor to " + path);
173 _rb_eim_eval(rb_eim_eval)
180 LOG_SCOPE(
"write_to_file()",
"RBEIMEvaluationSerialization");
184 capnp::MallocMessageBuilder message;
186 #ifndef LIBMESH_USE_COMPLEX_NUMBERS 187 RBData::RBEIMEvaluationReal::Builder rb_eim_eval_builder =
188 message.initRoot<RBData::RBEIMEvaluationReal>();
190 RBData::RBEIMEvaluationComplex::Builder rb_eim_eval_builder =
191 message.initRoot<RBData::RBEIMEvaluationComplex>();
195 rb_eim_eval_builder);
197 int fd = open(path.c_str(), O_WRONLY | O_CREAT | O_TRUNC, 0664);
198 libmesh_error_msg_if(!fd,
"Error opening a write-only file descriptor to " + path);
201 capnp::writePackedMessageToFd(fd, message);
203 capnp::writeMessageToFd(fd, message);
205 int error = close(fd);
206 libmesh_error_msg_if(error,
"Error closing a write-only file descriptor to " + path);
215 #if defined(LIBMESH_HAVE_SLEPC) && (LIBMESH_HAVE_GLPK) 219 _rb_scm_eval(rb_scm_eval)
227 LOG_SCOPE(
"write_to_file()",
"RBSCMEvaluationSerialization");
231 capnp::MallocMessageBuilder message;
233 RBData::RBSCMEvaluation::Builder rb_scm_eval_builder =
234 message.initRoot<RBData::RBSCMEvaluation>();
238 int fd = open(path.c_str(), O_WRONLY | O_CREAT | O_TRUNC, 0664);
239 libmesh_error_msg_if(!fd,
"Error opening a write-only file descriptor to " + path);
242 capnp::writePackedMessageToFd(fd, message);
244 capnp::writeMessageToFd(fd, message);
246 int error = close(fd);
247 libmesh_error_msg_if(error,
"Error closing a write-only file descriptor to " + path);
251 #endif // LIBMESH_HAVE_SLEPC && LIBMESH_HAVE_GLPK 259 RBData::ParameterRanges::Builder & parameter_ranges_list,
260 RBData::DiscreteParameterList::Builder & discrete_parameters_list)
265 auto names = parameter_ranges_list.initNames(n_continuous_parameters);
266 auto mins = parameter_ranges_list.initMinValues(n_continuous_parameters);
267 auto maxs = parameter_ranges_list.initMaxValues(n_continuous_parameters);
273 unsigned int count = 0;
274 for (
const auto & [key, val] : parameters_min)
277 names.set(count, key);
278 mins.set(count, parameters_min.get_value(key));
279 maxs.set(count, parameters_max.
get_value(key));
283 libmesh_error_msg_if(count != n_continuous_parameters,
"Mismatch in number of continuous parameters");
289 auto names = discrete_parameters_list.initNames(n_discrete_parameters);
290 auto values_outer = discrete_parameters_list.initValues(n_discrete_parameters);
292 const std::map<std::string, std::vector<Real>> & discrete_parameters =
295 unsigned int count = 0;
296 for (
const auto & discrete_parameter : discrete_parameters)
298 names.set(count, discrete_parameter.first);
300 const std::vector<Real> & values = discrete_parameter.second;
301 unsigned int n_values = values.size();
303 values_outer.init(count, n_values);
304 auto values_inner = values_outer[count];
305 for (
unsigned int i=0; i<n_values; ++i)
307 values_inner.set(i, values[i]);
313 libmesh_error_msg_if(count != n_discrete_parameters,
"Mismatch in number of discrete parameters");
317 template <
typename RBEvaluationBuilderNumber>
319 RBEvaluationBuilderNumber & rb_evaluation_builder)
328 rb_evaluation_builder.setNBfs(n_bfs);
332 unsigned int Q_f_hat = n_F_terms*(n_F_terms+1)/2;
334 auto fq_innerprods_list = rb_evaluation_builder.initFqInnerprods(Q_f_hat);
336 for (
unsigned int i=0; i<Q_f_hat; i++)
337 set_scalar_in_list(fq_innerprods_list,
344 auto fq_aq_innerprods_list =
345 rb_evaluation_builder.initFqAqInnerprods(n_F_terms*n_A_terms*n_bfs);
347 for (
unsigned int q_f=0; q_f < n_F_terms; ++q_f)
348 for (
unsigned int q_a=0; q_a < n_A_terms; ++q_a)
349 for (
unsigned int i=0; i < n_bfs; ++i)
351 unsigned int offset = q_f*n_A_terms*n_bfs + q_a*n_bfs + i;
353 fq_aq_innerprods_list, offset,
360 unsigned int Q_a_hat = n_A_terms*(n_A_terms+1)/2;
361 auto aq_aq_innerprods_list =
362 rb_evaluation_builder.initAqAqInnerprods(Q_a_hat*n_bfs*n_bfs);
364 for (
unsigned int i=0; i < Q_a_hat; ++i)
365 for (
unsigned int j=0; j < n_bfs; ++j)
366 for (
unsigned int l=0; l < n_bfs; ++l)
368 unsigned int offset = i*n_bfs*n_bfs + j*n_bfs + l;
370 aq_aq_innerprods_list,
379 auto output_innerprod_outer = rb_evaluation_builder.initOutputDualInnerprods(n_outputs);
380 auto output_vector_outer = rb_evaluation_builder.initOutputVectors(n_outputs);
382 for (
unsigned int output_id=0; output_id < n_outputs; ++output_id)
387 unsigned int Q_l_hat = n_output_terms*(n_output_terms+1)/2;
388 auto output_innerprod_inner = output_innerprod_outer.init(output_id, Q_l_hat);
389 for (
unsigned int q=0; q < Q_l_hat; ++q)
397 auto output_vector_middle = output_vector_outer.init(output_id, n_output_terms);
398 for (
unsigned int q_l=0; q_l<n_output_terms; ++q_l)
400 auto output_vector_inner = output_vector_middle.init(q_l, n_bfs);
401 for (
unsigned int j=0; j<n_bfs; ++j)
416 auto rb_fq_vectors_outer_list = rb_evaluation_builder.initRbFqVectors(n_F_terms);
417 for (
unsigned int q_f=0; q_f < n_F_terms; ++q_f)
419 auto rb_fq_vectors_inner_list = rb_fq_vectors_outer_list.init(q_f, n_bfs);
420 for (
unsigned int i=0; i<n_bfs; i++)
421 set_scalar_in_list(rb_fq_vectors_inner_list, i, rb_eval.
RB_Fq_vector[q_f](i));
424 auto rb_Aq_matrices_outer_list = rb_evaluation_builder.initRbAqMatrices(n_A_terms);
425 for (
unsigned int q_a=0; q_a < n_A_terms; ++q_a)
427 auto rb_Aq_matrices_inner_list = rb_Aq_matrices_outer_list.init(q_a, n_bfs*n_bfs);
428 for (
unsigned int i=0; i < n_bfs; ++i)
429 for (
unsigned int j=0; j < n_bfs; ++j)
431 unsigned int offset = i*n_bfs+j;
432 set_scalar_in_list(rb_Aq_matrices_inner_list, offset, rb_eval.
RB_Aq_vector[q_a](i,j));
440 auto rb_inner_product_matrix_list =
441 rb_evaluation_builder.initRbInnerProductMatrix(n_bfs*n_bfs);
443 for (
unsigned int i=0; i < n_bfs; ++i)
444 for (
unsigned int j=0; j < n_bfs; ++j)
446 unsigned int offset = i*n_bfs + j;
448 rb_inner_product_matrix_list,
454 auto parameter_ranges_list =
455 rb_evaluation_builder.initParameterRanges();
456 auto discrete_parameters_list =
457 rb_evaluation_builder.initDiscreteParameters();
459 parameter_ranges_list,
460 discrete_parameters_list);
463 template <
typename RBEvaluationBuilderNumber,
typename TransRBEvaluationBuilderNumber>
465 RBEvaluationBuilderNumber & rb_eval_builder,
466 TransRBEvaluationBuilderNumber & trans_rb_eval_builder)
470 trans_rb_eval_builder.setDeltaT(trans_rb_eval.
get_delta_t());
473 trans_rb_eval_builder.setTimeStep(trans_rb_eval.
get_time_step());
479 auto rb_L2_matrix_list =
480 trans_rb_eval_builder.initRbL2Matrix(n_bfs*n_bfs);
482 for (
unsigned int i=0; i<n_bfs; ++i)
483 for (
unsigned int j=0; j<n_bfs; ++j)
485 unsigned int offset = i*n_bfs + j;
486 set_scalar_in_list(rb_L2_matrix_list,
494 unsigned int n_M_terms = trans_theta_expansion.
get_n_M_terms();
497 auto rb_Mq_matrices_outer_list = trans_rb_eval_builder.initRbMqMatrices(n_M_terms);
498 for (
unsigned int q_m=0; q_m < n_M_terms; ++q_m)
500 auto rb_Mq_matrices_inner_list = rb_Mq_matrices_outer_list.init(q_m, n_bfs*n_bfs);
501 for (
unsigned int i=0; i < n_bfs; ++i)
502 for (
unsigned int j=0; j < n_bfs; ++j)
504 unsigned int offset = i*n_bfs+j;
505 set_scalar_in_list(rb_Mq_matrices_inner_list,
515 auto initial_l2_errors_builder =
516 trans_rb_eval_builder.initInitialL2Errors(n_bfs);
517 auto initial_conditions_outer_list =
518 trans_rb_eval_builder.initInitialConditions(n_bfs);
520 for (
unsigned int i=0; i<n_bfs; i++)
524 auto initial_conditions_inner_list =
525 initial_conditions_outer_list.init(i, i+1);
526 for (
unsigned int j=0; j<=i; j++)
528 set_scalar_in_list(initial_conditions_inner_list,
537 unsigned int n_F_terms = trans_theta_expansion.
get_n_F_terms();
538 auto fq_mq_innerprods_list =
539 trans_rb_eval_builder.initFqMqInnerprods(n_F_terms*n_M_terms*n_bfs);
541 for (
unsigned int q_f=0; q_f<n_F_terms; ++q_f)
542 for (
unsigned int q_m=0; q_m<n_M_terms; ++q_m)
543 for (
unsigned int i=0; i<n_bfs; ++i)
545 unsigned int offset = q_f*n_M_terms*n_bfs + q_m*n_bfs + i;
546 set_scalar_in_list(fq_mq_innerprods_list,
554 unsigned int Q_m_hat = n_M_terms*(n_M_terms+1)/2;
555 auto mq_mq_innerprods_list =
556 trans_rb_eval_builder.initMqMqInnerprods(Q_m_hat*n_bfs*n_bfs);
558 for (
unsigned int i=0; i < Q_m_hat; ++i)
559 for (
unsigned int j=0; j < n_bfs; ++j)
560 for (
unsigned int l=0; l < n_bfs; ++l)
562 unsigned int offset = i*n_bfs*n_bfs + j*n_bfs + l;
563 set_scalar_in_list(mq_mq_innerprods_list,
571 unsigned int n_A_terms = trans_theta_expansion.
get_n_A_terms();
573 auto aq_mq_innerprods_list =
574 trans_rb_eval_builder.initAqMqInnerprods(n_A_terms*n_M_terms*n_bfs*n_bfs);
576 for (
unsigned int q_a=0; q_a<n_A_terms; q_a++)
577 for (
unsigned int q_m=0; q_m<n_M_terms; q_m++)
578 for (
unsigned int i=0; i<n_bfs; i++)
579 for (
unsigned int j=0; j<n_bfs; j++)
581 unsigned int offset =
582 q_a*(n_M_terms*n_bfs*n_bfs) + q_m*(n_bfs*n_bfs) + i*n_bfs + j;
583 set_scalar_in_list(aq_mq_innerprods_list,
591 template <
typename RBEIMEvaluationBuilderNumber>
593 RBEIMEvaluationBuilderNumber & rb_eim_evaluation_builder)
597 rb_eim_evaluation_builder.setNBfs(n_bfs);
603 unsigned int half_matrix_size = n_bfs*(n_bfs+1)/2;
605 auto interpolation_matrix_list =
606 rb_eim_evaluation_builder.initInterpolationMatrix(half_matrix_size);
607 for (
unsigned int i=0; i < n_bfs; ++i)
608 for (
unsigned int j=0; j <= i; ++j)
610 unsigned int offset = i*(i+1)/2 + j;
611 set_scalar_in_list(interpolation_matrix_list,
623 bool use_error_indicator =
627 if (use_error_indicator)
629 auto error_indicator_data_list =
630 rb_eim_evaluation_builder.initEimErrorIndicatorInterpData(n_bfs);
631 const auto & error_indicator_row =
633 for (
unsigned int i=0; i < n_bfs; ++i)
635 set_scalar_in_list(error_indicator_data_list,
637 error_indicator_row(i));
648 if (use_error_indicator)
652 "Number of basis functions should match number of interpolation points");
654 auto parameter_ranges_list =
655 rb_eim_evaluation_builder.initParameterRanges();
656 auto discrete_parameters_list =
657 rb_eim_evaluation_builder.initDiscreteParameters();
659 parameter_ranges_list,
660 discrete_parameters_list);
664 auto interpolation_points_list =
665 rb_eim_evaluation_builder.initInterpolationXyz(n_bfs);
666 for (
unsigned int i=0; i < n_bfs; ++i)
668 interpolation_points_list[i]);
673 auto interpolation_points_comp_list =
674 rb_eim_evaluation_builder.initInterpolationComp(n_bfs);
675 for (
unsigned int i=0; i<n_bfs; ++i)
676 interpolation_points_comp_list.set(i,
682 auto interpolation_points_subdomain_id_list =
683 rb_eim_evaluation_builder.initInterpolationSubdomainId(n_bfs);
684 for (
unsigned int i=0; i<n_bfs; ++i)
685 interpolation_points_subdomain_id_list.set(i,
694 auto interpolation_points_boundary_id_list =
695 rb_eim_evaluation_builder.initInterpolationBoundaryId(n_bfs);
696 for (
unsigned int i=0; i<n_bfs; ++i)
697 interpolation_points_boundary_id_list.set(i,
703 auto interpolation_points_elem_id_list =
704 rb_eim_evaluation_builder.initInterpolationElemId(n_bfs);
705 for (
unsigned int i=0; i<n_bfs; ++i)
706 interpolation_points_elem_id_list.set(i,
713 auto interpolation_points_node_id_list =
714 rb_eim_evaluation_builder.initInterpolationNodeId(n_bfs);
715 for (
unsigned int i=0; i<n_bfs; ++i)
716 interpolation_points_node_id_list.set(i,
723 auto interpolation_points_side_index_list =
724 rb_eim_evaluation_builder.initInterpolationSideIndex(n_bfs);
725 for (
unsigned int i=0; i<n_bfs; ++i)
726 interpolation_points_side_index_list.set(i,
732 auto interpolation_points_qp_list =
733 rb_eim_evaluation_builder.initInterpolationQp(n_bfs);
734 for (
unsigned int i=0; i<n_bfs; ++i)
735 interpolation_points_qp_list.set(i,
741 auto interpolation_points_list_outer =
742 rb_eim_evaluation_builder.initInterpolationXyzPerturb(n_bfs);
743 for (
unsigned int i=0; i < n_bfs; ++i)
746 auto interpolation_points_list_inner = interpolation_points_list_outer.init(i, perturbs.size());
755 unsigned int n_elems = rb_eim_evaluation.
get_n_elems();
756 rb_eim_evaluation_builder.setNElems(n_elems);
759 auto elem_id_to_local_index_list =
760 rb_eim_evaluation_builder.initElemIdToLocalIndex(n_elems);
761 unsigned int counter = 0;
764 elem_id_to_local_index_list[counter].setFirst(elem_id);
765 elem_id_to_local_index_list[counter].setSecond(local_index);
772 auto interpolation_points_list_outer =
773 rb_eim_evaluation_builder.initInterpolationJxWAllQp(n_elems);
774 for (
unsigned int i=0; i < n_elems; ++i)
777 auto interpolation_points_list_inner = interpolation_points_list_outer.init(i, JxW.size());
783 interpolation_points_list_inner.set(j, JxW[j]);
790 auto interpolation_points_list_outer =
791 rb_eim_evaluation_builder.initInterpolationPhiValuesAllQp(n_elems);
793 for (
unsigned int i=0; i < n_elems; ++i)
796 auto interpolation_points_list_middle = interpolation_points_list_outer.init(i, phi_i_all_qp.size());
800 auto interpolation_points_list_inner = interpolation_points_list_middle.init(j, phi_i_all_qp[j].size());
803 interpolation_points_list_inner.set(k, phi_i_all_qp[j][k]);
810 auto dxyzdxi_elem_center =
811 rb_eim_evaluation_builder.initInterpolationDxyzDxiElem(n_elems);
820 auto dxyzdeta_elem_center =
821 rb_eim_evaluation_builder.initInterpolationDxyzDetaElem(n_elems);
830 auto interpolation_points_qrule_order_list =
831 rb_eim_evaluation_builder.initInterpolationQruleOrder(n_elems);
832 for (
unsigned int i=0; i<n_elems; ++i)
833 interpolation_points_qrule_order_list.set(i,
839 auto interpolation_points_elem_type_list =
840 rb_eim_evaluation_builder.initInterpolationElemType(n_bfs);
841 for (
unsigned int i=0; i<n_bfs; ++i)
842 interpolation_points_elem_type_list.set(i,
849 auto interpolation_points_property_list =
850 rb_eim_evaluation_builder.initPropertyMap(n_properties);
851 unsigned int property_counter = 0;
854 interpolation_points_property_list[property_counter].setName(property_name);
856 unsigned int entity_counter = 0;
857 auto property_entity_ids = interpolation_points_property_list[property_counter].initEntityIds(entity_ids.size());
858 for (
const auto entity_id : entity_ids)
860 property_entity_ids.set(entity_counter, entity_id);
872 auto eim_rhs_list_outer =
873 rb_eim_evaluation_builder.initEimSolutionsForTrainingSet(eim_solutions.size());
874 for (
auto i :
make_range(eim_solutions.size()))
877 auto eim_rhs_list_inner = eim_rhs_list_outer.init(i, values.
size());
881 set_scalar_in_list(eim_rhs_list_inner,
891 auto interpolation_points_list_outer =
892 rb_eim_evaluation_builder.initInterpolationPhiValues(n_bfs);
893 for (
unsigned int i=0; i < n_bfs; ++i)
896 auto interpolation_points_list_inner = interpolation_points_list_outer.init(i, phi_i_qp_vec.size());
902 interpolation_points_list_inner.set(j, phi_i_qp_vec[j]);
913 "Error: Number of spatial indices should match number of EIM basis functions");
915 auto interpolation_points_spatial_indices_builder =
916 rb_eim_evaluation_builder.initInterpolationSpatialIndices(n_bfs);
917 for (
unsigned int i=0; i<n_bfs; ++i)
919 const std::vector<unsigned int> & spatial_indices =
921 unsigned int n_spatial_indices = spatial_indices.size();
923 auto spatial_indices_builder =
924 interpolation_points_spatial_indices_builder.init(i, n_spatial_indices);
928 spatial_indices_builder.set(j, spatial_indices[j]);
934 #if defined(LIBMESH_HAVE_SLEPC) && (LIBMESH_HAVE_GLPK) 936 RBData::RBSCMEvaluation::Builder & rb_scm_eval_builder)
938 auto parameter_ranges_list =
939 rb_scm_eval_builder.initParameterRanges();
940 auto discrete_parameters_list =
941 rb_scm_eval_builder.initDiscreteParameters();
943 parameter_ranges_list,
944 discrete_parameters_list);
948 "Size error while writing B_min");
949 auto b_min_list = rb_scm_eval_builder.initBMin( rb_scm_eval.
B_min.size() );
951 b_min_list.set(i, rb_scm_eval.
get_B_min(i));
956 "Size error while writing B_max");
958 auto b_max_list = rb_scm_eval_builder.initBMax( rb_scm_eval.
B_max.size() );
960 b_max_list.set(i, rb_scm_eval.
get_B_max(i));
964 auto cj_stability_vector =
971 auto cj_parameters_outer =
972 rb_scm_eval_builder.initCJ( rb_scm_eval.
C_J.size() );
976 auto cj_parameters_inner =
977 cj_parameters_outer.init(i, rb_scm_eval.
C_J[i].n_parameters());
979 unsigned int count = 0;
980 for (
const auto & [key, val] : rb_scm_eval.
C_J[i])
982 cj_parameters_inner[count].setName(key);
983 cj_parameters_inner[count].setValue(val);
991 unsigned int n_C_J_values = rb_scm_eval.
C_J.size();
993 unsigned int n_values = n_C_J_values*n_A_terms;
994 auto scm_ub_vectors =
995 rb_scm_eval_builder.initScmUbVectors( n_values );
997 for (
unsigned int i=0; i<n_C_J_values; i++)
998 for (
unsigned int j=0; j<n_A_terms; j++)
1000 unsigned int offset = i*n_A_terms + j;
1005 #endif // LIBMESH_HAVE_SLEPC && LIBMESH_HAVE_GLPK 1009 point_builder.setX(point(0));
1011 if (LIBMESH_DIM >= 2)
1012 point_builder.setY(point(1));
1014 if (LIBMESH_DIM >= 3)
1015 point_builder.setZ(point(2));
1024 #endif // #if defined(LIBMESH_HAVE_CAPNPROTO) std::vector< Real > B_min
B_min, B_max define the bounding box.
Real get_value(const std::string ¶m_name) const
Get the value of the specified parameter, throw an error if it does not exist.
std::vector< std::vector< std::vector< Number > > > Mq_Mq_representor_innerprods
void write_to_file(const std::string &path, bool use_packing=false)
Write the Cap'n'Proto buffer to disk.
bool compute_RB_inner_product
Boolean flag to indicate whether we compute the RB_inner_product_matrix.
RBSCMEvaluationSerialization(RBSCMEvaluation &rb_eval)
Initialize a new buffer using the structure from the Cap'n'Proto schema described in rb_data...
virtual ~RBEvaluationSerialization()
std::vector< Number > Fq_representor_innerprods
Vectors storing the residual representor inner products to be used in computing the residuals online...
boost::multiprecision::float128 real(const boost::multiprecision::float128 in)
const DenseMatrix< Number > & get_interpolation_matrix() const
Get the EIM interpolation matrix.
unsigned int get_n_F_terms() const
Get Q_f, the number of terms in the affine expansion for the right-hand side.
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.
subdomain_id_type get_interpolation_points_subdomain_id(unsigned int index) const
Order get_interpolation_points_qrule_order(unsigned int index) const
std::vector< DenseMatrix< Number > > RB_M_q_vector
Dense matrices for the RB mass matrices.
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.
unsigned int get_n_basis_functions() const
Return the current number of EIM basis functions.
DenseMatrix< Number > RB_inner_product_matrix
The inner product matrix.
Real get_B_min(unsigned int i) const
Get B_min and B_max.
unsigned int get_n_interpolation_points_spatial_indices() const
_interpolation_points_spatial_indices is optional data, so we need to be able to check how many _inte...
const RBParameters & get_parameters_max() const
Get an RBParameters object that specifies the maximum allowable value for each parameter.
const std::vector< std::vector< Real > > & get_interpolation_points_phi_i_all_qp(unsigned int index) const
processor_id_type rank() const
unsigned int get_n_A_terms() const
Get Q_a, the number of terms in the affine expansion for the bilinear form.
const Parallel::Communicator & comm() const
This class is part of the rbOOmit framework.
This class stores the set of RBTheta functor objects that define the "parameter-dependent expansion" ...
void add_parameter_ranges_to_builder(const RBParametrized &rb_evaluation, RBData::ParameterRanges::Builder ¶meter_ranges, RBData::DiscreteParameterList::Builder &discrete_parameters_list)
Add parameter ranges for continuous and discrete parameters.
std::vector< std::vector< Number > > output_dual_innerprods
The vector storing the dual norm inner product terms for each output.
void add_rb_evaluation_data_to_builder(RBEvaluation &rb_eval, RBEvaluationBuilderNumber &rb_eval_builder)
Add data for an RBEvaluation to the builder.
const std::vector< Point > & get_interpolation_points_xyz_perturbations(unsigned int index) const
unsigned int get_interpolation_points_side_index(unsigned int index) const
unsigned int get_n_interpolation_points() const
Return the number of interpolation points.
RBEIMEvaluation & _rb_eim_eval
The RBEvaluation object that will be written to disk.
The libMesh namespace provides an interface to certain functionality in the library.
unsigned int get_interpolation_points_comp(unsigned int index) const
const std::vector< DenseVector< Number > > & get_eim_solutions_for_training_set() const
Return a const reference to the EIM solutions for the parameters in the training set.
RBEvaluationSerialization(RBEvaluation &rb_eval)
Initialize a new buffer using the structure from the Cap'n'Proto schema described in rb_data...
Real get_delta_t() const
Get/set delta_t, the time-step size.
RBThetaExpansion & get_rb_theta_expansion()
Get a reference to the rb_theta_expansion.
std::vector< Real > B_max
unsigned int get_n_outputs() const
Get n_outputs, the number output functionals.
This class stores the set of RBTheta functor objects that define the "parameter-dependent expansion" ...
Real get_euler_theta() const
Get/set euler_theta, parameter that determines the temporal discretization.
unsigned int get_n_discrete_params() const
Get the number of discrete parameters.
virtual ~RBSCMEvaluationSerialization()
std::vector< Real > C_J_stability_vector
Vector storing the (truth) stability values at the parameters in C_J.
const std::vector< Real > & get_interpolation_points_JxW_all_qp(unsigned int index) const
void write_to_file(const std::string &path, bool use_packing=false)
Write the Cap'n'Proto buffer to disk.
TransientRBEvaluationSerialization(TransientRBEvaluation &rb_eval)
Initialize a new buffer using the structure from the Cap'n'Proto schema described in rb_data...
const Point & get_elem_center_dxyzdeta(unsigned int index) const
std::vector< DenseVector< Number > > RB_Fq_vector
Dense vector for the RHS.
bool is_lookup_table
Boolean to indicate if this parametrized function is defined based on a lookup table or not...
std::vector< Real > initial_L2_error_all_N
Vector storing initial L2 error for all 1 <= N <= RB_size.
const RBParameters & get_parameters_min() const
Get an RBParameters object that specifies the minimum allowable value for each parameter.
boundary_id_type get_interpolation_points_boundary_id(unsigned int index) const
unsigned int get_n_properties() const
Return the number of properties stored in the rb_property_map.
std::vector< RBParameters > C_J
Vector storing the greedily selected parameters during SCM training.
const std::unordered_map< std::string, std::set< dof_id_type > > & get_rb_property_map() const
virtual unsigned int get_n_M_terms() const
Get Q_m, the number of terms in the affine expansion for the mass operator.
const std::vector< unsigned int > & get_interpolation_points_spatial_indices(unsigned int index) const
const std::map< std::string, std::vector< Real > > & get_discrete_parameter_values() const
Get a const reference to the discrete parameter values.
const DenseVector< Number > & get_error_indicator_interpolation_row() const
Get/set _extra_points_interpolation_matrix.
ElemType get_interpolation_points_elem_type(unsigned int index) const
std::vector< std::vector< std::vector< Number > > > Fq_Mq_representor_innerprods
Vectors storing the residual representor inner products to be used in computing the residuals online...
bool on_mesh_nodes() const
void write_to_file(const std::string &path, bool use_packing=false)
Write the Cap'n'Proto buffer to disk.
unsigned int get_time_step() const
Get/set the current time-step.
This class is part of the rbOOmit framework.
unsigned int get_n_continuous_params() const
Get the number of continuous parameters.
RBEIMEvaluationSerialization(RBEIMEvaluation &rb_eval)
Initialize a new buffer using the structure from the Cap'n'Proto schema described in rb_data...
void add_point_to_builder(const Point &point, RBData::Point3D::Builder point_builder)
Helper function that adds point data.
std::vector< std::vector< std::vector< Number > > > Fq_Aq_representor_innerprods
Vectors storing the residual representor inner products to be used in computing the residuals online...
std::vector< DenseVector< Number > > RB_initial_condition_all_N
The RB initial conditions (i.e.
RBSCMEvaluation & _rb_scm_eval
The RBEvaluation object that will be written to disk.
This class is part of the rbOOmit framework.
Point get_interpolation_points_xyz(unsigned int index) const
Get the data associated with EIM interpolation points.
dof_id_type get_interpolation_points_node_id(unsigned int index) const
virtual bool use_eim_error_indicator() const
Virtual function to indicate if we use the EIM error indicator in this case.
std::vector< DenseMatrix< Number > > RB_Aq_vector
Dense matrices for the RB computations.
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.
const std::map< dof_id_type, unsigned int > & get_elem_id_to_local_index_map() const
void write_to_file(const std::string &path, bool use_packing=false)
Write the Cap'n'Proto buffer to disk.
unsigned int get_n_output_terms(unsigned int output_index) const
Get the number of affine terms associated with the specified output.
Real get_C_J_stability_constraint(unsigned int j) const
Get stability constraints (i.e.
This class is part of the rbOOmit framework.
virtual unsigned int get_n_basis_functions() const
Get the current number of basis functions.
Real get_B_max(unsigned int i) const
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...
dof_id_type get_interpolation_points_elem_id(unsigned int index) const
TransientRBEvaluation & _trans_rb_eval
The RBEvaluation object that will be written to disk.
virtual unsigned int size() const override final
unsigned int get_n_time_steps() const
Get/set the total number of time-steps.
virtual ~RBEIMEvaluationSerialization()
unsigned int get_n_elems() const
Return the number of unique elements containing interpolation points.
std::vector< std::vector< std::vector< std::vector< Number > > > > Aq_Mq_representor_innerprods
bool on_mesh_sides() const
boost::multiprecision::float128 imag(const boost::multiprecision::float128)
virtual ~TransientRBEvaluationSerialization()
RBEvaluation & _rb_eval
The RBEvaluation object that will be written to disk.
RBParametrizedFunction & get_parametrized_function()
Get a reference to the parametrized function.
const std::vector< Real > & get_interpolation_points_phi_i_qp(unsigned int index) const
This class enables evaluation of an Empirical Interpolation Method (EIM) approximation.
const Point & get_elem_center_dxyzdxi(unsigned int index) const
void initialize_interpolation_points_spatial_indices()
Initialize _interpolation_points_spatial_indices.
A Point defines a location in LIBMESH_DIM dimensional Real space.
auto index_range(const T &sizable)
Helper function that returns an IntRange<std::size_t> representing all the indices of the passed-in v...
RBThetaExpansion & get_rb_theta_expansion()
Get a reference to the rb_theta_expansion.
unsigned int get_interpolation_points_qp(unsigned int index) const
DenseMatrix< Number > RB_L2_matrix
Dense RB L2 matrix.
Real get_SCM_UB_vector(unsigned int j, unsigned int q)
Get entries of SCM_UB_vector, which stores the vector y, corresponding to the minimizing eigenvectors...
bool is_discrete_parameter(const std::string &mu_name) const
Is parameter mu_name discrete?
This class is part of the rbOOmit framework.
std::vector< std::vector< std::vector< Number > > > Aq_Aq_representor_innerprods
std::vector< std::vector< DenseVector< Number > > > RB_output_vectors
The vectors storing the RB output vectors.