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/string_to_enum.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 
34 // Cap'n'Proto includes
35 #include <capnp/serialize.h>
36 
37 // C++ includes
38 #include <iostream>
39 #include <fstream>
40 #include <unistd.h>
41 #include <fcntl.h>
42 
43 namespace libMesh
44 {
45 
46 namespace
47 {
48 
53 template <typename T, typename U>
54 void set_scalar_in_list(T list, unsigned int i, U value)
55 {
56 #ifdef LIBMESH_USE_COMPLEX_NUMBERS
57  list[i].setReal(std::real(value));
58  list[i].setImag(std::imag(value));
59 #else
60  list.set(i, value);
61 #endif
62 }
63 
64 }
65 
66 namespace RBDataSerialization
67 {
68 
69 // ---- RBEvaluationSerialization (BEGIN) ----
70 
72  :
73  _rb_eval(rb_eval)
74 {
75 }
76 
78 {
79 }
80 
81 void RBEvaluationSerialization::write_to_file(const std::string & path)
82 {
83  LOG_SCOPE("write_to_file()", "RBEvaluationSerialization");
84 
85  if (_rb_eval.comm().rank() == 0)
86  {
87  capnp::MallocMessageBuilder message;
88 
89 #ifndef LIBMESH_USE_COMPLEX_NUMBERS
90  RBData::RBEvaluationReal::Builder rb_eval_builder =
91  message.initRoot<RBData::RBEvaluationReal>();
92 #else
93  RBData::RBEvaluationComplex::Builder rb_eval_builder =
94  message.initRoot<RBData::RBEvaluationComplex>();
95 #endif
96 
98 
99  int fd = open(path.c_str(), O_WRONLY | O_CREAT | O_TRUNC, 0664);
100  if (!fd)
101  libmesh_error_msg("Error opening a write-only file descriptor to " + path);
102 
103  capnp::writeMessageToFd(fd, message);
104 
105  int error = close(fd);
106  if (error)
107  libmesh_error_msg("Error closing a write-only file descriptor to " + path);
108  }
109 }
110 
111 // ---- RBEvaluationSerialization (END) ----
112 
113 
114 // ---- TransientRBEvaluationSerialization (BEGIN) ----
115 
118  _trans_rb_eval(trans_rb_eval)
119 {
120 }
121 
123 {
124 }
125 
127 {
128  LOG_SCOPE("write_to_file()", "TransientRBEvaluationSerialization");
129 
130  if (_trans_rb_eval.comm().rank() == 0)
131  {
132  capnp::MallocMessageBuilder message;
133 
134 #ifndef LIBMESH_USE_COMPLEX_NUMBERS
135  RBData::TransientRBEvaluationReal::Builder trans_rb_eval_builder =
136  message.initRoot<RBData::TransientRBEvaluationReal>();
137  RBData::RBEvaluationReal::Builder rb_eval_builder =
138  trans_rb_eval_builder.initRbEvaluation();
139 #else
140  RBData::TransientRBEvaluationComplex::Builder trans_rb_eval_builder =
141  message.initRoot<RBData::TransientRBEvaluationComplex>();
142  RBData::RBEvaluationComplex::Builder rb_eval_builder =
143  trans_rb_eval_builder.initRbEvaluation();
144 #endif
145 
147  rb_eval_builder,
148  trans_rb_eval_builder);
149 
150  int fd = open(path.c_str(), O_WRONLY | O_CREAT | O_TRUNC, 0664);
151  if (!fd)
152  libmesh_error_msg("Error opening a write-only file descriptor to " + path);
153 
154  capnp::writeMessageToFd(fd, message);
155 
156  int error = close(fd);
157  if (error)
158  libmesh_error_msg("Error closing a write-only file descriptor to " + path);
159  }
160 }
161 
162 // ---- TransientRBEvaluationSerialization (END) ----
163 
164 
165 // ---- RBEIMEvaluationSerialization (BEGIN) ----
166 
168  :
169  _rb_eim_eval(rb_eim_eval)
170 {
171 }
172 
174 {
175 }
176 
177 void RBEIMEvaluationSerialization::write_to_file(const std::string & path)
178 {
179  LOG_SCOPE("write_to_file()", "RBEIMEvaluationSerialization");
180 
181  if (_rb_eim_eval.comm().rank() == 0)
182  {
183  capnp::MallocMessageBuilder message;
184 
185 #ifndef LIBMESH_USE_COMPLEX_NUMBERS
186  RBData::RBEIMEvaluationReal::Builder rb_eim_eval_builder =
187  message.initRoot<RBData::RBEIMEvaluationReal>();
188  RBData::RBEvaluationReal::Builder rb_eval_builder =
189  rb_eim_eval_builder.initRbEvaluation();
190 #else
191  RBData::RBEIMEvaluationComplex::Builder rb_eim_eval_builder =
192  message.initRoot<RBData::RBEIMEvaluationComplex>();
193  RBData::RBEvaluationComplex::Builder rb_eval_builder =
194  rb_eim_eval_builder.initRbEvaluation();
195 #endif
196 
198  rb_eval_builder,
199  rb_eim_eval_builder);
200 
201  int fd = open(path.c_str(), O_WRONLY | O_CREAT | O_TRUNC, 0664);
202  if (!fd)
203  libmesh_error_msg("Error opening a write-only file descriptor to " + path);
204 
205  capnp::writeMessageToFd(fd, message);
206 
207  int error = close(fd);
208  if (error)
209  libmesh_error_msg("Error closing a write-only file descriptor to " + path);
210  }
211 }
212 
213 // ---- RBEIMEvaluationSerialization (END) ----
214 
215 
216 // ---- RBSCMEvaluationSerialization (BEGIN) ----
217 
218 #if defined(LIBMESH_HAVE_SLEPC) && (LIBMESH_HAVE_GLPK)
219 
221  :
222  _rb_scm_eval(rb_scm_eval)
223 {
224 }
225 
227 {
228 }
229 
230 void RBSCMEvaluationSerialization::write_to_file(const std::string & path)
231 {
232  LOG_SCOPE("write_to_file()", "RBSCMEvaluationSerialization");
233 
234  if (_rb_scm_eval.comm().rank() == 0)
235  {
236  capnp::MallocMessageBuilder message;
237 
238  RBData::RBSCMEvaluation::Builder rb_scm_eval_builder =
239  message.initRoot<RBData::RBSCMEvaluation>();
240 
242 
243  int fd = open(path.c_str(), O_WRONLY | O_CREAT | O_TRUNC, 0664);
244  if (!fd)
245  libmesh_error_msg("Error opening a write-only file descriptor to " + path);
246 
247  capnp::writeMessageToFd(fd, message);
248 
249  int error = close(fd);
250  if (error)
251  libmesh_error_msg("Error closing a write-only file descriptor to " + path);
252  }
253 }
254 
255 #endif // LIBMESH_HAVE_SLEPC && LIBMESH_HAVE_GLPK
256 
257 // ---- RBSCMEvaluationSerialization (END) ----
258 
259 
260 // ---- Helper functions for adding data to capnp Builders (BEGIN) ----
261 
263  RBData::ParameterRanges::Builder & parameter_ranges_list,
264  RBData::DiscreteParameterList::Builder & discrete_parameters_list)
265 {
266  // Continuous parameters
267  {
268  unsigned int n_continuous_parameters = rb_evaluation.get_n_continuous_params();
269  auto names = parameter_ranges_list.initNames(n_continuous_parameters);
270  auto mins = parameter_ranges_list.initMinValues(n_continuous_parameters);
271  auto maxs = parameter_ranges_list.initMaxValues(n_continuous_parameters);
272 
273  std::set<std::string> parameter_names = rb_evaluation.get_parameter_names();
274  const RBParameters & parameters_min = rb_evaluation.get_parameters_min();
275  const RBParameters & parameters_max = rb_evaluation.get_parameters_max();
276 
277  unsigned int count = 0;
278  for (const auto & parameter_name : parameter_names)
279  {
280  if (!rb_evaluation.is_discrete_parameter(parameter_name))
281  {
282  names.set(count, parameter_name);
283  mins.set(count, parameters_min.get_value(parameter_name));
284  maxs.set(count, parameters_max.get_value(parameter_name));
285 
286  ++count;
287  }
288  }
289 
290  if (count != n_continuous_parameters)
291  libmesh_error_msg("Mismatch in number of continuous parameters");
292  }
293 
294  // Discrete parameters
295  {
296  unsigned int n_discrete_parameters = rb_evaluation.get_n_discrete_params();
297  auto names = discrete_parameters_list.initNames(n_discrete_parameters);
298  auto values_outer = discrete_parameters_list.initValues(n_discrete_parameters);
299 
300  const std::map<std::string, std::vector<Real>> & discrete_parameters =
301  rb_evaluation.get_discrete_parameter_values();
302 
303  unsigned int count = 0;
304  for (const auto & discrete_parameter : discrete_parameters)
305  {
306  names.set(count, discrete_parameter.first);
307 
308  const std::vector<Real> & values = discrete_parameter.second;
309  unsigned int n_values = values.size();
310 
311  values_outer.init(count, n_values);
312  auto values_inner = values_outer[count];
313  for (unsigned int i=0; i<n_values; ++i)
314  {
315  values_inner.set(i, values[i]);
316  }
317 
318  ++count;
319  }
320 
321  if (count != n_discrete_parameters)
322  libmesh_error_msg("Mismatch in number of discrete parameters");
323  }
324 }
325 
326 template <typename RBEvaluationBuilderNumber>
328  RBEvaluationBuilderNumber & rb_evaluation_builder)
329 {
330  const RBThetaExpansion & rb_theta_expansion = rb_eval.get_rb_theta_expansion();
331 
332  unsigned int n_F_terms = rb_theta_expansion.get_n_F_terms();
333  unsigned int n_A_terms = rb_theta_expansion.get_n_A_terms();
334 
335  // Number of basis functions
336  unsigned int n_bfs = rb_eval.get_n_basis_functions();
337  rb_evaluation_builder.setNBfs(n_bfs);
338 
339  // Fq representor inner-product data
340  {
341  unsigned int Q_f_hat = n_F_terms*(n_F_terms+1)/2;
342 
343  auto fq_innerprods_list = rb_evaluation_builder.initFqInnerprods(Q_f_hat);
344 
345  for (unsigned int i=0; i<Q_f_hat; i++)
346  set_scalar_in_list(fq_innerprods_list,
347  i,
348  rb_eval.Fq_representor_innerprods[i]);
349  }
350 
351  // FqAq representor inner-product data
352  {
353  auto fq_aq_innerprods_list =
354  rb_evaluation_builder.initFqAqInnerprods(n_F_terms*n_A_terms*n_bfs);
355 
356  for (unsigned int q_f=0; q_f < n_F_terms; ++q_f)
357  for (unsigned int q_a=0; q_a < n_A_terms; ++q_a)
358  for (unsigned int i=0; i < n_bfs; ++i)
359  {
360  unsigned int offset = q_f*n_A_terms*n_bfs + q_a*n_bfs + i;
361  set_scalar_in_list(
362  fq_aq_innerprods_list, offset,
363  rb_eval.Fq_Aq_representor_innerprods[q_f][q_a][i]);
364  }
365  }
366 
367  // AqAq representor inner-product data
368  {
369  unsigned int Q_a_hat = n_A_terms*(n_A_terms+1)/2;
370  auto aq_aq_innerprods_list =
371  rb_evaluation_builder.initAqAqInnerprods(Q_a_hat*n_bfs*n_bfs);
372 
373  for (unsigned int i=0; i < Q_a_hat; ++i)
374  for (unsigned int j=0; j < n_bfs; ++j)
375  for (unsigned int l=0; l < n_bfs; ++l)
376  {
377  unsigned int offset = i*n_bfs*n_bfs + j*n_bfs + l;
378  set_scalar_in_list(
379  aq_aq_innerprods_list,
380  offset,
381  rb_eval.Aq_Aq_representor_innerprods[i][j][l]);
382  }
383  }
384 
385  // Output dual inner-product data, and output vectors
386  {
387  unsigned int n_outputs = rb_theta_expansion.get_n_outputs();
388  auto output_innerprod_outer = rb_evaluation_builder.initOutputDualInnerprods(n_outputs);
389  auto output_vector_outer = rb_evaluation_builder.initOutputVectors(n_outputs);
390 
391  for (unsigned int output_id=0; output_id < n_outputs; ++output_id)
392  {
393  unsigned int n_output_terms = rb_theta_expansion.get_n_output_terms(output_id);
394 
395  {
396  unsigned int Q_l_hat = n_output_terms*(n_output_terms+1)/2;
397  auto output_innerprod_inner = output_innerprod_outer.init(output_id, Q_l_hat);
398  for (unsigned int q=0; q < Q_l_hat; ++q)
399  {
400  set_scalar_in_list(
401  output_innerprod_inner, q, rb_eval.output_dual_innerprods[output_id][q]);
402  }
403  }
404 
405  {
406  auto output_vector_middle = output_vector_outer.init(output_id, n_output_terms);
407  for (unsigned int q_l=0; q_l<n_output_terms; ++q_l)
408  {
409  auto output_vector_inner = output_vector_middle.init(q_l, n_bfs);
410  for (unsigned int j=0; j<n_bfs; ++j)
411  {
412  set_scalar_in_list(
413  output_vector_inner, j, rb_eval.RB_output_vectors[output_id][q_l](j));
414  }
415  }
416  }
417  }
418  }
419 
420  // Fq vectors and Aq matrices
421  {
422  unsigned int n_F_terms = rb_theta_expansion.get_n_F_terms();
423  unsigned int n_A_terms = rb_theta_expansion.get_n_A_terms();
424 
425  auto rb_fq_vectors_outer_list = rb_evaluation_builder.initRbFqVectors(n_F_terms);
426  for (unsigned int q_f=0; q_f < n_F_terms; ++q_f)
427  {
428  auto rb_fq_vectors_inner_list = rb_fq_vectors_outer_list.init(q_f, n_bfs);
429  for (unsigned int i=0; i<n_bfs; i++)
430  set_scalar_in_list(rb_fq_vectors_inner_list, i, rb_eval.RB_Fq_vector[q_f](i));
431  }
432 
433  auto rb_Aq_matrices_outer_list = rb_evaluation_builder.initRbAqMatrices(n_A_terms);
434  for (unsigned int q_a=0; q_a < n_A_terms; ++q_a)
435  {
436  auto rb_Aq_matrices_inner_list = rb_Aq_matrices_outer_list.init(q_a, n_bfs*n_bfs);
437  for (unsigned int i=0; i < n_bfs; ++i)
438  for (unsigned int j=0; j < n_bfs; ++j)
439  {
440  unsigned int offset = i*n_bfs+j;
441  set_scalar_in_list(rb_Aq_matrices_inner_list, offset, rb_eval.RB_Aq_vector[q_a](i,j));
442  }
443  }
444  }
445 
446  // Inner-product matrix
447  if (rb_eval.compute_RB_inner_product)
448  {
449  auto rb_inner_product_matrix_list =
450  rb_evaluation_builder.initRbInnerProductMatrix(n_bfs*n_bfs);
451 
452  for (unsigned int i=0; i < n_bfs; ++i)
453  for (unsigned int j=0; j < n_bfs; ++j)
454  {
455  unsigned int offset = i*n_bfs + j;
456  set_scalar_in_list(
457  rb_inner_product_matrix_list,
458  offset,
459  rb_eval.RB_inner_product_matrix(i,j) );
460  }
461  }
462 
463  auto parameter_ranges_list =
464  rb_evaluation_builder.initParameterRanges();
465  auto discrete_parameters_list =
466  rb_evaluation_builder.initDiscreteParameters();
468  parameter_ranges_list,
469  discrete_parameters_list);
470 }
471 
472 template <typename RBEvaluationBuilderNumber, typename TransRBEvaluationBuilderNumber>
474  RBEvaluationBuilderNumber & rb_eval_builder,
475  TransRBEvaluationBuilderNumber & trans_rb_eval_builder)
476 {
477  add_rb_evaluation_data_to_builder(trans_rb_eval, rb_eval_builder);
478 
479  trans_rb_eval_builder.setDeltaT(trans_rb_eval.get_delta_t());
480  trans_rb_eval_builder.setEulerTheta(trans_rb_eval.get_euler_theta());
481  trans_rb_eval_builder.setNTimeSteps(trans_rb_eval.get_n_time_steps());
482  trans_rb_eval_builder.setTimeStep(trans_rb_eval.get_time_step());
483 
484  unsigned int n_bfs = trans_rb_eval.get_n_basis_functions();
485 
486  // L2-inner-product matrix
487  {
488  auto rb_L2_matrix_list =
489  trans_rb_eval_builder.initRbL2Matrix(n_bfs*n_bfs);
490 
491  for (unsigned int i=0; i<n_bfs; ++i)
492  for (unsigned int j=0; j<n_bfs; ++j)
493  {
494  unsigned int offset = i*n_bfs + j;
495  set_scalar_in_list(rb_L2_matrix_list,
496  offset,
497  trans_rb_eval.RB_L2_matrix(i,j));
498  }
499  }
500 
501  TransientRBThetaExpansion & trans_theta_expansion =
502  cast_ref<TransientRBThetaExpansion &>(trans_rb_eval.get_rb_theta_expansion());
503  unsigned int n_M_terms = trans_theta_expansion.get_n_M_terms();
504  // Mq matrices
505  {
506  auto rb_Mq_matrices_outer_list = trans_rb_eval_builder.initRbMqMatrices(n_M_terms);
507  for (unsigned int q_m=0; q_m < n_M_terms; ++q_m)
508  {
509  auto rb_Mq_matrices_inner_list = rb_Mq_matrices_outer_list.init(q_m, n_bfs*n_bfs);
510  for (unsigned int i=0; i < n_bfs; ++i)
511  for (unsigned int j=0; j < n_bfs; ++j)
512  {
513  unsigned int offset = i*n_bfs+j;
514  set_scalar_in_list(rb_Mq_matrices_inner_list,
515  offset,
516  trans_rb_eval.RB_M_q_vector[q_m](i,j));
517  }
518  }
519  }
520 
521  // The initial condition and L2 error at t=0.
522  // We store the values for each RB space of dimension (0,...,n_basis_functions).
523  {
524  auto initial_l2_errors_builder =
525  trans_rb_eval_builder.initInitialL2Errors(n_bfs);
526  auto initial_conditions_outer_list =
527  trans_rb_eval_builder.initInitialConditions(n_bfs);
528 
529  for (unsigned int i=0; i<n_bfs; i++)
530  {
531  initial_l2_errors_builder.set(i, trans_rb_eval.initial_L2_error_all_N[i]);
532 
533  auto initial_conditions_inner_list =
534  initial_conditions_outer_list.init(i, i+1);
535  for (unsigned int j=0; j<=i; j++)
536  {
537  set_scalar_in_list(initial_conditions_inner_list,
538  j,
539  trans_rb_eval.RB_initial_condition_all_N[i](j));
540  }
541  }
542  }
543 
544  // FqMq representor inner-product data
545  {
546  unsigned int n_F_terms = trans_theta_expansion.get_n_F_terms();
547  auto fq_mq_innerprods_list =
548  trans_rb_eval_builder.initFqMqInnerprods(n_F_terms*n_M_terms*n_bfs);
549 
550  for (unsigned int q_f=0; q_f<n_F_terms; ++q_f)
551  for (unsigned int q_m=0; q_m<n_M_terms; ++q_m)
552  for (unsigned int i=0; i<n_bfs; ++i)
553  {
554  unsigned int offset = q_f*n_M_terms*n_bfs + q_m*n_bfs + i;
555  set_scalar_in_list(fq_mq_innerprods_list,
556  offset,
557  trans_rb_eval.Fq_Mq_representor_innerprods[q_f][q_m][i]);
558  }
559  }
560 
561  // MqMq representor inner-product data
562  {
563  unsigned int Q_m_hat = n_M_terms*(n_M_terms+1)/2;
564  auto mq_mq_innerprods_list =
565  trans_rb_eval_builder.initMqMqInnerprods(Q_m_hat*n_bfs*n_bfs);
566 
567  for (unsigned int i=0; i < Q_m_hat; ++i)
568  for (unsigned int j=0; j < n_bfs; ++j)
569  for (unsigned int l=0; l < n_bfs; ++l)
570  {
571  unsigned int offset = i*n_bfs*n_bfs + j*n_bfs + l;
572  set_scalar_in_list(mq_mq_innerprods_list,
573  offset,
574  trans_rb_eval.Mq_Mq_representor_innerprods[i][j][l]);
575  }
576  }
577 
578  // AqMq representor inner-product data
579  {
580  unsigned int n_A_terms = trans_theta_expansion.get_n_A_terms();
581 
582  auto aq_mq_innerprods_list =
583  trans_rb_eval_builder.initAqMqInnerprods(n_A_terms*n_M_terms*n_bfs*n_bfs);
584 
585  for (unsigned int q_a=0; q_a<n_A_terms; q_a++)
586  for (unsigned int q_m=0; q_m<n_M_terms; q_m++)
587  for (unsigned int i=0; i<n_bfs; i++)
588  for (unsigned int j=0; j<n_bfs; j++)
589  {
590  unsigned int offset =
591  q_a*(n_M_terms*n_bfs*n_bfs) + q_m*(n_bfs*n_bfs) + i*n_bfs + j;
592  set_scalar_in_list(aq_mq_innerprods_list,
593  offset,
594  trans_rb_eval.Aq_Mq_representor_innerprods[q_a][q_m][i][j]);
595  }
596  }
597 
598 }
599 
600 template <typename RBEvaluationBuilderNumber, typename RBEIMEvaluationBuilderNumber>
602  RBEvaluationBuilderNumber & rb_evaluation_builder,
603  RBEIMEvaluationBuilderNumber & rb_eim_evaluation_builder)
604 {
605  add_rb_evaluation_data_to_builder(rb_eim_evaluation, rb_evaluation_builder);
606 
607  unsigned int n_bfs = rb_eim_evaluation.get_n_basis_functions();
608 
609  // EIM interpolation matrix
610  {
611  // We store the lower triangular part of an NxN matrix, the size of which is given by
612  // (N(N + 1))/2
613  unsigned int half_matrix_size = n_bfs*(n_bfs+1)/2;
614 
615  auto interpolation_matrix_list =
616  rb_eim_evaluation_builder.initInterpolationMatrix(half_matrix_size);
617  for (unsigned int i=0; i < n_bfs; ++i)
618  for (unsigned int j=0; j <= i; ++j)
619  {
620  unsigned int offset = i*(i+1)/2 + j;
621  set_scalar_in_list(interpolation_matrix_list,
622  offset,
623  rb_eim_evaluation.interpolation_matrix(i,j));
624  }
625  }
626 
627  // Interpolation points
628  {
629  auto interpolation_points_list =
630  rb_eim_evaluation_builder.initInterpolationPoints(n_bfs);
631  for (unsigned int i=0; i < n_bfs; ++i)
632  add_point_to_builder(rb_eim_evaluation.interpolation_points[i],
633  interpolation_points_list[i]);
634  }
635 
636  // Interpolation points variables
637  {
638  auto interpolation_points_var_list =
639  rb_eim_evaluation_builder.initInterpolationPointsVar(n_bfs);
640  for (unsigned int i=0; i<n_bfs; ++i)
641  interpolation_points_var_list.set(i,
642  rb_eim_evaluation.interpolation_points_var[i]);
643  }
644 
645  // Interpolation elements
646  {
647  unsigned int n_interpolation_elems =
648  rb_eim_evaluation.interpolation_points_elem.size();
649  auto interpolation_points_elem_list =
650  rb_eim_evaluation_builder.initInterpolationPointsElems(n_interpolation_elems);
651 
652  if (n_interpolation_elems != n_bfs)
653  libmesh_error_msg("The number of elements should match the number of basis functions");
654 
655  for (unsigned int i=0; i<n_interpolation_elems; ++i)
656  {
657  const libMesh::Elem & elem = *rb_eim_evaluation.interpolation_points_elem[i];
658  auto mesh_elem_builder = interpolation_points_elem_list[i];
659  add_elem_to_builder(elem, mesh_elem_builder);
660  }
661  }
662 }
663 
664 #if defined(LIBMESH_HAVE_SLEPC) && (LIBMESH_HAVE_GLPK)
666  RBData::RBSCMEvaluation::Builder & rb_scm_eval_builder)
667 {
668  auto parameter_ranges_list =
669  rb_scm_eval_builder.initParameterRanges();
670  auto discrete_parameters_list =
671  rb_scm_eval_builder.initDiscreteParameters();
673  parameter_ranges_list,
674  discrete_parameters_list);
675 
676  {
677  if (rb_scm_eval.B_min.size() != rb_scm_eval.get_rb_theta_expansion().get_n_A_terms())
678  libmesh_error_msg("Size error while writing B_min");
679  auto b_min_list = rb_scm_eval_builder.initBMin( rb_scm_eval.B_min.size() );
680  for (std::size_t i=0; i<rb_scm_eval.B_min.size(); i++)
681  b_min_list.set(i, rb_scm_eval.get_B_min(i));
682  }
683 
684  {
685  if (rb_scm_eval.B_max.size() != rb_scm_eval.get_rb_theta_expansion().get_n_A_terms())
686  libmesh_error_msg("Size error while writing B_max");
687 
688  auto b_max_list = rb_scm_eval_builder.initBMax( rb_scm_eval.B_max.size() );
689  for (std::size_t i=0; i<rb_scm_eval.B_max.size(); i++)
690  b_max_list.set(i, rb_scm_eval.get_B_max(i));
691  }
692 
693  {
694  auto cj_stability_vector =
695  rb_scm_eval_builder.initCJStabilityVector( rb_scm_eval.C_J_stability_vector.size() );
696  for (std::size_t i=0; i<rb_scm_eval.C_J_stability_vector.size(); i++)
697  cj_stability_vector.set(i, rb_scm_eval.get_C_J_stability_constraint(i));
698  }
699 
700  {
701  auto cj_parameters_outer =
702  rb_scm_eval_builder.initCJ( rb_scm_eval.C_J.size() );
703 
704  for (std::size_t i=0; i<rb_scm_eval.C_J.size(); i++)
705  {
706  auto cj_parameters_inner =
707  cj_parameters_outer.init(i, rb_scm_eval.C_J[i].n_parameters());
708 
709  RBParameters::const_iterator it = rb_scm_eval.C_J[i].begin();
710  RBParameters::const_iterator it_end = rb_scm_eval.C_J[i].end();
711 
712  unsigned int count = 0;
713  for ( ; it != it_end; ++it)
714  {
715  cj_parameters_inner[count].setName( it->first );
716  cj_parameters_inner[count].setValue( it->second );
717  count++;
718  }
719 
720  }
721  }
722 
723  {
724  unsigned int n_C_J_values = rb_scm_eval.C_J.size();
725  unsigned int n_A_terms = rb_scm_eval.get_rb_theta_expansion().get_n_A_terms();
726  unsigned int n_values = n_C_J_values*n_A_terms;
727  auto scm_ub_vectors =
728  rb_scm_eval_builder.initScmUbVectors( n_values );
729 
730  for (unsigned int i=0; i<n_C_J_values; i++)
731  for (unsigned int j=0; j<n_A_terms; j++)
732  {
733  unsigned int offset = i*n_A_terms + j;
734  scm_ub_vectors.set(offset, rb_scm_eval.get_SCM_UB_vector(i,j));
735  }
736  }
737 }
738 #endif // LIBMESH_HAVE_SLEPC && LIBMESH_HAVE_GLPK
739 
740 void add_point_to_builder(const Point & point, RBData::Point3D::Builder point_builder)
741 {
742  point_builder.setX(point(0));
743 
744  if (LIBMESH_DIM >= 2)
745  point_builder.setY(point(1));
746 
747  if (LIBMESH_DIM >= 3)
748  point_builder.setZ(point(2));
749 }
750 
751 void add_elem_to_builder(const libMesh::Elem & elem, RBData::MeshElem::Builder mesh_elem_builder)
752 {
753  std::string elem_type_string = libMesh::Utility::enum_to_string(elem.type());
754 
755  mesh_elem_builder.setType(elem_type_string.c_str());
756  mesh_elem_builder.setSubdomainId(elem.subdomain_id());
757 
758  unsigned int n_points = elem.n_nodes();
759  auto mesh_elem_point_list = mesh_elem_builder.initPoints(n_points);
760 
761  for (unsigned int j=0; j < n_points; ++j)
762  {
763  add_point_to_builder(elem.node_ref(j), mesh_elem_point_list[j]);
764  }
765 }
766 
767 // ---- Helper functions for adding data to capnp Builders (END) ----
768 
769 } // namespace RBDataSerialization
770 
771 } // namespace libMesh
772 
773 #endif // #if defined(LIBMESH_HAVE_CAPNPROTO)
std::vector< Real > B_min
B_min, B_max define the bounding box.
std::vector< std::vector< std::vector< Number > > > Mq_Mq_representor_innerprods
DenseMatrix< Number > interpolation_matrix
Dense matrix that stores the lower triangular interpolation matrix that can be used.
const RBParameters & get_parameters_max() const
Get an RBParameters object that specifies the maximum allowable value for each parameter.
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...
Real get_B_min(unsigned int i) const
Get B_min and B_max.
std::vector< Number > Fq_representor_innerprods
Vectors storing the residual representor inner products to be used in computing the residuals online...
virtual ElemType type() const =0
unsigned int get_n_A_terms() const
Get Q_a, the number of terms in the affine expansion for the bilinear form.
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.
Real get_value(const std::string &param_name) const
Get the value of the specific parameter.
Definition: rb_parameters.C:45
DenseMatrix< Number > RB_inner_product_matrix
The inner product matrix.
unsigned int get_n_output_terms(unsigned int output_index) const
Get the number of affine terms associated with the specified output.
This is the base class from which all geometric element types are derived.
Definition: elem.h:89
bool is_discrete_parameter(const std::string &mu_name) const
Is parameter mu_name discrete?
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.
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.
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_B_max(unsigned int i) const
RBThetaExpansion & get_rb_theta_expansion()
Get a reference to the rb_theta_expansion.
std::vector< Real > B_max
This class stores the set of RBTheta functor objects that define the "parameter-dependent expansion" ...
void add_rb_eim_evaluation_data_to_builder(RBEIMEvaluation &rb_eim_eval, RBEvaluationBuilderNumber &rb_eval_builder, RBEIMEvaluationBuilderNumber &rb_eim_eval_builder)
Add data for an RBEIMEvaluation to the builder.
std::vector< Real > C_J_stability_vector
Vector storing the (truth) stability values at the parameters in C_J.
virtual unsigned int n_nodes() const =0
TransientRBEvaluationSerialization(TransientRBEvaluation &rb_eval)
Initialize a new buffer using the structure from the Cap&#39;n&#39;Proto schema described in rb_data...
Real get_euler_theta() const
Get/set euler_theta, parameter that determines the temporal discretization.
std::vector< DenseVector< Number > > RB_Fq_vector
Dense vector for the RHS.
Real get_C_J_stability_constraint(unsigned int j) const
Get stability constraints (i.e.
std::vector< Real > initial_L2_error_all_N
Vector storing initial L2 error for all 1 <= N <= RB_size.
std::vector< RBParameters > C_J
Vector storing the greedily selected parameters during SCM training.
unsigned int get_n_F_terms() const
Get Q_f, the number of terms in the affine expansion for the right-hand side.
const Node & node_ref(const unsigned int i) const
Definition: elem.h:1896
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)
Write the Cap&#39;n&#39;Proto buffer to disk.
subdomain_id_type subdomain_id() const
Definition: elem.h:1951
void write_to_file(const std::string &path)
Write the Cap&#39;n&#39;Proto buffer to disk.
std::vector< unsigned int > interpolation_points_var
The corresponding list of variables indices at which the interpolation points were identified...
This class is part of the rbOOmit framework.
Definition: rb_parameters.h:42
std::vector< Elem * > interpolation_points_elem
The corresponding list of elements at which the interpolation points were identified.
virtual unsigned int get_n_M_terms()
Get Q_m, the number of terms in the affine expansion for the mass operator.
void write_to_file(const std::string &path)
Write the Cap&#39;n&#39;Proto buffer to disk.
void write_to_file(const std::string &path)
Write the Cap&#39;n&#39;Proto buffer to disk.
RBEIMEvaluationSerialization(RBEIMEvaluation &rb_eval)
Initialize a new buffer using the structure from the Cap&#39;n&#39;Proto schema described in rb_data...
const RBParameters & get_parameters_min() const
Get an RBParameters object that specifies the minimum allowable value for each parameter.
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::string enum_to_string(const T e)
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.
Real get_delta_t() const
Get/set delta_t, the time-step size.
unsigned int get_n_outputs() const
Get n_outputs, the number output functionals.
This class is part of the rbOOmit framework.
Definition: rb_evaluation.h:50
unsigned int get_n_discrete_params() const
Get the number of discrete parameters.
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.
unsigned int get_n_continuous_params() const
Get the number of continuous parameters.
void add_elem_to_builder(const Elem &elem, RBData::MeshElem::Builder mesh_elem_builder)
Helper function that adds element data.
unsigned int get_time_step() const
Get/set the current time-step.
This class is part of the rbOOmit framework.
const Parallel::Communicator & comm() const
static const bool value
Definition: xdr_io.C:108
TransientRBEvaluation & _trans_rb_eval
The RBEvaluation object that will be written to disk.
unsigned int rank() const
Definition: parallel.h:724
const std::map< std::string, std::vector< Real > > & get_discrete_parameter_values() const
Get a const reference to the discrete parameter values.
unsigned int get_n_time_steps() const
Get/set the total number of time-steps.
std::vector< std::vector< std::vector< std::vector< Number > > > > Aq_Mq_representor_innerprods
RBEvaluation & _rb_eval
The RBEvaluation object that will be written to disk.
virtual unsigned int get_n_basis_functions() const
Get the current number of basis functions.
std::map< std::string, Real >::const_iterator const_iterator
Definition: rb_parameters.h:57
std::vector< Point > interpolation_points
The list of interpolation points, i.e.
This class is part of the rbOOmit framework.
A Point defines a location in LIBMESH_DIM dimensional Real space.
Definition: point.h:38
RBThetaExpansion & get_rb_theta_expansion()
Get a reference to the rb_theta_expansion.
Definition: rb_evaluation.C:89
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...
This class is part of the rbOOmit framework.
std::vector< std::vector< std::vector< Number > > > Aq_Aq_representor_innerprods
std::set< std::string > get_parameter_names() const
Get a set that stores the parameter names.
std::vector< std::vector< DenseVector< Number > > > RB_output_vectors
The vectors storing the RB output vectors.