libMesh
rb_data_serialization.C
Go to the documentation of this file.
1 // rbOOmit: An implementation of the Certified Reduced Basis method.
2 // Copyright (C) 2009, 2010, 2015 David J. Knezevic
3 
4 // This file is part of rbOOmit.
5 
6 // rbOOmit is free software; you can redistribute it and/or
7 // modify it under the terms of the GNU Lesser General Public
8 // License as published by the Free Software Foundation; either
9 // version 2.1 of the License, or (at your option) any later version.
10 
11 // rbOOmit is distributed in the hope that it will be useful,
12 // but WITHOUT ANY WARRANTY; without even the implied warranty of
13 // MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
14 // Lesser General Public License for more details.
15 
16 // You should have received a copy of the GNU Lesser General Public
17 // License along with this library; if not, write to the Free Software
18 // Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA
19 
20 #include "libmesh/libmesh_config.h"
21 #if defined(LIBMESH_HAVE_CAPNPROTO)
22 
23 // libMesh includes
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"
36 
37 // Cap'n'Proto includes
38 #include <capnp/serialize.h>
39 #include <capnp/serialize-packed.h> // for writePackedMessageToFd()
40 
41 // C++ includes
42 #include <iostream>
43 #include <fstream>
44 #ifdef LIBMESH_HAVE_UNISTD_H
45 #include <unistd.h> // for close()
46 #endif
47 #include <fcntl.h>
48 
49 namespace libMesh
50 {
51 
52 namespace
53 {
54 
59 template <typename T, typename U>
60 void set_scalar_in_list(T list, unsigned int i, U value)
61 {
62 #ifdef LIBMESH_USE_COMPLEX_NUMBERS
63  list[i].setReal(std::real(value));
64  list[i].setImag(std::imag(value));
65 #else
66  list.set(i, value);
67 #endif
68 }
69 
70 }
71 
72 namespace RBDataSerialization
73 {
74 
75 // ---- RBEvaluationSerialization (BEGIN) ----
76 
78  :
79  _rb_eval(rb_eval)
80 {
81 }
82 
84 
85 void RBEvaluationSerialization::write_to_file(const std::string & path, bool use_packing)
86 {
87  LOG_SCOPE("write_to_file()", "RBEvaluationSerialization");
88 
89  if (_rb_eval.comm().rank() == 0)
90  {
91  capnp::MallocMessageBuilder message;
92 
93 #ifndef LIBMESH_USE_COMPLEX_NUMBERS
94  RBData::RBEvaluationReal::Builder rb_eval_builder =
95  message.initRoot<RBData::RBEvaluationReal>();
96 #else
97  RBData::RBEvaluationComplex::Builder rb_eval_builder =
98  message.initRoot<RBData::RBEvaluationComplex>();
99 #endif
100 
102 
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);
105 
106  if(use_packing)
107  capnp::writePackedMessageToFd(fd, message);
108  else
109  capnp::writeMessageToFd(fd, message);
110 
111  int error = close(fd);
112  libmesh_error_msg_if(error, "Error closing a write-only file descriptor to " + path);
113  }
114 }
115 
116 // ---- RBEvaluationSerialization (END) ----
117 
118 
119 // ---- TransientRBEvaluationSerialization (BEGIN) ----
120 
123  _trans_rb_eval(trans_rb_eval)
124 {
125 }
126 
128 
129 void TransientRBEvaluationSerialization::write_to_file(const std::string & path, bool use_packing)
130 {
131  LOG_SCOPE("write_to_file()", "TransientRBEvaluationSerialization");
132 
133  if (_trans_rb_eval.comm().rank() == 0)
134  {
135  capnp::MallocMessageBuilder message;
136 
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();
142 #else
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();
147 #endif
148 
150  rb_eval_builder,
151  trans_rb_eval_builder);
152 
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);
155 
156  if(use_packing)
157  capnp::writePackedMessageToFd(fd, message);
158  else
159  capnp::writeMessageToFd(fd, message);
160 
161  int error = close(fd);
162  libmesh_error_msg_if(error, "Error closing a write-only file descriptor to " + path);
163  }
164 }
165 
166 // ---- TransientRBEvaluationSerialization (END) ----
167 
168 
169 // ---- RBEIMEvaluationSerialization (BEGIN) ----
170 
172  :
173  _rb_eim_eval(rb_eim_eval)
174 {
175 }
176 
178 
179 void RBEIMEvaluationSerialization::write_to_file(const std::string & path, bool use_packing) {
180  LOG_SCOPE("write_to_file()", "RBEIMEvaluationSerialization");
181 
182  if (_rb_eim_eval.comm().rank() == 0)
183  {
184  capnp::MallocMessageBuilder message;
185 
186 #ifndef LIBMESH_USE_COMPLEX_NUMBERS
187  RBData::RBEIMEvaluationReal::Builder rb_eim_eval_builder =
188  message.initRoot<RBData::RBEIMEvaluationReal>();
189 #else
190  RBData::RBEIMEvaluationComplex::Builder rb_eim_eval_builder =
191  message.initRoot<RBData::RBEIMEvaluationComplex>();
192 #endif
193 
195  rb_eim_eval_builder);
196 
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);
199 
200  if(use_packing)
201  capnp::writePackedMessageToFd(fd, message);
202  else
203  capnp::writeMessageToFd(fd, message);
204 
205  int error = close(fd);
206  libmesh_error_msg_if(error, "Error closing a write-only file descriptor to " + path);
207  }
208 }
209 
210 // ---- RBEIMEvaluationSerialization (END) ----
211 
212 
213 // ---- RBSCMEvaluationSerialization (BEGIN) ----
214 
215 #if defined(LIBMESH_HAVE_SLEPC) && (LIBMESH_HAVE_GLPK)
216 
218  :
219  _rb_scm_eval(rb_scm_eval)
220 {
221 }
222 
224 
225 void RBSCMEvaluationSerialization::write_to_file(const std::string & path, bool use_packing)
226 {
227  LOG_SCOPE("write_to_file()", "RBSCMEvaluationSerialization");
228 
229  if (_rb_scm_eval.comm().rank() == 0)
230  {
231  capnp::MallocMessageBuilder message;
232 
233  RBData::RBSCMEvaluation::Builder rb_scm_eval_builder =
234  message.initRoot<RBData::RBSCMEvaluation>();
235 
237 
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);
240 
241  if(use_packing)
242  capnp::writePackedMessageToFd(fd, message);
243  else
244  capnp::writeMessageToFd(fd, message);
245 
246  int error = close(fd);
247  libmesh_error_msg_if(error, "Error closing a write-only file descriptor to " + path);
248  }
249 }
250 
251 #endif // LIBMESH_HAVE_SLEPC && LIBMESH_HAVE_GLPK
252 
253 // ---- RBSCMEvaluationSerialization (END) ----
254 
255 
256 // ---- Helper functions for adding data to capnp Builders (BEGIN) ----
257 
259  RBData::ParameterRanges::Builder & parameter_ranges_list,
260  RBData::DiscreteParameterList::Builder & discrete_parameters_list)
261 {
262  // Continuous parameters
263  {
264  unsigned int n_continuous_parameters = rb_evaluation.get_n_continuous_params();
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);
268 
269  const RBParameters & parameters_min = rb_evaluation.get_parameters_min();
270  const RBParameters & parameters_max = rb_evaluation.get_parameters_max();
271 
272  // We could loop over either parameters_min or parameters_max, they should have the same keys.
273  unsigned int count = 0;
274  for (const auto & [key, val] : parameters_min)
275  if (!rb_evaluation.is_discrete_parameter(key))
276  {
277  names.set(count, key);
278  mins.set(count, parameters_min.get_value(key));
279  maxs.set(count, parameters_max.get_value(key));
280  ++count;
281  }
282 
283  libmesh_error_msg_if(count != n_continuous_parameters, "Mismatch in number of continuous parameters");
284  }
285 
286  // Discrete parameters
287  {
288  unsigned int n_discrete_parameters = rb_evaluation.get_n_discrete_params();
289  auto names = discrete_parameters_list.initNames(n_discrete_parameters);
290  auto values_outer = discrete_parameters_list.initValues(n_discrete_parameters);
291 
292  const std::map<std::string, std::vector<Real>> & discrete_parameters =
293  rb_evaluation.get_discrete_parameter_values();
294 
295  unsigned int count = 0;
296  for (const auto & discrete_parameter : discrete_parameters)
297  {
298  names.set(count, discrete_parameter.first);
299 
300  const std::vector<Real> & values = discrete_parameter.second;
301  unsigned int n_values = values.size();
302 
303  values_outer.init(count, n_values);
304  auto values_inner = values_outer[count];
305  for (unsigned int i=0; i<n_values; ++i)
306  {
307  values_inner.set(i, values[i]);
308  }
309 
310  ++count;
311  }
312 
313  libmesh_error_msg_if(count != n_discrete_parameters, "Mismatch in number of discrete parameters");
314  }
315 }
316 
317 template <typename RBEvaluationBuilderNumber>
319  RBEvaluationBuilderNumber & rb_evaluation_builder)
320 {
321  const RBThetaExpansion & rb_theta_expansion = rb_eval.get_rb_theta_expansion();
322 
323  unsigned int n_F_terms = rb_theta_expansion.get_n_F_terms();
324  unsigned int n_A_terms = rb_theta_expansion.get_n_A_terms();
325 
326  // Number of basis functions
327  unsigned int n_bfs = rb_eval.get_n_basis_functions();
328  rb_evaluation_builder.setNBfs(n_bfs);
329 
330  // Fq representor inner-product data
331  {
332  unsigned int Q_f_hat = n_F_terms*(n_F_terms+1)/2;
333 
334  auto fq_innerprods_list = rb_evaluation_builder.initFqInnerprods(Q_f_hat);
335 
336  for (unsigned int i=0; i<Q_f_hat; i++)
337  set_scalar_in_list(fq_innerprods_list,
338  i,
339  rb_eval.Fq_representor_innerprods[i]);
340  }
341 
342  // FqAq representor inner-product data
343  {
344  auto fq_aq_innerprods_list =
345  rb_evaluation_builder.initFqAqInnerprods(n_F_terms*n_A_terms*n_bfs);
346 
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)
350  {
351  unsigned int offset = q_f*n_A_terms*n_bfs + q_a*n_bfs + i;
352  set_scalar_in_list(
353  fq_aq_innerprods_list, offset,
354  rb_eval.Fq_Aq_representor_innerprods[q_f][q_a][i]);
355  }
356  }
357 
358  // AqAq representor inner-product data
359  {
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);
363 
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)
367  {
368  unsigned int offset = i*n_bfs*n_bfs + j*n_bfs + l;
369  set_scalar_in_list(
370  aq_aq_innerprods_list,
371  offset,
372  rb_eval.Aq_Aq_representor_innerprods[i][j][l]);
373  }
374  }
375 
376  // Output dual inner-product data, and output vectors
377  {
378  unsigned int n_outputs = rb_theta_expansion.get_n_outputs();
379  auto output_innerprod_outer = rb_evaluation_builder.initOutputDualInnerprods(n_outputs);
380  auto output_vector_outer = rb_evaluation_builder.initOutputVectors(n_outputs);
381 
382  for (unsigned int output_id=0; output_id < n_outputs; ++output_id)
383  {
384  unsigned int n_output_terms = rb_theta_expansion.get_n_output_terms(output_id);
385 
386  {
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)
390  {
391  set_scalar_in_list(
392  output_innerprod_inner, q, rb_eval.output_dual_innerprods[output_id][q]);
393  }
394  }
395 
396  {
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)
399  {
400  auto output_vector_inner = output_vector_middle.init(q_l, n_bfs);
401  for (unsigned int j=0; j<n_bfs; ++j)
402  {
403  set_scalar_in_list(
404  output_vector_inner, j, rb_eval.RB_output_vectors[output_id][q_l](j));
405  }
406  }
407  }
408  }
409  }
410 
411  // Fq vectors and Aq matrices
412  {
413  unsigned int n_F_terms = rb_theta_expansion.get_n_F_terms();
414  unsigned int n_A_terms = rb_theta_expansion.get_n_A_terms();
415 
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)
418  {
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));
422  }
423 
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)
426  {
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)
430  {
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));
433  }
434  }
435  }
436 
437  // Inner-product matrix
438  if (rb_eval.compute_RB_inner_product)
439  {
440  auto rb_inner_product_matrix_list =
441  rb_evaluation_builder.initRbInnerProductMatrix(n_bfs*n_bfs);
442 
443  for (unsigned int i=0; i < n_bfs; ++i)
444  for (unsigned int j=0; j < n_bfs; ++j)
445  {
446  unsigned int offset = i*n_bfs + j;
447  set_scalar_in_list(
448  rb_inner_product_matrix_list,
449  offset,
450  rb_eval.RB_inner_product_matrix(i,j) );
451  }
452  }
453 
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);
461 }
462 
463 template <typename RBEvaluationBuilderNumber, typename TransRBEvaluationBuilderNumber>
465  RBEvaluationBuilderNumber & rb_eval_builder,
466  TransRBEvaluationBuilderNumber & trans_rb_eval_builder)
467 {
468  add_rb_evaluation_data_to_builder(trans_rb_eval, rb_eval_builder);
469 
470  trans_rb_eval_builder.setDeltaT(trans_rb_eval.get_delta_t());
471  trans_rb_eval_builder.setEulerTheta(trans_rb_eval.get_euler_theta());
472  trans_rb_eval_builder.setNTimeSteps(trans_rb_eval.get_n_time_steps());
473  trans_rb_eval_builder.setTimeStep(trans_rb_eval.get_time_step());
474 
475  unsigned int n_bfs = trans_rb_eval.get_n_basis_functions();
476 
477  // L2-inner-product matrix
478  {
479  auto rb_L2_matrix_list =
480  trans_rb_eval_builder.initRbL2Matrix(n_bfs*n_bfs);
481 
482  for (unsigned int i=0; i<n_bfs; ++i)
483  for (unsigned int j=0; j<n_bfs; ++j)
484  {
485  unsigned int offset = i*n_bfs + j;
486  set_scalar_in_list(rb_L2_matrix_list,
487  offset,
488  trans_rb_eval.RB_L2_matrix(i,j));
489  }
490  }
491 
492  TransientRBThetaExpansion & trans_theta_expansion =
493  cast_ref<TransientRBThetaExpansion &>(trans_rb_eval.get_rb_theta_expansion());
494  unsigned int n_M_terms = trans_theta_expansion.get_n_M_terms();
495  // Mq matrices
496  {
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)
499  {
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)
503  {
504  unsigned int offset = i*n_bfs+j;
505  set_scalar_in_list(rb_Mq_matrices_inner_list,
506  offset,
507  trans_rb_eval.RB_M_q_vector[q_m](i,j));
508  }
509  }
510  }
511 
512  // The initial condition and L2 error at t=0.
513  // We store the values for each RB space of dimension (0,...,n_basis_functions).
514  {
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);
519 
520  for (unsigned int i=0; i<n_bfs; i++)
521  {
522  initial_l2_errors_builder.set(i, trans_rb_eval.initial_L2_error_all_N[i]);
523 
524  auto initial_conditions_inner_list =
525  initial_conditions_outer_list.init(i, i+1);
526  for (unsigned int j=0; j<=i; j++)
527  {
528  set_scalar_in_list(initial_conditions_inner_list,
529  j,
530  trans_rb_eval.RB_initial_condition_all_N[i](j));
531  }
532  }
533  }
534 
535  // FqMq representor inner-product data
536  {
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);
540 
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)
544  {
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,
547  offset,
548  trans_rb_eval.Fq_Mq_representor_innerprods[q_f][q_m][i]);
549  }
550  }
551 
552  // MqMq representor inner-product data
553  {
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);
557 
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)
561  {
562  unsigned int offset = i*n_bfs*n_bfs + j*n_bfs + l;
563  set_scalar_in_list(mq_mq_innerprods_list,
564  offset,
565  trans_rb_eval.Mq_Mq_representor_innerprods[i][j][l]);
566  }
567  }
568 
569  // AqMq representor inner-product data
570  {
571  unsigned int n_A_terms = trans_theta_expansion.get_n_A_terms();
572 
573  auto aq_mq_innerprods_list =
574  trans_rb_eval_builder.initAqMqInnerprods(n_A_terms*n_M_terms*n_bfs*n_bfs);
575 
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++)
580  {
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,
584  offset,
585  trans_rb_eval.Aq_Mq_representor_innerprods[q_a][q_m][i][j]);
586  }
587  }
588 
589 }
590 
591 template <typename RBEIMEvaluationBuilderNumber>
593  RBEIMEvaluationBuilderNumber & rb_eim_evaluation_builder)
594 {
595  // Number of basis functions
596  unsigned int n_bfs = rb_eim_evaluation.get_n_basis_functions();
597  rb_eim_evaluation_builder.setNBfs(n_bfs);
598 
599  // EIM interpolation matrix
600  {
601  // We store the lower triangular part of an NxN matrix, the size of which is given by
602  // (N(N + 1))/2
603  unsigned int half_matrix_size = n_bfs*(n_bfs+1)/2;
604 
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)
609  {
610  unsigned int offset = i*(i+1)/2 + j;
611  set_scalar_in_list(interpolation_matrix_list,
612  offset,
613  rb_eim_evaluation.get_interpolation_matrix()(i,j));
614  }
615  }
616 
617  // We check use_eim_error_indicator() and the number of interpolation
618  // points in order to set use_error_indicator. This is because even if
619  // use_error_indicator() is true, we may not be using an error indicator,
620  // e.g. if there are no parameters in the model then we will not use an
621  // error indicator, and we can detect this by comparing the number of
622  // interpolation points with n_bfs.
623  bool use_error_indicator =
624  (rb_eim_evaluation.use_eim_error_indicator() &&
625  (rb_eim_evaluation.get_n_interpolation_points() > n_bfs));
626 
627  if (use_error_indicator)
628  {
629  auto error_indicator_data_list =
630  rb_eim_evaluation_builder.initEimErrorIndicatorInterpData(n_bfs);
631  const auto & error_indicator_row =
632  rb_eim_evaluation.get_error_indicator_interpolation_row();
633  for (unsigned int i=0; i < n_bfs; ++i)
634  {
635  set_scalar_in_list(error_indicator_data_list,
636  i,
637  error_indicator_row(i));
638  }
639  }
640 
641  // If we're using the EIM error indicator then we store one extra
642  // interpolation point and associated data, hence we increment n_bfs
643  // here so that we write out the extra data below.
644  //
645  // However the interpolation matrix and _error_indicator_interpolation_row
646  // does not include this extra point, which is why we write it out above
647  // before n_bfs is incremented.
648  if (use_error_indicator)
649  n_bfs++;
650 
651  libmesh_error_msg_if(n_bfs != rb_eim_evaluation.get_n_interpolation_points(),
652  "Number of basis functions should match number of interpolation points");
653 
654  auto parameter_ranges_list =
655  rb_eim_evaluation_builder.initParameterRanges();
656  auto discrete_parameters_list =
657  rb_eim_evaluation_builder.initDiscreteParameters();
658  add_parameter_ranges_to_builder(rb_eim_evaluation,
659  parameter_ranges_list,
660  discrete_parameters_list);
661 
662  // Interpolation points
663  {
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]);
669  }
670 
671  // Interpolation points comps
672  {
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,
677  rb_eim_evaluation.get_interpolation_points_comp(i));
678  }
679 
680  // Interpolation points subdomain IDs
681  {
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,
686  rb_eim_evaluation.get_interpolation_points_subdomain_id(i));
687  }
688 
689  // Interpolation points boundary IDs, relevant if the parametrized function is defined
690  // on mesh sides or nodesets
691  if (rb_eim_evaluation.get_parametrized_function().on_mesh_sides() ||
692  rb_eim_evaluation.get_parametrized_function().on_mesh_nodes())
693  {
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,
698  rb_eim_evaluation.get_interpolation_points_boundary_id(i));
699  }
700 
701  // Interpolation points element IDs
702  {
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,
707  rb_eim_evaluation.get_interpolation_points_elem_id(i));
708  }
709 
710  // Interpolation points node IDs, relevant if the parametrized function is defined on mesh sides
711  if (rb_eim_evaluation.get_parametrized_function().on_mesh_nodes())
712  {
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,
717  rb_eim_evaluation.get_interpolation_points_node_id(i));
718  }
719 
720  // Interpolation points side indices, relevant if the parametrized function is defined on mesh sides
721  if (rb_eim_evaluation.get_parametrized_function().on_mesh_sides())
722  {
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,
727  rb_eim_evaluation.get_interpolation_points_side_index(i));
728  }
729 
730  // Interpolation points quadrature point indices
731  {
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,
736  rb_eim_evaluation.get_interpolation_points_qp(i));
737  }
738 
739  // Interpolation points perturbations
740  {
741  auto interpolation_points_list_outer =
742  rb_eim_evaluation_builder.initInterpolationXyzPerturb(n_bfs);
743  for (unsigned int i=0; i < n_bfs; ++i)
744  {
745  const std::vector<Point> & perturbs = rb_eim_evaluation.get_interpolation_points_xyz_perturbations(i);
746  auto interpolation_points_list_inner = interpolation_points_list_outer.init(i, perturbs.size());
747 
748  for (unsigned int j : index_range(perturbs))
749  {
750  add_point_to_builder(perturbs[j], interpolation_points_list_inner[j]);
751  }
752  }
753  }
754 
755  unsigned int n_elems = rb_eim_evaluation.get_n_elems();
756  rb_eim_evaluation_builder.setNElems(n_elems);
757  // Elem id to local index map
758  {
759  auto elem_id_to_local_index_list =
760  rb_eim_evaluation_builder.initElemIdToLocalIndex(n_elems);
761  unsigned int counter = 0;
762  for (auto const & [elem_id, local_index] : rb_eim_evaluation.get_elem_id_to_local_index_map())
763  {
764  elem_id_to_local_index_list[counter].setFirst(elem_id);
765  elem_id_to_local_index_list[counter].setSecond(local_index);
766  counter++;
767  }
768  }
769 
770  // Interpolation points JxW values at each qp
771  {
772  auto interpolation_points_list_outer =
773  rb_eim_evaluation_builder.initInterpolationJxWAllQp(n_elems);
774  for (unsigned int i=0; i < n_elems; ++i)
775  {
776  const std::vector<Real> & JxW = rb_eim_evaluation.get_interpolation_points_JxW_all_qp(i);
777  auto interpolation_points_list_inner = interpolation_points_list_outer.init(i, JxW.size());
778 
779  for (unsigned int j : index_range(JxW))
780  {
781  // Here we can use set() instead of set_scalar_in_list() because
782  // phi stores real-valued data only.
783  interpolation_points_list_inner.set(j, JxW[j]);
784  }
785  }
786  }
787 
788  // Interpolation points phi values at each qp
789  {
790  auto interpolation_points_list_outer =
791  rb_eim_evaluation_builder.initInterpolationPhiValuesAllQp(n_elems);
792 
793  for (unsigned int i=0; i < n_elems; ++i)
794  {
795  const auto & phi_i_all_qp = rb_eim_evaluation.get_interpolation_points_phi_i_all_qp(i);
796  auto interpolation_points_list_middle = interpolation_points_list_outer.init(i, phi_i_all_qp.size());
797 
798  for (unsigned int j : index_range(phi_i_all_qp))
799  {
800  auto interpolation_points_list_inner = interpolation_points_list_middle.init(j, phi_i_all_qp[j].size());
801 
802  for (auto k : index_range(phi_i_all_qp[j]))
803  interpolation_points_list_inner.set(k, phi_i_all_qp[j][k]);
804  }
805  }
806  }
807 
808  // Dxyzdxi at the element center for the element that contains each interpolation point.
809  {
810  auto dxyzdxi_elem_center =
811  rb_eim_evaluation_builder.initInterpolationDxyzDxiElem(n_elems);
812  for (auto i : make_range(n_elems))
813  {
814  add_point_to_builder(rb_eim_evaluation.get_elem_center_dxyzdxi(i), dxyzdxi_elem_center[i]);
815  }
816  }
817 
818  // Dxyzdeta at the element center for the element that contains each interpolation point.
819  {
820  auto dxyzdeta_elem_center =
821  rb_eim_evaluation_builder.initInterpolationDxyzDetaElem(n_elems);
822  for (auto i : make_range(n_elems))
823  {
824  add_point_to_builder(rb_eim_evaluation.get_elem_center_dxyzdeta(i), dxyzdeta_elem_center[i]);
825  }
826  }
827 
828  // Quadrature rule order associated to the element that contains each interpolation point.
829  {
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,
834  rb_eim_evaluation.get_interpolation_points_qrule_order(i));
835  }
836 
837  // Element type for the element that contains each interpolation point
838  {
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,
843  rb_eim_evaluation.get_interpolation_points_elem_type(i));
844  }
845 
846  unsigned int n_properties = rb_eim_evaluation.get_n_properties();
847  // Property map used to store generic properties by flaging entites like elements, nodes etc...
848  {
849  auto interpolation_points_property_list =
850  rb_eim_evaluation_builder.initPropertyMap(n_properties);
851  unsigned int property_counter = 0;
852  for (const auto& [property_name, entity_ids] : rb_eim_evaluation.get_rb_property_map())
853  {
854  interpolation_points_property_list[property_counter].setName(property_name);
855 
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)
859  {
860  property_entity_ids.set(entity_counter, entity_id);
861  entity_counter++;
862  }
863  property_counter++;
864  }
865  }
866 
867  // Optionally store EIM solutions for the training set
868  if (rb_eim_evaluation.get_parametrized_function().is_lookup_table)
869  {
870  const std::vector<DenseVector<Number>> & eim_solutions = rb_eim_evaluation.get_eim_solutions_for_training_set();
871 
872  auto eim_rhs_list_outer =
873  rb_eim_evaluation_builder.initEimSolutionsForTrainingSet(eim_solutions.size());
874  for (auto i : make_range(eim_solutions.size()))
875  {
876  const DenseVector<Number> & values = eim_solutions[i];
877  auto eim_rhs_list_inner = eim_rhs_list_outer.init(i, values.size());
878 
879  for (auto j : index_range(values))
880  {
881  set_scalar_in_list(eim_rhs_list_inner,
882  j,
883  values(j));
884  }
885  }
886  }
887 
888  // The shape function values at the interpolation points. This can be used to evaluate nodal data
889  // at EIM interpolation points, which are at quadrature points.
890  {
891  auto interpolation_points_list_outer =
892  rb_eim_evaluation_builder.initInterpolationPhiValues(n_bfs);
893  for (unsigned int i=0; i < n_bfs; ++i)
894  {
895  const std::vector<Real> & phi_i_qp_vec = rb_eim_evaluation.get_interpolation_points_phi_i_qp(i);
896  auto interpolation_points_list_inner = interpolation_points_list_outer.init(i, phi_i_qp_vec.size());
897 
898  for (unsigned int j : index_range(phi_i_qp_vec))
899  {
900  // Here we can use set() instead of set_scalar_in_list() because
901  // phi stores real-valued data only.
902  interpolation_points_list_inner.set(j, phi_i_qp_vec[j]);
903  }
904  }
905  }
906 
907  // Initialize EIM spatial indices for the interpolation points, and if there are
908  // any spatial indices then we store them in the buffer.
910  if (rb_eim_evaluation.get_n_interpolation_points_spatial_indices() > 0)
911  {
912  libmesh_error_msg_if(n_bfs != rb_eim_evaluation.get_n_interpolation_points_spatial_indices(),
913  "Error: Number of spatial indices should match number of EIM basis functions");
914 
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)
918  {
919  const std::vector<unsigned int> & spatial_indices =
920  rb_eim_evaluation.get_interpolation_points_spatial_indices(i);
921  unsigned int n_spatial_indices = spatial_indices.size();
922 
923  auto spatial_indices_builder =
924  interpolation_points_spatial_indices_builder.init(i, n_spatial_indices);
925 
926  for (auto j : make_range(n_spatial_indices))
927  {
928  spatial_indices_builder.set(j, spatial_indices[j]);
929  }
930  }
931  }
932 }
933 
934 #if defined(LIBMESH_HAVE_SLEPC) && (LIBMESH_HAVE_GLPK)
936  RBData::RBSCMEvaluation::Builder & rb_scm_eval_builder)
937 {
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);
945 
946  {
947  libmesh_error_msg_if(rb_scm_eval.B_min.size() != rb_scm_eval.get_rb_theta_expansion().get_n_A_terms(),
948  "Size error while writing B_min");
949  auto b_min_list = rb_scm_eval_builder.initBMin( rb_scm_eval.B_min.size() );
950  for (auto i : index_range(rb_scm_eval.B_min))
951  b_min_list.set(i, rb_scm_eval.get_B_min(i));
952  }
953 
954  {
955  libmesh_error_msg_if(rb_scm_eval.B_max.size() != rb_scm_eval.get_rb_theta_expansion().get_n_A_terms(),
956  "Size error while writing B_max");
957 
958  auto b_max_list = rb_scm_eval_builder.initBMax( rb_scm_eval.B_max.size() );
959  for (auto i : index_range(rb_scm_eval.B_max))
960  b_max_list.set(i, rb_scm_eval.get_B_max(i));
961  }
962 
963  {
964  auto cj_stability_vector =
965  rb_scm_eval_builder.initCJStabilityVector( rb_scm_eval.C_J_stability_vector.size() );
966  for (auto i : index_range(rb_scm_eval.C_J_stability_vector))
967  cj_stability_vector.set(i, rb_scm_eval.get_C_J_stability_constraint(i));
968  }
969 
970  {
971  auto cj_parameters_outer =
972  rb_scm_eval_builder.initCJ( rb_scm_eval.C_J.size() );
973 
974  for (auto i : index_range(rb_scm_eval.C_J))
975  {
976  auto cj_parameters_inner =
977  cj_parameters_outer.init(i, rb_scm_eval.C_J[i].n_parameters());
978 
979  unsigned int count = 0;
980  for (const auto & [key, val] : rb_scm_eval.C_J[i])
981  {
982  cj_parameters_inner[count].setName(key);
983  cj_parameters_inner[count].setValue(val);
984  count++;
985  }
986 
987  }
988  }
989 
990  {
991  unsigned int n_C_J_values = rb_scm_eval.C_J.size();
992  unsigned int n_A_terms = rb_scm_eval.get_rb_theta_expansion().get_n_A_terms();
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 );
996 
997  for (unsigned int i=0; i<n_C_J_values; i++)
998  for (unsigned int j=0; j<n_A_terms; j++)
999  {
1000  unsigned int offset = i*n_A_terms + j;
1001  scm_ub_vectors.set(offset, rb_scm_eval.get_SCM_UB_vector(i,j));
1002  }
1003  }
1004 }
1005 #endif // LIBMESH_HAVE_SLEPC && LIBMESH_HAVE_GLPK
1006 
1007 void add_point_to_builder(const Point & point, RBData::Point3D::Builder point_builder)
1008 {
1009  point_builder.setX(point(0));
1010 
1011  if (LIBMESH_DIM >= 2)
1012  point_builder.setY(point(1));
1013 
1014  if (LIBMESH_DIM >= 3)
1015  point_builder.setZ(point(2));
1016 }
1017 
1018 // ---- Helper functions for adding data to capnp Builders (END) ----
1019 
1020 } // namespace RBDataSerialization
1021 
1022 } // namespace libMesh
1023 
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 &param_name) const
Get the value of the specified parameter, throw an error if it does not exist.
Definition: rb_parameters.C:65
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&#39;n&#39;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&#39;n&#39;Proto schema described in rb_data...
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 &parameter_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&#39;n&#39;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.
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&#39;n&#39;Proto buffer to disk.
TransientRBEvaluationSerialization(TransientRBEvaluation &rb_eval)
Initialize a new buffer using the structure from the Cap&#39;n&#39;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...
void write_to_file(const std::string &path, bool use_packing=false)
Write the Cap&#39;n&#39;Proto buffer to disk.
unsigned int get_time_step() const
Get/set the current time-step.
This class is part of the rbOOmit framework.
Definition: rb_parameters.h:52
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&#39;n&#39;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.
Definition: rb_evaluation.h:50
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&#39;n&#39;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.
static const bool value
Definition: xdr_io.C:54
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...
Definition: int_range.h:140
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
Definition: dense_vector.h:104
unsigned int get_n_time_steps() const
Get/set the total number of time-steps.
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
boost::multiprecision::float128 imag(const boost::multiprecision::float128)
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.
Definition: point.h:39
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:117
RBThetaExpansion & get_rb_theta_expansion()
Get a reference to the rb_theta_expansion.
Definition: rb_evaluation.C:85
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.