HyPar  1.0
Finite-Difference Hyperbolic-Parabolic PDE Solver on Cartesian Grids
 All Data Structures Namespaces Files Functions Variables Typedefs Enumerations Enumerator Macros Pages
Solve.cpp
Go to the documentation of this file.
1 
6 #include <stdio.h>
7 #include <math.h>
8 #include <string.h>
9 #include <vector>
10 #include <common_cpp.h>
11 #include <io_cpp.h>
12 #include <timeintegration_cpp.h>
13 #include <mpivars_cpp.h>
14 #include <simulation_object.h>
15 
16 #ifdef with_librom
17 #include <librom_interface.h>
18 #endif
19 
20 #ifdef compute_rhs_operators
21 extern "C" int ComputeRHSOperators(void*,void*,double);
22 #endif
23 extern "C" int CalculateError(void*,void*);
24 int OutputSolution(void*,int,double);
25 extern "C" void ResetFilenameIndex(char*, int);
26 #ifdef with_librom
27 extern "C" int CalculateROMDiff(void*,void*);
28 int OutputROMSolution(void*,int,double);
29 #endif
30 
37 int Solve( void *s,
38  int nsims,
39  int rank,
40  int nproc
41  )
42 {
44 
45  /* make sure none of the simulation objects sent in the array
46  * are "barebones" type */
47  for (int ns = 0; ns < nsims; ns++) {
48  if (sim[ns].is_barebones == 1) {
49  fprintf(stderr, "Error in Solve(): simulation object %d on rank %d is barebones!\n",
50  ns, rank );
51  return 1;
52  }
53  }
54 
55  /* write out iblank to file for visualization */
56  for (int ns = 0; ns < nsims; ns++) {
57  if (sim[ns].solver.flag_ib) {
58 
59  char fname_root[_MAX_STRING_SIZE_] = "iblank";
60  if (nsims > 1) {
61  char index[_MAX_STRING_SIZE_];
62  GetStringFromInteger(ns, index, (int)log10((nsims)+1));
63  strcat(fname_root, "_");
64  strcat(fname_root, index);
65  }
66 
67  WriteArray( sim[ns].solver.ndims,
68  1,
69  sim[ns].solver.dim_global,
70  sim[ns].solver.dim_local,
71  sim[ns].solver.ghosts,
72  sim[ns].solver.x,
73  sim[ns].solver.iblank,
74  &(sim[ns].solver),
75  &(sim[ns].mpi),
76  fname_root );
77  }
78  }
79 
80 #ifdef with_librom
81  if (!rank) printf("Setting up libROM interface.\n");
82  libROMInterface rom_interface( sim, nsims, rank, nproc, sim[0].solver.dt );
83  const std::string& rom_mode( rom_interface.mode() );
84  std::vector<double> op_times_arr(0);
85 #endif
86 
87 #ifdef with_librom
88  if ((rom_mode == _ROM_MODE_TRAIN_) || (rom_mode == _ROM_MODE_NONE_)) {
89 #endif
90  /* Define and initialize the time-integration object */
91  TimeIntegration TS;
92  if (!rank) printf("Setting up time integration.\n");
93  TimeInitialize(sim, nsims, rank, nproc, &TS);
94  double ti_runtime = 0.0;
95 
96  if (!rank) printf("Solving in time (from %d to %d iterations)\n",TS.restart_iter,TS.n_iter);
97  for (TS.iter = TS.restart_iter; TS.iter < TS.n_iter; TS.iter++) {
98 
99  /* Write initial solution to file if this is the first iteration */
100  if (!TS.iter) {
101  for (int ns = 0; ns < nsims; ns++) {
102  if (sim[ns].solver.PhysicsOutput) {
103  sim[ns].solver.PhysicsOutput( &(sim[ns].solver),
104  &(sim[ns].mpi),
105  TS.waqt );
106  }
107  }
108  OutputSolution(sim, nsims, TS.waqt);
109 #ifdef with_librom
110  op_times_arr.push_back(TS.waqt);
111 #endif
112  }
113 
114 #ifdef with_librom
115  if ((rom_mode == _ROM_MODE_TRAIN_) && (TS.iter%rom_interface.samplingFrequency() == 0)) {
116  rom_interface.takeSample( sim, TS.waqt );
117  }
118 #endif
119 
120  /* Call pre-step function */
121  TimePreStep (&TS);
122 #ifdef compute_rhs_operators
123  /* compute and write (to file) matrix operators representing the right-hand side */
124 // if (((TS.iter+1)%solver->file_op_iter == 0) || (!TS.iter))
125 // { ComputeRHSOperators(solver,mpi,TS.waqt);
126 #endif
127 
128  /* Step in time */
129  TimeStep (&TS);
130 
131  /* Call post-step function */
132  TimePostStep (&TS);
133 
134  ti_runtime += TS.iter_wctime;
135 
136  /* Print information to screen */
137  TimePrintStep(&TS);
138 
139  /* Write intermediate solution to file */
140  if ( ((TS.iter+1)%sim[0].solver.file_op_iter == 0)
141  && ((TS.iter+1) < TS.n_iter) ) {
142  for (int ns = 0; ns < nsims; ns++) {
143  if (sim[ns].solver.PhysicsOutput) {
144  sim[ns].solver.PhysicsOutput( &(sim[ns].solver),
145  &(sim[ns].mpi),
146  TS.waqt );
147  }
148  }
149  OutputSolution(sim, nsims, TS.waqt);
150 #ifdef with_librom
151  op_times_arr.push_back(TS.waqt);
152 #endif
153  }
154 
155  }
156 
157  double t_final = TS.waqt;
158  TimeCleanup(&TS);
159 
160  if (!rank) {
161  printf( "Completed time integration (Final time: %f), total wctime: %f (seconds).\n",
162  t_final, ti_runtime );
163  if (nsims > 1) printf("\n");
164  }
165 
166  /* calculate error if exact solution has been provided */
167  for (int ns = 0; ns < nsims; ns++) {
168  CalculateError(&(sim[ns].solver),
169  &(sim[ns].mpi) );
170  }
171 
172  /* write a final solution file */
173  for (int ns = 0; ns < nsims; ns++) {
174  if (sim[ns].solver.PhysicsOutput) {
175  sim[ns].solver.PhysicsOutput( &(sim[ns].solver),
176  &(sim[ns].mpi),
177  t_final );
178  }
179  }
180  OutputSolution(sim, nsims, t_final);
181 
182 #ifdef with_librom
183  op_times_arr.push_back(TS.waqt);
184 
185  for (int ns = 0; ns < nsims; ns++) {
186  ResetFilenameIndex( sim[ns].solver.filename_index,
187  sim[ns].solver.index_length );
188  }
189 
190  if (rom_interface.mode() == _ROM_MODE_TRAIN_) {
191 
192  rom_interface.train();
193  if (!rank) printf("libROM: total training wallclock time: %f (seconds).\n",
194  rom_interface.trainWallclockTime() );
195 
196  double total_rom_predict_time = 0;
197  for (int iter = 0; iter < op_times_arr.size(); iter++) {
198 
199  double waqt = op_times_arr[iter];
200 
201  rom_interface.predict(sim, waqt);
202  if (!rank) printf( "libROM: Predicted solution at time %1.4e using ROM, wallclock time: %f.\n",
203  waqt, rom_interface.predictWallclockTime() );
204  total_rom_predict_time += rom_interface.predictWallclockTime();
205 
206  /* calculate diff between ROM and PDE solutions */
207  if (iter == (op_times_arr.size()-1)) {
208  if (!rank) printf("libROM: Calculating diff between PDE and ROM solutions.\n");
209  for (int ns = 0; ns < nsims; ns++) {
210  CalculateROMDiff( &(sim[ns].solver),
211  &(sim[ns].mpi) );
212  }
213  }
214  /* write the ROM solution to file */
215  OutputROMSolution(sim, nsims,waqt);
216 
217  }
218 
219  if (!rank) {
220  printf( "libROM: total prediction/query wallclock time: %f (seconds).\n",
221  total_rom_predict_time );
222  }
223 
224  rom_interface.saveROM();
225 
226  } else {
227 
228  for (int ns = 0; ns < nsims; ns++) {
229  sim[ns].solver.rom_diff_norms[0]
230  = sim[ns].solver.rom_diff_norms[1]
231  = sim[ns].solver.rom_diff_norms[2]
232  = -1;
233  }
234 
235  }
236 
237  } else if (rom_mode == _ROM_MODE_PREDICT_) {
238 
239  for (int ns = 0; ns < nsims; ns++) {
240  sim[ns].solver.rom_diff_norms[0]
241  = sim[ns].solver.rom_diff_norms[1]
242  = sim[ns].solver.rom_diff_norms[2]
243  = -1;
244  strcpy(sim[ns].solver.ConservationCheck,"no");
245  }
246 
247  rom_interface.loadROM();
248  rom_interface.projectInitialSolution(sim);
249 
250  {
251  int start_iter = sim[0].solver.restart_iter;
252  int n_iter = sim[0].solver.n_iter;
253  double dt = sim[0].solver.dt;
254 
255  double cur_time = start_iter * dt;
256  op_times_arr.push_back(cur_time);
257 
258  for (int iter = start_iter; iter < n_iter; iter++) {
259  cur_time += dt;
260  if ( ( (iter+1)%sim[0].solver.file_op_iter == 0)
261  && ( (iter+1) < n_iter) ) {
262  op_times_arr.push_back(cur_time);
263  }
264  }
265 
266  double t_final = n_iter*dt;
267  op_times_arr.push_back(t_final);
268  }
269 
270  double total_rom_predict_time = 0;
271  for (int iter = 0; iter < op_times_arr.size(); iter++) {
272 
273  double waqt = op_times_arr[iter];
274 
275  rom_interface.predict(sim, waqt);
276  if (!rank) printf( "libROM: Predicted solution at time %1.4e using ROM, wallclock time: %f.\n",
277  waqt, rom_interface.predictWallclockTime() );
278  total_rom_predict_time += rom_interface.predictWallclockTime();
279 
280  /* write the solution to file */
281  for (int ns = 0; ns < nsims; ns++) {
282  if (sim[ns].solver.PhysicsOutput) {
283  sim[ns].solver.PhysicsOutput( &(sim[ns].solver),
284  &(sim[ns].mpi),
285  waqt );
286  }
287  }
288  OutputSolution(sim, nsims, waqt);
289 
290  }
291 
292  /* calculate error if exact solution has been provided */
293  for (int ns = 0; ns < nsims; ns++) {
294  CalculateError(&(sim[ns].solver),
295  &(sim[ns].mpi) );
296  }
297 
298  if (!rank) {
299  printf( "libROM: total prediction/query wallclock time: %f (seconds).\n",
300  total_rom_predict_time );
301  }
302 
303  }
304 #endif
305 
306  return 0;
307 }
int TimePostStep(void *)
Definition: TimePostStep.c:28
int TimePreStep(void *)
Definition: TimePreStep.c:22
int CalculateError(void *s, void *m)
int * dim_global
Definition: hypar.h:33
double * iblank
Definition: hypar.h:436
int TimeStep(void *)
Definition: TimeStep.c:13
Function declarations for file I/O functions.
void GetStringFromInteger(int, char *, int)
void projectInitialSolution(void *a_s)
void saveROM(const std::string &a_fname_root="") const
Structure of variables/parameters and function pointers for time integration.
int * dim_local
Definition: hypar.h:37
int(* PhysicsOutput)(void *, void *, double)
Definition: hypar.h:347
int TimeCleanup(void *)
Definition: TimeCleanup.c:18
void loadROM(const std::string &a_fname_root="")
libROM interface class
int OutputROMSolution(void *s, int nsims, double a_time)
Contains C++ function declarations for time integration.
int ghosts
Definition: hypar.h:52
double dt
Definition: hypar.h:67
void ResetFilenameIndex(char *f, int len)
int WriteArray(int, int, int *, int *, int, double *, double *, void *, void *, char *)
Definition: WriteArray.c:27
#define _MAX_STRING_SIZE_
Definition: basic.h:14
int samplingFrequency() const
double predictWallclockTime() const
Class implementing interface with libROM.
int CalculateROMDiff(void *s, void *m)
Simulation object.
Some common functions used here and there (C++ declarations)
int file_op_iter
Definition: hypar.h:171
int Solve(void *, int, int, int)
Definition: Solve.cpp:37
int OutputSolution(void *, int, double)
#define _ROM_MODE_PREDICT_
double trainWallclockTime() const
int n_iter
Definition: hypar.h:55
void predict(void *a_s, const double a_t) const
void takeSample(void *a_s, const double a_t)
Structure defining a simulation.
const std::string & mode() const
int TimePrintStep(void *)
Definition: TimePrintStep.c:16
#define _ROM_MODE_NONE_
int index_length
Definition: hypar.h:199
double * x
Definition: hypar.h:107
#define _ROM_MODE_TRAIN_
C++ declarations for MPI-related functions.
double rom_diff_norms[3]
Definition: hypar.h:405
int restart_iter
Definition: hypar.h:58
int TimeInitialize(void *, int, int, int, void *)
int ComputeRHSOperators(void *s, void *m, double t)