Skip to content

Commit

Permalink
Cleanup and port metgrid-related code from dirty old branches into a …
Browse files Browse the repository at this point in the history
…single clean branch. These changes include a new vertical interpolation and quality control algorithm for init_type="metgrid", similar to that in WRF's real.exe. Also includes updates to docs and cleanup of some code in touched files.
  • Loading branch information
wiersema1 committed May 24, 2024
1 parent 80f78f6 commit 40b1355
Show file tree
Hide file tree
Showing 8 changed files with 1,011 additions and 211 deletions.
75 changes: 73 additions & 2 deletions Docs/sphinx_doc/Inputs.rst
Original file line number Diff line number Diff line change
Expand Up @@ -954,7 +954,8 @@ List of Parameters
| **erf.init_type** | Initialization | “custom”, | “*custom*” |
| | type | “ideal”, | |
| | | "real", | |
| | |"input_sounding" | |
| | | "input_sounding", | |
| | | "metgrid", | |
+----------------------------------+-------------------+--------------------+------------------+
| **erf.input_sounding_file** | Path to WRF-style | String | "input_sounding" |
| | input sounding | | |
Expand All @@ -979,9 +980,70 @@ List of Parameters
| | mesoscale data at | | |
| | lateral boundaries| | |
+----------------------------------+-------------------+--------------------+------------------+
| **erf.project_initial_velocity** | project initial | Integer | 1 |
| **erf.project_initial_velocity** | Project initial | Integer | 1 |
| | velocity? | | |
+----------------------------------+-------------------+--------------------+------------------+
| **erf.real_width** | Lateral boundary | Integer | 0 |
| | total width if | | |
| | use_real_bcs is | | |
| | true | | |
+----------------------------------+-------------------+--------------------+------------------+
| **erf.real_set_width** | Lateral boundary | Integer | 0 |
| | specified zone | | |
| | width if | | |
| | use_real_bcs is | | |
| | true | | |
+----------------------------------+-------------------+--------------------+------------------+
| **erf.metgrid_basic_linear** | If init_type is | true or false | false |
| | metgrid, use | | |
| | linear vertical | | |
| | interpolation and | | |
| | no quality | | |
| | control? | | |
+----------------------------------+-------------------+--------------------+------------------+
| **erf.metgrid_use_below_sfc** | If init_type is | true or false | true |
| | metgrid, use the | | |
| | origin data levels| | |
| | below the surface?| | |
+----------------------------------+-------------------+--------------------+------------------+
| **erf.metgrid_use_sfc** | If init_type is | true or false | true |
| | metgrid, use the | | |
| | origin data level | | |
| | at the surface? | | |
+----------------------------------+-------------------+--------------------+------------------+
| **erf.metgrid_retain_sfc** | If init_type is | true or false | false |
| | metgrid, assign | | |
| | the lowest level | | |
| | directly using the| | |
| | surface value from| | |
| | the origin data? | | |
+----------------------------------+-------------------+--------------------+------------------+
| **erf.metgrid_proximity** | If init_type is | Real | 1000. |
| | metgrid, pressure | | |
| | differential for | | |
| | detecting origin | | |
| | levels that are | | |
| | problematically | | |
| | close together | | |
+----------------------------------+-------------------+--------------------+------------------+
| **erf.metgrid_order** | If init_type is | Integer | 2 |
| | metgrid, order of | | |
| | the Lagrange | | |
| | polynomial | | |
| | interpolation | | |
| | scheme for | | |
| | vertical | | |
| | interpolation | | |
+----------------------------------+-------------------+--------------------+------------------+
| **erf.metgrid_force_sfc_k** | If init_type is | Integer | 0 |
| | metgrid, force the| | |
| | origin data | | |
| | surface level to | | |
| | be included in the| | |
| | interpolation for | | |
| | this many ERF | | |
| | vertical levels | | |
+----------------------------------+-------------------+--------------------+------------------+

Notes
-----------------
Expand All @@ -992,6 +1054,15 @@ If **erf.init_type = real**, the problem is initialized with mesoscale data cont
provided via ``erf.nc_init_file``. The mesoscale data are realistic with variation in all three directions.
In addition, the lateral boundary conditions must be supplied in a NetCDF files specified by **erf.nc_bdy_file = wrfbdy_d01**

If **erf.init_type = metgrid**, the problem is initialized with data
contained in the first NetCDF file provided via ``erf.nc_init_file_0``.
Lateral boundary conditions are derived from the sequence of NetCDF
files provided via ``erf.nc_init_file_0``. The sequence of
``erf.nc_init_file_0`` should be output from the WRF Preprocessing
System (WPS) listed chronologically starting with the earliest
timestamp. A minimum of two files are required to derive lateral
boundary conditions.

If **erf.init_type = input_sounding**, a WRF-style input sounding is read from
``erf.input_sounding_file``. This text file includes any set of levels that
goes at least up to the model top height. The first line includes the surface
Expand Down
2 changes: 1 addition & 1 deletion Exec/DevTests/MetGrid/inputs
Original file line number Diff line number Diff line change
Expand Up @@ -44,7 +44,7 @@ erf.check_int = 100 # number of timesteps between checkpoints
# PLOTFILES
erf.plot_file_1 = plt # prefix of plotfile name
erf.plot_int_1 = 10 # number of timesteps between plotfiles
erf.plot_vars_1 = density dens_hse rhoadv_0 x_velocity y_velocity z_velocity pressure temp theta z_phys mapfac pres_hse pert_pres
erf.plot_vars_1 = lat_m lon_m density dens_hse rhoadv_0 x_velocity y_velocity z_velocity pressure temp theta z_phys mapfac pres_hse pert_pres

# SOLVER CHOICE
erf.alpha_T = 1.0
Expand Down
9 changes: 9 additions & 0 deletions Source/ERF.H
Original file line number Diff line number Diff line change
Expand Up @@ -882,6 +882,15 @@ private:
// used with init_type == "input_sounding"
static bool init_sounding_ideal;

// Options for vertical interpolation of met_em*.nc data.
bool metgrid_basic_linear{false};
bool metgrid_use_below_sfc{true};
bool metgrid_use_sfc{true};
bool metgrid_retain_sfc{false};
amrex::Real metgrid_proximity{1000.0};
int metgrid_order{2};
int metgrid_force_sfc_k{0};

// 1D CDF output (for ingestion in AMR-Wind)
static int output_1d_column;
static int column_interval;
Expand Down
9 changes: 9 additions & 0 deletions Source/ERF.cpp
Original file line number Diff line number Diff line change
Expand Up @@ -1358,6 +1358,15 @@ ERF::ReadParameters ()
// Flag to trigger initialization from input_sounding like WRF's ideal.exe
pp.query("init_sounding_ideal", init_sounding_ideal);

// Options for vertical interpolation of met_em*.nc data.
pp.query("metgrid_basic_linear", metgrid_basic_linear);
pp.query("metgrid_use_below_sfc", metgrid_use_below_sfc);
pp.query("metgrid_use_sfc", metgrid_use_sfc);
pp.query("metgrid_retain_sfc", metgrid_retain_sfc);
pp.query("metgrid_proximity", metgrid_proximity);
pp.query("metgrid_order", metgrid_order);
pp.query("metgrid_force_sfc_k", metgrid_force_sfc_k);

// Output format
pp.query("plotfile_type", plotfile_type);
if (plotfile_type != "amrex" &&
Expand Down
2 changes: 1 addition & 1 deletion Source/IO/NCWpsFile.H
Original file line number Diff line number Diff line change
Expand Up @@ -17,7 +17,7 @@ using PlaneVector = amrex::Vector<amrex::FArrayBox>;
NetCDF variables of dimensions Time_BT_SN_WE: "UU", "VV", "TT", "RH", "PRES", "GHT"
NetCDF variables of dimensions Time_SN_WE : "HGT", "MAPFAC_U", "MAPFAC_V", "MAPFAC_M", "PSFC"
NetCDF global attributes of type int : "FLAG_PSFC", "FLAG_MAPFAC_U", "FLAG_MAPFAC_V", "FLAG_MAPFAC_M",
"FLAG_HGT_M", "WEST-EAST_GRID_DIMENSION", "SOUTH-NORTH_GRID_DIMENSION"
"WEST-EAST_GRID_DIMENSION", "SOUTH-NORTH_GRID_DIMENSION"
NetCDF global attributes of type string : "SIMULATION_START_DATE"
NetCDF global attributes of type real : "DX", "DY"
Expand Down
63 changes: 34 additions & 29 deletions Source/IO/ReadFromMetgrid.cpp
Original file line number Diff line number Diff line change
Expand Up @@ -10,8 +10,8 @@ using namespace amrex;
void
read_from_metgrid (int lev, const Box& domain, const std::string& fname,
std::string& NC_dateTime, Real& NC_epochTime,
int& flag_psfc, int& flag_msfu, int& flag_msfv, int& flag_msfm,
int& flag_hgt, int& flag_sst, int& flag_lmask,
int& flag_psfc, int& flag_msf,
int& flag_sst, int& flag_lmask,
int& NC_nx, int& NC_ny,
Real& NC_dx, Real& NC_dy,
FArrayBox& NC_xvel_fab, FArrayBox& NC_yvel_fab,
Expand All @@ -32,24 +32,32 @@ read_from_metgrid (int lev, const Box& domain, const std::string& fname,
auto ncf = ncutils::NCFile::open(fname, NC_CLOBBER | NC_NETCDF4);
{ // Global Attributes (int)
std::vector<int> attr;
ncf.get_attr("FLAG_PSFC", attr); flag_psfc = attr[0];

/* UNCOMMENT FOR FLOWMAS
flag_msfu = 0; //ncf.get_attr("FLAG_MAPFAC_U", attr); flag_msfu = attr[0];
flag_msfv = 0; //ncf.get_attr("FLAG_MAPFAC_V", attr); flag_msfv = attr[0];
flag_msfm = 0; //ncf.get_attr("FLAG_MAPFAC_M", attr); flag_msfm = attr[0];
flag_hgt = 1; //ncf.get_attr("FLAG_HGT_M", attr); flag_hgt = attr[0];
flag_sst = 0; //ncf.get_attr("FLAG_SST", attr); flag_sst = attr[0];
flag_psfc = 0; //ncf.get_attr("FLAG_PSFC", attr); flag_psfc = attr[0];
flag_msf = 0; //ncf.get_attr("FLAG_MF_XY", attr); flag_msf = attr[0];
flag_sst = 0; //ncf.get_attr("FLAG_SST", attr); flag_sst = attr[0];
flag_lmask = 0; //ncf.get_attr("FLAG_LANDMASK", attr); flag_lmask = attr[0];
*/

ncf.get_attr("FLAG_MAPFAC_U", attr); flag_msfu = attr[0];
ncf.get_attr("FLAG_MAPFAC_V", attr); flag_msfv = attr[0];
ncf.get_attr("FLAG_MAPFAC_M", attr); flag_msfm = attr[0];
ncf.get_attr("FLAG_HGT_M", attr); flag_hgt = attr[0];
flag_sst = 0; //ncf.get_attr("FLAG_SST", attr); flag_sst = attr[0];
ncf.get_attr("FLAG_LANDMASK", attr); flag_lmask = attr[0];

if (ncf.has_attr("FLAG_PSFC")) {
ncf.get_attr("FLAG_PSFC", attr); flag_psfc = attr[0];
} else {
flag_psfc = 0;
}
if (ncf.has_attr("FLAG_MF_XY")) {
ncf.get_attr("FLAG_MF_XY", attr); flag_msf = attr[0];
} else {
flag_msf = 0;
}
if (ncf.has_attr("FLAG_SST")) {
ncf.get_attr("FLAG_SST", attr); flag_sst = attr[0];
} else {
flag_sst = 0;
}
if (ncf.has_attr("FLAG_LANDMASK")) {
ncf.get_attr("FLAG_LANDMASK", attr); flag_lmask = attr[0];
} else {
flag_lmask = 0;
}

ncf.get_attr("WEST-EAST_GRID_DIMENSION", attr); NC_nx = attr[0];
ncf.get_attr("SOUTH-NORTH_GRID_DIMENSION", attr); NC_ny = attr[0];
Expand Down Expand Up @@ -85,10 +93,7 @@ read_from_metgrid (int lev, const Box& domain, const std::string& fname,
}
int ioproc = ParallelDescriptor::IOProcessorNumber(); // I/O rank
ParallelDescriptor::Bcast(&flag_psfc, 1, ioproc);
ParallelDescriptor::Bcast(&flag_msfu, 1, ioproc);
ParallelDescriptor::Bcast(&flag_msfv, 1, ioproc);
ParallelDescriptor::Bcast(&flag_msfm, 1, ioproc);
ParallelDescriptor::Bcast(&flag_hgt, 1, ioproc);
ParallelDescriptor::Bcast(&flag_msf, 1, ioproc);
ParallelDescriptor::Bcast(&flag_sst, 1, ioproc);
ParallelDescriptor::Bcast(&flag_lmask, 1, ioproc);
ParallelDescriptor::Bcast(&NC_nx, 1, ioproc);
Expand All @@ -112,21 +117,21 @@ read_from_metgrid (int lev, const Box& domain, const std::string& fname,
NC_fabs.push_back(&NC_rhum_fab); NC_fnames.push_back("RH"); NC_fdim_types.push_back(NC_Data_Dims_Type::Time_BT_SN_WE);
NC_fabs.push_back(&NC_pres_fab); NC_fnames.push_back("PRES"); NC_fdim_types.push_back(NC_Data_Dims_Type::Time_BT_SN_WE);
NC_fabs.push_back(&NC_ght_fab); NC_fnames.push_back("GHT"); NC_fdim_types.push_back(NC_Data_Dims_Type::Time_BT_SN_WE);
NC_fabs.push_back(&NC_LAT_fab); NC_fnames.push_back("XLAT_V"); NC_fdim_types.push_back(NC_Data_Dims_Type::Time_SN_WE);
NC_fabs.push_back(&NC_LON_fab); NC_fnames.push_back("XLONG_U"); NC_fdim_types.push_back(NC_Data_Dims_Type::Time_SN_WE);
NC_fabs.push_back(&NC_LAT_fab); NC_fnames.push_back("XLAT_M"); NC_fdim_types.push_back(NC_Data_Dims_Type::Time_SN_WE);
NC_fabs.push_back(&NC_LON_fab); NC_fnames.push_back("XLONG_M"); NC_fdim_types.push_back(NC_Data_Dims_Type::Time_SN_WE);
NC_fabs.push_back(&NC_hgt_fab); NC_fnames.push_back("HGT_M"); NC_fdim_types.push_back(NC_Data_Dims_Type::Time_SN_WE);

if (flag_psfc) { NC_fabs.push_back(&NC_psfc_fab); NC_fnames.push_back("PSFC"); NC_fdim_types.push_back(NC_Data_Dims_Type::Time_SN_WE); }
if (flag_msfu) { NC_fabs.push_back(&NC_msfu_fab); NC_fnames.push_back("MAPFAC_U"); NC_fdim_types.push_back(NC_Data_Dims_Type::Time_SN_WE); }
if (flag_msfv) { NC_fabs.push_back(&NC_msfv_fab); NC_fnames.push_back("MAPFAC_V"); NC_fdim_types.push_back(NC_Data_Dims_Type::Time_SN_WE); }
if (flag_msfm) { NC_fabs.push_back(&NC_msfm_fab); NC_fnames.push_back("MAPFAC_M"); NC_fdim_types.push_back(NC_Data_Dims_Type::Time_SN_WE); }
if (flag_hgt) { NC_fabs.push_back(&NC_hgt_fab); NC_fnames.push_back("HGT_M"); NC_fdim_types.push_back(NC_Data_Dims_Type::Time_SN_WE); }
if (flag_msf) { NC_fabs.push_back(&NC_msfu_fab); NC_fnames.push_back("MAPFAC_U"); NC_fdim_types.push_back(NC_Data_Dims_Type::Time_SN_WE); }
if (flag_msf) { NC_fabs.push_back(&NC_msfv_fab); NC_fnames.push_back("MAPFAC_V"); NC_fdim_types.push_back(NC_Data_Dims_Type::Time_SN_WE); }
if (flag_msf) { NC_fabs.push_back(&NC_msfm_fab); NC_fnames.push_back("MAPFAC_M"); NC_fdim_types.push_back(NC_Data_Dims_Type::Time_SN_WE); }
if (flag_sst) { NC_fabs.push_back(&NC_sst_fab); NC_fnames.push_back("SST"); NC_fdim_types.push_back(NC_Data_Dims_Type::Time_SN_WE); }

if (flag_lmask) { NC_iabs.push_back(&NC_lmask_iab); NC_inames.push_back("LANDMASK"); NC_idim_types.push_back(NC_Data_Dims_Type::Time_SN_WE); }

// Read the netcdf file and fill these FABs
std::string Lat_var_name = "XLAT_V";
std::string Lon_var_name = "XLONG_U";
std::string Lat_var_name = "XLAT_M";
std::string Lon_var_name = "XLONG_M";
Print() << "Building initial FABS from file " << fname << std::endl;
BuildFABsFromNetCDFFile<FArrayBox,Real>(domain, Latitude, Longitude,
Lat_var_name, Lon_var_name,
Expand Down
Loading

0 comments on commit 40b1355

Please sign in to comment.