REMORA
Regional Modeling of Oceans Refined Adaptively
Loading...
Searching...
No Matches
REMORA.cpp
Go to the documentation of this file.
1/**
2 * \file REMORA.cpp
3 */
4
6#include <REMORA.H>
7
8#include <AMReX_buildInfo.H>
9
10using namespace amrex;
11
12amrex::Real REMORA::startCPUTime = 0.0_rt;
13amrex::Real REMORA::previousCPUTimeUsed = 0.0_rt;
14
15Vector<AMRErrorTag> REMORA::ref_tags;
16
18
19// Time step control
20amrex::Real REMORA::cfl = 0.8_rt;
21amrex::Real REMORA::fixed_dt = -1.0_rt;
22amrex::Real REMORA::fixed_fast_dt = -1.0_rt;
23amrex::Real REMORA::change_max = 1.1_rt;
24
26
27// Dictate verbosity in screen output
28int REMORA::verbose = 0;
29
30// Frequency of diagnostic output
32amrex::Real REMORA::sum_per = -1.0_rt;
33
34// Minimum number of digits in plotfile name
36
37// Native AMReX vs NetCDF
39
40#ifdef REMORA_USE_NETCDF
41
43
44// Do we write one file per timestep (false) or one file for all timesteps (true)
46
47// NetCDF initialization file
48amrex::Vector<std::string> REMORA::nc_bdry_file = {""}; // Must provide via input
49amrex::Vector<amrex::Vector<std::string>> REMORA::nc_init_file = {{""}}; // Must provide via input
50amrex::Vector<amrex::Vector<std::string>> REMORA::nc_grid_file = {{""}}; // Must provide via input
51#endif
52
53amrex::Vector<std::string> BCNames = {"xlo", "ylo", "zlo", "xhi", "yhi", "zhi"};
54
55/**
56 * constructor:
57 * - reads in parameters from inputs file
58 * - sizes multilevel arrays and data structures
59 * - initializes BCRec boundary condition object
60 */
62{
63 if (ParallelDescriptor::IOProcessor()) {
64 const char* remora_hash = amrex::buildInfoGetGitHash(1);
65 const char* amrex_hash = amrex::buildInfoGetGitHash(2);
66 const char* buildgithash = amrex::buildInfoGetBuildGitHash();
67 const char* buildgitname = amrex::buildInfoGetBuildGitName();
68
69 if (strlen(remora_hash) > 0) {
70 amrex::Print() << "\n"
71 << "REMORA git hash: " << remora_hash << "\n";
72 }
73 if (strlen(amrex_hash) > 0) {
74 amrex::Print() << "AMReX git hash: " << amrex_hash << "\n";
75 }
76 if (strlen(buildgithash) > 0) {
77 amrex::Print() << buildgitname << " git hash: " << buildgithash << "\n";
78 }
79
80 amrex::Print() << "\n";
81 }
82
84 const std::string& pv1 = "plot_vars"; setPlotVariables(pv1);
85
86 prob = amrex_probinit(geom[0].ProbLo(),geom[0].ProbHi());
87
88 // Geometry on all levels has been defined already.
89
90 // No valid BoxArray and DistributionMapping have been defined.
91 // But the arrays for them have been resized.
92
93 int nlevs_max = max_level + 1;
94
95 istep.resize(nlevs_max, 0);
96 nsubsteps.resize(nlevs_max, 1);
97 for (int lev = 1; lev <= max_level; ++lev) {
98 nsubsteps[lev] = do_substep ? MaxRefRatio(lev-1) : 1;
99 }
100
101 physbcs.resize(nlevs_max);
102
103 t_new.resize(nlevs_max, 0.0_rt);
104 t_old.resize(nlevs_max, -1.e100_rt);
105 dt.resize(nlevs_max, 1.e100_rt);
106
107 cons_new.resize(nlevs_max);
108 cons_old.resize(nlevs_max);
109 xvel_new.resize(nlevs_max);
110 xvel_old.resize(nlevs_max);
111 yvel_new.resize(nlevs_max);
112 yvel_old.resize(nlevs_max);
113 zvel_new.resize(nlevs_max);
114 zvel_old.resize(nlevs_max);
115
116 advflux_reg.resize(nlevs_max);
117
118 // Initialize tagging criteria for mesh refinement
120
121 // We have already read in the ref_Ratio (via amr.ref_ratio =) but we need to enforce
122 // that there is no refinement in the vertical so we test on that here.
123 for (int lev = 0; lev < max_level; ++lev)
124 {
125 amrex::Print() << "Refinement ratio at level " << lev << " set to be " <<
126 ref_ratio[lev][0] << " " << ref_ratio[lev][1] << " " << ref_ratio[lev][2] << std::endl;
127
128 if (ref_ratio[lev][2] != 1)
129 {
130 amrex::Error("We don't allow refinement in the vertical -- make sure to set ref_ratio = 1 in z");
131 }
132 }
133}
134
136{
137}
138
139void
141{
142 Real cur_time = t_new[0];
143
144 // Take one coarse timestep by calling timeStep -- which recursively calls timeStep
145 // for finer levels (with or without subcycling)
146 for (int step = istep[0]; step < max_step && cur_time < stop_time; ++step)
147 {
148 amrex::Print() << "\nCoarse STEP " << step+1 << " starts ..." << std::endl;
149
150 ComputeDt();
151
152 int lev = 0;
153 int iteration = 1;
154
155 if (max_level == 0) {
156 timeStep(lev, cur_time, iteration);
157 }
158 else {
159 timeStepML(cur_time, iteration);
160 }
161
162 cur_time += dt[0];
163
164 amrex::Print() << "Coarse STEP " << step+1 << " ends." << " TIME = " << cur_time
165 << " DT = " << dt[0] << std::endl;
166
167 if ((plot_int > 0 && (step+1 - last_plot_file_step) == plot_int)
168 || (plot_int_time > 0 && cur_time >= (last_plot_file_time + plot_int_time))) {
169 last_plot_file_step = step+1;
170 last_plot_file_time = cur_time;
172
174 }
175#ifdef REMORA_USE_NETCDF
177 WriteNCPlotFile(step+1);
179 }
180#endif
181 }
182
183 if ((check_int > 0 && (step+1 - last_check_file_step) == check_int)
184 || (check_int_time > 0 && cur_time >= (last_check_file_time + check_int_time))) {
185 last_check_file_step = step+1;
186 last_check_file_time = cur_time;
188 }
189
190 post_timestep(step, cur_time, dt[0]);
191
192#ifdef AMREX_MEM_PROFILING
193 {
194 std::ostringstream ss;
195 ss << "[STEP " << step+1 << "]";
196 MemProfiler::report(ss.str());
197 }
198#endif
199
200 if (cur_time >= stop_time - 1.e-6*dt[0]) break;
201 }
202
203 if ((plot_int > 0 || plot_int_time > 0.0) && istep[0] > last_plot_file_step) {
206 }
207#ifdef REMORA_USE_NETCDF
211 }
212#endif
213 }
214
215 if ((check_int > 0 || check_int_time > 0.0) && istep[0] > last_check_file_step) {
217 }
218}
219
220/**
221 * @param[in ] nstep which step we're on
222 * @param[in ] time current time
223 * @param[in ] dt_lev0 time step on level 0
224 */
225void
226REMORA::post_timestep (int nstep, Real time, Real dt_lev0)
227{
228 BL_PROFILE("REMORA::post_timestep()");
229
230#ifdef REMORA_USE_PARTICLES
231 particleData.Redistribute();
232#endif
233
235 {
236 for (int lev = finest_level-1; lev >= 0; lev--)
237 {
238 // This call refluxes from the lev/lev+1 interface onto lev
239 //getAdvFluxReg(lev+1)->Reflux(*cons_new[lev], 0, 0, NCONS);
240
241 // We need to do this before anything else because refluxing changes the
242 // values of coarse cells underneath fine grids with the assumption they'll
243 // be over-written by averaging down
244 //
245 AverageDownTo(lev);
246 }
247 }
248
249 if (is_it_time_for_action(nstep, time, dt_lev0, sum_interval, sum_per)) {
251 }
252}
253
254/**
255 * This is called from main.cpp and handles all initialization, whether from start or restart
256 */
257void
259{
260 // Initialize the start time for our CPU-time tracker
261 startCPUTime = Real(ParallelDescriptor::second());
262
263 // Map the words in the inputs file to BC types, then translate
264 // those types into what they mean for each variable
265 init_bcs();
266
269 last_plot_file_time = -1.0_rt;
270 last_check_file_time = -1.0_rt;
271
272 if (restart_chkfile == "") {
273 // start simulation from the beginning
274
275 InitFromScratch(start_time);
276
278 AverageDown();
279 }
280
281 } else { // Restart from a checkpoint
282
283 restart();
284
285 }
286#ifdef REMORA_USE_MOAB
287 InitMOABMesh();
288#endif
289 // Initialize flux registers (whether we start from scratch or restart)
291 advflux_reg[0] = nullptr;
292 for (int lev = 1; lev <= finest_level; lev++)
293 {
294 advflux_reg[lev] = new YAFluxRegister(grids[lev], grids[lev-1],
295 dmap[lev], dmap[lev-1],
296 geom[lev], geom[lev-1],
297 ref_ratio[lev-1], lev, NCONS);
298 }
299 }
300
301 // Fill ghost cells/faces
302 for (int lev = 0; lev <= finest_level; ++lev)
303 {
304 if (lev > 0 && cf_width >= 0) {
306 }
307
308 if (restart_chkfile == "") {
309 FillPatch(lev, t_new[lev], *cons_new[lev], cons_new, BCVars::cons_bc, BdyVars::t, 0, true, false,0,0,0.0,*cons_new[lev]);
310 FillPatch(lev, t_new[lev], *xvel_new[lev], xvel_new, BCVars::xvel_bc, BdyVars::u, 0, true, false,0,0,0.0,*xvel_new[lev]);
311 FillPatch(lev, t_new[lev], *yvel_new[lev], yvel_new, BCVars::yvel_bc, BdyVars::v, 0, true, false,0,0,0.0,*yvel_new[lev]);
312 FillPatch(lev, t_new[lev], *zvel_new[lev], zvel_new, BCVars::zvel_bc, BdyVars::null, 0, true, false);
313
314 // Copy from new into old just in case when initializing from scratch
315 int ngs = cons_new[lev]->nGrow();
316 int ngvel = xvel_new[lev]->nGrow();
317 MultiFab::Copy(*cons_old[lev],*cons_new[lev],0,0,NCONS,ngs);
318 MultiFab::Copy(*xvel_old[lev],*xvel_new[lev],0,0,1,ngvel);
319 MultiFab::Copy(*yvel_old[lev],*yvel_new[lev],0,0,1,ngvel);
320 MultiFab::Copy(*zvel_old[lev],*zvel_new[lev],0,0,1,IntVect(ngvel,ngvel,0));
321 }
322 } // lev
323
324 // Check for additional plotting variables that are available after
325 // particle containers are setup.
326 const std::string& pv1 = "plot_vars"; appendPlotVariables(pv1);
327
328 if (restart_chkfile == "" && (check_int > 0 || check_int_time > 0.0_rt))
329 {
332 }
333
334 if ( (restart_chkfile == "") ||
336 {
337 if (plot_int > 0 || plot_int_time > 0.0)
338 {
341#ifdef REMORA_USE_NETCDF
343 int step0 = 0;
344 WriteNCPlotFile(step0);
346 }
347#endif
349 }
350 }
351
354 }
355
356 ComputeDt();
357
358}
359
360/**
361 * @param[in ] lev level to operate on
362 */
363void
365{
366 amrex::Print() << ":::Construct_REMORAFillPatchers " << lev << std::endl;
367
368 auto& ba_fine = cons_new[lev ]->boxArray();
369 auto& ba_crse = cons_new[lev-1]->boxArray();
370 auto& dm_fine = cons_new[lev ]->DistributionMap();
371 auto& dm_crse = cons_new[lev-1]->DistributionMap();
372
373 BoxList bl2d_fine = ba_fine.boxList();
374 for (auto& b : bl2d_fine) {
375 b.setRange(2,0);
376 }
377 BoxArray ba2d_fine(std::move(bl2d_fine));
378
379 BoxList bl2d_crse = ba_crse.boxList();
380 for (auto& b : bl2d_crse) {
381 b.setRange(2,0);
382 }
383 BoxArray ba2d_crse(std::move(bl2d_crse));
384
385 int ncomp = cons_new[lev]->nComp();
386
387 FPr_c.emplace_back(ba_fine, dm_fine, geom[lev] ,
388 ba_crse, dm_crse, geom[lev-1],
389 -cf_width, -cf_set_width, ncomp, &cell_cons_interp);
390 FPr_u.emplace_back(convert(ba_fine, IntVect(1,0,0)), dm_fine, geom[lev] ,
391 convert(ba_crse, IntVect(1,0,0)), dm_crse, geom[lev-1],
392 -cf_width, -cf_set_width, 1, &face_cons_linear_interp);
393 FPr_v.emplace_back(convert(ba_fine, IntVect(0,1,0)), dm_fine, geom[lev] ,
394 convert(ba_crse, IntVect(0,1,0)), dm_crse, geom[lev-1],
395 -cf_width, -cf_set_width, 1, &face_cons_linear_interp);
396 FPr_w.emplace_back(convert(ba_fine, IntVect(0,0,1)), dm_fine, geom[lev] ,
397 convert(ba_crse, IntVect(0,0,1)), dm_crse, geom[lev-1],
398 -cf_width, -cf_set_width, 1, &face_cons_linear_interp);
399
400 FPr_ubar.emplace_back(convert(ba2d_fine, IntVect(1,0,0)), dm_fine, geom[lev] ,
401 convert(ba2d_crse, IntVect(1,0,0)), dm_crse, geom[lev-1],
402 -cf_width, -cf_set_width, 3, &face_cons_linear_interp);
403 FPr_vbar.emplace_back(convert(ba2d_fine, IntVect(0,1,0)), dm_fine, geom[lev] ,
404 convert(ba2d_crse, IntVect(0,1,0)), dm_crse, geom[lev-1],
405 -cf_width, -cf_set_width, 3, &face_cons_linear_interp);
406}
407
408/**
409 * @param[in ] lev level to operate on
410 */
411void
413{
414 amrex::Print() << ":::Define_REMORAFillPatchers " << lev << std::endl;
415
416 auto& ba_fine = cons_new[lev ]->boxArray();
417 auto& ba_crse = cons_new[lev-1]->boxArray();
418 auto& dm_fine = cons_new[lev ]->DistributionMap();
419 auto& dm_crse = cons_new[lev-1]->DistributionMap();
420
421 BoxList bl2d_fine = ba_fine.boxList();
422 for (auto& b : bl2d_fine) {
423 b.setRange(2,0);
424 }
425 BoxArray ba2d_fine(std::move(bl2d_fine));
426
427 BoxList bl2d_crse = ba_crse.boxList();
428 for (auto& b : bl2d_crse) {
429 b.setRange(2,0);
430 }
431 BoxArray ba2d_crse(std::move(bl2d_crse));
432
433
434 int ncomp = cons_new[lev]->nComp();
435
436 FPr_c[lev-1].Define(ba_fine, dm_fine, geom[lev] ,
437 ba_crse, dm_crse, geom[lev-1],
438 -cf_width, -cf_set_width, ncomp, &cell_cons_interp);
439 FPr_u[lev-1].Define(convert(ba_fine, IntVect(1,0,0)), dm_fine, geom[lev] ,
440 convert(ba_crse, IntVect(1,0,0)), dm_crse, geom[lev-1],
441 -cf_width, -cf_set_width, 1, &face_cons_linear_interp);
442 FPr_v[lev-1].Define(convert(ba_fine, IntVect(0,1,0)), dm_fine, geom[lev] ,
443 convert(ba_crse, IntVect(0,1,0)), dm_crse, geom[lev-1],
444 -cf_width, -cf_set_width, 1, &face_cons_linear_interp);
445 FPr_w[lev-1].Define(convert(ba_fine, IntVect(0,0,1)), dm_fine, geom[lev] ,
446 convert(ba_crse, IntVect(0,0,1)), dm_crse, geom[lev-1],
447 -cf_width, -cf_set_width, 1, &face_cons_linear_interp);
448
449 FPr_ubar[lev-1].Define(convert(ba2d_fine, IntVect(1,0,0)), dm_fine, geom[lev] ,
450 convert(ba2d_crse, IntVect(1,0,0)), dm_crse, geom[lev-1],
451 -cf_width, -cf_set_width, 3, &face_cons_linear_interp);
452 FPr_vbar[lev-1].Define(convert(ba2d_fine, IntVect(0,1,0)), dm_fine, geom[lev] ,
453 convert(ba2d_crse, IntVect(0,1,0)), dm_crse, geom[lev-1],
454 -cf_width, -cf_set_width, 3, &face_cons_linear_interp);
455}
456
457void
459{
461
462 // We set this here so that we don't over-write the checkpoint file we just started from
464}
465
466/**
467 * @param[in ] lev level to operate on
468 */
469void
471{
473 prob->init_analytic_zeta(lev, geom[lev], solverChoice, *this, *vec_zeta[lev]);
474
475#ifdef REMORA_USE_NETCDF
476 } else if (solverChoice.ic_type == IC_Type::netcdf) {
477 amrex::Print() << "Calling init_zeta_from_netcdf on level " << lev << std::endl;
479 amrex::Print() << "Sea surface height loaded from netcdf file \n " << std::endl;
480#endif
481 } else {
482 Abort("Don't know this ic_type!");
483 }
484 vec_zeta[lev]->FillBoundary(geom[lev].periodicity());
485 set_zeta_average(lev);
486}
487
488/**
489 * @param[in ] lev level to operate on
490 */
491void
493{
494 // Only set bathymetry on level 0, and interpolate for finer levels
496 if (lev==0) {
499 prob->init_analytic_bathymetry(lev, geom[lev], solverChoice, *this, *vec_hOfTheConfusingName[lev]);
500 } else {
502 }
503#ifdef REMORA_USE_NETCDF
504 } else if (solverChoice.ic_type == IC_Type::netcdf) {
506 amrex::Print() << "Calling init_bathymetry_from_netcdf " << std::endl;
508 amrex::Print() << "Bathymetry loaded from netcdf file \n " << std::endl;
509 } else {
511 }
512#endif
513 } else {
514 Abort("Don't know this ic_type!");
515 }
516 // Need FillBoundary to fill at grid-grid boundaries, and EnforcePeriodicity
517 // to make sure ghost cells in the domain corners are consistent.
518 vec_hOfTheConfusingName[lev]->FillBoundary(geom[lev].periodicity());
519 vec_hOfTheConfusingName[lev]->EnforcePeriodicity(geom[lev].periodicity());
520 } else {
521 Real dummy_time = 0.0_rt;
522 FillCoarsePatch(lev,dummy_time,vec_hOfTheConfusingName[lev].get(), vec_hOfTheConfusingName[lev-1].get());
523 }
524 } else if (solverChoice.init_ana_h) {
527 prob->init_analytic_bathymetry(lev, geom[lev], solverChoice, *this,*vec_hOfTheConfusingName[lev]);
528 } else {
530 }
531#ifdef REMORA_USE_NETCDF
532 } else if (solverChoice.ic_type == IC_Type::netcdf) {
533 amrex::Print() << "Calling init_bathymetry_from_netcdf level " << lev << std::endl;
535 amrex::Print() << "Bathymetry loaded from netcdf file \n " << std::endl;
536#endif
537 } else {
538 Abort("Don't know this ic_type!");
539 }
540 // Need FillBoundary to fill at grid-grid boundaries, and EnforcePeriodicity
541 // to make sure ghost cells in the domain corners are consistent.
542 vec_hOfTheConfusingName[lev]->FillBoundary(geom[lev].periodicity());
543 vec_hOfTheConfusingName[lev]->EnforcePeriodicity(geom[lev].periodicity());
544 } else if (solverChoice.init_l1ad_h) {
547 prob->init_analytic_bathymetry(lev, geom[lev], solverChoice, *this, *vec_hOfTheConfusingName[lev]);
548 } else {
550 }
551#ifdef REMORA_USE_NETCDF
552 } else if (solverChoice.ic_type == IC_Type::netcdf) {
554 amrex::Print() << "Calling init_bathymetry_from_netcdf " << std::endl;
556 amrex::Print() << "Bathymetry loaded from netcdf file \n " << std::endl;
557 } else {
559 }
560#endif
561 } else {
562 Abort("Don't know this ic_type!");
563 }
564 // Need FillBoundary to fill at grid-grid boundaries, and EnforcePeriodicity
565 // to make sure ghost cells in the domain corners are consistent.
566 vec_hOfTheConfusingName[lev]->FillBoundary(geom[lev].periodicity());
567 vec_hOfTheConfusingName[lev]->EnforcePeriodicity(geom[lev].periodicity());
568 } else {
569 amrex::Abort("Don't know this h init type");
570 }
571}
572
573/**
574 * @param[in ] lev level to operate on
575 */
576void
580 prob->init_analytic_coriolis(lev, geom[lev], solverChoice, *this, *vec_fcor[lev]);
583#ifdef REMORA_USE_NETCDF
585 amrex::Print() << "Calling init_coriolis_from_netcdf " << std::endl;
587 amrex::Print() << "Coriolis loaded from netcdf file \n" << std::endl;
588#endif
589 } else {
590 Abort("Don't know this coriolis_type!");
591 }
592
593 Real time = 0.0_rt;
594 FillPatch(lev, time, *vec_fcor[lev], GetVecOfPtrs(vec_fcor),BCVars::foextrap_bc);
595 vec_fcor[lev]->EnforcePeriodicity(geom[lev].periodicity());
596 }
597}
598
599void
605 // The GLS initialization just sets the multifab to a value, so there's
606 // no need to call FillPatch here
607 } else {
608 Abort("Don't know this vertical mixing type");
609 }
610}
611
612/**
613 * @param[in ] lev level to operate on
614 */
615void
617 Real time = 0.0_rt;
618 prob->init_analytic_vmix(lev, geom[lev], solverChoice, *this,*vec_Akv[lev], *vec_Akt[lev]);
619 FillPatch(lev, time, *vec_Akv[lev], GetVecOfPtrs(vec_Akv),BCVars::zvel_bc,BdyVars::null,0,true,false);
620 for (int n=0; n<NCONS;n++) {
621 FillPatch(lev, time, *vec_Akt[lev], GetVecOfPtrs(vec_Akt),BCVars::zvel_bc,BdyVars::null,0,false,false);
622 }
623}
624
625/**
626 * @param[in ] lev level to operate on
627 */
628void
630{
632 prob->init_analytic_hmix(lev, geom[lev], solverChoice, *this, *vec_visc2_p[lev], *vec_visc2_r[lev], *vec_diff2[lev]);
634 vec_visc2_p[lev]->setVal(solverChoice.visc2);
635 vec_visc2_r[lev]->setVal(solverChoice.visc2);
636 for (int n=0; n<NCONS; n++) {
637 vec_diff2[lev]->setVal(solverChoice.tnu2[n],n,1);
638 }
639 } else {
640 Abort("Don't know this horizontal mixing type");
641 }
642 Real time = 0.0_rt;
643 FillPatch(lev, time, *vec_visc2_p[lev], GetVecOfPtrs(vec_visc2_p),BCVars::foextrap_periodic_bc);
644 FillPatch(lev, time, *vec_visc2_r[lev], GetVecOfPtrs(vec_visc2_r),BCVars::foextrap_periodic_bc);
645 for (int n=0; n<NCONS; n++) {
646 FillPatch(lev, time, *vec_diff2[lev] , GetVecOfPtrs(vec_diff2),BCVars::foextrap_periodic_bc,BdyVars::null,n,false);
647 }
648}
649
650/**
651 * @param[in ] lev level to operate on
652 */
653void
655{
656 vec_hOfTheConfusingName[lev]->setVal(-geom[0].ProbLo()[2]);
657}
658
659/**
660 * @param[in ] lev level to operate on
661 */
662void
664{
666 prob->init_analytic_smflux(lev, geom[lev], solverChoice, *this,*vec_sustr[lev], *vec_svstr[lev]);
668#ifdef REMORA_USE_NETCDF
671 FillPatch(lev, t_old[lev], *vec_sustr[lev], GetVecOfPtrs(vec_sustr),BCVars::foextrap_periodic_bc,BdyVars::null,0,false);
672 FillPatch(lev, t_old[lev], *vec_svstr[lev], GetVecOfPtrs(vec_svstr),BCVars::foextrap_periodic_bc,BdyVars::null,0,false);
673#endif
674 }
675}
676
677/**
678 * @param[in ] lev level to operate on
679 */
680void
682{
684 prob->init_analytic_wind(lev,geom[lev], solverChoice, *this, *vec_uwind[lev], *vec_vwind[lev]);
686#ifdef REMORA_USE_NETCDF
689 FillPatch(lev, t_old[lev], *vec_uwind[lev], GetVecOfPtrs(vec_uwind),BCVars::foextrap_periodic_bc,BdyVars::null,0,false);
690 FillPatch(lev, t_old[lev], *vec_vwind[lev], GetVecOfPtrs(vec_vwind),BCVars::foextrap_periodic_bc,BdyVars::null,0,false);
691#endif
692 }
693}
694
695/**
696 * @param[in ] lev level to operate on
697 */
698void
700{
702 prob->init_analytic_masks(lev,geom[lev], solverChoice, *this, *vec_mskr[lev]);
705#ifdef REMORA_USE_NETCDF
706 amrex::Print() << "Calling init_masks_from_netcdf level " << lev << std::endl;
708 amrex::Print() << "Masks loaded from netcdf file \n " << std::endl;
709#endif
710 }
711}
712
713/**
714 * @param[in ] lev level to operate on
715 * @param[in ] time current time for initialization
716 */
717void
718REMORA::init_only (int lev, Real time)
719{
720 t_new[lev] = time;
721 t_old[lev] = time - 1.e200_rt;
722
723 cons_new[lev]->setVal(0.0_rt);
724 xvel_new[lev]->setVal(0.0_rt);
725 yvel_new[lev]->setVal(0.0_rt);
726 zvel_new[lev]->setVal(0.0_rt);
727
728 xvel_old[lev]->setVal(0.0_rt);
729 yvel_old[lev]->setVal(0.0_rt);
730 zvel_old[lev]->setVal(0.0_rt);
731
732 vec_ru[lev]->setVal(0.0_rt);
733 vec_rv[lev]->setVal(0.0_rt);
734
735 vec_ru2d[lev]->setVal(0.0_rt);
736 vec_rv2d[lev]->setVal(0.0_rt);
737
739 set_grid_scale(lev);
740 }
741 set_masks(lev);
742
743#ifdef REMORA_USE_NETCDF
746
748 if (nc_clim_his_file.empty()) {
749 amrex::Error("NetCDF climatology file name must be provided via input");
750 }
752 ubar_clim_data_from_file = new NCTimeSeries(nc_clim_his_file, "ubar", clim_ubar_time_varname, geom[lev].Domain(),vec_ubar[lev].get(),true,true);
753 vbar_clim_data_from_file = new NCTimeSeries(nc_clim_his_file, "vbar", clim_ubar_time_varname, geom[lev].Domain(),vec_vbar[lev].get(),true,true);
756 }
758 u_clim_data_from_file = new NCTimeSeries(nc_clim_his_file, "u", clim_u_time_varname, geom[lev].Domain(),xvel_new[lev],false,true);
759 v_clim_data_from_file = new NCTimeSeries(nc_clim_his_file, "v", clim_v_time_varname, geom[lev].Domain(),yvel_new[lev],false,true);
762 }
763 // Since the NCTimeSeries object isn't filling the cons_new MultiFab directly, we don't have to specify a component.
764 // It just needs to know the shape of the MultiFab
766 temp_clim_data_from_file = new NCTimeSeries(nc_clim_his_file, "temp", clim_temp_time_varname,geom[lev].Domain(),cons_new[lev],false,true);
768 }
770 salt_clim_data_from_file = new NCTimeSeries(nc_clim_his_file, "salt", clim_salt_time_varname,geom[lev].Domain(),cons_new[lev],false,true);
772 }
773 }
774 }
775
777 amrex::Print() << "Calling init_bdry_from_netcdf " << std::endl;
779 amrex::Print() << "Boundary data loaded from netcdf file \n " << std::endl;
780 }
781
782 // This will be a non-op if forcings specified analytically
784 if (nc_frc_file.empty()) {
785 amrex::Error("NetCDF forcing file name must be provided via input for winds");
786 }
787 Uwind_data_from_file = new NCTimeSeries(nc_frc_file, "Uwind", frc_time_varname, geom[lev].Domain(),vec_uwind[lev].get(), true, false);
788 Vwind_data_from_file = new NCTimeSeries(nc_frc_file, "Vwind", frc_time_varname, geom[lev].Domain(),vec_vwind[lev].get(), true, false);
792 if (nc_frc_file.empty()) {
793 amrex::Error("NetCDF forcing file name must be provided via input for surface momentum fluxes");
794 }
795 sustr_data_from_file = new NCTimeSeries(nc_frc_file, "sustr", frc_time_varname, geom[lev].Domain(),vec_sustr[lev].get(), true, false);
796 svstr_data_from_file = new NCTimeSeries(nc_frc_file, "svstr", frc_time_varname, geom[lev].Domain(),vec_svstr[lev].get(), true, false);
799 }
800
802 auto dom = geom[0].Domain();
803 int nz = dom.length(2);
804 river_source_cons.resize(NCONS);
805 Print() << solverChoice.do_rivers_cons[0] << std::endl;
808 river_source_cons[Salt_comp]->Initialize();
809 }
812 river_source_cons[Temp_comp]->Initialize();
813 }
816 river_source_cons[Scalar_comp]->Initialize();
817 }
823 }
824#else
826 Abort("Not compiled with NetCDF, but selected boundary conditions require NetCDF");
827 }
828#endif
829
830 set_bathymetry(lev);
831 set_zeta(lev);
833
835 if (lev==0) {
837 {
838 init_analytic(lev);
839#ifdef REMORA_USE_NETCDF
840 } else if (solverChoice.ic_type == IC_Type::netcdf) {
841 amrex::Print() << "Calling init_data_from_netcdf " << std::endl;
844 amrex::Print() << "Initial data loaded from netcdf file \n " << std::endl;
845
846#endif
847 } else {
848 Abort("Need to specify ic_type");
849 }
850 } else {
851 FillCoarsePatch(lev, time, cons_new[lev], cons_new[lev-1]);
852 FillCoarsePatch(lev, time, xvel_new[lev], xvel_new[lev-1]);
853 FillCoarsePatch(lev, time, yvel_new[lev], yvel_new[lev-1]);
854 FillCoarsePatch(lev, time, zvel_new[lev], zvel_new[lev-1]);
855 }
858 {
859 init_analytic(lev);
860#ifdef REMORA_USE_NETCDF
861 } else if (solverChoice.ic_type == IC_Type::netcdf) {
862 amrex::Print() << "Calling init_data_from_netcdf on level " << lev << std::endl;
865 amrex::Print() << "Initial data loaded from netcdf file on level " << lev << std::endl;
866
867#endif
868 } else {
869 Abort("Need to specify ic_type");
870 }
871 } else {
872 amrex::Abort("Need to specify T init procedure");
873 }
874
875 // Ensure that the face-based data are the same on both sides of a periodic domain.
876 // The data associated with the lower grid ID is considered the correct value.
877 xvel_new[lev]->OverrideSync(geom[lev].periodicity());
878 yvel_new[lev]->OverrideSync(geom[lev].periodicity());
879 zvel_new[lev]->OverrideSync(geom[lev].periodicity());
880
881 set_2darrays(lev);
882
883 init_set_vmix(lev);
884 set_hmixcoef(lev);
885 set_coriolis(lev);
886
887 // Previously set smflux here with OverrideSync:
888// set_smflux(lev);
889// prob->init_analytic_smflux(lev, geom[lev], solverChoice, *this, *vec_sustr[lev], *vec_svstr[lev]);
890// vec_sustr[lev]->OverrideSync(geom[lev].periodicity());
891// vec_svstr[lev]->OverrideSync(geom[lev].periodicity());
892
893}
894
895void
897{
898 {
899 ParmParse pp; // Traditionally, max_step and stop_time do not have prefix, so allow it for now.
900 bool noprefix_max_step = pp.query("max_step", max_step);
901 bool noprefix_stop_time = pp.query("stop_time", stop_time);
902 bool remora_max_step = pp.query("remora.max_step", max_step);
903 bool remora_stop_time = pp.query("remora.stop_time", stop_time);
904 if (remora_max_step and noprefix_max_step) {
905 Abort("remora.max_step and max_step are both specified. Please use only one!");
906 }
907 if (remora_stop_time and noprefix_stop_time) {
908 Abort("remora.stop_time and stop_time are both specified. Please use only one!");
909 }
910 }
911
912 ParmParse pp(pp_prefix);
913 ParmParse pp_amr("amr");
914 {
915 pp_amr.query("regrid_int", regrid_int);
916 pp.query("check_file", check_file);
917 pp.query("check_int", check_int);
918 pp_amr.query("check_int", check_int);
919 pp.query("check_int_time", check_int_time);
920 pp_amr.query("check_int_time", check_int_time);
921
922 pp.query("restart", restart_chkfile);
923 pp_amr.query("restart", restart_chkfile);
924 pp.query("start_time",start_time);
925
926 if (pp.contains("data_log"))
927 {
928 int num_datalogs = pp.countval("data_log");
929 datalog.resize(num_datalogs);
930 datalogname.resize(num_datalogs);
931 pp.queryarr("data_log",datalogname,0,num_datalogs);
932 for (int i = 0; i < num_datalogs; i++)
934 }
935
936 // Verbosity
937 pp.query("v", verbose);
938
939 // Frequency of diagnostic output
940 pp.query("sum_interval", sum_interval);
941 pp.query("sum_period" , sum_per);
942 pp.query("file_min_digits", file_min_digits);
943
944 if (file_min_digits < 0) {
945 amrex::Abort("remora.file_min_digits must be non-negative");
946 }
947
948 // Time step controls
949 pp.query("cfl", cfl);
950 pp.query("change_max", change_max);
951
952 pp.query("fixed_dt", fixed_dt);
953 pp.query("fixed_fast_dt", fixed_fast_dt);
954
955 pp.query("fixed_ndtfast_ratio", fixed_ndtfast_ratio);
956
957 // If all three are specified, they must be consistent
958 if (fixed_dt > 0. && fixed_fast_dt > 0. && fixed_ndtfast_ratio > 0)
959 {
961 {
962 amrex::Abort("Dt is over-specfied");
963 }
964 }
965 // If two are specified, initialize fixed_ndtfast_ratio
966 else if (fixed_dt > 0. && fixed_fast_dt > 0. && fixed_ndtfast_ratio <= 0)
967 {
968 fixed_ndtfast_ratio = static_cast<int>(fixed_dt / fixed_fast_dt);
969 }
970
971 AMREX_ASSERT(cfl > 0. || fixed_dt > 0.);
972
973 pp_amr.query("do_substep", do_substep);
974 if (do_substep) {
975 amrex::Abort("Time substepping is not yet implemented. amr.do_substep must be 0");
976 }
977
978 // We use this to keep track of how many boxes we read in from WRF initialization
979 num_files_at_level.resize(max_level+1,0);
980
981 // We use this to keep track of how many boxes are specified thru the refinement indicators
982 num_boxes_at_level.resize(max_level+1,0);
983 boxes_at_level.resize(max_level+1);
984
985 // We always have exactly one file at level 0
986 num_boxes_at_level[0] = 1;
987 boxes_at_level[0].resize(1);
988 boxes_at_level[0][0] = geom[0].Domain();
989
990 // Plotfile name and frequency
991 pp.query("plot_file", plot_file_name);
992 pp.query("plot_int", plot_int);
993 pp.query("plot_int_time", plot_int_time);
994 // Output format
995 std::string plotfile_type_str = "amrex";
996 pp.query("plotfile_type", plotfile_type_str);
997 if (plotfile_type_str == "amrex") {
999 } else if (plotfile_type_str == "netcdf" || plotfile_type_str == "NetCDF") {
1001#ifdef REMORA_USE_NETCDF
1002 pp.query("write_history_file",write_history_file);
1003 pp.query("chunk_history_file",chunk_history_file);
1004 pp.query("steps_per_history_file",steps_per_history_file);
1005 // Estimate size of domain for one timestep of netcdf
1006 auto dom = geom[0].Domain();
1007 int nx = dom.length(0) + 2;
1008 int ny = dom.length(1) + 2;
1009 int nz = dom.length(2);
1011 // Estimate number of steps that will fit into a 2GB file.
1012 steps_per_history_file = int((1.6e10 - NCH2D * nx * ny * 64.0_rt)
1013 / (nx * ny * 64.0_rt * (NC3D*nz + NC2D)));
1014 // If we calculate that a single step will exceed 2GB and the user has
1015 // requested automatic history file sizing, warn about a possible impending
1016 // error, and set steps_per_history_file = 1 to attempt output anyway.
1017 if (steps_per_history_file == 0) {
1018 amrex::Warning("NetCDF output for a single timestep appears to exceed 2GB. NetCDF output may not work. See Documentation for information about tested MPICH versions.");
1020 }
1021 } else if (write_history_file and !chunk_history_file) {
1022 // Estimate number of output steps we'll need
1023 int nt_out = int((max_step) / plot_int) + 1;
1024 Real est_hist_file_size = NCH2D * nx * ny * 64.0_rt + nt_out * nx * ny * 64.0_rt * (NC3D*nz + NC2D);
1025 if (est_hist_file_size > 1.6e10) {
1026 amrex::Warning("WARNING: NetCDF history file may be larger than 2GB limit. Consider setting remora.chunk_history_file=true");
1027 }
1028 }
1030 Print() << "NetCDF history files will have " << steps_per_history_file << " steps per file." << std::endl;
1031 }
1032#endif
1033 } else {
1034 amrex::Print() << "User selected plotfile_type = " << plotfile_type_str << std::endl;
1035 amrex::Abort("Dont know this plotfile_type");
1036 }
1037#ifndef REMORA_USE_NETCDF
1039 {
1040 amrex::Abort("Please compile with NetCDF in order to enable NetCDF plotfiles");
1041 }
1042
1043#endif
1044
1045#ifdef REMORA_USE_NETCDF
1046 nc_init_file.resize(max_level+1);
1047 nc_grid_file.resize(max_level+1);
1048
1049 // NetCDF initialization files -- possibly multiple files at each of multiple levels
1050 // but we always have exactly one file at level 0
1051 for (int lev = 0; lev <= max_level; lev++)
1052 {
1053 const std::string nc_file_names = amrex::Concatenate("nc_init_file_",lev,1);
1054 const std::string nc_bathy_file_names = amrex::Concatenate("nc_grid_file_",lev,1);
1055
1056 if (pp.contains(nc_file_names.c_str()))
1057 {
1058 int num_files = pp.countval(nc_file_names.c_str());
1059 int num_bathy_files = pp.countval(nc_bathy_file_names.c_str());
1060 if (num_files != num_bathy_files) {
1061 amrex::Error("Must have same number of netcdf files for grid info as for solution");
1062 }
1063
1064 num_files_at_level[lev] = num_files;
1065 nc_init_file[lev].resize(num_files);
1066 nc_grid_file[lev].resize(num_files);
1067
1068 pp.queryarr(nc_file_names.c_str() , nc_init_file[lev] ,0,num_files);
1069 pp.queryarr(nc_bathy_file_names.c_str(), nc_grid_file[lev],0,num_files);
1070 }
1071 }
1072 // We only read boundary data at level 0
1073 pp.queryarr("nc_bdry_file", nc_bdry_file);
1074
1075 // Also only read forcings at level 0 (for now)
1076 pp.query("nc_frc_file", nc_frc_file);
1077
1078 // Get river file
1079 pp.query("nc_river_file", nc_riv_file);
1080
1081 // Read in file names for climatology history and nudging weights
1082 pp.query("nc_clim_his_file", nc_clim_his_file);
1083 pp.query("nc_clim_coeff_file", nc_clim_coeff_file);
1084
1085 pp.query("bdy_time_varname",bdry_time_varname);
1086 pp.query("frc_time_varname",frc_time_varname);
1087 pp.query("riv_time_varname",riv_time_varname);
1088 pp.query("clim_ubar_time_varname",clim_ubar_time_varname);
1089 pp.query("clim_vbar_time_varname",clim_vbar_time_varname);
1090 pp.query("clim_u_time_varname",clim_u_time_varname);
1091 pp.query("clim_v_time_varname",clim_v_time_varname);
1092 pp.query("clim_salt_time_varname",clim_salt_time_varname);
1093 pp.query("clim_temp_time_varname",clim_temp_time_varname);
1094
1095#endif
1096
1097#ifdef REMORA_USE_PARTICLES
1098 readTracersParams();
1099#endif
1100 }
1101
1103}
1104
1105void
1107{
1108 for (int lev = finest_level-1; lev >= 0; --lev)
1109 {
1110 AverageDownTo(lev);
1111 }
1112}
1113
1114/**
1115 * @param[in ] crse_lev level to average down to
1116 */
1117void
1119{
1120 average_down(*cons_new[crse_lev+1], *cons_new[crse_lev],
1121 0, cons_new[crse_lev]->nComp(), refRatio(crse_lev));
1122 average_down(*vec_zeta[crse_lev+1].get(), *vec_zeta[crse_lev].get(),
1123 0, vec_zeta[crse_lev]->nComp(), refRatio(crse_lev));
1124
1125 Array<MultiFab*,AMREX_SPACEDIM> faces_crse;
1126 Array<MultiFab*,AMREX_SPACEDIM> faces_fine;
1127 faces_crse[0] = xvel_new[crse_lev];
1128 faces_crse[1] = yvel_new[crse_lev];
1129 faces_crse[2] = zvel_new[crse_lev];
1130
1131 faces_fine[0] = xvel_new[crse_lev+1];
1132 faces_fine[1] = yvel_new[crse_lev+1];
1133 faces_fine[2] = zvel_new[crse_lev+1];
1134
1135 average_down_faces(GetArrOfConstPtrs(faces_fine), faces_crse,
1136 refRatio(crse_lev),geom[crse_lev]);
1137}
1138
1139/**
1140 * @param[in ] lev level at which to get time
1141 */
1142amrex::Real REMORA::get_t_old(int lev) const
1143{
1144 return t_old[lev];
1145}
amrex::Vector< std::string > BCNames
Definition REMORA.cpp:53
PlotfileType
plotfile format
#define Temp_comp
#define Scalar_comp
#define NC3D
#define Salt_comp
#define NC2D
#define NCH2D
#define NCONS
std::unique_ptr< ProblemBase > amrex_probinit(const amrex_real *problo, const amrex_real *probhi) AMREX_ATTRIBUTE_WEAK
Function to init the physical bounds of the domain and instantiate a Problem derived from ProblemBase...
A class to hold and interpolate time series data read from a NetCDF file.
void Initialize()
Read in time array from file and allocate data arrays.
void update_interpolated_to_time(amrex::Real time)
Calculate interpolated values at time, reading in data as necessary.
static PlotfileType plotfile_type
Native or NetCDF plotfile output.
Definition REMORA.H:1275
NCTimeSeriesRiver * river_source_transportbar
Data container for vertically integrated momentum transport in rivers.
Definition REMORA.H:1049
std::string riv_time_varname
Name of time field for river time.
Definition REMORA.H:1308
amrex::Vector< std::unique_ptr< amrex::MultiFab > > vec_rv2d
v velocity RHS (2D, includes horizontal and vertical advection)
Definition REMORA.H:247
static amrex::Real fixed_dt
User specified fixed baroclinic time step.
Definition REMORA.H:1190
amrex::Real last_plot_file_time
Simulation time when we last output a plotfile.
Definition REMORA.H:1161
amrex::Vector< NCTimeSeriesRiver * > river_source_cons
Vector of data containers for scalar data in rivers.
Definition REMORA.H:1045
void init_bathymetry_from_netcdf(int lev)
Bathymetry data initialization from NetCDF file.
void init_bcs()
Read in boundary parameters from input file and set up data structures.
void calculate_nodal_masks(int lev)
Calculate u-, v-, and psi-point masks based on rho-point masks after analytic initialization.
static amrex::Real previousCPUTimeUsed
Accumulator variable for CPU time used thusfar.
Definition REMORA.H:1375
amrex::Vector< std::unique_ptr< amrex::MultiFab > > vec_fcor
coriolis factor (2D)
Definition REMORA.H:363
void init_gls_vmix(int lev, SolverChoice solver_choice)
Initialize GLS variables.
amrex::Vector< std::unique_ptr< amrex::MultiFab > > vec_hOfTheConfusingName
Bathymetry data (2D, positive valued, h in ROMS)
Definition REMORA.H:232
amrex::Vector< REMORAFillPatcher > FPr_v
Vector over levels of FillPatchers for v (3D)
Definition REMORA.H:1071
amrex::Vector< amrex::MultiFab * > cons_new
multilevel data container for current step's scalar data: temperature, salinity, passive scalar
Definition REMORA.H:220
static bool write_history_file
Whether to output NetCDF files as a single history file with several time steps.
Definition REMORA.H:1020
void stretch_transform(int lev)
Calculate vertical stretched coordinates.
amrex::Vector< std::unique_ptr< amrex::MultiFab > > vec_vwind
Wind in the v direction, defined at rho-points.
Definition REMORA.H:292
std::unique_ptr< ProblemBase > prob
Pointer to container of analytical functions for problem definition.
Definition REMORA.H:1100
amrex::Vector< std::unique_ptr< amrex::MultiFab > > vec_mskr
land/sea mask at cell centers (2D)
Definition REMORA.H:349
void Construct_REMORAFillPatchers(int lev)
Construct FillPatchers.
Definition REMORA.cpp:364
static int sum_interval
Diagnostic sum output interval in number of steps.
Definition REMORA.H:1267
int history_count
Counter for which time index we are writing to in the netcdf history file.
Definition REMORA.H:1226
amrex::Real stop_time
Time to stop.
Definition REMORA.H:1176
int do_substep
Whether to substep fine levels in time.
Definition REMORA.H:1199
NCTimeSeries * vbar_clim_data_from_file
Data container for vbar climatology data read from file.
Definition REMORA.H:1034
void Evolve()
Advance solution to final time.
Definition REMORA.cpp:140
NCTimeSeries * Uwind_data_from_file
Data container for u-direction wind read from file.
Definition REMORA.H:1027
std::string bdry_time_varname
Name of time field for boundary data.
Definition REMORA.H:1294
void ReadCheckpointFile()
read checkpoint file from disk
bool chunk_history_file
Whether to chunk netcdf history file.
Definition REMORA.H:1219
amrex::Vector< std::unique_ptr< amrex::MultiFab > > vec_sustr
Surface stress in the u direction.
Definition REMORA.H:285
amrex::Real get_t_old(int lev) const
Accessor method for t_old to expose to outside classes.
Definition REMORA.cpp:1142
NCTimeSeries * v_clim_data_from_file
Data container for v-velocity climatology data read from file.
Definition REMORA.H:1038
amrex::Vector< std::unique_ptr< amrex::MultiFab > > vec_ru2d
u velocity RHS (2D, includes horizontal and vertical advection)
Definition REMORA.H:245
amrex::Vector< std::string > datalogname
Definition REMORA.H:1407
amrex::Vector< amrex::MultiFab * > zvel_new
multilevel data container for current step's z velocities (largely unused; W stored separately)
Definition REMORA.H:226
void init_only(int lev, amrex::Real time)
Init (NOT restart or regrid)
Definition REMORA.cpp:718
void init_set_vmix(int lev)
Initialize vertical mixing coefficients from file or analytic.
Definition REMORA.cpp:600
std::string clim_u_time_varname
Name of time field for u climatology data.
Definition REMORA.H:1300
void set_grid_scale(int lev)
Set pm and pn arrays on level lev. Only works if using analytic initialization.
void set_coriolis(int lev)
Initialize Coriolis factor from file or analytic.
Definition REMORA.cpp:577
amrex::Vector< REMORAFillPatcher > FPr_u
Vector over levels of FillPatchers for u (3D)
Definition REMORA.H:1069
std::string clim_temp_time_varname
Name of time field for temperature climatology data.
Definition REMORA.H:1306
amrex::Vector< amrex::Vector< amrex::Box > > boxes_at_level
the boxes specified at each level by tagging criteria
Definition REMORA.H:1107
static amrex::Vector< amrex::AMRErrorTag > ref_tags
Holds info for dynamically generated tagging criteria.
Definition REMORA.H:1317
std::string nc_frc_file
NetCDF forcing file.
Definition REMORA.H:1285
amrex::Vector< std::unique_ptr< amrex::MultiFab > > vec_Akt
Vertical diffusion coefficient (3D)
Definition REMORA.H:255
std::string clim_ubar_time_varname
Name of time field for ubar climatology data.
Definition REMORA.H:1296
void WriteNCPlotFile(int istep)
Write plotfile using NetCDF (wrapper)
std::string check_file
Checkpoint file prefix.
Definition REMORA.H:1212
static amrex::Real startCPUTime
Variable for CPU timing.
Definition REMORA.H:1373
void setPlotVariables(const std::string &pp_plot_var_names)
NCTimeSeries * u_clim_data_from_file
Data container for u-velocity climatology data read from file.
Definition REMORA.H:1036
amrex::Vector< amrex::MultiFab * > xvel_old
multilevel data container for last step's x velocities (u in ROMS)
Definition REMORA.H:213
amrex::Real start_time
Time of the start of the simulation, in seconds.
Definition REMORA.H:1179
NCTimeSeries * ubar_clim_data_from_file
Data container for ubar climatology data read from file.
Definition REMORA.H:1032
void init_data_from_netcdf(int lev)
Problem initialization from NetCDF file.
void init_masks_from_netcdf(int lev)
Mask data initialization from NetCDF file.
amrex::Vector< amrex::MultiFab * > yvel_new
multilevel data container for current step's y velocities (v in ROMS)
Definition REMORA.H:224
amrex::Vector< std::unique_ptr< amrex::MultiFab > > vec_uwind
Wind in the u direction, defined at rho-points.
Definition REMORA.H:290
static amrex::Real fixed_fast_dt
User specified fixed barotropic time step.
Definition REMORA.H:1192
int regrid_int
how often each level regrids the higher levels of refinement (after a level advances that many time s...
Definition REMORA.H:1202
amrex::Real check_int_time
Checkpoint output interval in seconds.
Definition REMORA.H:1216
void init_bdry_from_netcdf()
Boundary data initialization from NetCDF file.
NCTimeSeries * sustr_data_from_file
Data container for u-component surface momentum flux read from file.
Definition REMORA.H:1023
void Define_REMORAFillPatchers(int lev)
Define FillPatchers.
Definition REMORA.cpp:412
amrex::Vector< std::unique_ptr< amrex::MultiFab > > vec_visc2_p
Harmonic viscosity defined on the psi points (corners of horizontal grid cells)
Definition REMORA.H:257
amrex::Real plot_int_time
Plotfile output interval in seconds.
Definition REMORA.H:1210
amrex::Vector< int > num_files_at_level
how many netcdf input files specified at each level
Definition REMORA.H:1105
amrex::Vector< REMORAFillPatcher > FPr_vbar
Vector over levels of FillPatchers for vbar (2D)
Definition REMORA.H:1079
void AverageDownTo(int crse_lev)
more flexible version of AverageDown() that lets you average down across multiple levels
Definition REMORA.cpp:1118
int steps_per_history_file
Number of time steps per netcdf history file.
Definition REMORA.H:1224
void post_timestep(int nstep, amrex::Real time, amrex::Real dt_lev)
Called after every level 0 timestep.
Definition REMORA.cpp:226
int max_step
maximum number of steps
Definition REMORA.H:1174
amrex::Vector< amrex::MultiFab * > zvel_old
multilevel data container for last step's z velocities (largely unused; W stored separately)
Definition REMORA.H:217
amrex::Vector< int > num_boxes_at_level
how many boxes specified at each level by tagging criteria
Definition REMORA.H:1103
amrex::Vector< amrex::MultiFab * > xvel_new
multilevel data container for current step's x velocities (u in ROMS)
Definition REMORA.H:222
void refinement_criteria_setup()
Set refinement criteria.
int last_check_file_step
Step when we last output a checkpoint file.
Definition REMORA.H:1164
void init_beta_plane_coriolis(int lev)
Calculate Coriolis parameters from beta plane parametrization.
std::string clim_vbar_time_varname
Name of time field for vbar climatology data.
Definition REMORA.H:1298
amrex::Vector< int > nsubsteps
How many substeps on each level?
Definition REMORA.H:1112
amrex::Vector< std::unique_ptr< REMORAPhysBCFunct > > physbcs
Vector (over level) of functors to apply physical boundary conditions.
Definition REMORA.H:1124
void ComputeDt()
a wrapper for estTimeStep()
int plot_int
Plotfile output interval in iterations.
Definition REMORA.H:1208
void FillPatch(int lev, amrex::Real time, amrex::MultiFab &mf_to_be_filled, amrex::Vector< amrex::MultiFab * > const &mfs, const int bccomp, const int bdy_var_type=BdyVars::null, const int icomp=0, const bool fill_all=true, const bool fill_set=true, const int n_not_fill=0, const int icomp_calc=0, const amrex::Real dt=amrex::Real(0.0), const amrex::MultiFab &mf_calc=amrex::MultiFab())
Fill a new MultiFab by copying in phi from valid region and filling ghost cells.
void InitData()
Initialize multilevel data.
Definition REMORA.cpp:258
amrex::Vector< int > istep
which step?
Definition REMORA.H:1110
void WriteCheckpointFile()
write checkpoint file to disk
std::string nc_clim_coeff_file
NetCDF climatology coefficient file.
Definition REMORA.H:1291
void setRecordDataInfo(int i, const std::string &filename)
Definition REMORA.H:1393
void set_analytic_vmix(int lev)
Set vertical mixing coefficients from analytic.
Definition REMORA.cpp:616
void init_flat_bathymetry(int lev)
Initialize flat bathymetry to value from problo.
Definition REMORA.cpp:654
static int file_min_digits
Minimum number of digits in plotfile name or chunked history file.
Definition REMORA.H:1272
void init_riv_pos_from_netcdf(int lev)
static amrex::Vector< std::string > nc_bdry_file
NetCDF boundary data.
Definition REMORA.H:48
amrex::Vector< std::unique_ptr< amrex::MultiFab > > vec_visc2_r
Harmonic viscosity defined on the rho points (centers)
Definition REMORA.H:259
std::string nc_riv_file
Definition REMORA.H:1286
amrex::Vector< std::unique_ptr< amrex::MultiFab > > vec_svstr
Surface stress in the v direction.
Definition REMORA.H:287
REMORA()
Definition REMORA.cpp:61
void set_zeta(int lev)
Initialize zeta from file or analytic.
Definition REMORA.cpp:470
static amrex::Real change_max
Fraction maximum change in subsequent time steps.
Definition REMORA.H:1188
void init_zeta_from_netcdf(int lev)
Sea-surface height data initialization from NetCDF file.
void set_zeta_average(int lev)
Set Zt_avg1 to zeta.
void init_coriolis_from_netcdf(int lev)
Coriolis parameter data initialization from NetCDF file.
std::string pp_prefix
default prefix for input file parameters
Definition REMORA.H:208
void set_bathymetry(int lev)
Initialize bathymetry from file or analytic.
Definition REMORA.cpp:492
amrex::Vector< amrex::MultiFab * > yvel_old
multilevel data container for last step's y velocities (v in ROMS)
Definition REMORA.H:215
amrex::Vector< REMORAFillPatcher > FPr_c
Vector over levels of FillPatchers for scalars.
Definition REMORA.H:1067
std::string clim_v_time_varname
Name of time field for v climatology data.
Definition REMORA.H:1302
amrex::Vector< REMORAFillPatcher > FPr_w
Vector over levels of FillPatchers for w.
Definition REMORA.H:1073
amrex::Vector< amrex::YAFluxRegister * > advflux_reg
array of flux registers for refluxing in multilevel
Definition REMORA.H:1127
std::string clim_salt_time_varname
Name of time field for salinity climatology data.
Definition REMORA.H:1304
NCTimeSeries * Vwind_data_from_file
Data container for v-direction wind read from file.
Definition REMORA.H:1029
amrex::Vector< amrex::Real > t_new
new time at each level
Definition REMORA.H:1114
NCTimeSeriesRiver * river_source_transport
Data container for momentum transport in rivers.
Definition REMORA.H:1047
NCTimeSeries * temp_clim_data_from_file
Data container for temperature climatology data read from file.
Definition REMORA.H:1040
static SolverChoice solverChoice
Container for algorithmic choices.
Definition REMORA.H:1237
void set_masks(int lev)
Initialize land-sea masks from file or analytic.
Definition REMORA.cpp:699
int cf_set_width
Width for fixing values at coarse-fine interface.
Definition REMORA.H:1064
NCTimeSeries * salt_clim_data_from_file
Data container for salinity climatology data read from file.
Definition REMORA.H:1042
void ReadParameters()
read in some parameters from inputs file
Definition REMORA.cpp:896
amrex::Vector< std::unique_ptr< amrex::MultiFab > > vec_ru
u velocity RHS (3D, includes horizontal and vertical advection)
Definition REMORA.H:241
NCTimeSeries * svstr_data_from_file
Data container for v-component surface momentum flux read from file.
Definition REMORA.H:1025
void sum_integrated_quantities(amrex::Real time)
Integrate conserved quantities for diagnostics.
static int total_nc_plot_file_step
Definition REMORA.H:926
static amrex::Vector< amrex::Vector< std::string > > nc_grid_file
NetCDF grid file.
Definition REMORA.H:50
amrex::Vector< std::unique_ptr< amrex::MultiFab > > vec_zeta
free surface height (2D)
Definition REMORA.H:346
amrex::Vector< std::unique_ptr< amrex::MultiFab > > vec_vbar
barotropic y velocity (2D)
Definition REMORA.H:344
void set_zeta_to_Ztavg(int lev)
Set zeta components to be equal to time-averaged Zt_avg1.
int plot_file_on_restart
Whether to output a plotfile on restart from checkpoint.
Definition REMORA.H:1168
void set_2darrays(int lev)
Set 2D momentum arrays from 3D momentum.
void init_analytic(int lev)
Initialize initial problem data from analytic functions.
amrex::Vector< std::unique_ptr< amrex::MultiFab > > vec_ubar
barotropic x velocity (2D)
Definition REMORA.H:342
amrex::Vector< amrex::MultiFab * > cons_old
multilevel data container for last step's scalar data: temperature, salinity, passive scalar
Definition REMORA.H:211
std::string frc_time_varname
Name of time field for forcing data.
Definition REMORA.H:1310
void set_wind(int lev)
Initialize or calculate wind speed from file or analytic.
Definition REMORA.cpp:681
amrex::Vector< REMORAFillPatcher > FPr_ubar
Vector over levels of FillPatchers for ubar (2D)
Definition REMORA.H:1077
bool is_it_time_for_action(int nstep, amrex::Real time, amrex::Real dt, int action_interval, amrex::Real action_per)
Decide if it is time to take an action.
void appendPlotVariables(const std::string &pp_plot_var_names)
void set_hmixcoef(int lev)
Initialize horizontal mixing coefficients.
Definition REMORA.cpp:629
std::string nc_clim_his_file
NetCDF climatology history file.
Definition REMORA.H:1289
void timeStep(int lev, amrex::Real time, int iteration)
advance a level by dt, includes a recursive call for finer levels
void AverageDown()
set covered coarse cells to be the average of overlying fine cells
Definition REMORA.cpp:1106
static int fixed_ndtfast_ratio
User specified, number of barotropic steps per baroclinic step.
Definition REMORA.H:1194
void timeStepML(amrex::Real time, int iteration)
advance all levels by dt, loops over finer levels
amrex::Vector< std::unique_ptr< std::fstream > > datalog
Definition REMORA.H:1406
amrex::Vector< std::unique_ptr< amrex::MultiFab > > vec_Akv
Vertical viscosity coefficient (3D)
Definition REMORA.H:253
static amrex::Real cfl
CFL condition.
Definition REMORA.H:1186
static int verbose
Verbosity level of output.
Definition REMORA.H:1264
void FillCoarsePatch(int lev, amrex::Real time, amrex::MultiFab *mf_fine, amrex::MultiFab *mf_crse, const int bdy_var_type=BdyVars::null, const int icomp=0, const bool fill_all=true, const int n_not_fill=0, const int icomp_calc=0, const amrex::Real dt=amrex::Real(0.0), const amrex::MultiFab &mf_calc=amrex::MultiFab())
fill an entire multifab by interpolating from the coarser level
std::string plot_file_name
Plotfile prefix.
Definition REMORA.H:1206
int check_int
Checkpoint output interval in iterations.
Definition REMORA.H:1214
void WritePlotFile()
main driver for writing AMReX plotfiles
void set_smflux(int lev)
Initialize or calculate surface momentum flux from file or analytic.
Definition REMORA.cpp:663
std::string restart_chkfile
If set, restart from this checkpoint file.
Definition REMORA.H:1182
void init_clim_nudg_coeff(int lev)
Wrapper to initialize climatology nudging coefficient.
amrex::Vector< std::unique_ptr< amrex::MultiFab > > vec_rv
v velocity RHS (3D, includes horizontal and vertical advection)
Definition REMORA.H:243
int cf_width
Nudging width at coarse-fine interface.
Definition REMORA.H:1062
static amrex::Vector< amrex::Vector< std::string > > nc_init_file
NetCDF initialization file.
Definition REMORA.H:49
int last_plot_file_step
Step when we last output a plotfile.
Definition REMORA.H:1159
amrex::Vector< amrex::Real > t_old
old time at each level
Definition REMORA.H:1116
amrex::Real last_check_file_time
Simulation time when we last output a checkpoint file.
Definition REMORA.H:1166
amrex::Vector< amrex::Real > dt
time step at each level
Definition REMORA.H:1118
static amrex::Real sum_per
Diagnostic sum output interval in time.
Definition REMORA.H:1269
virtual ~REMORA()
Definition REMORA.cpp:135
void restart()
Definition REMORA.cpp:458
amrex::Vector< std::unique_ptr< amrex::MultiFab > > vec_diff2
Harmonic diffusivity for temperature / salinity.
Definition REMORA.H:261
const char * buildInfoGetGitHash(int i)
HorizMixingType horiz_mixing_type
amrex::Vector< amrex::Real > tnu2
amrex::Vector< int > do_rivers_cons
void init_params()
read in and initialize parameters
amrex::Real visc2
SMFluxType smflux_type
VertMixingType vert_mixing_type
CouplingType coupling_type