REMORA
Regional Modeling of Oceans Refined Adaptively
Loading...
Searching...
No Matches
REMORA.H
Go to the documentation of this file.
1#ifndef REMORA_H_
2#define REMORA_H_
3
4#include <string>
5#include <limits>
6#include <memory>
7
8#ifdef _OPENMP
9#include <omp.h>
10#endif
11
12#include <AMReX_AmrCore.H>
13#include <AMReX_BCRec.H>
14#include <AMReX_InterpFaceRegister.H>
15
16#include <AMReX_ParallelDescriptor.H>
17#include <AMReX_ParmParse.H>
18#include <AMReX_MultiFabUtil.H>
19#include <AMReX_FillPatchUtil.H>
20#include <AMReX_VisMF.H>
21#include <AMReX_PhysBCFunct.H>
22#include <AMReX_YAFluxRegister.H>
23#include <AMReX_ErrorList.H>
24
25#ifdef AMREX_MEM_PROFILING
26#include <AMReX_MemProfiler.H>
27#endif
28
29#include <REMORA_Math.H>
30#include <REMORA_PhysBCFunct.H>
31#include <REMORA_FillPatcher.H>
32
33#include <REMORA_IndexDefines.H>
35#include <REMORA_DataStruct.H>
36#include <REMORA_Derive.H>
37#include "REMORA_prob_common.H"
38
39#ifdef REMORA_USE_PARTICLES
40#include "REMORA_ParticleData.H"
41#endif
42
43#ifdef REMORA_USE_NETCDF
44#include "REMORA_NCInterface.H"
45#include "REMORA_NCTimeSeries.H"
47#endif
48
49#ifdef REMORA_USE_MOAB
50#include "REMORA_MOAB.H"
51#endif
53#include <iostream>
54
55#ifdef AMREX_LAZY
56#include <AMReX_Lazy.H>
57#endif
58
59#ifndef AMREX_USE_MPI
60using amrex::MPI_COMM_WORLD;
61using amrex::MPI_Comm;
62#endif
63
75
76// Forward declarations to resolve circular dependency
77class ProblemBase;
78#ifdef REMORA_USE_NETCDF
79class NCTimeSeries;
81#endif
82
83/*!
84 * \brief Class that stores all relevant simulation state data with methods for time stepping
85 */
86class REMORA
87 : public amrex::AmrCore
88{
89public:
90
91 REMORA ();
92 virtual ~REMORA ();
93
94 /** \brief Advance solution to final time */
95 void Evolve ();
96
97 /** \brief Tag cells for refinement */
98 virtual void ErrorEst (int lev, amrex::TagBoxArray& tags, amrex::Real time, int ngrow) override;
99
100 /** \brief Initialize multilevel data */
101 void InitData ();
102
103 /** \brief Init (NOT restart or regrid) */
104 void init_only (int lev, amrex::Real time);
105
106 /** \briefRestart */
107 void restart ();
108
109 /** \brief Called after every level 0 timestep */
110 void post_timestep (int nstep, amrex::Real time, amrex::Real dt_lev);
111
112 // Diagnostics
113
114 /** \brief Integrate conserved quantities for diagnostics */
115 void sum_integrated_quantities (amrex::Real time);
116
117 /** \brief Perform the volume-weighted sum */
118 amrex::Real
119 volWgtSumMF (int lev,
120 const amrex::MultiFab& mf, int comp, bool local, bool finemask);
121
122 /** \brief Decide if it is time to take an action */
123 bool is_it_time_for_action (int nstep, amrex::Real time, amrex::Real dt,
124 int action_interval, amrex::Real action_per);
125
126 /** \brief Make a new level using provided BoxArray and DistributionMapping and
127 * fill with interpolated coarse level data.
128 * Overrides the pure virtual function in AmrCore
129 */
130 virtual void MakeNewLevelFromCoarse (int lev, amrex::Real time, const amrex::BoxArray& ba,
131 const amrex::DistributionMapping& dm) override;
132
133 /** \brief Remake an existing level using provided BoxArray and DistributionMapping and
134 * fill with existing fine and coarse data.
135 * Overrides the pure virtual function in AmrCore
136 */
137 virtual void RemakeLevel (int lev, amrex::Real time, const amrex::BoxArray& ba,
138 const amrex::DistributionMapping& dm) override;
139
140 /** \brief Delete level data
141 * Overrides the pure virtual function in AmrCore
142 */
143 virtual void ClearLevel (int lev) override;
144
145 /** \brief Make a new level from scratch using provided BoxArray and DistributionMapping.
146 * Only used during initialization.
147 * Overrides the pure virtual function in AmrCore
148 */
149 virtual void MakeNewLevelFromScratch (int lev, amrex::Real time, const amrex::BoxArray& ba,
150 const amrex::DistributionMapping& dm) override;
151
152 /** \brief Set pm and pn arrays on level lev. Only works if using analytic initialization */
153 void set_grid_scale (int lev);
154
155 /** \brief Set zeta components to be equal to time-averaged Zt_avg1 */
156 void set_zeta_to_Ztavg (int lev);
157
158 /** \brief Calculate u-, v-, and psi-point masks based on rho-point masks after analytic initialization */
159 void calculate_nodal_masks (int lev);
160
161 /** \brief Set psi-point mask to be consistent with rho-point mask */
162 void update_mskp (int lev);
163
164 /** \brief compute dt from CFL considerations */
165 amrex::Real estTimeStep (int lev) const;
166
167 /** \brief Interface for advancing the data at one level by one "slow" timestep */
168 void remora_advance(int level,
169 amrex::MultiFab& cons_old, amrex::MultiFab& cons_new,
170 amrex::MultiFab& xvel_old, amrex::MultiFab& yvel_old, amrex::MultiFab& zvel_old,
171 amrex::MultiFab& xvel_new, amrex::MultiFab& yvel_new, amrex::MultiFab& zvel_new,
172 amrex::MultiFab& source,
173 const amrex::Geometry fine_geom,
174 const amrex::Real dt, const amrex::Real time);
175
176 /** \brief Make mask to zero out covered cells (for mesh refinement) */
177 amrex::MultiFab& build_fine_mask(int lev);
178
179 /** \brief main driver for writing AMReX plotfiles */
180 void WritePlotFile ();
181
182 /** \brief write out particular data to an AMReX plotfile */
183 void WriteMultiLevelPlotfileWithBathymetry (const std::string &plotfilename,
184 int nlevels,
185 const amrex::Vector<const amrex::MultiFab*> &mf,
186 const amrex::Vector<const amrex::MultiFab*> &mf_nd,
187 const amrex::Vector<const amrex::MultiFab*> &mf_u,
188 const amrex::Vector<const amrex::MultiFab*> &mf_v,
189 const amrex::Vector<const amrex::MultiFab*> &mf_w,
190 const amrex::Vector<std::string> &varnames,
191 amrex::Real time,
192 const amrex::Vector<int> &level_steps,
193 const std::string &versionName = "HyperCLaw-V1.1",
194 const std::string &levelPrefix = "Level_",
195 const std::string &mfPrefix = "Cell",
196 const amrex::Vector<std::string>& extra_dirs = amrex::Vector<std::string>()) const;
197
198
199 /** \brief write out header data for an AMReX plotfile */
200 void WriteGenericPlotfileHeaderWithBathymetry (std::ostream &HeaderFile,
201 int nlevels,
202 const amrex::Vector<amrex::BoxArray> &bArray,
203 const amrex::Vector<std::string> &varnames,
204 amrex::Real time,
205 const amrex::Vector<int> &level_steps,
206 const std::string &versionName,
207 const std::string &levelPrefix,
208 const std::string &mfPrefix) const;
209
210 /** \brief default prefix for input file parameters */
211 std::string pp_prefix {"remora"};
212
213 /** \brief multilevel data container for last step's scalar data: temperature, salinity, passive scalar */
214 amrex::Vector<amrex::MultiFab*> cons_old;
215 /** \brief multilevel data container for last step's x velocities (u in ROMS) */
216 amrex::Vector<amrex::MultiFab*> xvel_old;
217 /** \brief multilevel data container for last step's y velocities (v in ROMS) */
218 amrex::Vector<amrex::MultiFab*> yvel_old;
219 /** \brief multilevel data container for last step's z velocities (largely unused; W stored separately) */
220 amrex::Vector<amrex::MultiFab*> zvel_old;
221
222 /** \brief multilevel data container for current step's scalar data: temperature, salinity, passive scalar */
223 amrex::Vector<amrex::MultiFab*> cons_new;
224 /** \brief multilevel data container for current step's x velocities (u in ROMS) */
225 amrex::Vector<amrex::MultiFab*> xvel_new;
226 /** \brief multilevel data container for current step's y velocities (v in ROMS) */
227 amrex::Vector<amrex::MultiFab*> yvel_new;
228 /** \brief multilevel data container for current step's z velocities (largely unused; W stored separately) */
229 amrex::Vector<amrex::MultiFab*> zvel_new;
230
231 // Program state data is represented by vectors of pointers to AMReX Multifabs.
232 // There is one pointer per level
233
234 /** \brief Bathymetry data (2D, positive valued, h in ROMS) **/
235 amrex::Vector<std::unique_ptr<amrex::MultiFab>> vec_hOfTheConfusingName;
236
237 /** \brief Width of cells in the vertical (z-) direction (3D, Hz in ROMS) */
238 amrex::Vector<std::unique_ptr<amrex::MultiFab>> vec_Hz;
239 /** \brief u-volume flux (3D) */
240 amrex::Vector<std::unique_ptr<amrex::MultiFab>> vec_Huon;
241 /** \brief v-volume flux (3D) */
242 amrex::Vector<std::unique_ptr<amrex::MultiFab>> vec_Hvom;
243 /** \brief u velocity RHS (3D, includes horizontal and vertical advection) */
244 amrex::Vector<std::unique_ptr<amrex::MultiFab>> vec_ru;
245 /** \brief v velocity RHS (3D, includes horizontal and vertical advection) */
246 amrex::Vector<std::unique_ptr<amrex::MultiFab>> vec_rv;
247 /** \brief u velocity RHS (2D, includes horizontal and vertical advection) */
248 amrex::Vector<std::unique_ptr<amrex::MultiFab>> vec_ru2d;
249 /** \brief v velocity RHS (2D, includes horizontal and vertical advection) */
250 amrex::Vector<std::unique_ptr<amrex::MultiFab>> vec_rv2d;
251 /** \brief u velocity RHS, integrated, including advection and bottom/surface stresses (2D) */
252 amrex::Vector<std::unique_ptr<amrex::MultiFab>> vec_rufrc;
253 /** \brief v velocity RHS, integrated, including advection and bottom/surface stresses (2D) */
254 amrex::Vector<std::unique_ptr<amrex::MultiFab>> vec_rvfrc;
255 /** \brief Vertical viscosity coefficient (3D) */
256 amrex::Vector<std::unique_ptr<amrex::MultiFab>> vec_Akv;
257 /** \brief Vertical diffusion coefficient (3D) */
258 amrex::Vector<std::unique_ptr<amrex::MultiFab>> vec_Akt;
259 /** \brief Harmonic viscosity defined on the psi points (corners of horizontal grid cells) */
260 amrex::Vector<std::unique_ptr<amrex::MultiFab>> vec_visc2_p;
261 /** \brief Harmonic viscosity defined on the rho points (centers) */
262 amrex::Vector<std::unique_ptr<amrex::MultiFab>> vec_visc2_r;
263 /** \brief Harmonic diffusivity for temperature / salinity */
264 amrex::Vector<std::unique_ptr<amrex::MultiFab>> vec_diff2;
265
266 /** \brief z coordinates at rho points (cell centers) */
267 amrex::Vector<std::unique_ptr<amrex::MultiFab>> vec_z_r;
268
269 /** \brief z coordinates at w points (faces between z-cells) */
270 amrex::Vector<std::unique_ptr<amrex::MultiFab>> vec_z_w;
271 /** \brief Scaled vertical coordinate (range [0,1]) that transforms to z, defined at rho points (cell centers) */
272 amrex::Gpu::DeviceVector<amrex::Real> s_r;
273 /** \brief Scaled vertical coordinate (range [0,1]) that transforms to z, defined at w-points (cell faces) */
274 amrex::Gpu::DeviceVector<amrex::Real> s_w;
275
276 //amrex::Vector<std::unique_ptr<amrex::MultiFab>> vec_s_r;
277 //amrex::Vector<std::unique_ptr<amrex::MultiFab>> vec_s_w;
278 //amrex::Vector<std::unique_ptr<amrex::MultiFab>> vec_Cs_r;
279 //amrex::Vector<std::unique_ptr<amrex::MultiFab>> vec_Cs_w;
280
281 /** \brief Stretching coefficients at rho points */
282 amrex::Gpu::DeviceVector<amrex::Real> Cs_r;
283 /** \brief Stretching coefficients at w points */
284 amrex::Gpu::DeviceVector<amrex::Real> Cs_w;
285
286 /** \brief z coordinates at psi points (cell nodes) **/
287 amrex::Vector<std::unique_ptr<amrex::MultiFab>> vec_z_phys_nd;
288
289 /** \brief Average of the free surface, zeta (2D) */
290 amrex::Vector<std::unique_ptr<amrex::MultiFab>> vec_Zt_avg1;
291
292 /** \brief Surface stress in the u direction */
293 amrex::Vector<std::unique_ptr<amrex::MultiFab>> vec_sustr;
294 /** \brief Surface stress in the v direction */
295 amrex::Vector<std::unique_ptr<amrex::MultiFab>> vec_svstr;
296
297 /** \brief Wind in the u direction, defined at rho-points */
298 amrex::Vector<std::unique_ptr<amrex::MultiFab>> vec_uwind;
299 /** \brief Wind in the v direction, defined at rho-points */
300 amrex::Vector<std::unique_ptr<amrex::MultiFab>> vec_vwind;
301
302 /** \brief longwave radiation */
303 amrex::Vector<std::unique_ptr<amrex::MultiFab>> vec_lrflx;
304 /** \brief latent heat flux */
305 amrex::Vector<std::unique_ptr<amrex::MultiFab>> vec_lhflx;
306 /** \brief sensible heat flux */
307 amrex::Vector<std::unique_ptr<amrex::MultiFab>> vec_shflx;
308
309 /** \brief Surface tracer flux; working arrays */
310 amrex::Vector<std::unique_ptr<amrex::MultiFab>> vec_stflx;
311 /** \brief Surface tracer flux; input arrays */
312 amrex::Vector<std::unique_ptr<amrex::MultiFab>> vec_stflux;
313 /** \brief Bottom tracer flux; working arrays */
314 amrex::Vector<std::unique_ptr<amrex::MultiFab>> vec_btflx;
315 /** \brief Bottom tracer flux; input arrays */
316 amrex::Vector<std::unique_ptr<amrex::MultiFab>> vec_btflux;
317
318 /** \brief precipitation rate [kg/m^2/s] */
319 amrex::Vector<std::unique_ptr<amrex::MultiFab>> vec_rain;
320 /** \brief evaporation rate [kg/m^2/s] */
321 amrex::Vector<std::unique_ptr<amrex::MultiFab>> vec_evap;
322
323 /** \brief Linear drag coefficient [m/s], defined at rho points */
324 amrex::Vector<std::unique_ptr<amrex::MultiFab>> vec_rdrag;
325 /** \brief Quadratic drag coefficient [unitless], defined at rho points */
326 amrex::Vector<std::unique_ptr<amrex::MultiFab>> vec_rdrag2;
327 /** \brief Bottom roughness length [m], defined at rho points */
328 amrex::Vector<std::unique_ptr<amrex::MultiFab>> vec_ZoBot;
329
330 /** \brief Bottom stress in the u direction */
331 amrex::Vector<std::unique_ptr<amrex::MultiFab>> vec_bustr;
332 /** \brief Bottom stress in the v direction */
333 amrex::Vector<std::unique_ptr<amrex::MultiFab>> vec_bvstr;
334
335 /** \brief time average of barotropic x velocity flux (2D) */
336 amrex::Vector<std::unique_ptr<amrex::MultiFab>> vec_DU_avg1;
337 /** \brief correct time average of barotropic x velocity flux for coupling (2D) */
338 amrex::Vector<std::unique_ptr<amrex::MultiFab>> vec_DU_avg2;
339 /** \brief time average of barotropic y velocity flux */
340 amrex::Vector<std::unique_ptr<amrex::MultiFab>> vec_DV_avg1;
341 /** \brief correct time average of barotropic y velocity flux for coupling (2D) */
342 amrex::Vector<std::unique_ptr<amrex::MultiFab>> vec_DV_avg2;
343 /** \brief barotropic x velocity for the RHS (2D) */
344 amrex::Vector<std::unique_ptr<amrex::MultiFab>> vec_rubar;
345 /** \brief barotropic y velocity for the RHS (2D) */
346 amrex::Vector<std::unique_ptr<amrex::MultiFab>> vec_rvbar;
347 /** \brief free surface height for the RHS (2D) */
348 amrex::Vector<std::unique_ptr<amrex::MultiFab>> vec_rzeta;
349 /** \brief barotropic x velocity (2D) */
350 amrex::Vector<std::unique_ptr<amrex::MultiFab>> vec_ubar;
351 /** \brief barotropic y velocity (2D) */
352 amrex::Vector<std::unique_ptr<amrex::MultiFab>> vec_vbar;
353 /** \brief free surface height (2D) */
354 amrex::Vector<std::unique_ptr<amrex::MultiFab>> vec_zeta;
355
356 /** \brief land/sea mask at cell centers (2D) */
357 amrex::Vector<std::unique_ptr<amrex::MultiFab>> vec_mskr;
358 /** \brief land/sea mask at x-faces (2D) */
359 amrex::Vector<std::unique_ptr<amrex::MultiFab>> vec_msku;
360 /** \brief land/sea mask at y-faces (2D) */
361 amrex::Vector<std::unique_ptr<amrex::MultiFab>> vec_mskv;
362 /** \brief land/sea mask at cell corners (2D) */
363 amrex::Vector<std::unique_ptr<amrex::MultiFab>> vec_mskp;
364
365 /** \brief horizontal scaling factor: 1 / dx (2D) */
366 amrex::Vector<std::unique_ptr<amrex::MultiFab>> vec_pm;
367 /** \brief horizontal scaling factor: 1 / dy (2D) */
368 amrex::Vector<std::unique_ptr<amrex::MultiFab>> vec_pn;
369
370 /** \brief coriolis factor (2D) */
371 amrex::Vector<std::unique_ptr<amrex::MultiFab>> vec_fcor;
372
373 /** \brief x_grid on rho points (2D) */
374 amrex::Vector<std::unique_ptr<amrex::MultiFab>> vec_xr;
375 /** \brief y_grid on rho points (2D) */
376 amrex::Vector<std::unique_ptr<amrex::MultiFab>> vec_yr;
377
378 /** \brief x_grid on u-points (2D) */
379 amrex::Vector<std::unique_ptr<amrex::MultiFab>> vec_xu;
380 /** \brief y_grid on u-points (2D) */
381 amrex::Vector<std::unique_ptr<amrex::MultiFab>> vec_yu;
382
383 /** \brief x_grid on v-points (2D) */
384 amrex::Vector<std::unique_ptr<amrex::MultiFab>> vec_xv;
385 /** \brief y_grid on v-points (2D) */
386 amrex::Vector<std::unique_ptr<amrex::MultiFab>> vec_yv;
387
388 /** \brief x_grid on psi-points (2D) */
389 amrex::Vector<std::unique_ptr<amrex::MultiFab>> vec_xp;
390 /** \brief y_grid on psi-points (2D) */
391 amrex::Vector<std::unique_ptr<amrex::MultiFab>> vec_yp;
392
393 /** \brief additional scratch space for calculations on temp, salt, etc */
394 amrex::Vector<std::unique_ptr<amrex::MultiFab>> vec_sstore;
395
396 /** \brief density perturbation */
397 amrex::Vector<std::unique_ptr<amrex::MultiFab>> vec_rhoS;
398 /** \brief vertically-averaged density */
399 amrex::Vector<std::unique_ptr<amrex::MultiFab>> vec_rhoA;
400 /** \brief Brunt-Vaisala frequency (3D) */
401 amrex::Vector<std::unique_ptr<amrex::MultiFab>> vec_bvf;
402 /** \brief Thermal expansion coefficient (3D) */
403 amrex::Vector<std::unique_ptr<amrex::MultiFab>> vec_alpha;
404 /** \brief Saline contraction coefficient (3D) */
405 amrex::Vector<std::unique_ptr<amrex::MultiFab>> vec_beta;
406
407 /** \brief Weights for calculating avg1 in 2D advance */
408 amrex::Vector<amrex::Real> vec_weight1;
409 /** \brief Weights for calculating avg2 in 2D advance */
410 amrex::Vector<amrex::Real> vec_weight2;
411
412 // GLS
413 /** \brief Turbulent kinetic energy */
414 amrex::Vector<std::unique_ptr<amrex::MultiFab>> vec_tke;
415 /** \brief Turbulent generic length scale */
416 amrex::Vector<std::unique_ptr<amrex::MultiFab>> vec_gls;
417 /** \brief Vertical mixing turbulent length scale */
418 amrex::Vector<std::unique_ptr<amrex::MultiFab>> vec_Lscale;
419 /** \brief Turbulent kinetic energy vertical diffusion coefficient */
420 amrex::Vector<std::unique_ptr<amrex::MultiFab>> vec_Akk;
421 /** \brief Turbulent length scale vertical diffusion coefficient */
422 amrex::Vector<std::unique_ptr<amrex::MultiFab>> vec_Akp;
423
424 /** \brief Climatology nudging coefficients */
425 amrex::Vector<amrex::Vector<std::unique_ptr<amrex::MultiFab>>> vec_nudg_coeff;
426
427 /** \brief advance a single level for a single time step */
428 void Advance (int lev, amrex::Real time, amrex::Real dt_lev, int iteration, int ncycle);
429
430 /** \brief Set everything up for a step on a level */
431 void setup_step(int lev, amrex::Real time, amrex::Real dt_lev);
432
433 /** \brief 3D advance on a single level */
434 void advance_3d_ml (int lev, amrex::Real dt_lev);
435
436 /** \brief 2D advance, one predictor/corrector step */
437 void advance_2d_onestep (int lev, amrex::Real dt_lev, amrex::Real dtfast_lev, int my_iif, int nfast_counter);
438
439 /** \brief Perform a 2D predictor (predictor_2d_step=True) or corrector (predictor_2d_step=False) step */
440 void advance_2d (int lev,
441 amrex::MultiFab const* mf_rhoS,
442 amrex::MultiFab const* mf_rhoA,
443 amrex::MultiFab * mf_ru2d,
444 amrex::MultiFab * mf_rv2d,
445 amrex::MultiFab * mf_rufrc,
446 amrex::MultiFab * mf_rvfrc,
447 amrex::MultiFab * mf_Zt_avg1,
448 std::unique_ptr<amrex::MultiFab>& mf_DU_avg1,
449 std::unique_ptr<amrex::MultiFab>& mf_DU_avg2,
450 std::unique_ptr<amrex::MultiFab>& mf_DV_avg1,
451 std::unique_ptr<amrex::MultiFab>& mf_DV_avg2,
452 std::unique_ptr<amrex::MultiFab>& mf_rubar,
453 std::unique_ptr<amrex::MultiFab>& mf_rvbar,
454 std::unique_ptr<amrex::MultiFab>& mf_rzeta,
455 std::unique_ptr<amrex::MultiFab>& mf_ubar,
456 std::unique_ptr<amrex::MultiFab>& mf_vbar,
457 amrex::MultiFab * mf_zeta,
458 amrex::MultiFab const* mf_h,
459 amrex::MultiFab const* mf_pm,
460 amrex::MultiFab const* mf_pn,
461 amrex::MultiFab const* mf_fcor,
462 amrex::MultiFab const* mf_visc2_p,
463 amrex::MultiFab const* mf_visc2_r,
464 amrex::MultiFab const* mf_mskr,
465 amrex::MultiFab const* mf_msku,
466 amrex::MultiFab const* mf_mskv,
467 amrex::MultiFab const* mf_mskp,
468 amrex::Real dtfast_lev,
469 bool predictor_2d_step,
470 bool first_2d_step, int my_iif,
471 int & next_indx1);
472
473 /** \brief Advance the 3D variables */
474 void advance_3d (int lev,
475 amrex::MultiFab& mf_cons,
476 amrex::MultiFab& mf_u , amrex::MultiFab& mf_v,
477 amrex::MultiFab* mf_sstore,
478 amrex::MultiFab* mf_ru , amrex::MultiFab* mf_rv,
479 std::unique_ptr<amrex::MultiFab>& mf_DU_avg1,
480 std::unique_ptr<amrex::MultiFab>& mf_DU_avg2,
481 std::unique_ptr<amrex::MultiFab>& mf_DV_avg1,
482 std::unique_ptr<amrex::MultiFab>& mf_DV_avg2,
483 std::unique_ptr<amrex::MultiFab>& mf_ubar,
484 std::unique_ptr<amrex::MultiFab>& mf_vbar,
485 std::unique_ptr<amrex::MultiFab>& mf_Akv,
486 std::unique_ptr<amrex::MultiFab>& mf_Akt,
487 std::unique_ptr<amrex::MultiFab>& mf_Hz,
488 std::unique_ptr<amrex::MultiFab>& mf_Huon,
489 std::unique_ptr<amrex::MultiFab>& mf_Hvom,
490 std::unique_ptr<amrex::MultiFab>& mf_z_w,
491 amrex::MultiFab const* mf_h,
492 amrex::MultiFab const* mf_pm,
493 amrex::MultiFab const* mf_pn,
494 amrex::MultiFab const* mf_mskr,
495 amrex::MultiFab const* mf_msku,
496 amrex::MultiFab const* mf_mskv,
497 const int N,
498 const amrex::Real dt_lev);
499
500 /** \brief Calculate bulk temperature, salinity, wind fluxes */
501 void bulk_fluxes (int lev,
502 amrex::MultiFab* mf_cons,
503 amrex::MultiFab* mf_uwind,
504 amrex::MultiFab* mf_vwind,
505 amrex::MultiFab* mf_evap,
506 amrex::MultiFab* mf_sustr,
507 amrex::MultiFab* mf_svstr,
508 amrex::MultiFab* mf_stflux,
509 amrex::MultiFab* mf_lrflx,
510 amrex::MultiFab* mf_lhflx,
511 amrex::MultiFab* mf_shflx,
512 const int N);
513
514 /** \brief Wrapper function for prestep */
515 void prestep (int lev,
516 amrex::MultiFab& mf_uold, amrex::MultiFab& mf_vold,
517 amrex::MultiFab& mf_u, amrex::MultiFab& mf_v,
518 amrex::MultiFab* mf_ru,
519 amrex::MultiFab* mf_rv,
520 amrex::MultiFab& S_old,
521 amrex::MultiFab& S_new,
522 amrex::MultiFab& mf_W, amrex::MultiFab& mf_DC,
523 /* MF mf_FC? */
524 const amrex::MultiFab* mf_z_r,
525 const amrex::MultiFab* mf_z_w,
526 const amrex::MultiFab* mf_h,
527 const amrex::MultiFab* mf_pm,
528 const amrex::MultiFab* mf_pn,
529 const amrex::MultiFab* mf_sustr,
530 const amrex::MultiFab* mf_svstr,
531 const amrex::MultiFab* mf_bustr,
532 const amrex::MultiFab* mf_bvstr,
533 const amrex::MultiFab* mf_msku,
534 const amrex::MultiFab* mf_mskv,
535 const int iic, const int nfirst,
536 const int nnew, int nstp, int nrhs,
537 int N, const amrex::Real dt_lev);
538
539 /** \brief Prestep advection calculations for the tracers */
540 void prestep_t_advection (int lev, const amrex::Box& tbx,
541 const amrex::Box& gbx,
542 const amrex::Array4<amrex::Real >& tempold,
543 const amrex::Array4<amrex::Real >& tempcache,
544 const amrex::Array4<amrex::Real >& Hz,
545 const amrex::Array4<amrex::Real >& Huon,
546 const amrex::Array4<amrex::Real >& Hvom,
547 const amrex::Array4<amrex::Real >& W,
548 const amrex::Array4<amrex::Real >& DC,
549 const amrex::Array4<amrex::Real >& FC,
550 const amrex::Array4<amrex::Real >& sstore,
551 const amrex::Array4<amrex::Real const>& z_w,
552 const amrex::Array4<amrex::Real const>& h,
553 const amrex::Array4<amrex::Real const>& pm,
554 const amrex::Array4<amrex::Real const>& pn,
555 const amrex::Array4<amrex::Real const>& msku,
556 const amrex::Array4<amrex::Real const>& mskv,
557 const amrex::Array4<int const>& river_pos,
558 const amrex::Array4<amrex::Real const>& river_source,
559 int iic, int ntfirst, int nrhs, int N,
560 const amrex::Real dt_lev);
561
562 /** \brief RHS terms for tracer */
563 void rhs_t_3d (int lev,
564 const amrex::Box& bx,
565 const amrex::Array4<amrex::Real >& t,
566 const amrex::Array4<amrex::Real const>& tempstore,
567 const amrex::Array4<amrex::Real const>& Huon,
568 const amrex::Array4<amrex::Real const>& Hvom,
569 const amrex::Array4<amrex::Real const>& Hz,
570 const amrex::Array4<amrex::Real const>& pn,
571 const amrex::Array4<amrex::Real const>& pm,
572 const amrex::Array4<amrex::Real const>& W,
573 const amrex::Array4<amrex::Real >& FC,
574 const amrex::Array4<amrex::Real const>& mskr,
575 const amrex::Array4<amrex::Real const>& msku,
576 const amrex::Array4<amrex::Real const>& mskv,
577 const amrex::Array4<int const>& river_pos,
578 const amrex::Array4<amrex::Real const>& river_source,
579 int nrhs, int nnew, int N, const amrex::Real dt_lev);
580
581 /** \brief RHS terms for 3D momentum */
582 void rhs_uv_3d (int lev,
583 const amrex::Box& xbx,
584 const amrex::Box& ybx,
585 const amrex::Array4<amrex::Real const>& uold,
586 const amrex::Array4<amrex::Real const>& vold,
587 const amrex::Array4<amrex::Real >& ru,
588 const amrex::Array4<amrex::Real >& rv,
589 const amrex::Array4<amrex::Real >& rufrc,
590 const amrex::Array4<amrex::Real >& rvfrc,
591 const amrex::Array4<amrex::Real const>& sustr,
592 const amrex::Array4<amrex::Real const>& svstr,
593 const amrex::Array4<amrex::Real const>& bustr,
594 const amrex::Array4<amrex::Real const>& bvstr,
595 const amrex::Array4<amrex::Real const>& Huon,
596 const amrex::Array4<amrex::Real const>& Hvom,
597 const amrex::Array4<amrex::Real const>& pm,
598 const amrex::Array4<amrex::Real const>& pn,
599 const amrex::Array4<amrex::Real const>& W,
600 const amrex::Array4<amrex::Real >& FC,
601 int nrhs, int N);
602
603 /** \brief RHS terms for 2D momentum */
604 void rhs_uv_2d (int lev,
605 const amrex::Box& xbx,
606 const amrex::Box& ybx,
607 const amrex::Array4<amrex::Real const>& uold,
608 const amrex::Array4<amrex::Real const>& vold,
609 const amrex::Array4<amrex::Real >& ru,
610 const amrex::Array4<amrex::Real >& rv,
611 const amrex::Array4<amrex::Real const>& Duon,
612 const amrex::Array4<amrex::Real const>& Dvom,
613 const int nrhs);
614
615 /** \brief Wrapper around equation of state calculation */
616 void rho_eos (const amrex::Box& bx,
617 const amrex::Array4<amrex::Real const>& state,
618 const amrex::Array4<amrex::Real >& rho,
619 const amrex::Array4<amrex::Real >& rhoA,
620 const amrex::Array4<amrex::Real >& rhoS,
621 const amrex::Array4<amrex::Real >& bvf,
622 const amrex::Array4<amrex::Real >& alpha,
623 const amrex::Array4<amrex::Real >& beta,
624 const amrex::Array4<amrex::Real const>& Hz,
625 const amrex::Array4<amrex::Real const>& z_w,
626 const amrex::Array4<amrex::Real const>& z_r,
627 const amrex::Array4<amrex::Real const>& h,
628 const amrex::Array4<amrex::Real const>& mskr,
629 const int N);
630
631 /** \brief Calculate density and related quantities from linear equation of state */
632 void lin_eos (const amrex::Box& bx,
633 const amrex::Array4<amrex::Real const>& state,
634 const amrex::Array4<amrex::Real >& rho,
635 const amrex::Array4<amrex::Real >& rhoA,
636 const amrex::Array4<amrex::Real >& rhoS,
637 const amrex::Array4<amrex::Real >& bvf,
638 const amrex::Array4<amrex::Real const>& Hz,
639 const amrex::Array4<amrex::Real const>& z_w,
640 const amrex::Array4<amrex::Real const>& z_r,
641 const amrex::Array4<amrex::Real const>& h,
642 const amrex::Array4<amrex::Real const>& mskr,
643 const int N);
644
645 /** \brief Calculate density and related quantities from nonlinear equation of state */
646 void nonlin_eos (const amrex::Box& bx,
647 const amrex::Array4<amrex::Real const>& state,
648 const amrex::Array4<amrex::Real >& rho,
649 const amrex::Array4<amrex::Real >& rhoA,
650 const amrex::Array4<amrex::Real >& rhoS,
651 const amrex::Array4<amrex::Real >& bvf,
652 const amrex::Array4<amrex::Real >& alpha,
653 const amrex::Array4<amrex::Real >& beta,
654 const amrex::Array4<amrex::Real const>& Hz,
655 const amrex::Array4<amrex::Real const>& z_w,
656 const amrex::Array4<amrex::Real const>& z_r,
657 const amrex::Array4<amrex::Real const>& h,
658 const amrex::Array4<amrex::Real const>& mskr,
659 const int N);
660
661 /** \brief Calculate pressure gradient */
662 void prsgrd (const amrex::Box& bx,
663 const amrex::Box& gbx,
664 const amrex::Box& utbx,
665 const amrex::Box& vtbx,
666 const amrex::Array4<amrex::Real >& ru,
667 const amrex::Array4<amrex::Real >& rv,
668 const amrex::Array4<amrex::Real const>& pn,
669 const amrex::Array4<amrex::Real const>& pm,
670 const amrex::Array4<amrex::Real const>& rho,
671 const amrex::Array4<amrex::Real >& FC,
672 const amrex::Array4<amrex::Real const>& Hz,
673 const amrex::Array4<amrex::Real const>& z_r,
674 const amrex::Array4<amrex::Real const>& z_w,
675 const amrex::Array4<amrex::Real const>& msku,
676 const amrex::Array4<amrex::Real const>& mskv,
677 const int nrhs, const int N);
678
679 /** \brief Update velocities or tracers with diffusion/viscosity
680 * as the last part of the prestep */
681 void prestep_diffusion (const amrex::Box& bx,
682 const amrex::Box& gbx,
683 const int ioff, const int joff,
684 const amrex::Array4<amrex::Real >& vel,
685 const amrex::Array4<amrex::Real const>& vel_old,
686 const amrex::Array4<amrex::Real >& rvel,
687 const amrex::Array4<amrex::Real const>& Hz,
688 const amrex::Array4<amrex::Real const>& Akv,
689 const amrex::Array4<amrex::Real >& FC,
690 const amrex::Array4<amrex::Real const>& sstr,
691 const amrex::Array4<amrex::Real const>& bstr,
692 const amrex::Array4<amrex::Real const>& z_r,
693 const amrex::Array4<amrex::Real const>& pm,
694 const amrex::Array4<amrex::Real const>& pn,
695 const int iic, const int ntfirst, const int nnew, int nstp, int nrhs, int N,
696 const amrex::Real lambda, const amrex::Real dt_lev);
697
698 /** \brief Calculate effects of vertical viscosity or diffusivity */
699 void vert_visc_3d (const amrex::Box& bx,
700 const int ioff, const int joff,
701 const amrex::Array4<amrex::Real >& phi,
702 const amrex::Array4<amrex::Real const>& Hz,
703 const amrex::Array4<amrex::Real >& Hzk,
704 const amrex::Array4<amrex::Real >& AK,
705 const amrex::Array4<amrex::Real const>& Akv,
706 const amrex::Array4<amrex::Real >& BC,
707 const amrex::Array4<amrex::Real >& DC,
708 const amrex::Array4<amrex::Real >& FC,
709 const amrex::Array4<amrex::Real >& CF,
710 const int nnew, const int N,
711 const amrex::Real dt_lev);
712
713 /** \brief Correct mass flux */
714 void update_massflux_3d (int lev, const amrex::Box& bx,
715 const int ioff, const int joff,
716 const amrex::Array4<amrex::Real >& phi,
717 const amrex::Array4<amrex::Real >& phibar,
718 const amrex::Array4<amrex::Real >& Hphi,
719 const amrex::Array4<amrex::Real const>& Hz,
720 const amrex::Array4<amrex::Real const>& pm_or_pn,
721 const amrex::Array4<amrex::Real const>& Dphi1,
722 const amrex::Array4<amrex::Real const>& Dphi2,
723 const amrex::Array4<amrex::Real >& DC,
724 const amrex::Array4<amrex::Real >& FC,
725 const amrex::Array4<amrex::Real const>& msk,
726 const int nnew);
727
728 /** \brief Adjust 3D momentum variables based on vertical mean momentum */
729 void vert_mean_3d (const amrex::Box& bx,
730 const int ioff, const int joff,
731 const amrex::Array4<amrex::Real >& phi,
732 const amrex::Array4<amrex::Real const>& Hz,
733 const amrex::Array4<amrex::Real const>& Dphi_avg1,
734 const amrex::Array4<amrex::Real >& DC,
735 const amrex::Array4<amrex::Real >& CF,
736 const amrex::Array4<amrex::Real const>& pm_or_pn,
737 const amrex::Array4<amrex::Real const>& msk,
738 const int nnew, const int N);
739
740 /** \brief Harmonic viscosity */
741 void uv3dmix (const amrex::Box& xbx,
742 const amrex::Box& ybx,
743 const amrex::Array4<amrex::Real >& u,
744 const amrex::Array4<amrex::Real >& v,
745 const amrex::Array4<amrex::Real const>& uold,
746 const amrex::Array4<amrex::Real const>& vold,
747 const amrex::Array4<amrex::Real >& rufrc,
748 const amrex::Array4<amrex::Real >& rvfrc,
749 const amrex::Array4<amrex::Real const>& visc2_p,
750 const amrex::Array4<amrex::Real const>& visc2_r,
751 const amrex::Array4<amrex::Real const>& Hz,
752 const amrex::Array4<amrex::Real const>& pm,
753 const amrex::Array4<amrex::Real const>& pn,
754 const amrex::Array4<amrex::Real const>& mskp,
755 int nrhs, int nnew,
756 const amrex::Real dt_lev);
757
758 /** \brief Harmonic diffusivity for tracers */
759 void t3dmix (const amrex::Box& bx,
760 const amrex::Array4<amrex::Real >& state,
761 const amrex::Array4<amrex::Real >& state_rhs,
762 const amrex::Array4<amrex::Real const>& diff2,
763 const amrex::Array4<amrex::Real const>& Hz,
764 const amrex::Array4<amrex::Real const>& pm,
765 const amrex::Array4<amrex::Real const>& pn,
766 const amrex::Array4<amrex::Real const>& msku,
767 const amrex::Array4<amrex::Real const>& mskv,
768 const amrex::Real dt_lev,
769 const int ncomp);
770
771 /** \brief Calculate Coriolis terms */
772 void coriolis (const amrex::Box& xbx,
773 const amrex::Box& ybx,
774 const amrex::Array4<amrex::Real const>& uold,
775 const amrex::Array4<amrex::Real const>& vold,
776 const amrex::Array4<amrex::Real >& ru,
777 const amrex::Array4<amrex::Real >& rv,
778 const amrex::Array4<amrex::Real const>& Hz,
779 const amrex::Array4<amrex::Real const>& fomn,
780 int nrhs, int nr);
781
782 /** \brief Apply climatology nudging */
783 void apply_clim_nudg (const amrex::Box& bx,
784 int ioff, int joff,
785 const amrex::Array4<amrex::Real >& var,
786 const amrex::Array4<amrex::Real const>& var_old,
787 const amrex::Array4<amrex::Real const>& var_clim,
788 const amrex::Array4<amrex::Real const>& clim_coeff,
789 const amrex::Array4<amrex::Real const>& Hz,
790 const amrex::Array4<amrex::Real const>& pm,
791 const amrex::Array4<amrex::Real const>& pn,
792 const amrex::Real dt_lev = amrex::Real(0.0));
793
794 /** \brief Set 2D momentum arrays from 3D momentum */
795 void set_2darrays (int lev);
796
797 /** \brief Set Zt_avg1 to zeta */
798 void set_zeta_average (int lev);
799
800 /** \brief Initialize zeta from file or analytic */
801 void set_zeta (int lev);
802
803 /** \brief Initialize bathymetry from file or analytic */
804 void set_bathymetry (int lev);
805
806 /** \brief Initialize Coriolis factor from file or analytic */
807 void set_coriolis (int lev);
808
809 /** \brief Initialize land-sea masks from file or analytic */
810 void set_masks (int lev);
811
812 /** \brief Calculate vertical stretched coordinates */
813 void stretch_transform (int lev);
814
815 /** \brief Initialize vertical mixing coefficients from file or analytic */
816 void init_set_vmix (int lev);
817 /** \brief Set vertical mixing coefficients from analytic */
818 void set_analytic_vmix (int lev);
819 /** \brief Initialize GLS variables */
820 void init_gls_vmix (int lev, SolverChoice solver_choice);
821 /** \brief Prestep for GLS calculation */
822 void gls_prestep (int lev, amrex::MultiFab* mf_gls, amrex::MultiFab* mf_tke,
823 amrex::MultiFab& mf_W,
824 amrex::MultiFab* mf_msku, amrex::MultiFab* mf_mskv,
825 const int nstp, const int nnew, const int iic, const int ntfirst,
826 const int N, const amrex::Real dt_lev);
827 /** \brief Corrector step for GLS calculation */
828 void gls_corrector (int lev, amrex::MultiFab* mf_gls, amrex::MultiFab* mf_tke,
829 amrex::MultiFab& mf_W, amrex::MultiFab* mf_Akv, amrex::MultiFab* mf_Akt,
830 amrex::MultiFab* mf_Akk, amrex::MultiFab* mf_Akp,
831 amrex::MultiFab* mf_mskr,
832 amrex::MultiFab* mf_msku, amrex::MultiFab* mf_mskv,
833 const int nstp, const int nnew,
834 const int N, const amrex::Real dt_lev);
835
836 /** \brief Scale RHS momentum variables by 1/cell area, needed before FillPatch to different levels*/
837 void scale_rhs_vars ();
838 /** \brief Scale RHS momentum variables by cell area, needed after FillPatch to different levels*/
839 void scale_rhs_vars_inv ();
840
841 /** \brief Initialize or calculate surface momentum flux from file or analytic */
842 void set_smflux (int lev);
843
844 /** \brief Initialize or calculate wind speed from file or analytic */
845 void set_wind (int lev);
846
847 /** \brief Initialize horizontal mixing coefficients */
848 void set_hmixcoef (int lev);
849
850 /** \brief Initialize flat bathymetry to value from problo */
851 void init_flat_bathymetry (int lev);
852
853 /** \brief Initialize or calculate bottom drag */
854 void set_drag (int lev);
855
856 /** \brief Set weights for averaging 3D variables to 2D */
857 void set_weights (int lev);
858
859 /** \brief Fill the physical boundary conditions for cell-centered velocity (diagnostic only) */
860 void FillBdyCCVels (int lev, amrex::MultiFab& mf_cc_vel);
861
862 /** \brief Fill a new MultiFab by copying in phi from valid region and filling ghost cells */
863 void FillPatch (int lev, amrex::Real time,
864 amrex::MultiFab& mf_to_be_filled,
865 amrex::Vector<amrex::MultiFab*> const& mfs,
866 const int bccomp,
867 const int bdy_var_type = BdyVars::null,
868 const int icomp=0,
869 const bool fill_all=true,
870 const bool fill_set=true,
871 const int n_not_fill=0,
872 const int icomp_calc=0,
873 const amrex::Real dt = amrex::Real(0.0),
874 const amrex::MultiFab& mf_calc = amrex::MultiFab());
875
876 /** \brief Fill a new MultiFab by copying in phi from valid region and filling ghost cells without applying boundary conditions */
877 void FillPatchNoBC (int lev, amrex::Real time,
878 amrex::MultiFab& mf_to_be_filled,
879 amrex::Vector<amrex::MultiFab*> const& mfs,
880 const int bdy_var_type = BdyVars::null,
881 const int icomp=0,
882 const bool fill_all=true,
883 const bool fill_set=true);
884 /** \brief Fill boundary data from netcdf file */
885 void fill_from_bdyfiles (amrex::MultiFab& mf_to_fill,
886 const amrex::MultiFab& mf_mask,
887 const amrex::Real time,
888 const int bccomp,
889 const int bdy_var_type,
890 const int icomp_to_fill,
891 const int icomp_calc = 0,
892 const amrex::MultiFab& mf_calc = amrex::MultiFab(),
893 const amrex::Real = amrex::Real(0.0));
894
895 /** \brief fill an entire multifab by interpolating from the coarser level
896 * using the piecewise constant interpolater */
897 void FillCoarsePatchPC (int lev, amrex::Real time, amrex::MultiFab* mf_fine,
898 amrex::MultiFab* mf_crse,
899 const int bccomp,
900 const int bdy_var_type = BdyVars::null,
901 const int icomp = 0,
902 const bool fill_all = true,
903 const int n_not_fill=0,
904 const int icomp_calc=0,
905 const amrex::Real dt = amrex::Real(0.0),
906 const amrex::MultiFab& mf_calc = amrex::MultiFab());
907
908 /** \brief fill an entire multifab by interpolating from the coarser level */
909 void FillCoarsePatch (int lev, amrex::Real time, amrex::MultiFab* mf_fine,
910 amrex::MultiFab* mf_crse,
911 const int bccomp,
912 const int bdy_var_type = BdyVars::null,
913 const int icomp = 0,
914 const bool fill_all = true,
915 const int n_not_fill=0,
916 const int icomp_calc=0,
917 const amrex::Real dt = amrex::Real(0.0),
918 const amrex::MultiFab& mf_calc = amrex::MultiFab());
919
920 /** \brief fill an entire multifab by interpolating from the coarser level,
921 * explicitly specifying interpolator to use */
922 void FillCoarsePatchMap (int lev, amrex::Real time, amrex::MultiFab* mf_fine,
923 amrex::MultiFab* mf_crse,
924 const int bccomp,
925 const int bdy_var_type = BdyVars::null,
926 const int icomp = 0,
927 const bool fill_all = true,
928 const int n_not_fill=0,
929 const int icomp_calc=0,
930 const amrex::Real dt = amrex::Real(0.0),
931 const amrex::MultiFab& mf_calc = amrex::MultiFab(),
932 amrex::Interpolater* mapper = nullptr);
933
934 /** \brief initialize and calculate stretch coefficients */
935 void init_stretch_coeffs ();
936 /** \brief calculate vertical stretch coefficients */
937 void calc_stretch_coeffs ();
938
939 /** \brief Calculate Coriolis parameters from beta plane parametrization */
940 void init_beta_plane_coriolis (int lev);
941
942 /** \brief Problem initialization from NetCDF file */
943 void init_data_from_netcdf (int lev);
944 /** \brief Boundary data initialization from NetCDF file */
945 void init_bdry_from_netcdf ();
946 /** \brief Mask data initialization from NetCDF file */
947 void init_masks_from_netcdf (int lev);
948 /** \brief Bathymetry data initialization from NetCDF file */
949 void init_bathymetry_from_netcdf (int lev);
950 /** \brief Sea-surface height data initialization from NetCDF file */
951 void init_zeta_from_netcdf (int lev);
952 /** \brief Coriolis parameter data initialization from NetCDF file */
953 void init_coriolis_from_netcdf (int lev);
954 /** \brief Climatology nudging coefficient initialization from NetCDF file */
956 /** \brief Wrapper to initialize climatology nudging coefficient */
957 void init_clim_nudg_coeff (int lev);
958 void init_riv_pos_from_netcdf (int lev);
959
960 /** \brief Convert data in a multifab from inverse days to inverse seconds */
961 void convert_inv_days_to_inv_s (amrex::MultiFab*);
962
963 /** \brief Mask data arrays before writing output */
964 void mask_arrays_for_write (int lev, amrex::Real fill_value, amrex::Real fill_where);
965
966
967#ifdef REMORA_USE_NETCDF
969#endif
970
971private:
972
973 ///////////////////////////
974 // private member functions
975 ///////////////////////////
976
977 /** \brief read in some parameters from inputs file */
978 void ReadParameters();
979
980 /** \brief set covered coarse cells to be the average of overlying fine cells */
981 void AverageDown ();
982
983 /** \brief Read in boundary parameters from input file and set up data structures */
984 void init_bcs();
985
986 /** \brief Initialize initial problem data from analytic functions */
987 void init_analytic(int lev);
988
989 /** \brief Allocate MultiFabs for state and evolution variables */
990 void init_stuff (int lev, const amrex::BoxArray& ba, const amrex::DistributionMapping& dm);
991
992 /** \brief Allocate MultiFabs for masks */
993 void init_masks (int lev, const amrex::BoxArray& ba, const amrex::DistributionMapping& dm);
994
995 /** \brief Resize variable containers to accommodate data on levels 0 to max_lev */
996 void resize_stuff (int lev);
997
998 /** \brief more flexible version of AverageDown() that lets you average down across multiple levels */
999 void AverageDownTo (int crse_lev);
1000
1001 /** \brief Construct FillPatchers */
1002 void Construct_REMORAFillPatchers (int lev);
1003
1004 /** \brief Define FillPatchers */
1005 void Define_REMORAFillPatchers (int lev);
1006
1007 /** \brief utility to copy in data from old and/or new state into another multifab */
1008 TimeInterpolatedData GetDataAtTime (int lev, amrex::Real time);
1009
1010 /** \brief advance a level by dt, includes a recursive call for finer levels */
1011 void timeStep (int lev, amrex::Real time, int iteration);
1012
1013 /** \brief advance all levels by dt, loops over finer levels */
1014 void timeStepML (amrex::Real time, int iteration);
1015
1016 /** \brief a wrapper for estTimeStep() */
1017 void ComputeDt ();
1018
1019 /** \brief get plotfile name */
1020 std::string PlotFileName (int lev) const;
1021
1022 /* set which variables and derived quantities go into plotfiles */
1023 void setPlotVariables (const std::string& pp_plot_var_names);
1024
1025 /* append variables to plot */
1026 void appendPlotVariables (const std::string& pp_plot_var_names);
1027
1028#ifdef REMORA_USE_NETCDF
1029 /** \brief Write plotfile using NetCDF (wrapper) */
1030 void WriteNCPlotFile (int istep);
1031 /** \brief Write a particular NetCDF plotfile */
1032 void WriteNCPlotFile_which (int lev, int which_subdomain,
1033 bool write_header, ncutils::NCFile& ncf,
1034 bool is_history);
1035
1036 /** \brief Write MultiFab in NetCDF format */
1037 void WriteNCMultiFab (const amrex::FabArray<amrex::FArrayBox> &fab,
1038 const std::string& name,
1039 bool set_ghost = false) const;
1040
1041 /** \brief Read MultiFab in NetCDF format */
1042 void ReadNCMultiFab (amrex::FabArray<amrex::FArrayBox> &fab,
1043 const std::string &name,
1044 int coordinatorProc = amrex::ParallelDescriptor::IOProcessorNumber(),
1045 int allow_empty_mf = 0);
1046
1047 /** \brief Vectors (over time) of Vector (over variables) of FArrayBoxs for holding Western boundary data from file */
1048 amrex::Vector<amrex::Vector<amrex::FArrayBox>> bdy_data_xlo;
1049 /** \brief Vectors (over time) of Vector (over variables) of FArrayBoxs for holding Eastern boundary data from file */
1050 amrex::Vector<amrex::Vector<amrex::FArrayBox>> bdy_data_xhi;
1051 /** \brief Vectors (over time) of Vector (over variables) of FArrayBoxs for holding Southern boundary data from file */
1052 amrex::Vector<amrex::Vector<amrex::FArrayBox>> bdy_data_ylo;
1053 /** \brief Vectors (over time) of Vector (over variables) of FArrayBoxs for holding Northern boundary data from file */
1054 amrex::Vector<amrex::Vector<amrex::FArrayBox>> bdy_data_yhi;
1055
1056 /** \brief Start time in the time series of boundary data */
1057 amrex::Real start_bdy_time;
1058 /** \brief Interval between boundary data times */
1060
1061 /** \brief Whether to output NetCDF files as a single history file with several time steps */
1063
1064 /** \brief Data container for u-component surface momentum flux read from file */
1066 /** \brief Data container for v-component surface momentum flux read from file */
1068 /** \brief Data container for u-direction wind read from file */
1070 /** \brief Data container for v-direction wind read from file */
1072
1073 /** \brief Data container for ubar climatology data read from file */
1075 /** \brief Data container for vbar climatology data read from file */
1077 /** \brief Data container for u-velocity climatology data read from file */
1079 /** \brief Data container for v-velocity climatology data read from file */
1081 /** \brief Data container for temperature climatology data read from file */
1083 /** \brief Data container for salinity climatology data read from file */
1085
1086 /** \brief Vector of data containers for scalar data in rivers */
1087 amrex::Vector<NCTimeSeriesRiver*> river_source_cons;
1088 /** \brief Data container for momentum transport in rivers */
1090 /** \brief Data container for vertically integrated momentum transport in rivers */
1092#endif // REMORA_USE_NETCDF
1093 /** \brief iMultiFab for river positions; contents are indices of rivers */
1094 amrex::Vector<std::unique_ptr<amrex::iMultiFab>> vec_river_position;
1095 /** \brief Vector over rivers of river direction: 0: u-face; 1: v-face; 2: w-face */
1096 amrex::Gpu::DeviceVector<int> river_direction;
1097
1098
1099#ifdef REMORA_USE_MOAB
1100 /** \brief Initialize connection to MOAB */
1101 void InitMOABMesh();
1102#endif
1103 /** \brief Nudging width at coarse-fine interface */
1104 int cf_width{0};
1105 /** \brief Width for fixing values at coarse-fine interface */
1107 // Fillpatcher classes for coarse-fine boundaries
1108 /** \brief Vector over levels of FillPatchers for scalars */
1109 amrex::Vector<REMORAFillPatcher> FPr_c;
1110 /** \brief Vector over levels of FillPatchers for u (3D) */
1111 amrex::Vector<REMORAFillPatcher> FPr_u;
1112 /** \brief Vector over levels of FillPatchers for v (3D) */
1113 amrex::Vector<REMORAFillPatcher> FPr_v;
1114 /** \brief Vector over levels of FillPatchers for w */
1115 amrex::Vector<REMORAFillPatcher> FPr_w;
1116
1117 // Fillpatcher classes for 2d velocity variables
1118 /** \brief Vector over levels of FillPatchers for ubar (2D) */
1119 amrex::Vector<REMORAFillPatcher> FPr_ubar;
1120 /** \brief Vector over levels of FillPatchers for vbar (2D) */
1121 amrex::Vector<REMORAFillPatcher> FPr_vbar;
1122
1123 /** \brief write checkpoint file to disk */
1124 void WriteCheckpointFile ();
1125
1126 /** \brief read checkpoint file from disk */
1127 void ReadCheckpointFile ();
1128
1129 /** \brief Read the file passed to remora.restart and use it as an initial condition for the current simulation */
1131
1132 /** \brief Initialize the new-time data at a level from the initial_data MultiFab */
1133 void InitializeLevelFromData (int lev, const amrex::MultiFab& initial_data);
1134
1135 /** \brief utility to skip to next line in Header */
1136 static void GotoNextLine (std::istream& is);
1137
1138 ////////////////
1139 // private data members
1140
1141 /** \brief Pointer to container of analytical functions for problem definition */
1142 std::unique_ptr<ProblemBase> prob = nullptr;
1143
1144 /** \brief how many boxes specified at each level by tagging criteria */
1145 amrex::Vector<int> num_boxes_at_level;
1146 /** \brief how many netcdf input files specified at each level */
1147 amrex::Vector<int> num_files_at_level;
1148 /** \brief the boxes specified at each level by tagging criteria */
1149 amrex::Vector<amrex::Vector<amrex::Box>> boxes_at_level;
1150
1151 /** \brief which step? */
1152 amrex::Vector<int> istep;
1153 /** \brief How many substeps on each level? */
1154 amrex::Vector<int> nsubsteps;
1155 /** \brief new time at each level */
1156 amrex::Vector<amrex::Real> t_new;
1157 /** \brief old time at each level */
1158 amrex::Vector<amrex::Real> t_old;
1159 /** \brief time step at each level */
1160 amrex::Vector<amrex::Real> dt;
1161
1162 /** \brief whether to set boundary conditions by variable rather than just by side */
1164
1165 /** \brief Vector (over level) of functors to apply physical boundary conditions */
1166 amrex::Vector<std::unique_ptr<REMORAPhysBCFunct>> physbcs;
1167
1168 /** \brief array of flux registers for refluxing in multilevel */
1169 amrex::Vector<amrex::YAFluxRegister*> advflux_reg;
1170
1171 ////////////////
1172 // boundary data members
1173
1174 // A BCRec is essentially a 2*DIM integer array storing the boundary
1175 // condition type at each lo/hi walls in each direction. We have one BCRec
1176 // for each component of the cell-centered variables and each velocity component.
1177 /** \brief vector (over BCVars) of BCRecs */
1178 amrex::Vector <amrex::BCRec> domain_bcs_type;
1179 /** \brief GPU vector (over BCVars) of BCRecs */
1180 amrex::Gpu::DeviceVector<amrex::BCRec> domain_bcs_type_d;
1181
1182 /** \brief Array of strings describing domain boundary conditions */
1183 amrex::Array<std::string,2*AMREX_SPACEDIM> domain_bc_type;
1184
1185 /** \brief Array holding the Dirichlet values at walls which need them */
1186 amrex::Array<amrex::Array<amrex::Real, AMREX_SPACEDIM*2>,AMREX_SPACEDIM+NCONS+8> m_bc_extdir_vals;
1187
1188 /** \brief Array holding the "physical" boundary condition types (e.g. "inflow") */
1189 amrex::GpuArray<amrex::GpuArray<REMORA_BC, AMREX_SPACEDIM*2>,BCVars::NumTypes> phys_bc_type;
1190
1191 /** \brief These are flags that indicate whether we need to read in boundary data from file */
1192 amrex::GpuArray<amrex::GpuArray<bool, AMREX_SPACEDIM*2>,BdyVars::NumTypes+1> phys_bc_need_data;
1193
1194 /** \brief Container to connect boundary data being read in boundary condition containers.
1195 *
1196 * Needed for phys_bc_need_data to work correctly
1197 */
1198 amrex::Vector<int> bdy_index;
1199
1200 /** \brief Step when we last output a plotfile */
1202 /** \brief Simulation time when we last output a plotfile */
1204
1205 /** \brief Step when we last output a checkpoint file */
1207 /** \brief Simulation time when we last output a checkpoint file */
1209 /** \brief Whether to output a plotfile on restart from checkpoint */
1211
1212 ////////////////
1213 // runtime parameters
1214
1215 /** \brief maximum number of steps */
1216 int max_step = std::numeric_limits<int>::max();
1217 /** \brief Time to stop */
1218 amrex::Real stop_time = std::numeric_limits<amrex::Real>::max();
1219
1220 /** \brief Time of the start of the simulation, in seconds */
1221 amrex::Real start_time = 0.0;
1222
1223 /** \brief If set, restart from this checkpoint file */
1224 std::string restart_chkfile = "";
1225
1226 // Time step controls
1227 /** \brief CFL condition */
1228 static amrex::Real cfl;
1229 /** \brief Fraction maximum change in subsequent time steps */
1230 static amrex::Real change_max;
1231 /** \brief User specified fixed baroclinic time step */
1232 static amrex::Real fixed_dt;
1233 /** \brief User specified fixed barotropic time step */
1234 static amrex::Real fixed_fast_dt;
1235 /** \brief User specified, number of barotropic steps per baroclinic step */
1237 /** \brief Number of fast steps to take */
1239
1240 /** \brief Whether to substep fine levels in time */
1241 int do_substep = 0;
1242
1243 /** \brief how often each level regrids the higher levels of refinement (after a level advances that many time steps) */
1244 int regrid_int = 2;
1245
1246 // I/O controls
1247 /** \brief Plotfile prefix */
1248 std::string plot_file_name {"plt_"};
1249 /** \brief Plotfile output interval in iterations */
1250 int plot_int = -1;
1251 /** \brief Plotfile output interval in seconds */
1252 amrex::Real plot_int_time = amrex::Real(-1.0);
1253 /** \brief Checkpoint file prefix */
1254 std::string check_file {"chk"};
1255 /** \brief Checkpoint output interval in iterations */
1256 int check_int = -1;
1257 /** \brief Checkpoint output interval in seconds */
1258 amrex::Real check_int_time = amrex::Real(-1.0);
1259
1260 /** \brief Whether to chunk netcdf history file */
1262 /** \brief Number of time steps per netcdf history file.
1263 *
1264 -1 is the default and means code will automatically compute
1265 number of steps per file based on grid size to keep files < 2 GB */
1267 /** \brief Counter for which time index we are writing to in the netcdf history file */
1269
1270 /** \brief Names of variables to output to AMReX plotfile */
1271 amrex::Vector<std::string> plot_var_names;
1272 /** \brief Names of scalars for plotfile output */
1273 const amrex::Vector<std::string> cons_names {"temp", "salt", "scalar"};
1274 // Note that the order of variable names here must match the order in PlotFile.cpp
1275 /** \brief Names of derived fields for plotfiles */
1276 const amrex::Vector<std::string> derived_names {"vorticity"};
1277
1278 /** \brief Container for algorithmic choices */
1280
1281#ifdef REMORA_USE_PARTICLES
1282 /** \brief Particle container with all particle species */
1283 ParticleData particleData;
1284
1285 // variables and functions for tracers particles
1286 /** \brief tracer particles that advect with flow */
1287 bool m_use_tracer_particles;
1288 /** \brief tracer particles that fall with gravity */
1289 bool m_use_hydro_particles;
1290
1291 /** \brief Read tracer and hydro particles parameters */
1292 void readTracersParams();
1293
1294 /** \brief Initialize tracer and hydro particles */
1295 void initializeTracers ( amrex::ParGDBBase*,
1296 const amrex::Vector<std::unique_ptr<amrex::MultiFab>>& );
1297
1298 /** \brief Evolve tracers and hydro particles */
1299 void evolveTracers( int lev,
1300 amrex::Real dt,
1301 amrex::Vector<amrex::MultiFab const*>& flow_vel,
1302 const amrex::Vector<std::unique_ptr<amrex::MultiFab>>& );
1303
1304#endif
1305 /** \brief Verbosity level of output */
1306 static int verbose;
1307
1308 /** \brief Diagnostic sum output interval in number of steps */
1309 static int sum_interval;
1310 /** \brief Diagnostic sum output interval in time */
1311 static amrex::Real sum_per;
1312
1313 /** \brief Minimum number of digits in plotfile name or chunked history file */
1315
1316 /** \brief Whether to write the staggered velocities (not averaged to cell centers) */
1318
1319 /** \brief Native or NetCDF plotfile output */
1321
1322 /** \brief NetCDF initialization file */
1323 static amrex::Vector<amrex::Vector<std::string>> nc_init_file;
1324 /** \brief NetCDF grid file */
1325 static amrex::Vector<amrex::Vector<std::string>> nc_grid_file;
1326 /** \brief NetCDF boundary data */
1327 static amrex::Vector<std::string> nc_bdry_file;
1328
1329 /** \brief NetCDF forcing file */
1330 std::string nc_frc_file;
1331 std::string nc_riv_file;
1332
1333 /** \brief NetCDF climatology history file */
1334 std::string nc_clim_his_file;
1335 /** \brief NetCDF climatology coefficient file */
1337
1338 /** \brief Name of time field for boundary data */
1339 std::string bdry_time_varname = "ocean_time";
1340 /** \brief Name of time field for ubar climatology data */
1341 std::string clim_ubar_time_varname = "ocean_time";
1342 /** \brief Name of time field for vbar climatology data */
1343 std::string clim_vbar_time_varname = "ocean_time";
1344 /** \brief Name of time field for u climatology data */
1345 std::string clim_u_time_varname = "ocean_time";
1346 /** \brief Name of time field for v climatology data */
1347 std::string clim_v_time_varname = "ocean_time";
1348 /** \brief Name of time field for salinity climatology data */
1349 std::string clim_salt_time_varname = "ocean_time";
1350 /** \brief Name of time field for temperature climatology data */
1351 std::string clim_temp_time_varname = "ocean_time";
1352 /** \brief Name of time field for river time */
1353 std::string riv_time_varname = "river_time";
1354 /** \brief Name of time field for forcing data */
1355 std::string frc_time_varname = "";
1356
1357
1358 /** \brief Set refinement criteria */
1360
1361 /** \brief Holds info for dynamically generated tagging criteria */
1362 static amrex::Vector<amrex::AMRErrorTag> ref_tags;
1363
1364 /** \brief Mask that zeroes out values on a coarse level underlying grids on the next finest level */
1365 amrex::MultiFab fine_mask;
1366
1367 /** \brief Helper function to determine number of ghost cells */
1368 AMREX_FORCE_INLINE
1369 int
1370 ComputeGhostCells(const int& spatial_order) {
1371 int nGhostCells = 0;
1372 switch (spatial_order) {
1373 case 2:
1374 nGhostCells = 2; // We need this many to compute the eddy viscosity in the ghost cells
1375 break;
1376 case 3:
1377 nGhostCells = 2;
1378 break;
1379 case 4:
1380 nGhostCells = 2;
1381 break;
1382 case 5:
1383 nGhostCells = 3;
1384 break;
1385 case 6:
1386 nGhostCells = 3;
1387 break;
1388 default:
1389 amrex::Error("Must specify spatial order to be 2,3,4,5 or 6");
1390 }
1391
1392 return nGhostCells;
1393 }
1394
1395 /** \brief Get flux register */
1396 AMREX_FORCE_INLINE
1397 amrex::YAFluxRegister* getAdvFluxReg (int lev)
1398 {
1399 return advflux_reg[lev];
1400 }
1401
1402 /** \brief Helper function for IO stream */
1403 AMREX_FORCE_INLINE
1404 std::ostream&
1405 DataLog (int i)
1406 {
1407 return *datalog[i];
1408 }
1409
1410 AMREX_FORCE_INLINE
1411 int
1412 NumDataLogs () noexcept
1413 {
1414 return datalog.size();
1415 }
1416
1417 /** \brief Variable for CPU timing */
1418 static amrex::Real startCPUTime;
1419 /** \brief Accumulator variable for CPU time used thusfar */
1420 static amrex::Real previousCPUTimeUsed;
1421
1422 /** \brief Get CPU time used */
1423 amrex::Real
1425 {
1426 int numCores = amrex::ParallelDescriptor::NProcs();
1427#ifdef _OPENMP
1428 numCores = numCores * omp_get_max_threads();
1429#endif
1430
1431 amrex::Real T =
1432 numCores * (amrex::Real(amrex::ParallelDescriptor::second()) - startCPUTime) +
1434
1435 return T;
1436 }
1437
1438 void setRecordDataInfo (int i, const std::string& filename)
1439 {
1440 if (amrex::ParallelDescriptor::IOProcessor())
1441 {
1442 datalog[i] = std::make_unique<std::fstream>();
1443 datalog[i]->open(filename.c_str(),std::ios::out|std::ios::app);
1444 if (!datalog[i]->good()) {
1445 amrex::FileOpenFailed(filename);
1446 }
1447 }
1448 amrex::ParallelDescriptor::Barrier("REMORA::setRecordDataInfo");
1449 }
1450
1451 amrex::Vector<std::unique_ptr<std::fstream> > datalog;
1452 amrex::Vector<std::string> datalogname;
1453
1454 /** \brief The filename of the ith datalog file. */
1455 const std::string DataLogName (int i) const noexcept { return datalogname[i]; }
1456
1457public:
1458 /** \brief Write job info to stdout */
1459 void writeJobInfo (const std::string& dir) const;
1460 /** \brief Write build info to os */
1461 static void writeBuildInfo (std::ostream& os);
1462
1463 static void print_banner(MPI_Comm /*comm*/, std::ostream& /*out*/);
1464 static void print_usage(MPI_Comm /*comm*/, std::ostream& /*out*/);
1465 static void print_error(MPI_Comm /*comm*/, const std::string& msg);
1466 static void print_summary(std::ostream&);
1467 static void print_tpls(std::ostream& /*out*/);
1468
1469 /** \brief Accessor method for t_old to expose to outside classes */
1470 amrex::Real get_t_old(int lev) const;
1471
1472 /** \brief Evaluate stability function psi for wind speed */
1473 AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE
1474 static amrex::Real bulk_psiu(amrex::Real ZoL) {
1475 // Compute stability function, psi
1476 using namespace amrex;
1477
1478 Real psiu;
1479 if (ZoL < 0.0_rt) {
1480 // Unstable conditions.
1481 Real x = std::pow(1.0_rt-15.0_rt*ZoL,0.25_rt);
1482 Real psik = 2.0_rt*std::log(0.5_rt*(1.0_rt+x))+
1483 std::log(0.5_rt*(1.0_rt+x*x))-
1484 2.0_rt*std::atan(x)+0.5_rt*PI;
1485 // For very unstable conditions, use free-convection (Fairall).
1486 Real sqrt3 = std::sqrt(3.0_rt);
1487 Real y = std::pow(1.0_rt-10.15_rt*ZoL,1.0_rt/3.0_rt);
1488 Real psic = 1.5_rt*std::log(1.0_rt/3.0_rt*(1.0_rt+y+y*y))-
1489 sqrt3*std::atan((1.0_rt+2.0_rt*y)/sqrt3)+PI/sqrt3;
1490 // Match Kansas and free-convection forms with weighting Fw.
1491 Real Zol2 = ZoL*ZoL;
1492 Real Fw = Zol2/(1.0_rt+Zol2);
1493
1494 psiu = (1.0_rt-Fw)*psik+Fw*psic;
1495 } else {
1496 // Stable conditions
1497 Real cff=std::min(50.0_rt,0.35_rt*ZoL);
1498 psiu = -((1.0_rt+ZoL)+0.6667_rt*(ZoL-14.28_rt)/
1499 std::exp(cff)+8.525_rt);
1500 }
1501 return psiu;
1502 }
1503
1504 /** \brief Evaluate stability function psi for moisture and heat */
1505 AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE
1506 static amrex::Real bulk_psit(amrex::Real ZoL) {
1507 // Compute stability function, psi
1508 using namespace amrex;
1509
1510 Real psit;
1511 if (ZoL < 0.0_rt) {
1512 // Unstable conditions.
1513 Real x = std::pow(1.0_rt-15.0_rt*ZoL,0.5_rt);
1514 Real psik = 2.0_rt*std::log(0.5_rt*(1.0_rt+x));
1515 // For very unstable conditions, use free-convection (Fairall).
1516 Real sqrt3=std::sqrt(3.0_rt);
1517 Real y=std::pow(1.0_rt-34.15_rt*ZoL,1.0_rt/3.0_rt);
1518 Real psic=1.5_rt*std::log(1.0_rt/3.0_rt*(1.0_rt+y+y*y))-
1519 sqrt3*std::atan((1.0_rt+2.0_rt*y)/sqrt3)+PI/sqrt3;
1520
1521 // Match Kansas and free-convection forms with weighting Fw.
1522 Real ZoL2=ZoL*ZoL;
1523 Real Fw=ZoL2/(1.0_rt+ZoL2);
1524 psit = (1.0_rt-Fw)*psik+Fw*psic;
1525 } else {
1526 // Stable conditions.
1527 Real cff=std::min(50.0_rt,0.35_rt*ZoL);
1528 psit=-(std::pow(1.0_rt+2.0_rt*ZoL,1.5_rt)+
1529 0.6667_rt*(ZoL-14.28_rt)/std::exp(cff)+8.525_rt);
1530 }
1531 return psit;
1532 }
1533};
1534
1535#endif
constexpr amrex::Real PI
PlotfileType
plotfile format
#define NCONS
A class to hold and interpolate time series data read from a NetCDF file.
Class that stores all relevant simulation state data with methods for time stepping.
Definition REMORA.H:88
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
const amrex::Vector< std::string > cons_names
Names of scalars for plotfile output.
Definition REMORA.H:1273
void FillCoarsePatchPC(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 using the piecewise constant interpol...
std::string riv_time_varname
Name of time field for river time.
Definition REMORA.H:1353
static void GotoNextLine(std::istream &is)
utility to skip to next line in Header
int nfast
Number of fast steps to take.
Definition REMORA.H:1238
void WriteNCMultiFab(const amrex::FabArray< amrex::FArrayBox > &fab, const std::string &name, bool set_ghost=false) const
Write MultiFab in NetCDF format.
amrex::Vector< std::unique_ptr< amrex::MultiFab > > vec_rv2d
v velocity RHS (2D, includes horizontal and vertical advection)
Definition REMORA.H:250
void prsgrd(const amrex::Box &bx, const amrex::Box &gbx, const amrex::Box &utbx, const amrex::Box &vtbx, const amrex::Array4< amrex::Real > &ru, const amrex::Array4< amrex::Real > &rv, const amrex::Array4< amrex::Real const > &pn, const amrex::Array4< amrex::Real const > &pm, const amrex::Array4< amrex::Real const > &rho, const amrex::Array4< amrex::Real > &FC, const amrex::Array4< amrex::Real const > &Hz, const amrex::Array4< amrex::Real const > &z_r, const amrex::Array4< amrex::Real const > &z_w, const amrex::Array4< amrex::Real const > &msku, const amrex::Array4< amrex::Real const > &mskv, const int nrhs, const int N)
Calculate pressure gradient.
static amrex::Real fixed_dt
User specified fixed baroclinic time step.
Definition REMORA.H:1232
amrex::Vector< std::unique_ptr< amrex::MultiFab > > vec_evap
evaporation rate [kg/m^2/s]
Definition REMORA.H:321
amrex::Real last_plot_file_time
Simulation time when we last output a plotfile.
Definition REMORA.H:1203
void update_massflux_3d(int lev, const amrex::Box &bx, const int ioff, const int joff, const amrex::Array4< amrex::Real > &phi, const amrex::Array4< amrex::Real > &phibar, const amrex::Array4< amrex::Real > &Hphi, const amrex::Array4< amrex::Real const > &Hz, const amrex::Array4< amrex::Real const > &pm_or_pn, const amrex::Array4< amrex::Real const > &Dphi1, const amrex::Array4< amrex::Real const > &Dphi2, const amrex::Array4< amrex::Real > &DC, const amrex::Array4< amrex::Real > &FC, const amrex::Array4< amrex::Real const > &msk, const int nnew)
Correct mass flux.
void scale_rhs_vars()
Scale RHS momentum variables by 1/cell area, needed before FillPatch to different levels.
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.
amrex::Vector< amrex::BCRec > domain_bcs_type
vector (over BCVars) of BCRecs
Definition REMORA.H:1178
bool set_bcs_by_var
whether to set boundary conditions by variable rather than just by side
Definition REMORA.H:1163
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
amrex::Vector< std::unique_ptr< amrex::MultiFab > > vec_btflux
Bottom tracer flux; input arrays.
Definition REMORA.H:316
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< std::unique_ptr< amrex::MultiFab > > vec_rubar
barotropic x velocity for the RHS (2D)
Definition REMORA.H:344
amrex::Vector< std::unique_ptr< amrex::MultiFab > > vec_pm
horizontal scaling factor: 1 / dx (2D)
Definition REMORA.H:366
amrex::MultiFab fine_mask
Mask that zeroes out values on a coarse level underlying grids on the next finest level.
Definition REMORA.H:1365
amrex::Vector< REMORAFillPatcher > FPr_v
Vector over levels of FillPatchers for v (3D)
Definition REMORA.H:1113
amrex::Vector< std::string > plot_var_names
Names of variables to output to AMReX plotfile.
Definition REMORA.H:1271
amrex::Vector< std::unique_ptr< amrex::MultiFab > > vec_ZoBot
Bottom roughness length [m], defined at rho points.
Definition REMORA.H:328
amrex::Vector< std::unique_ptr< amrex::MultiFab > > vec_DU_avg2
correct time average of barotropic x velocity flux for coupling (2D)
Definition REMORA.H:338
amrex::Vector< std::unique_ptr< amrex::MultiFab > > vec_lrflx
longwave radiation
Definition REMORA.H:303
amrex::Vector< std::unique_ptr< amrex::MultiFab > > vec_yv
y_grid on v-points (2D)
Definition REMORA.H:386
static void print_error(MPI_Comm, const std::string &msg)
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
virtual void MakeNewLevelFromCoarse(int lev, amrex::Real time, const amrex::BoxArray &ba, const amrex::DistributionMapping &dm) override
Make a new level using provided BoxArray and DistributionMapping and fill with interpolated coarse le...
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::Gpu::DeviceVector< amrex::Real > s_w
Scaled vertical coordinate (range [0,1]) that transforms to z, defined at w-points (cell faces)
Definition REMORA.H:274
amrex::GpuArray< amrex::GpuArray< REMORA_BC, AMREX_SPACEDIM *2 >, BCVars::NumTypes > phys_bc_type
Array holding the "physical" boundary condition types (e.g. "inflow")
Definition REMORA.H:1189
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
void advance_3d(int lev, amrex::MultiFab &mf_cons, amrex::MultiFab &mf_u, amrex::MultiFab &mf_v, amrex::MultiFab *mf_sstore, amrex::MultiFab *mf_ru, amrex::MultiFab *mf_rv, std::unique_ptr< amrex::MultiFab > &mf_DU_avg1, std::unique_ptr< amrex::MultiFab > &mf_DU_avg2, std::unique_ptr< amrex::MultiFab > &mf_DV_avg1, std::unique_ptr< amrex::MultiFab > &mf_DV_avg2, std::unique_ptr< amrex::MultiFab > &mf_ubar, std::unique_ptr< amrex::MultiFab > &mf_vbar, std::unique_ptr< amrex::MultiFab > &mf_Akv, std::unique_ptr< amrex::MultiFab > &mf_Akt, std::unique_ptr< amrex::MultiFab > &mf_Hz, std::unique_ptr< amrex::MultiFab > &mf_Huon, std::unique_ptr< amrex::MultiFab > &mf_Hvom, std::unique_ptr< amrex::MultiFab > &mf_z_w, amrex::MultiFab const *mf_h, amrex::MultiFab const *mf_pm, amrex::MultiFab const *mf_pn, amrex::MultiFab const *mf_mskr, amrex::MultiFab const *mf_msku, amrex::MultiFab const *mf_mskv, const int N, const amrex::Real dt_lev)
Advance the 3D variables.
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
amrex::Vector< std::unique_ptr< amrex::MultiFab > > vec_rain
precipitation rate [kg/m^2/s]
Definition REMORA.H:319
amrex::Vector< std::unique_ptr< amrex::MultiFab > > vec_tke
Turbulent kinetic energy.
Definition REMORA.H:414
void rho_eos(const amrex::Box &bx, const amrex::Array4< amrex::Real const > &state, const amrex::Array4< amrex::Real > &rho, const amrex::Array4< amrex::Real > &rhoA, const amrex::Array4< amrex::Real > &rhoS, const amrex::Array4< amrex::Real > &bvf, const amrex::Array4< amrex::Real > &alpha, const amrex::Array4< amrex::Real > &beta, const amrex::Array4< amrex::Real const > &Hz, const amrex::Array4< amrex::Real const > &z_w, const amrex::Array4< amrex::Real const > &z_r, const amrex::Array4< amrex::Real const > &h, const amrex::Array4< amrex::Real const > &mskr, const int N)
Wrapper around equation of state calculation.
static void print_summary(std::ostream &)
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
amrex::Vector< std::unique_ptr< amrex::MultiFab > > vec_stflx
Surface tracer flux; working arrays.
Definition REMORA.H:310
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
amrex::Vector< std::unique_ptr< amrex::MultiFab > > vec_gls
Turbulent generic length scale.
Definition REMORA.H:416
bool chunk_history_file
Whether to chunk netcdf history file.
Definition REMORA.H:1261
void writeJobInfo(const std::string &dir) const
Write job info to stdout.
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
amrex::Vector< std::unique_ptr< amrex::MultiFab > > vec_yp
y_grid on psi-points (2D)
Definition REMORA.H:391
amrex::Vector< std::unique_ptr< amrex::MultiFab > > vec_xr
x_grid on rho points (2D)
Definition REMORA.H:374
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
virtual void ClearLevel(int lev) override
Delete level data Overrides the pure virtual function in AmrCore.
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
AMREX_FORCE_INLINE int ComputeGhostCells(const int &spatial_order)
Helper function to determine number of ghost cells.
Definition REMORA.H:1370
amrex::Vector< std::unique_ptr< amrex::MultiFab > > vec_sstore
additional scratch space for calculations on temp, salt, etc
Definition REMORA.H:394
amrex::Vector< std::unique_ptr< amrex::MultiFab > > vec_xv
x_grid on v-points (2D)
Definition REMORA.H:384
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 rhs_uv_3d(int lev, const amrex::Box &xbx, const amrex::Box &ybx, const amrex::Array4< amrex::Real const > &uold, const amrex::Array4< amrex::Real const > &vold, const amrex::Array4< amrex::Real > &ru, const amrex::Array4< amrex::Real > &rv, const amrex::Array4< amrex::Real > &rufrc, const amrex::Array4< amrex::Real > &rvfrc, const amrex::Array4< amrex::Real const > &sustr, const amrex::Array4< amrex::Real const > &svstr, const amrex::Array4< amrex::Real const > &bustr, const amrex::Array4< amrex::Real const > &bvstr, const amrex::Array4< amrex::Real const > &Huon, const amrex::Array4< amrex::Real const > &Hvom, const amrex::Array4< amrex::Real const > &pm, const amrex::Array4< amrex::Real const > &pn, const amrex::Array4< amrex::Real const > &W, const amrex::Array4< amrex::Real > &FC, int nrhs, int N)
RHS terms for 3D momentum.
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< std::unique_ptr< amrex::MultiFab > > vec_Lscale
Vertical mixing turbulent length scale.
Definition REMORA.H:418
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
void prestep_diffusion(const amrex::Box &bx, const amrex::Box &gbx, const int ioff, const int joff, const amrex::Array4< amrex::Real > &vel, const amrex::Array4< amrex::Real const > &vel_old, const amrex::Array4< amrex::Real > &rvel, const amrex::Array4< amrex::Real const > &Hz, const amrex::Array4< amrex::Real const > &Akv, const amrex::Array4< amrex::Real > &FC, const amrex::Array4< amrex::Real const > &sstr, const amrex::Array4< amrex::Real const > &bstr, const amrex::Array4< amrex::Real const > &z_r, const amrex::Array4< amrex::Real const > &pm, const amrex::Array4< amrex::Real const > &pn, const int iic, const int ntfirst, const int nnew, int nstp, int nrhs, int N, const amrex::Real lambda, const amrex::Real dt_lev)
Update velocities or tracers with diffusion/viscosity as the last part of the prestep.
void vert_mean_3d(const amrex::Box &bx, const int ioff, const int joff, const amrex::Array4< amrex::Real > &phi, const amrex::Array4< amrex::Real const > &Hz, const amrex::Array4< amrex::Real const > &Dphi_avg1, const amrex::Array4< amrex::Real > &DC, const amrex::Array4< amrex::Real > &CF, const amrex::Array4< amrex::Real const > &pm_or_pn, const amrex::Array4< amrex::Real const > &msk, const int nnew, const int N)
Adjust 3D momentum variables based on vertical mean momentum.
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
void t3dmix(const amrex::Box &bx, const amrex::Array4< amrex::Real > &state, const amrex::Array4< amrex::Real > &state_rhs, const amrex::Array4< amrex::Real const > &diff2, const amrex::Array4< amrex::Real const > &Hz, const amrex::Array4< amrex::Real const > &pm, const amrex::Array4< amrex::Real const > &pn, const amrex::Array4< amrex::Real const > &msku, const amrex::Array4< amrex::Real const > &mskv, const amrex::Real dt_lev, const int ncomp)
Harmonic diffusivity for tracers.
amrex::Vector< std::unique_ptr< amrex::MultiFab > > vec_Hz
Width of cells in the vertical (z-) direction (3D, Hz in ROMS)
Definition REMORA.H:238
amrex::Vector< std::unique_ptr< amrex::MultiFab > > vec_Akt
Vertical diffusion coefficient (3D)
Definition REMORA.H:258
amrex::Vector< std::unique_ptr< amrex::MultiFab > > vec_msku
land/sea mask at x-faces (2D)
Definition REMORA.H:359
amrex::Vector< std::unique_ptr< amrex::MultiFab > > vec_rvfrc
v velocity RHS, integrated, including advection and bottom/surface stresses (2D)
Definition REMORA.H:254
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
void FillPatchNoBC(int lev, amrex::Real time, amrex::MultiFab &mf_to_be_filled, amrex::Vector< amrex::MultiFab * > const &mfs, const int bdy_var_type=BdyVars::null, const int icomp=0, const bool fill_all=true, const bool fill_set=true)
Fill a new MultiFab by copying in phi from valid region and filling ghost cells without applying boun...
static amrex::Real startCPUTime
Variable for CPU timing.
Definition REMORA.H:1418
void setPlotVariables(const std::string &pp_plot_var_names)
amrex::Vector< amrex::Vector< amrex::FArrayBox > > bdy_data_yhi
Vectors (over time) of Vector (over variables) of FArrayBoxs for holding Northern boundary data from ...
Definition REMORA.H:1054
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 WriteNCPlotFile_which(int lev, int which_subdomain, bool write_header, ncutils::NCFile &ncf, bool is_history)
Write a particular NetCDF plotfile.
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
amrex::Vector< std::unique_ptr< amrex::MultiFab > > vec_rufrc
u velocity RHS, integrated, including advection and bottom/surface stresses (2D)
Definition REMORA.H:252
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
amrex::Real bdy_time_interval
Interval between boundary data times.
Definition REMORA.H:1059
void init_bdry_from_netcdf()
Boundary data initialization from NetCDF file.
void gls_prestep(int lev, amrex::MultiFab *mf_gls, amrex::MultiFab *mf_tke, amrex::MultiFab &mf_W, amrex::MultiFab *mf_msku, amrex::MultiFab *mf_mskv, const int nstp, const int nnew, const int iic, const int ntfirst, const int N, const amrex::Real dt_lev)
Prestep for GLS calculation.
amrex::Real start_bdy_time
Start time in the time series of boundary data.
Definition REMORA.H:1057
NCTimeSeries * sustr_data_from_file
Data container for u-component surface momentum flux read from file.
Definition REMORA.H:1065
amrex::Vector< amrex::Real > vec_weight2
Weights for calculating avg2 in 2D advance.
Definition REMORA.H:410
amrex::Vector< int > bdy_index
Container to connect boundary data being read in boundary condition containers.
Definition REMORA.H:1198
amrex::Gpu::DeviceVector< amrex::Real > s_r
Scaled vertical coordinate (range [0,1]) that transforms to z, defined at rho points (cell centers)
Definition REMORA.H:272
void Define_REMORAFillPatchers(int lev)
Define FillPatchers.
Definition REMORA.cpp:422
amrex::Vector< std::unique_ptr< amrex::MultiFab > > vec_shflx
sensible heat flux
Definition REMORA.H:307
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
void set_drag(int lev)
Initialize or calculate bottom drag.
amrex::MultiFab & build_fine_mask(int lev)
Make mask to zero out covered cells (for mesh refinement)
amrex::Real plot_int_time
Plotfile output interval in seconds.
Definition REMORA.H:1252
amrex::Vector< std::unique_ptr< amrex::MultiFab > > vec_rvbar
barotropic y velocity for the RHS (2D)
Definition REMORA.H:346
amrex::Vector< int > num_files_at_level
how many netcdf input files specified at each level
Definition REMORA.H:1147
void advance_2d_onestep(int lev, amrex::Real dt_lev, amrex::Real dtfast_lev, int my_iif, int nfast_counter)
2D advance, one predictor/corrector step
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< std::unique_ptr< amrex::MultiFab > > vec_z_r
z coordinates at rho points (cell centers)
Definition REMORA.H:267
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
amrex::Vector< std::unique_ptr< amrex::MultiFab > > vec_lhflx
latent heat flux
Definition REMORA.H:305
void refinement_criteria_setup()
Set refinement criteria.
amrex::Vector< std::unique_ptr< amrex::MultiFab > > vec_mskp
land/sea mask at cell corners (2D)
Definition REMORA.H:363
void remora_advance(int level, amrex::MultiFab &cons_old, amrex::MultiFab &cons_new, amrex::MultiFab &xvel_old, amrex::MultiFab &yvel_old, amrex::MultiFab &zvel_old, amrex::MultiFab &xvel_new, amrex::MultiFab &yvel_new, amrex::MultiFab &zvel_new, amrex::MultiFab &source, const amrex::Geometry fine_geom, const amrex::Real dt, const amrex::Real time)
Interface for advancing the data at one level by one "slow" timestep.
int last_check_file_step
Step when we last output a checkpoint file.
Definition REMORA.H:1206
void prestep(int lev, amrex::MultiFab &mf_uold, amrex::MultiFab &mf_vold, amrex::MultiFab &mf_u, amrex::MultiFab &mf_v, amrex::MultiFab *mf_ru, amrex::MultiFab *mf_rv, amrex::MultiFab &S_old, amrex::MultiFab &S_new, amrex::MultiFab &mf_W, amrex::MultiFab &mf_DC, const amrex::MultiFab *mf_z_r, const amrex::MultiFab *mf_z_w, const amrex::MultiFab *mf_h, const amrex::MultiFab *mf_pm, const amrex::MultiFab *mf_pn, const amrex::MultiFab *mf_sustr, const amrex::MultiFab *mf_svstr, const amrex::MultiFab *mf_bustr, const amrex::MultiFab *mf_bvstr, const amrex::MultiFab *mf_msku, const amrex::MultiFab *mf_mskv, const int iic, const int nfirst, const int nnew, int nstp, int nrhs, int N, const amrex::Real dt_lev)
Wrapper function for prestep.
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< std::unique_ptr< amrex::MultiFab > > vec_bvf
Brunt-Vaisala frequency (3D)
Definition REMORA.H:401
virtual void RemakeLevel(int lev, amrex::Real time, const amrex::BoxArray &ba, const amrex::DistributionMapping &dm) override
Remake an existing level using provided BoxArray and DistributionMapping and fill with existing fine ...
amrex::Vector< std::unique_ptr< amrex::MultiFab > > vec_mskv
land/sea mask at y-faces (2D)
Definition REMORA.H:361
void init_masks(int lev, const amrex::BoxArray &ba, const amrex::DistributionMapping &dm)
Allocate MultiFabs for masks.
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()
void WriteGenericPlotfileHeaderWithBathymetry(std::ostream &HeaderFile, int nlevels, const amrex::Vector< amrex::BoxArray > &bArray, const amrex::Vector< std::string > &varnames, amrex::Real time, const amrex::Vector< int > &level_steps, const std::string &versionName, const std::string &levelPrefix, const std::string &mfPrefix) const
write out header data for an AMReX plotfile
static void writeBuildInfo(std::ostream &os)
Write build info to os.
virtual void ErrorEst(int lev, amrex::TagBoxArray &tags, amrex::Real time, int ngrow) override
Tag cells for refinement.
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.
amrex::Real getCPUTime() const
Get CPU time used.
Definition REMORA.H:1424
amrex::Vector< amrex::Vector< amrex::FArrayBox > > bdy_data_ylo
Vectors (over time) of Vector (over variables) of FArrayBoxs for holding Southern boundary data from ...
Definition REMORA.H:1052
void InitData()
Initialize multilevel data.
Definition REMORA.cpp:263
amrex::Vector< int > istep
which step?
Definition REMORA.H:1152
void apply_clim_nudg(const amrex::Box &bx, int ioff, int joff, const amrex::Array4< amrex::Real > &var, const amrex::Array4< amrex::Real const > &var_old, const amrex::Array4< amrex::Real const > &var_clim, const amrex::Array4< amrex::Real const > &clim_coeff, const amrex::Array4< amrex::Real const > &Hz, const amrex::Array4< amrex::Real const > &pm, const amrex::Array4< amrex::Real const > &pn, const amrex::Real dt_lev=amrex::Real(0.0))
Apply climatology nudging.
void uv3dmix(const amrex::Box &xbx, const amrex::Box &ybx, const amrex::Array4< amrex::Real > &u, const amrex::Array4< amrex::Real > &v, const amrex::Array4< amrex::Real const > &uold, const amrex::Array4< amrex::Real const > &vold, const amrex::Array4< amrex::Real > &rufrc, const amrex::Array4< amrex::Real > &rvfrc, const amrex::Array4< amrex::Real const > &visc2_p, const amrex::Array4< amrex::Real const > &visc2_r, const amrex::Array4< amrex::Real const > &Hz, const amrex::Array4< amrex::Real const > &pm, const amrex::Array4< amrex::Real const > &pn, const amrex::Array4< amrex::Real const > &mskp, int nrhs, int nnew, const amrex::Real dt_lev)
Harmonic viscosity.
void WriteCheckpointFile()
write checkpoint file to disk
std::string nc_clim_coeff_file
NetCDF climatology coefficient file.
Definition REMORA.H:1336
void lin_eos(const amrex::Box &bx, const amrex::Array4< amrex::Real const > &state, const amrex::Array4< amrex::Real > &rho, const amrex::Array4< amrex::Real > &rhoA, const amrex::Array4< amrex::Real > &rhoS, const amrex::Array4< amrex::Real > &bvf, const amrex::Array4< amrex::Real const > &Hz, const amrex::Array4< amrex::Real const > &z_w, const amrex::Array4< amrex::Real const > &z_r, const amrex::Array4< amrex::Real const > &h, const amrex::Array4< amrex::Real const > &mskr, const int N)
Calculate density and related quantities from linear equation of state.
void setRecordDataInfo(int i, const std::string &filename)
Definition REMORA.H:1438
void advance_3d_ml(int lev, amrex::Real dt_lev)
3D advance on a single level
const std::string DataLogName(int i) const noexcept
The filename of the ith datalog file.
Definition REMORA.H:1455
void set_analytic_vmix(int lev)
Set vertical mixing coefficients from analytic.
Definition REMORA.cpp:632
void mask_arrays_for_write(int lev, amrex::Real fill_value, amrex::Real fill_where)
Mask data arrays before writing output.
void init_flat_bathymetry(int lev)
Initialize flat bathymetry to value from problo.
Definition REMORA.cpp:672
amrex::Vector< std::unique_ptr< amrex::MultiFab > > vec_rhoS
density perturbation
Definition REMORA.H:397
void fill_from_bdyfiles(amrex::MultiFab &mf_to_fill, const amrex::MultiFab &mf_mask, const amrex::Real time, const int bccomp, const int bdy_var_type, const int icomp_to_fill, const int icomp_calc=0, const amrex::MultiFab &mf_calc=amrex::MultiFab(), const amrex::Real=amrex::Real(0.0))
Fill boundary data from netcdf file.
void rhs_uv_2d(int lev, const amrex::Box &xbx, const amrex::Box &ybx, const amrex::Array4< amrex::Real const > &uold, const amrex::Array4< amrex::Real const > &vold, const amrex::Array4< amrex::Real > &ru, const amrex::Array4< amrex::Real > &rv, const amrex::Array4< amrex::Real const > &Duon, const amrex::Array4< amrex::Real const > &Dvom, const int nrhs)
RHS terms for 2D momentum.
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
void update_mskp(int lev)
Set psi-point mask to be consistent with rho-point mask.
amrex::Vector< std::unique_ptr< amrex::MultiFab > > vec_svstr
Surface stress in the v direction.
Definition REMORA.H:295
amrex::Vector< std::unique_ptr< amrex::MultiFab > > vec_Huon
u-volume flux (3D)
Definition REMORA.H:240
REMORA()
Definition REMORA.cpp:64
void FillBdyCCVels(int lev, amrex::MultiFab &mf_cc_vel)
Fill the physical boundary conditions for cell-centered velocity (diagnostic only)
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 calc_stretch_coeffs()
calculate vertical stretch coefficients
void bulk_fluxes(int lev, amrex::MultiFab *mf_cons, amrex::MultiFab *mf_uwind, amrex::MultiFab *mf_vwind, amrex::MultiFab *mf_evap, amrex::MultiFab *mf_sustr, amrex::MultiFab *mf_svstr, amrex::MultiFab *mf_stflux, amrex::MultiFab *mf_lrflx, amrex::MultiFab *mf_lhflx, amrex::MultiFab *mf_shflx, const int N)
Calculate bulk temperature, salinity, wind fluxes.
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.
static void print_usage(MPI_Comm, std::ostream &)
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
void init_stuff(int lev, const amrex::BoxArray &ba, const amrex::DistributionMapping &dm)
Allocate MultiFabs for state and evolution variables.
TimeInterpolatedData GetDataAtTime(int lev, amrex::Real time)
utility to copy in data from old and/or new state into another multifab
void Advance(int lev, amrex::Real time, amrex::Real dt_lev, int iteration, int ncycle)
advance a single level for a single time step
amrex::Vector< std::unique_ptr< amrex::MultiFab > > vec_rhoA
vertically-averaged density
Definition REMORA.H:399
amrex::Vector< std::unique_ptr< amrex::MultiFab > > vec_DV_avg1
time average of barotropic y velocity flux
Definition REMORA.H:340
amrex::Vector< REMORAFillPatcher > FPr_c
Vector over levels of FillPatchers for scalars.
Definition REMORA.H:1109
void rhs_t_3d(int lev, const amrex::Box &bx, const amrex::Array4< amrex::Real > &t, const amrex::Array4< amrex::Real const > &tempstore, const amrex::Array4< amrex::Real const > &Huon, const amrex::Array4< amrex::Real const > &Hvom, const amrex::Array4< amrex::Real const > &Hz, const amrex::Array4< amrex::Real const > &pn, const amrex::Array4< amrex::Real const > &pm, const amrex::Array4< amrex::Real const > &W, const amrex::Array4< amrex::Real > &FC, const amrex::Array4< amrex::Real const > &mskr, const amrex::Array4< amrex::Real const > &msku, const amrex::Array4< amrex::Real const > &mskv, const amrex::Array4< int const > &river_pos, const amrex::Array4< amrex::Real const > &river_source, int nrhs, int nnew, int N, const amrex::Real dt_lev)
RHS terms for tracer.
AMREX_GPU_HOST_DEVICE static AMREX_FORCE_INLINE amrex::Real bulk_psiu(amrex::Real ZoL)
Evaluate stability function psi for wind speed.
Definition REMORA.H:1474
static void print_banner(MPI_Comm, std::ostream &)
std::string clim_v_time_varname
Name of time field for v climatology data.
Definition REMORA.H:1347
void scale_rhs_vars_inv()
Scale RHS momentum variables by cell area, needed after FillPatch to different levels.
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
amrex::Gpu::DeviceVector< amrex::Real > Cs_r
Stretching coefficients at rho points.
Definition REMORA.H:282
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 resize_stuff(int lev)
Resize variable containers to accommodate data on levels 0 to max_lev.
amrex::Vector< std::unique_ptr< amrex::iMultiFab > > vec_river_position
iMultiFab for river positions; contents are indices of rivers
Definition REMORA.H:1094
amrex::Vector< std::unique_ptr< amrex::MultiFab > > vec_Akk
Turbulent kinetic energy vertical diffusion coefficient.
Definition REMORA.H:420
void set_masks(int lev)
Initialize land-sea masks from file or analytic.
Definition REMORA.cpp:720
amrex::Vector< std::unique_ptr< amrex::MultiFab > > vec_rdrag2
Quadratic drag coefficient [unitless], defined at rho points.
Definition REMORA.H:326
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
static void print_tpls(std::ostream &)
void vert_visc_3d(const amrex::Box &bx, const int ioff, const int joff, const amrex::Array4< amrex::Real > &phi, const amrex::Array4< amrex::Real const > &Hz, const amrex::Array4< amrex::Real > &Hzk, const amrex::Array4< amrex::Real > &AK, const amrex::Array4< amrex::Real const > &Akv, const amrex::Array4< amrex::Real > &BC, const amrex::Array4< amrex::Real > &DC, const amrex::Array4< amrex::Real > &FC, const amrex::Array4< amrex::Real > &CF, const int nnew, const int N, const amrex::Real dt_lev)
Calculate effects of vertical viscosity or diffusivity.
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
amrex::Vector< amrex::Real > vec_weight1
Weights for calculating avg1 in 2D advance.
Definition REMORA.H:408
amrex::Vector< std::unique_ptr< amrex::MultiFab > > vec_xp
x_grid on psi-points (2D)
Definition REMORA.H:389
amrex::Vector< amrex::Vector< amrex::FArrayBox > > bdy_data_xhi
Vectors (over time) of Vector (over variables) of FArrayBoxs for holding Eastern boundary data from f...
Definition REMORA.H:1050
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::Gpu::DeviceVector< int > river_direction
Vector over rivers of river direction: 0: u-face; 1: v-face; 2: w-face.
Definition REMORA.H:1096
amrex::Vector< std::unique_ptr< amrex::MultiFab > > vec_vbar
barotropic y velocity (2D)
Definition REMORA.H:352
void ReadNCMultiFab(amrex::FabArray< amrex::FArrayBox > &fab, const std::string &name, int coordinatorProc=amrex::ParallelDescriptor::IOProcessorNumber(), int allow_empty_mf=0)
Read MultiFab in NetCDF format.
amrex::Vector< std::unique_ptr< amrex::MultiFab > > vec_DU_avg1
time average of barotropic x velocity flux (2D)
Definition REMORA.H:336
amrex::Gpu::DeviceVector< amrex::Real > Cs_w
Stretching coefficients at w points.
Definition REMORA.H:284
void set_zeta_to_Ztavg(int lev)
Set zeta components to be equal to time-averaged Zt_avg1.
void advance_2d(int lev, amrex::MultiFab const *mf_rhoS, amrex::MultiFab const *mf_rhoA, amrex::MultiFab *mf_ru2d, amrex::MultiFab *mf_rv2d, amrex::MultiFab *mf_rufrc, amrex::MultiFab *mf_rvfrc, amrex::MultiFab *mf_Zt_avg1, std::unique_ptr< amrex::MultiFab > &mf_DU_avg1, std::unique_ptr< amrex::MultiFab > &mf_DU_avg2, std::unique_ptr< amrex::MultiFab > &mf_DV_avg1, std::unique_ptr< amrex::MultiFab > &mf_DV_avg2, std::unique_ptr< amrex::MultiFab > &mf_rubar, std::unique_ptr< amrex::MultiFab > &mf_rvbar, std::unique_ptr< amrex::MultiFab > &mf_rzeta, std::unique_ptr< amrex::MultiFab > &mf_ubar, std::unique_ptr< amrex::MultiFab > &mf_vbar, amrex::MultiFab *mf_zeta, amrex::MultiFab const *mf_h, amrex::MultiFab const *mf_pm, amrex::MultiFab const *mf_pn, amrex::MultiFab const *mf_fcor, amrex::MultiFab const *mf_visc2_p, amrex::MultiFab const *mf_visc2_r, amrex::MultiFab const *mf_mskr, amrex::MultiFab const *mf_msku, amrex::MultiFab const *mf_mskv, amrex::MultiFab const *mf_mskp, amrex::Real dtfast_lev, bool predictor_2d_step, bool first_2d_step, int my_iif, int &next_indx1)
Perform a 2D predictor (predictor_2d_step=True) or corrector (predictor_2d_step=False) step.
int plot_file_on_restart
Whether to output a plotfile on restart from checkpoint.
Definition REMORA.H:1210
amrex::Vector< std::unique_ptr< amrex::MultiFab > > vec_alpha
Thermal expansion coefficient (3D)
Definition REMORA.H:403
void set_2darrays(int lev)
Set 2D momentum arrays from 3D momentum.
void FillCoarsePatchMap(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(), amrex::Interpolater *mapper=nullptr)
fill an entire multifab by interpolating from the coarser level, explicitly specifying interpolator t...
amrex::Array< amrex::Array< amrex::Real, AMREX_SPACEDIM *2 >, AMREX_SPACEDIM+NCONS+8 > m_bc_extdir_vals
Array holding the Dirichlet values at walls which need them.
Definition REMORA.H:1186
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
void WriteMultiLevelPlotfileWithBathymetry(const std::string &plotfilename, int nlevels, const amrex::Vector< const amrex::MultiFab * > &mf, const amrex::Vector< const amrex::MultiFab * > &mf_nd, const amrex::Vector< const amrex::MultiFab * > &mf_u, const amrex::Vector< const amrex::MultiFab * > &mf_v, const amrex::Vector< const amrex::MultiFab * > &mf_w, const amrex::Vector< std::string > &varnames, amrex::Real time, const amrex::Vector< int > &level_steps, const std::string &versionName="HyperCLaw-V1.1", const std::string &levelPrefix="Level_", const std::string &mfPrefix="Cell", const amrex::Vector< std::string > &extra_dirs=amrex::Vector< std::string >()) const
write out particular data to an AMReX plotfile
void InitializeLevelFromData(int lev, const amrex::MultiFab &initial_data)
Initialize the new-time data at a level from the initial_data MultiFab.
amrex::Vector< amrex::MultiFab * > cons_old
multilevel data container for last step's scalar data: temperature, salinity, passive scalar
Definition REMORA.H:214
amrex::Vector< std::unique_ptr< amrex::MultiFab > > vec_bustr
Bottom stress in the u direction.
Definition REMORA.H:331
AMREX_FORCE_INLINE int NumDataLogs() noexcept
Definition REMORA.H:1412
amrex::Real volWgtSumMF(int lev, const amrex::MultiFab &mf, int comp, bool local, bool finemask)
Perform the volume-weighted sum.
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
void init_clim_nudg_coeff_from_netcdf(int lev)
Climatology nudging coefficient initialization from NetCDF file.
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)
std::string PlotFileName(int lev) const
get plotfile name
amrex::Vector< std::unique_ptr< amrex::MultiFab > > vec_DV_avg2
correct time average of barotropic y velocity flux for coupling (2D)
Definition REMORA.H:342
void set_hmixcoef(int lev)
Initialize horizontal mixing coefficients.
Definition REMORA.cpp:646
amrex::Array< std::string, 2 *AMREX_SPACEDIM > domain_bc_type
Array of strings describing domain boundary conditions.
Definition REMORA.H:1183
amrex::Real estTimeStep(int lev) const
compute dt from CFL considerations
amrex::Vector< std::unique_ptr< amrex::MultiFab > > vec_yu
y_grid on u-points (2D)
Definition REMORA.H:381
amrex::Vector< std::unique_ptr< amrex::MultiFab > > vec_bvstr
Bottom stress in the v direction.
Definition REMORA.H:333
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
amrex::Vector< std::unique_ptr< amrex::MultiFab > > vec_rzeta
free surface height for the RHS (2D)
Definition REMORA.H:348
amrex::Vector< std::unique_ptr< amrex::MultiFab > > vec_z_phys_nd
z coordinates at psi points (cell nodes)
Definition REMORA.H:287
AMREX_FORCE_INLINE amrex::YAFluxRegister * getAdvFluxReg(int lev)
Get flux register.
Definition REMORA.H:1397
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
amrex::Vector< std::unique_ptr< amrex::MultiFab > > vec_xu
x_grid on u-points (2D)
Definition REMORA.H:379
void timeStepML(amrex::Real time, int iteration)
advance all levels by dt, loops over finer levels
amrex::GpuArray< amrex::GpuArray< bool, AMREX_SPACEDIM *2 >, BdyVars::NumTypes+1 > phys_bc_need_data
These are flags that indicate whether we need to read in boundary data from file.
Definition REMORA.H:1192
amrex::Vector< amrex::Vector< std::unique_ptr< amrex::MultiFab > > > vec_nudg_coeff
Climatology nudging coefficients.
Definition REMORA.H:425
void set_weights(int lev)
Set weights for averaging 3D variables to 2D.
void prestep_t_advection(int lev, const amrex::Box &tbx, const amrex::Box &gbx, const amrex::Array4< amrex::Real > &tempold, const amrex::Array4< amrex::Real > &tempcache, const amrex::Array4< amrex::Real > &Hz, const amrex::Array4< amrex::Real > &Huon, const amrex::Array4< amrex::Real > &Hvom, const amrex::Array4< amrex::Real > &W, const amrex::Array4< amrex::Real > &DC, const amrex::Array4< amrex::Real > &FC, const amrex::Array4< amrex::Real > &sstore, const amrex::Array4< amrex::Real const > &z_w, const amrex::Array4< amrex::Real const > &h, const amrex::Array4< amrex::Real const > &pm, const amrex::Array4< amrex::Real const > &pn, const amrex::Array4< amrex::Real const > &msku, const amrex::Array4< amrex::Real const > &mskv, const amrex::Array4< int const > &river_pos, const amrex::Array4< amrex::Real const > &river_source, int iic, int ntfirst, int nrhs, int N, const amrex::Real dt_lev)
Prestep advection calculations for the tracers.
AMREX_FORCE_INLINE std::ostream & DataLog(int i)
Helper function for IO stream.
Definition REMORA.H:1405
amrex::Vector< std::unique_ptr< amrex::MultiFab > > vec_pn
horizontal scaling factor: 1 / dy (2D)
Definition REMORA.H:368
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
virtual void MakeNewLevelFromScratch(int lev, amrex::Real time, const amrex::BoxArray &ba, const amrex::DistributionMapping &dm) override
Make a new level from scratch using provided BoxArray and DistributionMapping. Only used during initi...
amrex::Vector< std::unique_ptr< amrex::MultiFab > > vec_stflux
Surface tracer flux; input arrays.
Definition REMORA.H:312
amrex::Vector< amrex::Vector< amrex::FArrayBox > > bdy_data_xlo
Vectors (over time) of Vector (over variables) of FArrayBoxs for holding Western boundary data from f...
Definition REMORA.H:1048
static int verbose
Verbosity level of output.
Definition REMORA.H:1306
void setup_step(int lev, amrex::Real time, amrex::Real dt_lev)
Set everything up for a step on a level.
amrex::Vector< std::unique_ptr< amrex::MultiFab > > vec_rdrag
Linear drag coefficient [m/s], defined at rho points.
Definition REMORA.H:324
void nonlin_eos(const amrex::Box &bx, const amrex::Array4< amrex::Real const > &state, const amrex::Array4< amrex::Real > &rho, const amrex::Array4< amrex::Real > &rhoA, const amrex::Array4< amrex::Real > &rhoS, const amrex::Array4< amrex::Real > &bvf, const amrex::Array4< amrex::Real > &alpha, const amrex::Array4< amrex::Real > &beta, const amrex::Array4< amrex::Real const > &Hz, const amrex::Array4< amrex::Real const > &z_w, const amrex::Array4< amrex::Real const > &z_r, const amrex::Array4< amrex::Real const > &h, const amrex::Array4< amrex::Real const > &mskr, const int N)
Calculate density and related quantities from nonlinear equation of state.
std::string plot_file_name
Plotfile prefix.
Definition REMORA.H:1248
amrex::Vector< std::unique_ptr< amrex::MultiFab > > vec_Zt_avg1
Average of the free surface, zeta (2D)
Definition REMORA.H:290
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
amrex::Vector< std::unique_ptr< amrex::MultiFab > > vec_btflx
Bottom tracer flux; working arrays.
Definition REMORA.H:314
int cf_width
Nudging width at coarse-fine interface.
Definition REMORA.H:1104
AMREX_GPU_HOST_DEVICE static AMREX_FORCE_INLINE amrex::Real bulk_psit(amrex::Real ZoL)
Evaluate stability function psi for moisture and heat.
Definition REMORA.H:1506
static amrex::Vector< amrex::Vector< std::string > > nc_init_file
NetCDF initialization file.
Definition REMORA.H:52
amrex::Vector< std::unique_ptr< amrex::MultiFab > > vec_beta
Saline contraction coefficient (3D)
Definition REMORA.H:405
int last_plot_file_step
Step when we last output a plotfile.
Definition REMORA.H:1201
const amrex::Vector< std::string > derived_names
Names of derived fields for plotfiles.
Definition REMORA.H:1276
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
void InitializeFromFile()
Read the file passed to remora.restart and use it as an initial condition for the current simulation.
void coriolis(const amrex::Box &xbx, const amrex::Box &ybx, const amrex::Array4< amrex::Real const > &uold, const amrex::Array4< amrex::Real const > &vold, const amrex::Array4< amrex::Real > &ru, const amrex::Array4< amrex::Real > &rv, const amrex::Array4< amrex::Real const > &Hz, const amrex::Array4< amrex::Real const > &fomn, int nrhs, int nr)
Calculate Coriolis terms.
amrex::Gpu::DeviceVector< amrex::BCRec > domain_bcs_type_d
GPU vector (over BCVars) of BCRecs.
Definition REMORA.H:1180
virtual ~REMORA()
Definition REMORA.cpp:139
amrex::Vector< std::unique_ptr< amrex::MultiFab > > vec_z_w
z coordinates at w points (faces between z-cells)
Definition REMORA.H:270
void restart()
Definition REMORA.cpp:469
void convert_inv_days_to_inv_s(amrex::MultiFab *)
Convert data in a multifab from inverse days to inverse seconds.
amrex::Vector< std::unique_ptr< amrex::MultiFab > > vec_yr
y_grid on rho points (2D)
Definition REMORA.H:376
amrex::Vector< std::unique_ptr< amrex::MultiFab > > vec_Hvom
v-volume flux (3D)
Definition REMORA.H:242
amrex::Vector< std::unique_ptr< amrex::MultiFab > > vec_Akp
Turbulent length scale vertical diffusion coefficient.
Definition REMORA.H:422
void gls_corrector(int lev, amrex::MultiFab *mf_gls, amrex::MultiFab *mf_tke, amrex::MultiFab &mf_W, amrex::MultiFab *mf_Akv, amrex::MultiFab *mf_Akt, amrex::MultiFab *mf_Akk, amrex::MultiFab *mf_Akp, amrex::MultiFab *mf_mskr, amrex::MultiFab *mf_msku, amrex::MultiFab *mf_mskv, const int nstp, const int nnew, const int N, const amrex::Real dt_lev)
Corrector step for GLS calculation.
amrex::Vector< std::unique_ptr< amrex::MultiFab > > vec_diff2
Harmonic diffusivity for temperature / salinity.
Definition REMORA.H:264
@ CellConservativeQuartic
Definition REMORA.H:72
@ CellBilinear
Definition REMORA.H:69
@ PCInterp
Definition REMORA.H:66
@ NodeBilinear
Definition REMORA.H:67
@ CellConservativeLinear
Definition REMORA.H:68
@ CellQuadratic
Definition REMORA.H:70
@ CellConservativeProtected
Definition REMORA.H:71