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

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

Here is the call graph for this function:
Here is the caller graph for this function: