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// Do we include staggered velocities in the plotfile?
39
40// Native AMReX vs NetCDF
42
43#ifdef REMORA_USE_NETCDF
44
46
47// Do we write one file per timestep (false) or one file for all timesteps (true)
49
50// NetCDF initialization file
51amrex::Vector<std::string> REMORA::nc_bdry_file = {""}; // Must provide via input
52amrex::Vector<amrex::Vector<std::string>> REMORA::nc_init_file = {{""}}; // Must provide via input
53amrex::Vector<amrex::Vector<std::string>> REMORA::nc_grid_file = {{""}}; // Must provide via input
54#endif
55
56amrex::Vector<std::string> BCNames = {"xlo", "ylo", "zlo", "xhi", "yhi", "zhi"};
57
58/**
59 * constructor:
60 * - reads in parameters from inputs file
61 * - sizes multilevel arrays and data structures
62 * - initializes BCRec boundary condition object
63 */
65{
66 BL_PROFILE("REMORA::REMORA()");
67 if (ParallelDescriptor::IOProcessor()) {
68 const char* remora_hash = amrex::buildInfoGetGitHash(1);
69 const char* amrex_hash = amrex::buildInfoGetGitHash(2);
70 const char* buildgithash = amrex::buildInfoGetBuildGitHash();
71 const char* buildgitname = amrex::buildInfoGetBuildGitName();
72
73 if (strlen(remora_hash) > 0) {
74 amrex::Print() << "\n"
75 << "REMORA git hash: " << remora_hash << "\n";
76 }
77 if (strlen(amrex_hash) > 0) {
78 amrex::Print() << "AMReX git hash: " << amrex_hash << "\n";
79 }
80 if (strlen(buildgithash) > 0) {
81 amrex::Print() << buildgitname << " git hash: " << buildgithash << "\n";
82 }
83
84 amrex::Print() << "\n";
85 }
86
88 const std::string& pv1 = "plot_vars"; setPlotVariables(pv1);
89
90 prob = amrex_probinit(geom[0].ProbLo(),geom[0].ProbHi());
91
92 // Geometry on all levels has been defined already.
93
94 // No valid BoxArray and DistributionMapping have been defined.
95 // But the arrays for them have been resized.
96
97 int nlevs_max = max_level + 1;
98
99 istep.resize(nlevs_max, 0);
100 nsubsteps.resize(nlevs_max, 1);
101 for (int lev = 1; lev <= max_level; ++lev) {
102 nsubsteps[lev] = do_substep ? MaxRefRatio(lev-1) : 1;
103 }
104
105 physbcs.resize(nlevs_max);
106
107 t_new.resize(nlevs_max, 0.0_rt);
108 t_old.resize(nlevs_max, -1.e100_rt);
109 dt.resize(nlevs_max, 1.e100_rt);
110
111 cons_new.resize(nlevs_max);
112 cons_old.resize(nlevs_max);
113 xvel_new.resize(nlevs_max);
114 xvel_old.resize(nlevs_max);
115 yvel_new.resize(nlevs_max);
116 yvel_old.resize(nlevs_max);
117 zvel_new.resize(nlevs_max);
118 zvel_old.resize(nlevs_max);
119
120 advflux_reg.resize(nlevs_max);
121
122 // Initialize tagging criteria for mesh refinement
124
125 // We have already read in the ref_Ratio (via amr.ref_ratio =) but we need to enforce
126 // that there is no refinement in the vertical so we test on that here.
127 for (int lev = 0; lev < max_level; ++lev)
128 {
129 amrex::Print() << "Refinement ratio at level " << lev << " set to be " <<
130 ref_ratio[lev][0] << " " << ref_ratio[lev][1] << " " << ref_ratio[lev][2] << std::endl;
131
132 if (ref_ratio[lev][2] != 1)
133 {
134 amrex::Error("We don't allow refinement in the vertical -- make sure to set ref_ratio = 1 in z");
135 }
136 }
137}
138
140{
141}
142
143void
145{
146 BL_PROFILE("REMORA::Evolve()");
147 Real cur_time = t_new[0];
148
149 // Take one coarse timestep by calling timeStep -- which recursively calls timeStep
150 // for finer levels (with or without subcycling)
151 for (int step = istep[0]; step < max_step && cur_time < stop_time; ++step)
152 {
153 amrex::Print() << "\nCoarse STEP " << step+1 << " starts ..." << std::endl;
154
155 ComputeDt();
156
157 int lev = 0;
158 int iteration = 1;
159
160 if (max_level == 0) {
161 timeStep(lev, cur_time, iteration);
162 }
163 else {
164 timeStepML(cur_time, iteration);
165 }
166
167 cur_time += dt[0];
168
169 amrex::Print() << "Coarse STEP " << step+1 << " ends." << " TIME = " << cur_time
170 << " DT = " << dt[0] << std::endl;
171
172 if ((plot_int > 0 && (step+1 - last_plot_file_step) == plot_int)
173 || (plot_int_time > 0 && cur_time >= (last_plot_file_time + plot_int_time))) {
174 last_plot_file_step = step+1;
175 last_plot_file_time = cur_time;
177
179 }
180#ifdef REMORA_USE_NETCDF
182 WriteNCPlotFile(step+1);
184 }
185#endif
186 }
187
188 if ((check_int > 0 && (step+1 - last_check_file_step) == check_int)
189 || (check_int_time > 0 && cur_time >= (last_check_file_time + check_int_time))) {
190 last_check_file_step = step+1;
191 last_check_file_time = cur_time;
193 }
194
195 post_timestep(step, cur_time, dt[0]);
196
197#ifdef AMREX_MEM_PROFILING
198 {
199 std::ostringstream ss;
200 ss << "[STEP " << step+1 << "]";
201 MemProfiler::report(ss.str());
202 }
203#endif
204
205 if (cur_time >= stop_time - 1.e-6*dt[0]) break;
206 }
207
208 if ((plot_int > 0 || plot_int_time > 0.0) && istep[0] > last_plot_file_step) {
211 }
212#ifdef REMORA_USE_NETCDF
216 }
217#endif
218 }
219
220 if ((check_int > 0 || check_int_time > 0.0) && istep[0] > last_check_file_step) {
222 }
223}
224
225/**
226 * @param[in ] nstep which step we're on
227 * @param[in ] time current time
228 * @param[in ] dt_lev0 time step on level 0
229 */
230void
231REMORA::post_timestep (int nstep, Real time, Real dt_lev0)
232{
233 BL_PROFILE("REMORA::post_timestep()");
234
235#ifdef REMORA_USE_PARTICLES
236 particleData.Redistribute();
237#endif
238
240 {
241 for (int lev = finest_level-1; lev >= 0; lev--)
242 {
243 // This call refluxes from the lev/lev+1 interface onto lev
244 //getAdvFluxReg(lev+1)->Reflux(*cons_new[lev], 0, 0, NCONS);
245
246 // We need to do this before anything else because refluxing changes the
247 // values of coarse cells underneath fine grids with the assumption they'll
248 // be over-written by averaging down
249 //
250 AverageDownTo(lev);
251 }
252 }
253
254 if (is_it_time_for_action(nstep, time, dt_lev0, sum_interval, sum_per)) {
256 }
257}
258
259/**
260 * This is called from main.cpp and handles all initialization, whether from start or restart
261 */
262void
264{
265 BL_PROFILE("REMORA::InitData()");
266 // Initialize the start time for our CPU-time tracker
267 startCPUTime = Real(ParallelDescriptor::second());
268
269 // Map the words in the inputs file to BC types, then translate
270 // those types into what they mean for each variable
271 init_bcs();
272
273 // Init vertical stretching coeffs
275
278 last_plot_file_time = -1.0_rt;
279 last_check_file_time = -1.0_rt;
280
281 if (restart_chkfile == "") {
282 // start simulation from the beginning
283
284 InitFromScratch(start_time);
285
287 AverageDown();
288 }
289
290 } else { // Restart from a checkpoint
291
292 restart();
293
294 }
295#ifdef REMORA_USE_MOAB
296 InitMOABMesh();
297#endif
298 // Initialize flux registers (whether we start from scratch or restart)
300 advflux_reg[0] = nullptr;
301 for (int lev = 1; lev <= finest_level; lev++)
302 {
303 advflux_reg[lev] = new YAFluxRegister(grids[lev], grids[lev-1],
304 dmap[lev], dmap[lev-1],
305 geom[lev], geom[lev-1],
306 ref_ratio[lev-1], lev, NCONS);
307 }
308 }
309
310 // Fill ghost cells/faces
311 for (int lev = 0; lev <= finest_level; ++lev)
312 {
313 if (lev > 0 && cf_width >= 0) {
315 }
316
317 if (restart_chkfile == "") {
318 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]);
319 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]);
320 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]);
321 FillPatch(lev, t_new[lev], *zvel_new[lev], zvel_new, BCVars::zvel_bc, BdyVars::null, 0, true, false);
322
323 // Copy from new into old just in case when initializing from scratch
324 int ngs = cons_new[lev]->nGrow();
325 int ngvel = xvel_new[lev]->nGrow();
326 MultiFab::Copy(*cons_old[lev],*cons_new[lev],0,0,NCONS,ngs);
327 MultiFab::Copy(*xvel_old[lev],*xvel_new[lev],0,0,1,ngvel);
328 MultiFab::Copy(*yvel_old[lev],*yvel_new[lev],0,0,1,ngvel);
329 MultiFab::Copy(*zvel_old[lev],*zvel_new[lev],0,0,1,IntVect(ngvel,ngvel,0));
330 }
331 } // lev
332
333 // Check for additional plotting variables that are available after
334 // particle containers are setup.
335 const std::string& pv1 = "plot_vars"; appendPlotVariables(pv1);
336
337 if (restart_chkfile == "" && (check_int > 0 || check_int_time > 0.0_rt))
338 {
341 }
342
343 if ( (restart_chkfile == "") ||
345 {
346 if (plot_int > 0 || plot_int_time > 0.0)
347 {
350#ifdef REMORA_USE_NETCDF
352 int step0 = 0;
353 WriteNCPlotFile(step0);
355 }
356#endif
358 }
359 }
360
363 }
364
365 ComputeDt();
366
367}
368
369/**
370 * @param[in ] lev level to operate on
371 */
372void
374{
375 BL_PROFILE("REMORA::Construct_REMORAFillPatchers()");
376 amrex::Print() << ":::Construct_REMORAFillPatchers " << lev << std::endl;
377
378 auto& ba_fine = cons_new[lev ]->boxArray();
379 auto& ba_crse = cons_new[lev-1]->boxArray();
380 auto& dm_fine = cons_new[lev ]->DistributionMap();
381 auto& dm_crse = cons_new[lev-1]->DistributionMap();
382
383 BoxList bl2d_fine = ba_fine.boxList();
384 for (auto& b : bl2d_fine) {
385 b.setRange(2,0);
386 }
387 BoxArray ba2d_fine(std::move(bl2d_fine));
388
389 BoxList bl2d_crse = ba_crse.boxList();
390 for (auto& b : bl2d_crse) {
391 b.setRange(2,0);
392 }
393 BoxArray ba2d_crse(std::move(bl2d_crse));
394
395 int ncomp = cons_new[lev]->nComp();
396
397 FPr_c.emplace_back(ba_fine, dm_fine, geom[lev] ,
398 ba_crse, dm_crse, geom[lev-1],
399 -cf_width, -cf_set_width, ncomp, &cell_cons_interp);
400 FPr_u.emplace_back(convert(ba_fine, IntVect(1,0,0)), dm_fine, geom[lev] ,
401 convert(ba_crse, IntVect(1,0,0)), dm_crse, geom[lev-1],
402 -cf_width, -cf_set_width, 1, &face_cons_linear_interp);
403 FPr_v.emplace_back(convert(ba_fine, IntVect(0,1,0)), dm_fine, geom[lev] ,
404 convert(ba_crse, IntVect(0,1,0)), dm_crse, geom[lev-1],
405 -cf_width, -cf_set_width, 1, &face_cons_linear_interp);
406 FPr_w.emplace_back(convert(ba_fine, IntVect(0,0,1)), dm_fine, geom[lev] ,
407 convert(ba_crse, IntVect(0,0,1)), dm_crse, geom[lev-1],
408 -cf_width, -cf_set_width, 1, &face_cons_linear_interp);
409
410 FPr_ubar.emplace_back(convert(ba2d_fine, IntVect(1,0,0)), dm_fine, geom[lev] ,
411 convert(ba2d_crse, IntVect(1,0,0)), dm_crse, geom[lev-1],
412 -cf_width, -cf_set_width, 3, &face_cons_linear_interp);
413 FPr_vbar.emplace_back(convert(ba2d_fine, IntVect(0,1,0)), dm_fine, geom[lev] ,
414 convert(ba2d_crse, IntVect(0,1,0)), dm_crse, geom[lev-1],
415 -cf_width, -cf_set_width, 3, &face_cons_linear_interp);
416}
417
418/**
419 * @param[in ] lev level to operate on
420 */
421void
423{
424 BL_PROFILE("REMORA::Define_REMORAFillPatchers()");
425 amrex::Print() << ":::Define_REMORAFillPatchers " << lev << std::endl;
426
427 auto& ba_fine = cons_new[lev ]->boxArray();
428 auto& ba_crse = cons_new[lev-1]->boxArray();
429 auto& dm_fine = cons_new[lev ]->DistributionMap();
430 auto& dm_crse = cons_new[lev-1]->DistributionMap();
431
432 BoxList bl2d_fine = ba_fine.boxList();
433 for (auto& b : bl2d_fine) {
434 b.setRange(2,0);
435 }
436 BoxArray ba2d_fine(std::move(bl2d_fine));
437
438 BoxList bl2d_crse = ba_crse.boxList();
439 for (auto& b : bl2d_crse) {
440 b.setRange(2,0);
441 }
442 BoxArray ba2d_crse(std::move(bl2d_crse));
443
444
445 int ncomp = cons_new[lev]->nComp();
446
447 FPr_c[lev-1].Define(ba_fine, dm_fine, geom[lev] ,
448 ba_crse, dm_crse, geom[lev-1],
449 -cf_width, -cf_set_width, ncomp, &cell_cons_interp);
450 FPr_u[lev-1].Define(convert(ba_fine, IntVect(1,0,0)), dm_fine, geom[lev] ,
451 convert(ba_crse, IntVect(1,0,0)), dm_crse, geom[lev-1],
452 -cf_width, -cf_set_width, 1, &face_cons_linear_interp);
453 FPr_v[lev-1].Define(convert(ba_fine, IntVect(0,1,0)), dm_fine, geom[lev] ,
454 convert(ba_crse, IntVect(0,1,0)), dm_crse, geom[lev-1],
455 -cf_width, -cf_set_width, 1, &face_cons_linear_interp);
456 FPr_w[lev-1].Define(convert(ba_fine, IntVect(0,0,1)), dm_fine, geom[lev] ,
457 convert(ba_crse, IntVect(0,0,1)), dm_crse, geom[lev-1],
458 -cf_width, -cf_set_width, 1, &face_cons_linear_interp);
459
460 FPr_ubar[lev-1].Define(convert(ba2d_fine, IntVect(1,0,0)), dm_fine, geom[lev] ,
461 convert(ba2d_crse, IntVect(1,0,0)), dm_crse, geom[lev-1],
462 -cf_width, -cf_set_width, 3, &face_cons_linear_interp);
463 FPr_vbar[lev-1].Define(convert(ba2d_fine, IntVect(0,1,0)), dm_fine, geom[lev] ,
464 convert(ba2d_crse, IntVect(0,1,0)), dm_crse, geom[lev-1],
465 -cf_width, -cf_set_width, 3, &face_cons_linear_interp);
466}
467
468void
470{
471 BL_PROFILE("REMORA::restart()");
473
474 // We set this here so that we don't over-write the checkpoint file we just started from
476}
477
478/**
479 * @param[in ] lev level to operate on
480 */
481void
483{
484 BL_PROFILE("REMORA::set_zeta()");
486 prob->init_analytic_zeta(lev, geom[lev], solverChoice, *this, *vec_zeta[lev]);
487
488#ifdef REMORA_USE_NETCDF
489 } else if (solverChoice.ic_type == IC_Type::netcdf) {
490 amrex::Print() << "Calling init_zeta_from_netcdf on level " << lev << std::endl;
492 amrex::Print() << "Sea surface height loaded from netcdf file \n " << std::endl;
493#endif
494 } else {
495 Abort("Don't know this ic_type!");
496 }
497 vec_zeta[lev]->FillBoundary(geom[lev].periodicity());
498 set_zeta_average(lev);
499}
500
501/**
502 * @param[in ] lev level to operate on
503 */
504void
506{
507 BL_PROFILE("REMORA::bathymetry()");
508 // Only set bathymetry on level 0, and interpolate for finer levels
510 if (lev==0) {
513 prob->init_analytic_bathymetry(lev, geom[lev], solverChoice, *this, *vec_hOfTheConfusingName[lev]);
514 } else {
516 }
517#ifdef REMORA_USE_NETCDF
518 } else if (solverChoice.ic_type == IC_Type::netcdf) {
520 amrex::Print() << "Calling init_bathymetry_from_netcdf " << std::endl;
522 amrex::Print() << "Bathymetry loaded from netcdf file \n " << std::endl;
523 } else {
525 }
526#endif
527 } else {
528 Abort("Don't know this ic_type!");
529 }
530 // Need FillBoundary to fill at grid-grid boundaries, and EnforcePeriodicity
531 // to make sure ghost cells in the domain corners are consistent.
532 vec_hOfTheConfusingName[lev]->FillBoundary(geom[lev].periodicity());
533 vec_hOfTheConfusingName[lev]->EnforcePeriodicity(geom[lev].periodicity());
534 } else {
535 Real dummy_time = 0.0_rt;
537 }
538 } else if (solverChoice.init_ana_h) {
541 prob->init_analytic_bathymetry(lev, geom[lev], solverChoice, *this,*vec_hOfTheConfusingName[lev]);
542 } else {
544 }
545#ifdef REMORA_USE_NETCDF
546 } else if (solverChoice.ic_type == IC_Type::netcdf) {
547 amrex::Print() << "Calling init_bathymetry_from_netcdf level " << lev << std::endl;
549 amrex::Print() << "Bathymetry loaded from netcdf file \n " << std::endl;
550#endif
551 } else {
552 Abort("Don't know this ic_type!");
553 }
554 // Need FillBoundary to fill at grid-grid boundaries, and EnforcePeriodicity
555 // to make sure ghost cells in the domain corners are consistent.
556 vec_hOfTheConfusingName[lev]->FillBoundary(geom[lev].periodicity());
557 vec_hOfTheConfusingName[lev]->EnforcePeriodicity(geom[lev].periodicity());
558 } else if (solverChoice.init_l1ad_h) {
561 prob->init_analytic_bathymetry(lev, geom[lev], solverChoice, *this, *vec_hOfTheConfusingName[lev]);
562 } else {
564 }
565#ifdef REMORA_USE_NETCDF
566 } else if (solverChoice.ic_type == IC_Type::netcdf) {
568 amrex::Print() << "Calling init_bathymetry_from_netcdf " << std::endl;
570 amrex::Print() << "Bathymetry loaded from netcdf file \n " << std::endl;
571 } else {
573 }
574#endif
575 } else {
576 Abort("Don't know this ic_type!");
577 }
578 // Need FillBoundary to fill at grid-grid boundaries, and EnforcePeriodicity
579 // to make sure ghost cells in the domain corners are consistent.
580 vec_hOfTheConfusingName[lev]->FillBoundary(geom[lev].periodicity());
581 vec_hOfTheConfusingName[lev]->EnforcePeriodicity(geom[lev].periodicity());
582 } else {
583 amrex::Abort("Don't know this h init type");
584 }
585}
586
587/**
588 * @param[in ] lev level to operate on
589 */
590void
592 BL_PROFILE("REMORA::set_coriolis()");
595 prob->init_analytic_coriolis(lev, geom[lev], solverChoice, *this, *vec_fcor[lev]);
598#ifdef REMORA_USE_NETCDF
600 amrex::Print() << "Calling init_coriolis_from_netcdf " << std::endl;
602 amrex::Print() << "Coriolis loaded from netcdf file \n" << std::endl;
603#endif
604 } else {
605 Abort("Don't know this coriolis_type!");
606 }
607
608 Real time = 0.0_rt;
609 FillPatch(lev, time, *vec_fcor[lev], GetVecOfPtrs(vec_fcor),BCVars::foextrap_bc);
610 vec_fcor[lev]->EnforcePeriodicity(geom[lev].periodicity());
611 }
612}
613
614void
616 BL_PROFILE("REMORA::init_set_vmix()");
621 // The GLS initialization just sets the multifab to a value, so there's
622 // no need to call FillPatch here
623 } else {
624 Abort("Don't know this vertical mixing type");
625 }
626}
627
628/**
629 * @param[in ] lev level to operate on
630 */
631void
633 BL_PROFILE("REMORA::set_analytic_vmix()");
634 Real time = 0.0_rt;
635 prob->init_analytic_vmix(lev, geom[lev], solverChoice, *this,*vec_Akv[lev], *vec_Akt[lev]);
636 FillPatch(lev, time, *vec_Akv[lev], GetVecOfPtrs(vec_Akv),BCVars::zvel_bc,BdyVars::null,0,true,false);
637 for (int n=0; n<NCONS;n++) {
638 FillPatch(lev, time, *vec_Akt[lev], GetVecOfPtrs(vec_Akt),BCVars::zvel_bc,BdyVars::null,0,false,false);
639 }
640}
641
642/**
643 * @param[in ] lev level to operate on
644 */
645void
647{
648 BL_PROFILE("REMORA::set_hmixcoef()");
650 prob->init_analytic_hmix(lev, geom[lev], solverChoice, *this, *vec_visc2_p[lev], *vec_visc2_r[lev], *vec_diff2[lev]);
652 vec_visc2_p[lev]->setVal(solverChoice.visc2);
653 vec_visc2_r[lev]->setVal(solverChoice.visc2);
654 for (int n=0; n<NCONS; n++) {
655 vec_diff2[lev]->setVal(solverChoice.tnu2[n],n,1);
656 }
657 } else {
658 Abort("Don't know this horizontal mixing type");
659 }
660 Real time = 0.0_rt;
661 FillPatch(lev, time, *vec_visc2_p[lev], GetVecOfPtrs(vec_visc2_p),BCVars::foextrap_periodic_bc);
662 FillPatch(lev, time, *vec_visc2_r[lev], GetVecOfPtrs(vec_visc2_r),BCVars::foextrap_periodic_bc);
663 for (int n=0; n<NCONS; n++) {
664 FillPatch(lev, time, *vec_diff2[lev] , GetVecOfPtrs(vec_diff2),BCVars::foextrap_periodic_bc,BdyVars::null,n,false);
665 }
666}
667
668/**
669 * @param[in ] lev level to operate on
670 */
671void
673{
674 BL_PROFILE("REMORA::init_flat_bathymetry()");
675 vec_hOfTheConfusingName[lev]->setVal(-geom[0].ProbLo()[2]);
676}
677
678/**
679 * @param[in ] lev level to operate on
680 */
681void
683{
684 BL_PROFILE("REMORA::set_smflux()");
686 prob->init_analytic_smflux(lev, geom[lev], solverChoice, *this,*vec_sustr[lev], *vec_svstr[lev]);
688#ifdef REMORA_USE_NETCDF
691 FillPatch(lev, t_old[lev], *vec_sustr[lev], GetVecOfPtrs(vec_sustr),BCVars::foextrap_periodic_bc,BdyVars::null,0,false);
692 FillPatch(lev, t_old[lev], *vec_svstr[lev], GetVecOfPtrs(vec_svstr),BCVars::foextrap_periodic_bc,BdyVars::null,0,false);
693#endif
694 }
695}
696
697/**
698 * @param[in ] lev level to operate on
699 */
700void
702{
703 BL_PROFILE("REMORA::set_wind()");
705 prob->init_analytic_wind(lev,geom[lev], solverChoice, *this, *vec_uwind[lev], *vec_vwind[lev]);
707#ifdef REMORA_USE_NETCDF
710 FillPatch(lev, t_old[lev], *vec_uwind[lev], GetVecOfPtrs(vec_uwind),BCVars::foextrap_periodic_bc,BdyVars::null,0,false);
711 FillPatch(lev, t_old[lev], *vec_vwind[lev], GetVecOfPtrs(vec_vwind),BCVars::foextrap_periodic_bc,BdyVars::null,0,false);
712#endif
713 }
714}
715
716/**
717 * @param[in ] lev level to operate on
718 */
719void
721{
723 prob->init_analytic_masks(lev,geom[lev], solverChoice, *this, *vec_mskr[lev]);
726#ifdef REMORA_USE_NETCDF
727 amrex::Print() << "Calling init_masks_from_netcdf level " << lev << std::endl;
729 amrex::Print() << "Masks loaded from netcdf file \n " << std::endl;
730#endif
731 }
732}
733
734/**
735 * @param[in ] lev level to operate on
736 * @param[in ] time current time for initialization
737 */
738void
739REMORA::init_only (int lev, Real time)
740{
741 BL_PROFILE("REMORA::init_only()");
742 t_new[lev] = time;
743 t_old[lev] = time - 1.e200_rt;
744
745 cons_new[lev]->setVal(0.0_rt);
746 xvel_new[lev]->setVal(0.0_rt);
747 yvel_new[lev]->setVal(0.0_rt);
748 zvel_new[lev]->setVal(0.0_rt);
749
750 xvel_old[lev]->setVal(0.0_rt);
751 yvel_old[lev]->setVal(0.0_rt);
752 zvel_old[lev]->setVal(0.0_rt);
753
754 vec_ru[lev]->setVal(0.0_rt);
755 vec_rv[lev]->setVal(0.0_rt);
756
757 vec_ru2d[lev]->setVal(0.0_rt);
758 vec_rv2d[lev]->setVal(0.0_rt);
759
761 set_grid_scale(lev);
762 }
763 set_masks(lev);
764
765#ifdef REMORA_USE_NETCDF
768
770 if (nc_clim_his_file.empty()) {
771 amrex::Error("NetCDF climatology file name must be provided via input");
772 }
774 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);
775 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);
778 }
780 u_clim_data_from_file = new NCTimeSeries(nc_clim_his_file, "u", clim_u_time_varname, geom[lev].Domain(),xvel_new[lev],false,true);
781 v_clim_data_from_file = new NCTimeSeries(nc_clim_his_file, "v", clim_v_time_varname, geom[lev].Domain(),yvel_new[lev],false,true);
784 }
785 // Since the NCTimeSeries object isn't filling the cons_new MultiFab directly, we don't have to specify a component.
786 // It just needs to know the shape of the MultiFab
788 temp_clim_data_from_file = new NCTimeSeries(nc_clim_his_file, "temp", clim_temp_time_varname,geom[lev].Domain(),cons_new[lev],false,true);
790 }
792 salt_clim_data_from_file = new NCTimeSeries(nc_clim_his_file, "salt", clim_salt_time_varname,geom[lev].Domain(),cons_new[lev],false,true);
794 }
795 }
796 }
797
799 amrex::Print() << "Calling init_bdry_from_netcdf " << std::endl;
801 amrex::Print() << "Boundary data loaded from netcdf file \n " << std::endl;
802 }
803
804 // This will be a non-op if forcings specified analytically
806 if (nc_frc_file.empty()) {
807 amrex::Error("NetCDF forcing file name must be provided via input for winds");
808 }
809 Uwind_data_from_file = new NCTimeSeries(nc_frc_file, "Uwind", frc_time_varname, geom[lev].Domain(),vec_uwind[lev].get(), true, false);
810 Vwind_data_from_file = new NCTimeSeries(nc_frc_file, "Vwind", frc_time_varname, geom[lev].Domain(),vec_vwind[lev].get(), true, false);
814 if (nc_frc_file.empty()) {
815 amrex::Error("NetCDF forcing file name must be provided via input for surface momentum fluxes");
816 }
817 sustr_data_from_file = new NCTimeSeries(nc_frc_file, "sustr", frc_time_varname, geom[lev].Domain(),vec_sustr[lev].get(), true, false);
818 svstr_data_from_file = new NCTimeSeries(nc_frc_file, "svstr", frc_time_varname, geom[lev].Domain(),vec_svstr[lev].get(), true, false);
821 }
822
824 auto dom = geom[0].Domain();
825 int nz = dom.length(2);
826 river_source_cons.resize(NCONS);
827 Print() << solverChoice.do_rivers_cons[0] << std::endl;
830 river_source_cons[Salt_comp]->Initialize();
831 }
834 river_source_cons[Temp_comp]->Initialize();
835 }
838 river_source_cons[Scalar_comp]->Initialize();
839 }
845 }
846#else
848 Abort("Not compiled with NetCDF, but selected boundary conditions require NetCDF");
849 }
850#endif
851
852 set_bathymetry(lev);
853 set_zeta(lev);
855
857 if (lev==0) {
859 {
860 init_analytic(lev);
861#ifdef REMORA_USE_NETCDF
862 } else if (solverChoice.ic_type == IC_Type::netcdf) {
863 amrex::Print() << "Calling init_data_from_netcdf " << std::endl;
866 amrex::Print() << "Initial data loaded from netcdf file \n " << std::endl;
867
868#endif
869 } else {
870 Abort("Need to specify ic_type");
871 }
872 } else {
877 }
880 {
881 init_analytic(lev);
882#ifdef REMORA_USE_NETCDF
883 } else if (solverChoice.ic_type == IC_Type::netcdf) {
884 amrex::Print() << "Calling init_data_from_netcdf on level " << lev << std::endl;
887 amrex::Print() << "Initial data loaded from netcdf file on level " << lev << std::endl;
888
889#endif
890 } else {
891 Abort("Need to specify ic_type");
892 }
893 } else {
894 amrex::Abort("Need to specify T init procedure");
895 }
896
897 // Ensure that the face-based data are the same on both sides of a periodic domain.
898 // The data associated with the lower grid ID is considered the correct value.
899 xvel_new[lev]->OverrideSync(geom[lev].periodicity());
900 yvel_new[lev]->OverrideSync(geom[lev].periodicity());
901 zvel_new[lev]->OverrideSync(geom[lev].periodicity());
902
903 set_2darrays(lev);
904
905 init_set_vmix(lev);
906 set_hmixcoef(lev);
907 set_coriolis(lev);
908
909 // Previously set smflux here with OverrideSync:
910// set_smflux(lev);
911// prob->init_analytic_smflux(lev, geom[lev], solverChoice, *this, *vec_sustr[lev], *vec_svstr[lev]);
912// vec_sustr[lev]->OverrideSync(geom[lev].periodicity());
913// vec_svstr[lev]->OverrideSync(geom[lev].periodicity());
914
915}
916
917void
919{
920 BL_PROFILE("REMORA::ReadParameters()");
921 {
922 ParmParse pp; // Traditionally, max_step and stop_time do not have prefix, so allow it for now.
923 bool noprefix_max_step = pp.query("max_step", max_step);
924 bool noprefix_stop_time = pp.query("stop_time", stop_time);
925 bool remora_max_step = pp.query("remora.max_step", max_step);
926 bool remora_stop_time = pp.query("remora.stop_time", stop_time);
927 if (remora_max_step and noprefix_max_step) {
928 Abort("remora.max_step and max_step are both specified. Please use only one!");
929 }
930 if (remora_stop_time and noprefix_stop_time) {
931 Abort("remora.stop_time and stop_time are both specified. Please use only one!");
932 }
933 }
934
935 ParmParse pp(pp_prefix);
936 ParmParse pp_amr("amr");
937 {
938 pp_amr.query("regrid_int", regrid_int);
939 pp.query("check_file", check_file);
940 pp.query("check_int", check_int);
941 pp_amr.query("check_int", check_int);
942 pp.query("check_int_time", check_int_time);
943 pp_amr.query("check_int_time", check_int_time);
944
945 pp.query("restart", restart_chkfile);
946 pp_amr.query("restart", restart_chkfile);
947 pp.query("start_time",start_time);
948
949 if (pp.contains("data_log"))
950 {
951 int num_datalogs = pp.countval("data_log");
952 datalog.resize(num_datalogs);
953 datalogname.resize(num_datalogs);
954 pp.queryarr("data_log",datalogname,0,num_datalogs);
955 for (int i = 0; i < num_datalogs; i++)
957 }
958
959 // Verbosity
960 pp.query("v", verbose);
961
962 // Frequency of diagnostic output
963 pp.query("sum_interval", sum_interval);
964 pp.query("sum_period" , sum_per);
965 pp.query("file_min_digits", file_min_digits);
966
967 if (file_min_digits < 0) {
968 amrex::Abort("remora.file_min_digits must be non-negative");
969 }
970
971 // Time step controls
972 pp.query("cfl", cfl);
973 pp.query("change_max", change_max);
974
975 pp.query("fixed_dt", fixed_dt);
976 pp.query("fixed_fast_dt", fixed_fast_dt);
977
978 pp.query("fixed_ndtfast_ratio", fixed_ndtfast_ratio);
979
980 // If all three are specified, they must be consistent
981 if (fixed_dt > 0. && fixed_fast_dt > 0. && fixed_ndtfast_ratio > 0)
982 {
984 {
985 amrex::Abort("Dt is over-specfied");
986 }
987 }
988 // If two are specified, initialize fixed_ndtfast_ratio
989 else if (fixed_dt > 0. && fixed_fast_dt > 0. && fixed_ndtfast_ratio <= 0)
990 {
991 fixed_ndtfast_ratio = static_cast<int>(fixed_dt / fixed_fast_dt);
992 }
993
994 AMREX_ASSERT(cfl > 0. || fixed_dt > 0.);
995
996 pp_amr.query("do_substep", do_substep);
997 if (do_substep) {
998 amrex::Abort("Time substepping is not yet implemented. amr.do_substep must be 0");
999 }
1000
1001 // We use this to keep track of how many boxes we read in from WRF initialization
1002 num_files_at_level.resize(max_level+1,0);
1003
1004 // We use this to keep track of how many boxes are specified thru the refinement indicators
1005 num_boxes_at_level.resize(max_level+1,0);
1006 boxes_at_level.resize(max_level+1);
1007
1008 // We always have exactly one file at level 0
1009 num_boxes_at_level[0] = 1;
1010 boxes_at_level[0].resize(1);
1011 boxes_at_level[0][0] = geom[0].Domain();
1012
1013 // Plotfile name and frequency
1014 pp.query("plot_file", plot_file_name);
1015 pp.query("plot_int", plot_int);
1016 pp.query("plot_int_time", plot_int_time);
1017
1018 // Should we plot the staggered face velocities (without averaging to cell centers)
1019 pp.query("plot_staggered_vels", plot_staggered_vels);
1020
1021 // Output format
1022 std::string plotfile_type_str = "amrex";
1023 pp.query("plotfile_type", plotfile_type_str);
1024 if (plotfile_type_str == "amrex") {
1026 } else if (plotfile_type_str == "netcdf" || plotfile_type_str == "NetCDF") {
1028#ifdef REMORA_USE_NETCDF
1029 pp.query("write_history_file",write_history_file);
1030 pp.query("chunk_history_file",chunk_history_file);
1031 pp.query("steps_per_history_file",steps_per_history_file);
1032 // Estimate size of domain for one timestep of netcdf
1033 auto dom = geom[0].Domain();
1034 int nx = dom.length(0) + 2;
1035 int ny = dom.length(1) + 2;
1036 int nz = dom.length(2);
1038 // Estimate number of steps that will fit into a 2GB file.
1039 steps_per_history_file = int((1.6e10 - NCH2D * nx * ny * 64.0_rt)
1040 / (nx * ny * 64.0_rt * (NC3D*nz + NC2D)));
1041 // If we calculate that a single step will exceed 2GB and the user has
1042 // requested automatic history file sizing, warn about a possible impending
1043 // error, and set steps_per_history_file = 1 to attempt output anyway.
1044 if (steps_per_history_file == 0) {
1045 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.");
1047 }
1048 } else if (write_history_file and !chunk_history_file) {
1049 // Estimate number of output steps we'll need
1050 int nt_out = int((max_step) / plot_int) + 1;
1051 Real est_hist_file_size = NCH2D * nx * ny * 64.0_rt + nt_out * nx * ny * 64.0_rt * (NC3D*nz + NC2D);
1052 if (est_hist_file_size > 1.6e10) {
1053 amrex::Warning("WARNING: NetCDF history file may be larger than 2GB limit. Consider setting remora.chunk_history_file=true");
1054 }
1055 }
1057 Print() << "NetCDF history files will have " << steps_per_history_file << " steps per file." << std::endl;
1058 }
1059#endif
1060 } else {
1061 amrex::Print() << "User selected plotfile_type = " << plotfile_type_str << std::endl;
1062 amrex::Abort("Dont know this plotfile_type");
1063 }
1064#ifndef REMORA_USE_NETCDF
1066 {
1067 amrex::Abort("Please compile with NetCDF in order to enable NetCDF plotfiles");
1068 }
1069
1070#endif
1071
1072#ifdef REMORA_USE_NETCDF
1073 nc_init_file.resize(max_level+1);
1074 nc_grid_file.resize(max_level+1);
1075
1076 // NetCDF initialization files -- possibly multiple files at each of multiple levels
1077 // but we always have exactly one file at level 0
1078 for (int lev = 0; lev <= max_level; lev++)
1079 {
1080 const std::string nc_file_names = amrex::Concatenate("nc_init_file_",lev,1);
1081 const std::string nc_bathy_file_names = amrex::Concatenate("nc_grid_file_",lev,1);
1082
1083 if (pp.contains(nc_file_names.c_str()))
1084 {
1085 int num_files = pp.countval(nc_file_names.c_str());
1086 int num_bathy_files = pp.countval(nc_bathy_file_names.c_str());
1087 if (num_files != num_bathy_files) {
1088 amrex::Error("Must have same number of netcdf files for grid info as for solution");
1089 }
1090
1091 num_files_at_level[lev] = num_files;
1092 nc_init_file[lev].resize(num_files);
1093 nc_grid_file[lev].resize(num_files);
1094
1095 pp.queryarr(nc_file_names.c_str() , nc_init_file[lev] ,0,num_files);
1096 pp.queryarr(nc_bathy_file_names.c_str(), nc_grid_file[lev],0,num_files);
1097 }
1098 }
1099 // We only read boundary data at level 0
1100 pp.queryarr("nc_bdry_file", nc_bdry_file);
1101
1102 // Also only read forcings at level 0 (for now)
1103 pp.query("nc_frc_file", nc_frc_file);
1104
1105 // Get river file
1106 pp.query("nc_river_file", nc_riv_file);
1107
1108 // Read in file names for climatology history and nudging weights
1109 pp.query("nc_clim_his_file", nc_clim_his_file);
1110 pp.query("nc_clim_coeff_file", nc_clim_coeff_file);
1111
1112 pp.query("bdy_time_varname",bdry_time_varname);
1113 pp.query("frc_time_varname",frc_time_varname);
1114 pp.query("riv_time_varname",riv_time_varname);
1115 pp.query("clim_ubar_time_varname",clim_ubar_time_varname);
1116 pp.query("clim_vbar_time_varname",clim_vbar_time_varname);
1117 pp.query("clim_u_time_varname",clim_u_time_varname);
1118 pp.query("clim_v_time_varname",clim_v_time_varname);
1119 pp.query("clim_salt_time_varname",clim_salt_time_varname);
1120 pp.query("clim_temp_time_varname",clim_temp_time_varname);
1121
1122#endif
1123
1124#ifdef REMORA_USE_PARTICLES
1125 readTracersParams();
1126#endif
1127 }
1128
1130}
1131
1132void
1134{
1135 BL_PROFILE("REMORA::AverageDown()");
1136 for (int lev = finest_level-1; lev >= 0; --lev)
1137 {
1138 AverageDownTo(lev);
1139 }
1140}
1141
1142/**
1143 * @param[in ] crse_lev level to average down to
1144 */
1145void
1147{
1148 BL_PROFILE("REMORA::AverageDownTo()");
1149 average_down(*cons_new[crse_lev+1], *cons_new[crse_lev],
1150 0, cons_new[crse_lev]->nComp(), refRatio(crse_lev));
1151 average_down(*vec_zeta[crse_lev+1].get(), *vec_zeta[crse_lev].get(),
1152 0, vec_zeta[crse_lev]->nComp(), refRatio(crse_lev));
1153
1154 Array<MultiFab*,AMREX_SPACEDIM> faces_crse;
1155 Array<MultiFab*,AMREX_SPACEDIM> faces_fine;
1156 faces_crse[0] = xvel_new[crse_lev];
1157 faces_crse[1] = yvel_new[crse_lev];
1158 faces_crse[2] = zvel_new[crse_lev];
1159
1160 faces_fine[0] = xvel_new[crse_lev+1];
1161 faces_fine[1] = yvel_new[crse_lev+1];
1162 faces_fine[2] = zvel_new[crse_lev+1];
1163
1164 average_down_faces(GetArrOfConstPtrs(faces_fine), faces_crse,
1165 refRatio(crse_lev),geom[crse_lev]);
1166}
1167
1168/**
1169 * @param[in ] lev level at which to get time
1170 */
1171amrex::Real REMORA::get_t_old(int lev) const
1172{
1173 return t_old[lev];
1174}
amrex::Vector< std::string > BCNames
Definition REMORA.cpp:56
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:1320
NCTimeSeriesRiver * river_source_transportbar
Data container for vertically integrated momentum transport in rivers.
Definition REMORA.H:1091
std::string riv_time_varname
Name of time field for river time.
Definition REMORA.H:1353
amrex::Vector< std::unique_ptr< amrex::MultiFab > > vec_rv2d
v velocity RHS (2D, includes horizontal and vertical advection)
Definition REMORA.H:250
static amrex::Real fixed_dt
User specified fixed baroclinic time step.
Definition REMORA.H:1232
amrex::Real last_plot_file_time
Simulation time when we last output a plotfile.
Definition REMORA.H:1203
amrex::Vector< NCTimeSeriesRiver * > river_source_cons
Vector of data containers for scalar data in rivers.
Definition REMORA.H:1087
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:1420
amrex::Vector< std::unique_ptr< amrex::MultiFab > > vec_fcor
coriolis factor (2D)
Definition REMORA.H:371
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:235
amrex::Vector< REMORAFillPatcher > FPr_v
Vector over levels of FillPatchers for v (3D)
Definition REMORA.H:1113
amrex::Vector< amrex::MultiFab * > cons_new
multilevel data container for current step's scalar data: temperature, salinity, passive scalar
Definition REMORA.H:223
static bool write_history_file
Whether to output NetCDF files as a single history file with several time steps.
Definition REMORA.H:1062
void FillCoarsePatch(int lev, amrex::Real time, amrex::MultiFab *mf_fine, amrex::MultiFab *mf_crse, const int bccomp, 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
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:300
std::unique_ptr< ProblemBase > prob
Pointer to container of analytical functions for problem definition.
Definition REMORA.H:1142
amrex::Vector< std::unique_ptr< amrex::MultiFab > > vec_mskr
land/sea mask at cell centers (2D)
Definition REMORA.H:357
void Construct_REMORAFillPatchers(int lev)
Construct FillPatchers.
Definition REMORA.cpp:373
static int sum_interval
Diagnostic sum output interval in number of steps.
Definition REMORA.H:1309
int history_count
Counter for which time index we are writing to in the netcdf history file.
Definition REMORA.H:1268
amrex::Real stop_time
Time to stop.
Definition REMORA.H:1218
int do_substep
Whether to substep fine levels in time.
Definition REMORA.H:1241
NCTimeSeries * vbar_clim_data_from_file
Data container for vbar climatology data read from file.
Definition REMORA.H:1076
void Evolve()
Advance solution to final time.
Definition REMORA.cpp:144
NCTimeSeries * Uwind_data_from_file
Data container for u-direction wind read from file.
Definition REMORA.H:1069
std::string bdry_time_varname
Name of time field for boundary data.
Definition REMORA.H:1339
void ReadCheckpointFile()
read checkpoint file from disk
bool chunk_history_file
Whether to chunk netcdf history file.
Definition REMORA.H:1261
amrex::Vector< std::unique_ptr< amrex::MultiFab > > vec_sustr
Surface stress in the u direction.
Definition REMORA.H:293
amrex::Real get_t_old(int lev) const
Accessor method for t_old to expose to outside classes.
Definition REMORA.cpp:1171
NCTimeSeries * v_clim_data_from_file
Data container for v-velocity climatology data read from file.
Definition REMORA.H:1080
amrex::Vector< std::unique_ptr< amrex::MultiFab > > vec_ru2d
u velocity RHS (2D, includes horizontal and vertical advection)
Definition REMORA.H:248
amrex::Vector< std::string > datalogname
Definition REMORA.H:1452
amrex::Vector< amrex::MultiFab * > zvel_new
multilevel data container for current step's z velocities (largely unused; W stored separately)
Definition REMORA.H:229
void init_only(int lev, amrex::Real time)
Init (NOT restart or regrid)
Definition REMORA.cpp:739
void init_set_vmix(int lev)
Initialize vertical mixing coefficients from file or analytic.
Definition REMORA.cpp:615
std::string clim_u_time_varname
Name of time field for u climatology data.
Definition REMORA.H:1345
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:591
amrex::Vector< REMORAFillPatcher > FPr_u
Vector over levels of FillPatchers for u (3D)
Definition REMORA.H:1111
std::string clim_temp_time_varname
Name of time field for temperature climatology data.
Definition REMORA.H:1351
amrex::Vector< amrex::Vector< amrex::Box > > boxes_at_level
the boxes specified at each level by tagging criteria
Definition REMORA.H:1149
static amrex::Vector< amrex::AMRErrorTag > ref_tags
Holds info for dynamically generated tagging criteria.
Definition REMORA.H:1362
std::string nc_frc_file
NetCDF forcing file.
Definition REMORA.H:1330
amrex::Vector< std::unique_ptr< amrex::MultiFab > > vec_Akt
Vertical diffusion coefficient (3D)
Definition REMORA.H:258
std::string clim_ubar_time_varname
Name of time field for ubar climatology data.
Definition REMORA.H:1341
void WriteNCPlotFile(int istep)
Write plotfile using NetCDF (wrapper)
std::string check_file
Checkpoint file prefix.
Definition REMORA.H:1254
static amrex::Real startCPUTime
Variable for CPU timing.
Definition REMORA.H:1418
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:1078
amrex::Vector< amrex::MultiFab * > xvel_old
multilevel data container for last step's x velocities (u in ROMS)
Definition REMORA.H:216
amrex::Real start_time
Time of the start of the simulation, in seconds.
Definition REMORA.H:1221
NCTimeSeries * ubar_clim_data_from_file
Data container for ubar climatology data read from file.
Definition REMORA.H:1074
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:227
amrex::Vector< std::unique_ptr< amrex::MultiFab > > vec_uwind
Wind in the u direction, defined at rho-points.
Definition REMORA.H:298
static amrex::Real fixed_fast_dt
User specified fixed barotropic time step.
Definition REMORA.H:1234
int regrid_int
how often each level regrids the higher levels of refinement (after a level advances that many time s...
Definition REMORA.H:1244
amrex::Real check_int_time
Checkpoint output interval in seconds.
Definition REMORA.H:1258
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:1065
void Define_REMORAFillPatchers(int lev)
Define FillPatchers.
Definition REMORA.cpp:422
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:260
amrex::Real plot_int_time
Plotfile output interval in seconds.
Definition REMORA.H:1252
amrex::Vector< int > num_files_at_level
how many netcdf input files specified at each level
Definition REMORA.H:1147
amrex::Vector< REMORAFillPatcher > FPr_vbar
Vector over levels of FillPatchers for vbar (2D)
Definition REMORA.H:1121
void AverageDownTo(int crse_lev)
more flexible version of AverageDown() that lets you average down across multiple levels
Definition REMORA.cpp:1146
int steps_per_history_file
Number of time steps per netcdf history file.
Definition REMORA.H:1266
void post_timestep(int nstep, amrex::Real time, amrex::Real dt_lev)
Called after every level 0 timestep.
Definition REMORA.cpp:231
int max_step
maximum number of steps
Definition REMORA.H:1216
amrex::Vector< amrex::MultiFab * > zvel_old
multilevel data container for last step's z velocities (largely unused; W stored separately)
Definition REMORA.H:220
amrex::Vector< int > num_boxes_at_level
how many boxes specified at each level by tagging criteria
Definition REMORA.H:1145
amrex::Vector< amrex::MultiFab * > xvel_new
multilevel data container for current step's x velocities (u in ROMS)
Definition REMORA.H:225
void refinement_criteria_setup()
Set refinement criteria.
int last_check_file_step
Step when we last output a checkpoint file.
Definition REMORA.H:1206
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:1343
amrex::Vector< int > nsubsteps
How many substeps on each level?
Definition REMORA.H:1154
amrex::Vector< std::unique_ptr< REMORAPhysBCFunct > > physbcs
Vector (over level) of functors to apply physical boundary conditions.
Definition REMORA.H:1166
void ComputeDt()
a wrapper for estTimeStep()
int plot_int
Plotfile output interval in iterations.
Definition REMORA.H:1250
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:263
amrex::Vector< int > istep
which step?
Definition REMORA.H:1152
void WriteCheckpointFile()
write checkpoint file to disk
std::string nc_clim_coeff_file
NetCDF climatology coefficient file.
Definition REMORA.H:1336
void setRecordDataInfo(int i, const std::string &filename)
Definition REMORA.H:1438
void set_analytic_vmix(int lev)
Set vertical mixing coefficients from analytic.
Definition REMORA.cpp:632
void init_flat_bathymetry(int lev)
Initialize flat bathymetry to value from problo.
Definition REMORA.cpp:672
static int file_min_digits
Minimum number of digits in plotfile name or chunked history file.
Definition REMORA.H:1314
void init_riv_pos_from_netcdf(int lev)
static amrex::Vector< std::string > nc_bdry_file
NetCDF boundary data.
Definition REMORA.H:51
amrex::Vector< std::unique_ptr< amrex::MultiFab > > vec_visc2_r
Harmonic viscosity defined on the rho points (centers)
Definition REMORA.H:262
std::string nc_riv_file
Definition REMORA.H:1331
amrex::Vector< std::unique_ptr< amrex::MultiFab > > vec_svstr
Surface stress in the v direction.
Definition REMORA.H:295
REMORA()
Definition REMORA.cpp:64
void set_zeta(int lev)
Initialize zeta from file or analytic.
Definition REMORA.cpp:482
static amrex::Real change_max
Fraction maximum change in subsequent time steps.
Definition REMORA.H:1230
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:211
void set_bathymetry(int lev)
Initialize bathymetry from file or analytic.
Definition REMORA.cpp:505
amrex::Vector< amrex::MultiFab * > yvel_old
multilevel data container for last step's y velocities (v in ROMS)
Definition REMORA.H:218
amrex::Vector< REMORAFillPatcher > FPr_c
Vector over levels of FillPatchers for scalars.
Definition REMORA.H:1109
std::string clim_v_time_varname
Name of time field for v climatology data.
Definition REMORA.H:1347
amrex::Vector< REMORAFillPatcher > FPr_w
Vector over levels of FillPatchers for w.
Definition REMORA.H:1115
amrex::Vector< amrex::YAFluxRegister * > advflux_reg
array of flux registers for refluxing in multilevel
Definition REMORA.H:1169
std::string clim_salt_time_varname
Name of time field for salinity climatology data.
Definition REMORA.H:1349
NCTimeSeries * Vwind_data_from_file
Data container for v-direction wind read from file.
Definition REMORA.H:1071
amrex::Vector< amrex::Real > t_new
new time at each level
Definition REMORA.H:1156
NCTimeSeriesRiver * river_source_transport
Data container for momentum transport in rivers.
Definition REMORA.H:1089
void init_stretch_coeffs()
initialize and calculate stretch coefficients
NCTimeSeries * temp_clim_data_from_file
Data container for temperature climatology data read from file.
Definition REMORA.H:1082
static SolverChoice solverChoice
Container for algorithmic choices.
Definition REMORA.H:1279
void set_masks(int lev)
Initialize land-sea masks from file or analytic.
Definition REMORA.cpp:720
int cf_set_width
Width for fixing values at coarse-fine interface.
Definition REMORA.H:1106
NCTimeSeries * salt_clim_data_from_file
Data container for salinity climatology data read from file.
Definition REMORA.H:1084
void ReadParameters()
read in some parameters from inputs file
Definition REMORA.cpp:918
amrex::Vector< std::unique_ptr< amrex::MultiFab > > vec_ru
u velocity RHS (3D, includes horizontal and vertical advection)
Definition REMORA.H:244
NCTimeSeries * svstr_data_from_file
Data container for v-component surface momentum flux read from file.
Definition REMORA.H:1067
void sum_integrated_quantities(amrex::Real time)
Integrate conserved quantities for diagnostics.
static int total_nc_plot_file_step
Definition REMORA.H:968
static int plot_staggered_vels
Whether to write the staggered velocities (not averaged to cell centers)
Definition REMORA.H:1317
static amrex::Vector< amrex::Vector< std::string > > nc_grid_file
NetCDF grid file.
Definition REMORA.H:53
amrex::Vector< std::unique_ptr< amrex::MultiFab > > vec_zeta
free surface height (2D)
Definition REMORA.H:354
amrex::Vector< std::unique_ptr< amrex::MultiFab > > vec_vbar
barotropic y velocity (2D)
Definition REMORA.H:352
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:1210
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:350
amrex::Vector< amrex::MultiFab * > cons_old
multilevel data container for last step's scalar data: temperature, salinity, passive scalar
Definition REMORA.H:214
std::string frc_time_varname
Name of time field for forcing data.
Definition REMORA.H:1355
void set_wind(int lev)
Initialize or calculate wind speed from file or analytic.
Definition REMORA.cpp:701
amrex::Vector< REMORAFillPatcher > FPr_ubar
Vector over levels of FillPatchers for ubar (2D)
Definition REMORA.H:1119
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:646
std::string nc_clim_his_file
NetCDF climatology history file.
Definition REMORA.H:1334
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:1133
static int fixed_ndtfast_ratio
User specified, number of barotropic steps per baroclinic step.
Definition REMORA.H:1236
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:1451
amrex::Vector< std::unique_ptr< amrex::MultiFab > > vec_Akv
Vertical viscosity coefficient (3D)
Definition REMORA.H:256
static amrex::Real cfl
CFL condition.
Definition REMORA.H:1228
static int verbose
Verbosity level of output.
Definition REMORA.H:1306
std::string plot_file_name
Plotfile prefix.
Definition REMORA.H:1248
int check_int
Checkpoint output interval in iterations.
Definition REMORA.H:1256
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:682
std::string restart_chkfile
If set, restart from this checkpoint file.
Definition REMORA.H:1224
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:246
int cf_width
Nudging width at coarse-fine interface.
Definition REMORA.H:1104
static amrex::Vector< amrex::Vector< std::string > > nc_init_file
NetCDF initialization file.
Definition REMORA.H:52
int last_plot_file_step
Step when we last output a plotfile.
Definition REMORA.H:1201
amrex::Vector< amrex::Real > t_old
old time at each level
Definition REMORA.H:1158
amrex::Real last_check_file_time
Simulation time when we last output a checkpoint file.
Definition REMORA.H:1208
amrex::Vector< amrex::Real > dt
time step at each level
Definition REMORA.H:1160
static amrex::Real sum_per
Diagnostic sum output interval in time.
Definition REMORA.H:1311
virtual ~REMORA()
Definition REMORA.cpp:139
void restart()
Definition REMORA.cpp:469
amrex::Vector< std::unique_ptr< amrex::MultiFab > > vec_diff2
Harmonic diffusivity for temperature / salinity.
Definition REMORA.H:264
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