REMORA
Regional Modeling of Oceans Refined Adaptively
Loading...
Searching...
No Matches
REMORA_DataStruct.H
Go to the documentation of this file.
1#ifndef _REMORA_DATA_STRUCT_H_
2#define _REMORA_DATA_STRUCT_H_
3
4#include <string>
5#include <iostream>
6
7#include <AMReX_ParmParse.H>
8#include <AMReX_Print.H>
9#include <AMReX_Gpu.H>
10
11#include <REMORA_Constants.H>
12#include "REMORA_IndexDefines.H"
13
14/** \brief Type of coupling between levels in AMR */
15enum struct CouplingType {
17};
18
19/** \brief Coordinates */
20enum class Coord {
21 x, y, z
22};
23
24/** \brief Horizontal advection schemes */
28
29/** \brief Type of initial condition type. Analytic reads from prob.cpp. Netcdf is from file */
30enum class IC_Type {
32};
33
34/** \brief Coriolis factor */
35enum class Cor_Type {
37};
38
39/** \brief plotfile format */
40enum class PlotfileType {
42};
43
44/** \brief vertical mixing type */
45enum class VertMixingType {
47};
48
49/** \brief horizontal viscosity/diffusion type */
50enum class HorizMixingType {
52};
53
54/** \brief stability function for GLS */
58
59/** \brief equation of state */
60enum class EOSType {
62};
63
64/** \brief bottom stress formulation */
68
69/** \brief initialization for pm and pn */
70enum class GridScaleType {
72};
73
74/** \brief surface momentum flux */
75enum class SMFluxType {
77};
78
79/** \brief surface wind */
80enum class WindType {
82};
83
84/** \brief masks */
85enum class MaskType {
87};
88
90 public:
91 /** \brief read in and initialize parameters */
93 {
94 amrex::ParmParse pp(pp_prefix);
95
96 pp.query("flat_bathymetry", flat_bathymetry);
97
98 // Which horizontal advection scheme for tracers
99 static std::string tracer_hadv_string = "upstream3";
100 pp.query("tracer_horizontal_advection_scheme",tracer_hadv_string);
101 if (tracer_hadv_string == "centered4")
103 else if (tracer_hadv_string == "upstream3")
105 else
106 amrex::Error("Advection scheme unknown.");
107
108 // Which horizontal advection scheme
109 static std::string uv_hadv_string = "upstream3";
110 pp.query("uv_horizontal_advection_scheme",uv_hadv_string);
111 if (uv_hadv_string == "upstream3")
113 else if (uv_hadv_string == "centered2")
115 else
116 amrex::Error("UV advection scheme unknown.");
117
118 pp.query("rdrag", rdrag);
119 pp.query("rdrag2", rdrag2);
120 pp.query("Zos", Zos);
121 pp.query("Zob", Zob);
122 pp.query("Cdb_max", Cdb_max);
123 pp.query("Cdb_min", Cdb_min);
124
125 // Include salinity?
126 pp.query("use_salt", use_salt);
127
128 // Include Coriolis forcing?
129 pp.query("use_coriolis", use_coriolis);
130
131 // Include prestep / lagged predictor / corrections
132 pp.query("use_prestep", use_prestep);
133
134 //This affect forcing and some of the mixing terms for velocity
135 pp.query("use_uv3dmix", use_uv3dmix);
136
137 //This accounts for the main 2d loop but may not include coupling and copying properly
138 pp.query("use_barotropic", use_barotropic);
139
140 // Whether to do rivers. By default, rivers are temp and salt sources. Rivers
141 // have to be momentum sources.
142 pp.query("do_rivers", do_rivers);
143
144 pp.query("do_rivers_temp", do_rivers_temp);
145 pp.query("do_rivers_salt", do_rivers_salt);
146 pp.query("do_rivers_scalar", do_rivers_scalar);
147
148 do_rivers_cons.resize(NCONS);
149 // If we aren't doing rivers, set all of these flags to false
153
154 pp.query("init_l1ad_h", init_l1ad_h);
155 pp.query("init_l1ad_T", init_l1ad_T);
156
157 pp.query("init_ana_h", init_ana_h);
158 pp.query("init_ana_T", init_ana_T);
159
160 pp.query("init_l0int_h", init_l0int_h);
161 pp.query("init_l0int_T", init_l0int_T);
162
163 static std::string eos_type_string = "linear";
164 pp.query("eos_type",eos_type_string);
165 if (eos_type_string == "linear" || eos_type_string == "Linear" ||
166 eos_type_string == "lin" || eos_type_string == "Lin") {
168 pp.query("Tcoef",Tcoef);
169 pp.query("Scoef",Scoef);
170 } else if (eos_type_string == "nonlinear" || eos_type_string == "Nonlinear" ||
171 eos_type_string == "non-linear" || eos_type_string == "Non-linear" ||
172 eos_type_string == "nonlin" || eos_type_string == "Nonlin") {
174 } else {
175 amrex::Abort("Dont know this eos_type");
176 }
177 pp.query("R0",R0);
178 pp.query("S0",S0);
179 pp.query("T0",T0);
180 pp.query("rho0", rho0);
181
182 pp.query("bulk_fluxes",bulk_fluxes);
183 if (bulk_fluxes) {
184 do_salt_flux = true;
185 do_temp_flux = true;
186 }
187 pp.query("air_pressure",Pair);
188 pp.query("air_temperature",Tair);
189 pp.query("air_humidity",Hair);
190 pp.query("surface_radiation_flux",srflux);
191 pp.query("cloud",cloud);
192 pp.query("rain",rain);
193 pp.query("blk_ZQ",blk_ZQ);
194 pp.query("blk_ZT",blk_ZT);
195 pp.query("blk_ZW",blk_ZW);
196 pp.query("eminusp",eminusp);
197 pp.query("eminusp_correct_ssh",eminusp_correct_ssh);
198
199 if (eminusp_correct_ssh and !eminusp) {
200 amrex::Abort("If evaporation minus precipitation (E-P) sea surface height correct in is on, E-P must be on as well (remora.eminusp=true)");
201 }
202 if (eminusp and !bulk_fluxes) {
203 amrex::Abort("Evaporation minus precipitation (E-P) requires bulk flux parametrizations (remora.bulk_fluxes=true)");
204 }
205
206 //Read in linear eos parameters
207 //Grid stretching
208 pp.query("theta_s",theta_s);
209 pp.query("theta_b",theta_b);
210 pp.query("tcline",tcline);
211
212 //coriolis factor
213 pp.query("coriolis_f0",coriolis_f0);
214 pp.query("coriolis_beta",coriolis_beta);
215
216 pp.query("Akv_bak",Akv_bak);
217 pp.query("Akt_bak",Akt_bak);
218
219
220 static std::string grid_scale_type_string = "constant";
221 pp.query("grid_scale_type",grid_scale_type_string);
222
223 if (amrex::toLower(grid_scale_type_string) == "constant") {
225 } else if (amrex::toLower(grid_scale_type_string) == "custom") {
226 amrex::Warning("Initialization of grid scale from prob.cpp is now called 'analytic'. 'custom' will be deprecated");
228 } else if (amrex::toLower(grid_scale_type_string) == "analytic") {
230 } else {
231 amrex::Error("Don't know this grid_scale_type");
232 }
233
234 static std::string ic_type_string = "analytic";
235 bool found_ic_bc = pp.query("ic_bc_type", ic_type_string);
236 pp.query("ic_type", ic_type_string);
237
238 if (found_ic_bc) {
239 amrex::Warning("remora.ic_bc_type is now called remora.ic_type, and will eventually be deprecated");
240 }
241
242 if ( amrex::toLower(ic_type_string) == "custom") {
243 amrex::Warning("Problem initialization from prob.cpp is now called 'analytic'. 'custom' will be deprecated");
245 } else if ( amrex::toLower(ic_type_string) == "analytic") {
247 } else if ( amrex::toLower(ic_type_string) == "netcdf") {
249 } else if ( amrex::toLower(ic_type_string) == "real") {
250 amrex::Warning("Problem initialization from NetCDF (remora.ic_type) is now called 'netcdf'. 'real' will be deprecated");
252 } else {
253 amrex::Error("Don't know this ic_type");
254 }
255
256 // Which type of refinement
257 static std::string coupling_type_string = "TwoWay";
258 pp.query("coupling_type",coupling_type_string);
259 if (amrex::toLower(coupling_type_string) == "twoway" ||
260 amrex::toLower(coupling_type_string) == "two_way") {
262 } else if (amrex::toLower(coupling_type_string) == "oneway" ||
263 amrex::toLower(coupling_type_string) == "one_way") {
265 } else {
266 amrex::Abort("Dont know this coupling_type");
267 }
268
269 // Which type of coriolis forcing
270 if (use_coriolis) {
271 static std::string coriolis_type_string = "beta_plane";
272 pp.query("coriolis_type",coriolis_type_string);
273 if ( amrex::toLower(coriolis_type_string) == "custom") {
274 amrex::Warning("Coriolis initialization from prob.cpp is now called 'analytic'. 'custom' will be deprecated");
276 } else if ( amrex::toLower(coriolis_type_string) == "analytic") {
278 } else if ((amrex::toLower(coriolis_type_string) == "beta_plane") ||
279 (amrex::toLower(coriolis_type_string) == "betaplane")) {
281 } else if ( (amrex::toLower(coriolis_type_string) == "netcdf")) {
283 } else if ( (amrex::toLower(coriolis_type_string) == "real")) {
284 amrex::Warning("Coriolis initialization from NetCDF is now called 'netcdf'. 'real' will be deprecated");
286 } else {
287 amrex::Abort("Don't know this coriolis_type");
288 }
289 }
290
291 static std::string smflux_type_string = "analytic";
292 int smflux_specified = pp.query("smflux_type",smflux_type_string);
293 if ( amrex::toLower(smflux_type_string) == "custom") {
294 amrex::Warning("Surface momentum flux initialization from prob.cpp is now called 'analytic'. 'custom' will be deprecated");
296 } else if ( amrex::toLower(smflux_type_string) == "analytic") {
298 } else if ( amrex::toLower(smflux_type_string) == "netcdf") {
300 } else {
301 amrex::Abort("Don't know this smflux_type");
302 }
303
304 static std::string wind_type_string = "analytic";
305 int wind_specified = pp.query("wind_type",wind_type_string);
306 if ( amrex::toLower(wind_type_string) == "custom") {
307 amrex::Warning("Surface wind initialization from prob.cpp is now called 'analytic'. 'custom' will be deprecated");
309 } else if ( amrex::toLower(wind_type_string) == "analytic") {
311 } else if ( amrex::toLower(wind_type_string) == "netcdf") {
313 } else {
314 amrex::Abort("Don't know this smflux_type");
315 }
316
317 if (wind_specified && smflux_specified) {
318 amrex::Abort("Cannot specify both wind and surface momentum flux");
319 }
320
321 static std::string mask_type_string = "none";
322 // If initial condition type is netcdf, default to netcdf masks
323 if (ic_type == IC_Type::netcdf) {
324 mask_type_string = "netcdf";
325 }
326 pp.query("mask_type", mask_type_string);
327 if (amrex::toLower(mask_type_string) == "none") {
329 } else if (amrex::toLower(mask_type_string) == "analytic") {
331 } else if (amrex::toLower(mask_type_string) == "netcdf") {
333 } else {
334 amrex::Abort("Don't know this mask_type");
335 }
336
337 static std::string bottom_stress_type_string = "linear";
338 pp.query("bottom_stress_type", bottom_stress_type_string);
339 if (amrex::toLower(bottom_stress_type_string) == "linear") {
341 } else if (amrex::toLower(bottom_stress_type_string) == "quadratic") {
343 } else if (amrex::toLower(bottom_stress_type_string) == "logarithmic") {
345 } else {
346 amrex::Abort("Don't know this bottom_stress_type");
347 }
348
349 amrex::Real tnu2_salt = amrex::Real(0.0);
350 amrex::Real tnu2_temp = amrex::Real(0.0);
351 amrex::Real tnu2_scalar = amrex::Real(0.0);
352 static std::string horiz_mixing_type_string = "analytic";
353 pp.query("horizontal_mixing_type", horiz_mixing_type_string);
354 if (amrex::toLower(horiz_mixing_type_string) == "analytical" ||
355 amrex::toLower(horiz_mixing_type_string) == "analytic") {
357 } else if (amrex::toLower(horiz_mixing_type_string) == "constant") {
359 } else {
360 amrex::Abort("Don't know this horizontal mixing type");
361 }
362 pp.query("visc2",visc2);
363 pp.query("tnu2_salt",tnu2_salt);
364 pp.query("tnu2_temp",tnu2_temp);
365 pp.query("tnu2_scalar",tnu2_scalar);
366 tnu2.resize(NCONS);
367 if (use_salt) {
368 tnu2[0] = tnu2_temp;
369 tnu2[1] = tnu2_salt;
370 tnu2[2] = tnu2_scalar;
371 } else {
372 tnu2[0] = tnu2_temp;
373 tnu2[1] = tnu2_scalar;
374 }
375
376 pp.query("Akk_bak", Akk_bak);
377 pp.query("Akp_bak", Akp_bak);
378 static std::string vert_mixing_type_string = "analytic";
379 static std::string gls_stability_type_string = "Canuto_A";
380 pp.query("vertical_mixing_type", vert_mixing_type_string);
381 pp.query("gls_stability_type", gls_stability_type_string);
382 if (amrex::toLower(vert_mixing_type_string) == "analytical" ||
383 amrex::toLower(vert_mixing_type_string) == "analytic") {
385 } else if (amrex::toLower(vert_mixing_type_string) == "gls") {
387 if (amrex::toLower(gls_stability_type_string) == "canuto_a") {
389 }
390 else if (amrex::toLower(gls_stability_type_string) == "canuto_b") {
392 }
393 else if (amrex::toLower(gls_stability_type_string) == "galperin") {
395 }
396 else {
397 amrex::Abort("Don't know this GLS stability type");
398 }
399 } else {
400 amrex::Abort("Don't know this vertical mixing type");
401 }
402 // Read in GLS params
404 pp.query("gls_P", gls_p);
405 pp.query("gls_M", gls_m);
406 pp.query("gls_N", gls_n);
407 pp.query("gls_Kmin", gls_Kmin);
408 pp.query("gls_Pmin", gls_Pmin);
409
410 pp.query("gls_cmu0", gls_cmu0);
411 pp.query("gls_c1", gls_c1);
412 pp.query("gls_c2", gls_c2);
413 pp.query("gls_c3m", gls_c3m);
414 pp.query("gls_c3p", gls_c3p);
415 pp.query("gls_sigk", gls_sigk);
416 pp.query("gls_sigp", gls_sigp);
418 gls_Gh0 = amrex::Real(0.0329); // 0.0329 GOTM, 0.0673 Burchard
419 gls_Ghcri = amrex::Real(0.03);
420 gls_L1 = amrex::Real(0.107);
421 gls_L2 = amrex::Real(0.0032);
422 gls_L3 = amrex::Real(0.0864);
423 gls_L4 = amrex::Real(0.12);
424 gls_L5 = amrex::Real(11.9);
425 gls_L6 = amrex::Real(0.4);
426 gls_L7 = amrex::Real(0.0);
427 gls_L8 = amrex::Real(0.48);
429 gls_Gh0 = amrex::Real(0.0444); // 0.044 GOTM, 0.0673 Burchard
430 gls_Ghcri = amrex::Real(0.0414);
431 gls_L1 = amrex::Real(0.127);
432 gls_L2 = amrex::Real(0.00336);
433 gls_L3 = amrex::Real(0.0906);
434 gls_L4 = amrex::Real(0.101);
435 gls_L5 = amrex::Real(11.2);
436 gls_L6 = amrex::Real(0.4);
437 gls_L7 = amrex::Real(0.0);
438 gls_L8 = amrex::Real(0.318);
439 } else {
440 gls_Gh0 = amrex::Real(0.028);
441 gls_Ghcri = amrex::Real(0.02);
442 }
443 }
444
445 // Read and compute inverse nudging coeffs from inputs given in days,
446 // and store in a vector corresponding to BdyVars enum
447 amrex::Real tnudg = amrex::Real(0.0);
448 amrex::Real znudg = amrex::Real(0.0);
449 amrex::Real m2nudg = amrex::Real(0.0);
450 amrex::Real m3nudg = amrex::Real(0.0);
451 pp.query("tnudg",tnudg);
452 pp.query("znudg",znudg);
453 pp.query("m2nudg",m2nudg);
454 pp.query("m3nudg",m3nudg);
455 pp.query("obcfac",obcfac);
456
458 nudg_coeff[BdyVars::u ] = (m3nudg > 0.0) ? amrex::Real(1.0) / (m3nudg * amrex::Real(86400.0)) : 0.0;//BdyVars::u
459 nudg_coeff[BdyVars::v ] = (m3nudg > 0.0) ? amrex::Real(1.0) / (m3nudg * amrex::Real(86400.0)) : 0.0;//BdyVars::v
460 nudg_coeff[BdyVars::t ] = ( tnudg > 0.0) ? amrex::Real(1.0) / ( tnudg * amrex::Real(86400.0)) : 0.0;//BdyVars::t
461 nudg_coeff[BdyVars::s ] = ( tnudg > 0.0) ? amrex::Real(1.0) / ( tnudg * amrex::Real(86400.0)) : 0.0;//BdyVars::s
462 nudg_coeff[BdyVars::ubar ] = (m2nudg > 0.0) ? amrex::Real(1.0) / (m2nudg * amrex::Real(86400.0)) : 0.0;//BdyVars::ubar
463 nudg_coeff[BdyVars::vbar ] = (m2nudg > 0.0) ? amrex::Real(1.0) / (m2nudg * amrex::Real(86400.0)) : 0.0;//BdyVars::vbar
464 nudg_coeff[BdyVars::zeta ] = ( znudg > 0.0) ? amrex::Real(1.0) / ( znudg * amrex::Real(86400.0)) : 0.0;//BdyVars::zeta
465
466 pp.query("do_m3_clim_nudg", do_m3_clim_nudg);
467 pp.query("do_m2_clim_nudg", do_m2_clim_nudg);
468 pp.query("do_temp_clim_nudg", do_temp_clim_nudg);
469 pp.query("do_salt_clim_nudg", do_salt_clim_nudg);
470
472 do_any_clim_nudg = true;
473 }
474#ifndef REMORA_USE_NETCDF
475 if (do_any_clim_nudg) {
476 amrex::Abort("Climatology nudging requires building with NetCDF");
477 }
478#endif
479 }
480
481 void display()
482 {
483 amrex::Print() << "SOLVER CHOICE: " << std::endl;
484 amrex::Print() << "use_salt : " << use_salt << std::endl;
485 amrex::Print() << "use_coriolis : " << use_coriolis << std::endl;
486 amrex::Print() << "use_prestep : " << use_prestep << std::endl;
487 amrex::Print() << "use_uv3dmix : " << use_uv3dmix << std::endl;
488 amrex::Print() << "use_barotropic : " << use_barotropic << std::endl;
489 amrex::Print() << "flat_bathymetry : " << flat_bathymetry << std::endl;
490 amrex::Print() << "spatial_order : " << spatial_order << std::endl;
491
492 if (ic_type == IC_Type::analytic) {
493 amrex::Print() << "Using analytic initial onditions" << std::endl;
494 }
495 else if (ic_type == IC_Type::netcdf) {
496 amrex::Print() << "Using NetCDF initial conditions" << std::endl;
497 }
498
500 amrex::Print() << "Horizontal advection scheme for tracers: " << "Centered 4" << std::endl;
501 }
503 amrex::Print() << "Horizontal advection scheme for tracers: " << "Upstream 3" << std::endl;
504 }
505 else {
506 amrex::Error("Invalid horizontal advection scheme for tracers.");
507 }
508
510 amrex::Print() << "Horizontal advection scheme for momenta: " << "Centered 2" << std::endl;
511 }
513 amrex::Print() << "Horizontal advection scheme for momenta: " << "Upstream 3" << std::endl;
514 }
515 else {
516 amrex::Error("Invalid horizontal advection scheme for momenta.");
517 }
518
520 amrex::Print() << "Using two-way coupling " << std::endl;
521 } else if (coupling_type == CouplingType::one_way) {
522 amrex::Print() << "Using one-way coupling " << std::endl;
523 }
524
525 if (use_coriolis) {
527 amrex::Print() << "Using analytic coriolis forcing " << std::endl;
528 } else if (coriolis_type == Cor_Type::beta_plane) {
529 amrex::Print() << "Using beta plane coriolis forcing " << std::endl;
530 } else if (coriolis_type == Cor_Type::netcdf) {
531 amrex::Print() << "Using coriolis forcing loaded from file " << std::endl;
532 }
533 }
534 }
535
536 // Default prefix
537 std::string pp_prefix {"remora"};
538
539 bool flat_bathymetry = false;
540
541 bool use_salt = true;
542
543 // Specify what additional physics/forcing modules we use
544 bool use_coriolis = false;
545
546 // Specify whether terms are used for debugging purposes
547 bool use_prestep = true;
548 bool use_uv3dmix = true;
549 bool use_baroclinic = true;
550 bool use_barotropic = true;
551
552 bool bulk_fluxes = false;
553 bool do_temp_flux = false;
554 bool do_salt_flux = false;
555
556 bool do_rivers = false;
557 bool do_rivers_temp = true;
558 bool do_rivers_salt = true;
559 bool do_rivers_scalar = false;
560 amrex::Vector<int> do_rivers_cons;
561
562 bool init_l1ad_h = false;
563 bool init_l1ad_T = false;
564
565 bool init_ana_h = false;
566 bool init_ana_T = false;
567
568 bool init_l0int_h = true;
569 bool init_l0int_T = true;
570
573
574 // Coupling options are "OneWay" or "TwoWay"
576
577 // IC and BC Type: "analytic" or "netcdf"
579
580 // Coriolis forcing type
582
583 // Surface momentum flux type
585
586 // Surface wind speed type
588
589 // EOS type
591
592 // Bottom stress type
594
595 // Land/sea mask type
597
598 // Mixing type and parameters
602
603 // Type for grid scale (pm and pn)
605
606 // Stretching and depth parameters which may need to be read from inputs
607 amrex::Real theta_s = amrex::Real(3.0);
608 amrex::Real theta_b = amrex::Real(0.0);
609 amrex::Real tcline = amrex::Real(150.0);
610
611 // Linear drag coefficient [m/s]
612 amrex::Real rdrag = amrex::Real(3e-4);
613 // Quadratic drag coefficient [dimensionless]
614 amrex::Real rdrag2 = amrex::Real(3e-3);
615
616 // Momentum stress scales [m]
617 amrex::Real Zob = amrex::Real(2e-2);
618 amrex::Real Zos = amrex::Real(2e-2);
619
620 amrex::Real Cdb_max = amrex::Real(0.5);
621 amrex::Real Cdb_min = amrex::Real(1e-6);
622
623 // Linear equation of state parameters
624 amrex::Real R0 = amrex::Real(1028); // background density value (Kg/m3) used in Linear Equation of State
625 amrex::Real S0 = amrex::Real(35.0); // background salinity (nondimensional) constant
626 amrex::Real T0 = amrex::Real(5.0); // background potential temperature (Celsius) constant
627 amrex::Real Tcoef = amrex::Real(1.7e-4); // linear equation of state parameter (1/Celsius)
628 amrex::Real Scoef = amrex::Real(0.0); // linear equation of state parameter (nondimensional)
629 amrex::Real rho0 = amrex::Real(1025.0); // Mean density (Kg/m3) used when Boussinesq approx is inferred
630
631 // Coriolis forcing
632 amrex::Real coriolis_f0 = amrex::Real(0.0); // f-plane constant (1/s)
633 amrex::Real coriolis_beta = amrex::Real(0.0); // beta-plane constant (1/s/m)
634
635 // Air pressure
636 amrex::Real Pair = amrex::Real(1013.48);
637 // Air temperature
638 amrex::Real Tair = amrex::Real(23.567);
639 // Relative humidity (air)
640 amrex::Real Hair = amrex::Real(0.776);
641 // Cloud cover fraction (0=clear sky, 1=overcast)
642 amrex::Real cloud = amrex::Real(0.0);
643 // Precipitation rate (kg/m2/s)
644 amrex::Real rain = amrex::Real(0.0);
645 // Height (m) of atmospheric measurements for Bulk fluxes parametrization
646 amrex::Real blk_ZQ = amrex::Real(10.0); // air humidity
647 amrex::Real blk_ZT = amrex::Real(10.0); // air temperature
648 amrex::Real blk_ZW = amrex::Real(10.0); // winds
649
650 bool eminusp = false;
652
653 // Surface radiation flux
654 amrex::Real srflux = amrex::Real(0.0);
655
656 // Spatial discretization
658
659 // Horizontal mixing parameters
660 amrex::Real visc2 = amrex::Real(0.0);
661 amrex::Vector<amrex::Real> tnu2;
662
663 // GLS params
664 amrex::Real gls_p = amrex::Real(3.0);
665 amrex::Real gls_m = amrex::Real(1.5);
666 amrex::Real gls_n = amrex::Real(-1.0);
667 amrex::Real gls_Kmin = amrex::Real(7.6e-6);
668 amrex::Real gls_Pmin = amrex::Real(1.0e-12);
669
670 amrex::Real gls_cmu0 = amrex::Real(0.5477);
671 amrex::Real gls_c1 = amrex::Real(1.44);
672 amrex::Real gls_c2 = amrex::Real(1.92);
673 amrex::Real gls_c3m = amrex::Real(-0.4);
674 amrex::Real gls_c3p = amrex::Real(1.0);
675 amrex::Real gls_sigk = amrex::Real(1.0);
676 amrex::Real gls_sigp = amrex::Real(1.3);
677
678 // Turbulence closure
679 amrex::Real Akk_bak = amrex::Real(5.0e-6);
680 amrex::Real Akp_bak = amrex::Real(5.0e-6);
681 amrex::Real Akv_bak = amrex::Real(5.0e-6);
682 amrex::Real Akt_bak = amrex::Real(1.0e-6); // Note: this is a vector with one component per tracer in ROMS
683
684 // Params for stability functions.
685 amrex::Real gls_Gh0;
686 amrex::Real gls_Ghcri;
687 amrex::Real gls_Ghmin = amrex::Real(-0.28);
688 amrex::Real gls_E2 = amrex::Real(1.33);
689 // Params only for Canuto stability
690 amrex::Real gls_L1;
691 amrex::Real gls_L2;
692 amrex::Real gls_L3;
693 amrex::Real gls_L4;
694 amrex::Real gls_L5;
695 amrex::Real gls_L6;
696 amrex::Real gls_L7;
697 amrex::Real gls_L8;
698
699 // Params for some GLS and also Mellor-Yamada
700 amrex::Real my_A1 = amrex::Real(0.92);
701 amrex::Real my_A2 = amrex::Real(0.74);
702 amrex::Real my_B1 = amrex::Real(16.6);
703 amrex::Real my_B2 = amrex::Real(10.1);
704 amrex::Real my_C1 = amrex::Real(0.08);
705 amrex::Real my_C2 = amrex::Real(0.7);
706 amrex::Real my_C3 = amrex::Real(0.2);
707 amrex::Real my_E1 = amrex::Real(1.8);
708 amrex::Real my_E2 = amrex::Real(1.33);
709 amrex::Real my_Gh0 = amrex::Real(0.0233);
710 amrex::Real my_Sq = amrex::Real(0.2);
711 amrex::Real my_dtfac = amrex::Real(0.05);
712 amrex::Real my_lmax = amrex::Real(0.53);
713 amrex::Real my_qmin = amrex::Real(1.0E-8);
714
715 // Nudging time scales in 1/s
716 amrex::Vector<amrex::Real> nudg_coeff;
717
718 // Factor between passive (outflow) and active (inflow) open boundary
719 // conditions.
720 amrex::Real obcfac = amrex::Real(0.0);
721
722 // Whether to do climatoogy nudging
723 bool do_m2_clim_nudg = false;
724 bool do_m3_clim_nudg = false;
725 bool do_temp_clim_nudg = false;
726 bool do_salt_clim_nudg = false;
727 bool do_any_clim_nudg = false;
728
730};
731
732#endif
BottomStressType
bottom stress formulation
GridScaleType
initialization for pm and pn
SMFluxType
surface momentum flux
AdvectionScheme
Horizontal advection schemes.
PlotfileType
plotfile format
HorizMixingType
horizontal viscosity/diffusion type
Coord
Coordinates.
MaskType
masks
Cor_Type
Coriolis factor.
WindType
surface wind
GLS_StabilityType
stability function for GLS
CouplingType
Type of coupling between levels in AMR.
VertMixingType
vertical mixing type
IC_Type
Type of initial condition type. Analytic reads from prob.cpp. Netcdf is from file.
EOSType
equation of state
#define Temp_comp
#define Scalar_comp
#define Salt_comp
#define NCONS
GLS_StabilityType gls_stability_type
HorizMixingType horiz_mixing_type
amrex::Real Cdb_min
amrex::Real Akv_bak
amrex::Real blk_ZT
amrex::Real cloud
amrex::Real coriolis_beta
amrex::Real gls_sigp
amrex::Vector< amrex::Real > nudg_coeff
amrex::Real my_E2
amrex::Real rdrag2
amrex::Real gls_m
amrex::Real my_lmax
amrex::Real my_Gh0
amrex::Real my_B2
amrex::Real my_A1
amrex::Vector< amrex::Real > tnu2
amrex::Real gls_sigk
amrex::Real gls_L3
amrex::Real coriolis_f0
AdvectionScheme uv_Hadv_scheme
amrex::Vector< int > do_rivers_cons
amrex::Real Tcoef
amrex::Real gls_cmu0
void init_params()
read in and initialize parameters
std::string pp_prefix
amrex::Real gls_L6
amrex::Real gls_L2
AdvectionScheme tracer_Hadv_scheme
amrex::Real Akk_bak
amrex::Real gls_L1
amrex::Real blk_ZW
amrex::Real my_E1
amrex::Real my_C3
amrex::Real gls_n
amrex::Real theta_b
amrex::Real theta_s
amrex::Real Akt_bak
amrex::Real visc2
amrex::Real tcline
amrex::Real gls_c3m
BottomStressType bottom_stress_type
amrex::Real my_A2
amrex::Real gls_Gh0
amrex::Real gls_Ghmin
amrex::Real my_C1
amrex::Real gls_c1
amrex::Real gls_L5
amrex::Real gls_L8
amrex::Real my_dtfac
amrex::Real srflux
SMFluxType smflux_type
amrex::Real my_C2
amrex::Real gls_Kmin
amrex::Real rdrag
VertMixingType vert_mixing_type
amrex::Real gls_c3p
amrex::Real gls_p
amrex::Real gls_L4
amrex::Real my_B1
amrex::Real blk_ZQ
amrex::Real gls_L7
amrex::Real Cdb_max
GridScaleType grid_scale_type
amrex::Real my_qmin
amrex::Real gls_Ghcri
amrex::Real obcfac
amrex::Real gls_Pmin
amrex::Real gls_c2
amrex::Real Scoef
amrex::Real Akp_bak
amrex::Real my_Sq
CouplingType coupling_type
amrex::Real gls_E2