ERF
Energy Research and Forecasting: An Atmospheric Modeling Code
ERF_NCPlotFile.cpp File Reference
#include <iomanip>
#include <iostream>
#include <string>
#include <ctime>
#include <AMReX_Utility.H>
#include <AMReX_MultiFab.H>
#include "ERF_Constants.H"
#include "ERF_DataStruct.H"
#include "ERF_NCPlotFile.H"
#include "ERF_NCInterface.H"
Include dependency graph for ERF_NCPlotFile.cpp:

Functions

void writeNCPlotFile (int lev, int which_subdomain, const std::string &dir, const Vector< const MultiFab * > &plotMF, const Vector< std::string > &plot_var_names, const Vector< int > &, Array< Real, AMREX_SPACEDIM > prob_lo, Array< Real, AMREX_SPACEDIM > prob_hi, Array< Real, AMREX_SPACEDIM > dx_in, const Box &subdomain, const double &time, const double &start_bdy_time, const SolverChoice &solverChoice, const Vector< Real > &zlevels_stag)
 

Function Documentation

◆ writeNCPlotFile()

void writeNCPlotFile ( int  lev,
int  which_subdomain,
const std::string &  dir,
const Vector< const MultiFab * > &  plotMF,
const Vector< std::string > &  plot_var_names,
const Vector< int > &  ,
Array< Real, AMREX_SPACEDIM >  prob_lo,
Array< Real, AMREX_SPACEDIM >  prob_hi,
Array< Real, AMREX_SPACEDIM >  dx_in,
const Box &  subdomain,
const double &  time,
const double &  start_bdy_time,
const SolverChoice solverChoice,
const Vector< Real > &  zlevels_stag 
)
32 {
33  //
34  // Set the full IO path for NetCDF output
35  //
36  std::string FullPath = dir;
37  if (lev == 0) {
38  const std::string& extension = amrex::Concatenate("_d",lev+1,2);
39  FullPath += extension + ".nc";
40  } else {
41  const std::string& extension = amrex::Concatenate("_d",lev+1+which_subdomain,2);
42  FullPath += extension + ".nc";
43  }
44 
45  Print() << "Writing level " << lev << " NetCDF plot file " << FullPath << std::endl;
46 
47  //
48  // Open netcdf file to write data
49  //
50  auto ncf = ncutils::NCFile::create_par(FullPath, NC_NETCDF4 | NC_MPIIO,
51  amrex::ParallelContext::CommunicatorSub(), MPI_INFO_NULL);
52 
53  auto ba = plotMF[lev]->boxArray();
54  auto dm = plotMF[lev]->DistributionMap();
55 
56  int nblocks = ba.size();
57 
58  int nx = subdomain.length(0);
59  int ny = subdomain.length(1);
60  int nz = subdomain.length(2);
61 
62  int num_pts = nx*ny*nz;
63 
64  int n_data_items = plotMF[lev]->nComp();
65 
66  //
67  // Start Define stuff
68  //
69  ncf.enter_def_mode();
70 
71  const std::string nt_name = "num_time_steps";
72  const std::string nb_name = "num_blocks";
73  const std::string np_name = "num_pts";
74  const std::string nx_name = "nx";
75  const std::string ny_name = "ny";
76  const std::string nz_name = "nz";
77 
78  const std::string ndim_name = "num_geo_dims";
79 
80  ncf.put_attr("title", "ERF NetCDF Plot data output");
81 
82  ncf.def_dim(ndim_name, AMREX_SPACEDIM);
83  ncf.def_dim(np_name , num_pts);
84  ncf.def_dim(nb_name , nblocks);
85 
86  ncf.def_dim(nt_name, NC_UNLIMITED);
87  ncf.def_dim(nx_name, nx);
88  ncf.def_dim(ny_name, ny);
89  ncf.def_dim(nz_name, nz);
90 
91  ncf.def_var("probLo" , NC_FLOAT, {ndim_name});
92  ncf.def_var("probHi" , NC_FLOAT, {ndim_name});
93 
94  ncf.def_var("Geom.smallend", NC_INT, {ndim_name});
95  ncf.def_var("Geom.bigend" , NC_INT, {ndim_name});
96  ncf.def_var("CellSize" , NC_FLOAT, {ndim_name});
97 
98  ncf.def_var("x_grid", NC_DOUBLE, {np_name});
99  ncf.def_var("y_grid", NC_DOUBLE, {np_name});
100  ncf.def_var("z_grid", NC_DOUBLE, {np_name});
101 
102  for (int i = 0; i < plot_var_names.size(); i++) {
103  ncf.def_var(plot_var_names[i], NC_DOUBLE, {nz_name, ny_name, nx_name});
104  }
105 
106  ncf.exit_def_mode();
107  //
108  // End Define stuff
109  //
110 
111  // We are doing single-level writes but it doesn't have to be level 0
112  //
113  // Write out the netcdf plotfile head information.
114  //
115  if (n_data_items == 0) {
116  amrex::Error("Must specify at least one valid data item to plot");
117  }
118 
119  ncf.put_attr("number_variables", std::vector<int>{n_data_items});
120  ncf.put_attr("space_dimension", std::vector<int>{AMREX_SPACEDIM});
121  ncf.put_attr("current_time", std::vector<double>{time});
122  ncf.put_attr("start_time", std::vector<double>{start_bdy_time});
123  ncf.put_attr("CurrentLevel", std::vector<int>{lev});
124 
125  Real dx[AMREX_SPACEDIM];
126  for (int i = 0; i < AMREX_SPACEDIM; i++) {
127  dx[i] = dx_in[i];
128  }
129 
130  amrex::Vector<Real> probLo;
131  amrex::Vector<Real> probHi;
132  for (int i = 0; i < AMREX_SPACEDIM; i++) {
133  probLo.push_back(prob_lo[i]);
134  probHi.push_back(prob_hi[i]);
135  }
136 
137  auto nc_probLo = ncf.var("probLo");
138  nc_probLo.par_access(NC_COLLECTIVE);
139  nc_probLo.put(probLo.data(), {0}, {AMREX_SPACEDIM});
140 
141  auto nc_probHi = ncf.var("probHi");
142  nc_probHi.par_access(NC_COLLECTIVE);
143  nc_probHi.put(probHi.data(), {0}, {AMREX_SPACEDIM});
144 
145  amrex::Vector<int> smallend;
146  amrex::Vector<int> bigend;
147  smallend.clear(); bigend.clear();
148  for (int j = 0; j < AMREX_SPACEDIM; j++) {
149  smallend.push_back(subdomain.smallEnd(j));
150  bigend.push_back(subdomain.bigEnd(j));
151  }
152 
153  auto nc_Geom_smallend = ncf.var("Geom.smallend");
154  nc_Geom_smallend.par_access(NC_COLLECTIVE);
155  nc_Geom_smallend.put(smallend.data(), {0}, {AMREX_SPACEDIM});
156 
157  auto nc_Geom_bigend = ncf.var("Geom.bigend");
158  nc_Geom_bigend.par_access(NC_COLLECTIVE);
159  nc_Geom_bigend.put(bigend.data(), {0}, {AMREX_SPACEDIM});
160 
161  amrex::Vector<Real> CellSize;
162  CellSize.clear();
163  for (Real& j : dx) {
164  CellSize.push_back(j);
165  }
166  auto nc_CellSize = ncf.var("CellSize");
167  nc_CellSize.par_access(NC_COLLECTIVE);
168  nc_CellSize.put(CellSize.data(), {0}, {AMREX_SPACEDIM});
169 
170  ncf.put_attr("DefaultGeometry", std::vector<int>{amrex::DefaultGeometry().Coord()});
171 
172  std::vector<Real> x_grid;
173  std::vector<Real> y_grid;
174  std::vector<Real> z_grid;
175  long unsigned goffset = 0;
176  long unsigned glen = 0;
177 
178  // *******************************************************************************
179  // NOTE: the (x,y,z) output here are for a mesh withOUT terrain-fitted coordinates
180  // *******************************************************************************
181  if (solverChoice.mesh_type == MeshType::ConstantDz) {
182  for (int i = 0; i < ba.size(); ++i) {
183  auto bx = ba[i];
184  if (subdomain.contains(bx)) {
185  // The loop indices below are box-local, so they must be offset by the
186  // low corner of the box to give the global location of each cell
187  const auto bx_lo = bx.smallEnd();
188  x_grid.clear(); y_grid.clear(); z_grid.clear();
189  for (auto k3 = 0; k3 < bx.length(2); ++k3) {
190  for (auto k2 = 0; k2 < bx.length(1); ++k2) {
191  for (auto k1 = 0; k1 < bx.length(0); ++k1) {
192  x_grid.push_back(prob_lo[0]+dx[0]*(static_cast<Real>(bx_lo[0]+k1)+myhalf));
193  y_grid.push_back(prob_lo[1]+dx[1]*(static_cast<Real>(bx_lo[1]+k2)+myhalf));
194  z_grid.push_back(prob_lo[2]+dx[2]*(static_cast<Real>(bx_lo[2]+k3)+myhalf));
195  }
196  }
197  }
198 
199  goffset += glen;
200  glen = bx.numPts();
201 
202  auto nc_x_grid = ncf.var("x_grid");
203  auto nc_y_grid = ncf.var("y_grid");
204  auto nc_z_grid = ncf.var("z_grid");
205 
206  nc_x_grid.par_access(NC_COLLECTIVE);
207  nc_y_grid.par_access(NC_COLLECTIVE);
208  nc_z_grid.par_access(NC_COLLECTIVE);
209 
210  nc_x_grid.put(x_grid.data(), {goffset}, {glen});
211  nc_y_grid.put(y_grid.data(), {goffset}, {glen});
212  nc_z_grid.put(z_grid.data(), {goffset}, {glen});
213  }
214  }
215  }
216  else if (solverChoice.mesh_type == MeshType::StretchedDz)
217  {
218  for (int i = 0; i < ba.size(); ++i) {
219  auto bx = ba[i];
220  if (subdomain.contains(bx)) {
221  // The loop indices below are box-local, so they must be offset by the
222  // low corner of the box to give the global location of each cell;
223  // zlevels_stag is indexed globally as well
224  const auto bx_lo = bx.smallEnd();
225  x_grid.clear(); y_grid.clear(); z_grid.clear();
226  for (auto k3 = 0; k3 < bx.length(2); ++k3) {
227  for (auto k2 = 0; k2 < bx.length(1); ++k2) {
228  for (auto k1 = 0; k1 < bx.length(0); ++k1) {
229  x_grid.push_back(prob_lo[0]+dx[0]*(static_cast<Real>(bx_lo[0]+k1)+myhalf));
230  y_grid.push_back(prob_lo[1]+dx[1]*(static_cast<Real>(bx_lo[1]+k2)+myhalf));
231  z_grid.push_back(myhalf * (zlevels_stag[bx_lo[2]+k3] + zlevels_stag[bx_lo[2]+k3+1]));
232  }
233  }
234  }
235 
236  goffset += glen;
237  glen = bx.numPts();
238 
239  auto nc_x_grid = ncf.var("x_grid");
240  auto nc_y_grid = ncf.var("y_grid");
241  auto nc_z_grid = ncf.var("z_grid");
242 
243  nc_x_grid.par_access(NC_COLLECTIVE);
244  nc_y_grid.par_access(NC_COLLECTIVE);
245  nc_z_grid.par_access(NC_COLLECTIVE);
246 
247  nc_x_grid.put(x_grid.data(), {goffset}, {glen});
248  nc_y_grid.put(y_grid.data(), {goffset}, {glen});
249  nc_z_grid.put(z_grid.data(), {goffset}, {glen});
250  }
251  }
252  }
253 
254  const int ncomp = plotMF[lev]->nComp();
255 
256  for (MFIter mfi(*plotMF[lev]); mfi.isValid(); ++mfi)
257  {
258  auto bx = mfi.validbox();
259 
260  if (subdomain.contains(bx))
261 
262  {
263  //
264  // These are the dimensions of the data we write for only this box
265  //
266  long unsigned local_nx = bx.length()[0];
267  long unsigned local_ny = bx.length()[1];
268  long unsigned local_nz = bx.length()[2];
269 
270  long unsigned local_start_x = static_cast<long unsigned>(bx.smallEnd()[0]-subdomain.smallEnd()[0]);
271  long unsigned local_start_y = static_cast<long unsigned>(bx.smallEnd()[1]-subdomain.smallEnd()[1]);
272  long unsigned local_start_z = static_cast<long unsigned>(bx.smallEnd()[2]-subdomain.smallEnd()[2]);
273 
274  for (int k(0); k < ncomp; ++k) {
275  FArrayBox tmp;
276  tmp.resize(bx, 1, amrex::The_Pinned_Arena());
277  tmp.template copy<RunOn::Device>((*plotMF[lev])[mfi.index()], k, 0, 1);
278  Gpu::streamSynchronize();
279 
280  auto nc_plot_var = ncf.var(plot_var_names[k]);
281  nc_plot_var.par_access(NC_COLLECTIVE);
282  nc_plot_var.put(tmp.dataPtr(), {local_start_z,local_start_y,local_start_x},
283  {local_nz, local_ny, local_nx});
284  }
285  }
286  }
287  ncf.close();
288 }
constexpr amrex::Real myhalf
Definition: ERF_Constants.H:13
const Real dx
Definition: ERF_InitCustomPert_ABL.H:23
const GpuArray< Real, AMREX_SPACEDIM > prob_lo
Definition: ERF_InitCustomPert_CloudChamber.H:33
const amrex::Real * prob_hi
Definition: ERF_InitCustomPert_DataAssimilation_ISV.H:17
amrex::Real Real
Definition: ERF_ShocInterface.H:19
@ tmp
Definition: ERF_AdvanceWSM6.cpp:114
real(c_double), private k1
Definition: ERF_module_mp_morr_two_moment.F90:213
static MeshType mesh_type
Vertical mesh representation.
Definition: ERF_DataStruct.H:1377

Referenced by main_main(), ERF::Write2DPlotFile(), and ERF::Write3DPlotFile().

Here is the caller graph for this function: