ERF
Energy Research and Forecasting: An Atmospheric Modeling Code
All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Macros Pages
ERF_InitForecastData.cpp File Reference
Include dependency graph for ERF_InitForecastData.cpp:

Enumerations

enum class  BoundType { Lo , Hi }
 
enum class  MultiFabType { CC , NC }
 

Functions

void fill_weather_data_multifab (MultiFab &mf, const Geometry &geom_weather, const int nx, const int ny, const int nz, const Vector< Real > &latvec_h, const Vector< Real > &lonvec_h, const Vector< Real > &zvec_h, const Vector< Real > &rho_h, const Vector< Real > &uvel_h, const Vector< Real > &vvel_h, const Vector< Real > &wvel_h, const Vector< Real > &theta_h, const Vector< Real > &qv_h, const Vector< Real > &qc_h, const Vector< Real > &qr_h)
 
void PlotMultiFab (const MultiFab &mf, const Geometry &geom_mf, const std::string plotfilename, MultiFabType mftype)
 
IntVect find_bound_idx (const Real &x, const Real &y, const Real &z, const BoxList &bl_weather, const Geometry &geom_weather, BoundType bound_type)
 

Enumeration Type Documentation

◆ BoundType

enum BoundType
strong
Enumerator
Lo 
Hi 
99 { Lo, Hi };

◆ MultiFabType

enum MultiFabType
strong
Enumerator
CC 
NC 
100 { CC, NC };

Function Documentation

◆ fill_weather_data_multifab()

void fill_weather_data_multifab ( MultiFab &  mf,
const Geometry &  geom_weather,
const int  nx,
const int  ny,
const int  nz,
const Vector< Real > &  latvec_h,
const Vector< Real > &  lonvec_h,
const Vector< Real > &  zvec_h,
const Vector< Real > &  rho_h,
const Vector< Real > &  uvel_h,
const Vector< Real > &  vvel_h,
const Vector< Real > &  wvel_h,
const Vector< Real > &  theta_h,
const Vector< Real > &  qv_h,
const Vector< Real > &  qc_h,
const Vector< Real > &  qr_h 
)
26 {
27  const int nx_d = nx;
28  const int ny_d = ny;
29  const int nz_d = nz;
30  const int tot_size = nx*ny*nz;
31 
32  amrex::Gpu::DeviceVector<Real> latvec_d(nx*ny), lonvec_d(nx*ny), zvec_d(nz);
33  amrex::Gpu::DeviceVector<Real> rho_d(tot_size), uvel_d(tot_size),
34  vvel_d(tot_size), wvel_d(tot_size), theta_d(tot_size),
35  qv_d(tot_size), qc_d(tot_size), qr_d(tot_size);
36 
37  amrex::Gpu::copyAsync(amrex::Gpu::hostToDevice, latvec_h.begin(), latvec_h.end(), latvec_d.begin());
38  amrex::Gpu::copyAsync(amrex::Gpu::hostToDevice, lonvec_h.begin(), lonvec_h.end(), lonvec_d.begin());
39  amrex::Gpu::copyAsync(amrex::Gpu::hostToDevice, zvec_h.begin(), zvec_h.end(), zvec_d.begin());
40  amrex::Gpu::copyAsync(amrex::Gpu::hostToDevice, rho_h.begin(), rho_h.end(), rho_d.begin());
41  amrex::Gpu::copyAsync(amrex::Gpu::hostToDevice, theta_h.begin(), theta_h.end(), theta_d.begin());
42  amrex::Gpu::copyAsync(amrex::Gpu::hostToDevice, uvel_h.begin(), uvel_h.end(), uvel_d.begin());
43  amrex::Gpu::copyAsync(amrex::Gpu::hostToDevice, vvel_h.begin(), vvel_h.end(), vvel_d.begin());
44  amrex::Gpu::copyAsync(amrex::Gpu::hostToDevice, wvel_h.begin(), wvel_h.end(), wvel_d.begin());
45  amrex::Gpu::copyAsync(amrex::Gpu::hostToDevice, qv_h.begin(), qv_h.end(), qv_d.begin());
46  amrex::Gpu::copyAsync(amrex::Gpu::hostToDevice, qc_h.begin(), qc_h.end(), qc_d.begin());
47  amrex::Gpu::copyAsync(amrex::Gpu::hostToDevice, qr_h.begin(), qr_h.end(), qr_d.begin());
48 
49  Real* latvec_d_ptr = latvec_d.data();
50  Real* lonvec_d_ptr = lonvec_d.data();
51  Real* zvec_d_ptr = zvec_d.data();
52  Real* rho_d_ptr = rho_d.data();
53  Real* uvel_d_ptr = uvel_d.data();
54  Real* vvel_d_ptr = vvel_d.data();
55  Real* wvel_d_ptr = wvel_d.data();
56  Real* theta_d_ptr = theta_d.data();
57  Real* qv_d_ptr = qv_d.data();
58  Real* qc_d_ptr = qc_d.data();
59  Real* qr_d_ptr = qr_d.data();
60 
61  const auto prob_lo = geom_weather.ProbLoArray();
62  const auto dx = geom_weather.CellSizeArray();
63 
64  for (MFIter mfi(mf, TilingIfNotGPU()); mfi.isValid(); ++mfi)
65  {
66  const Box& bx = mfi.nodaltilebox();
67  Array4<Real> const& arr = mf.array(mfi);
68 
69  ParallelFor(bx, [=] AMREX_GPU_DEVICE (int i, int j, int k)
70  {
71  const Real z = prob_lo[2] + k * dx[2];
72  int kloc = -1;
73  for (int kk = 0; kk < nz_d; kk++) {
74  if (zvec_d_ptr[kk] > z) {
75  kloc = kk;
76  break;
77  }
78  }
79  // Replace this with logic using your coarse input data
80  int idx1 = get_single_index(i,j,kloc-1,nx_d,ny_d);
81  int idx2 = get_single_index(i,j,kloc,nx_d,ny_d);
82  Real dz = zvec_d_ptr[kloc] - zvec_d_ptr[kloc-1];
83  Real fac = (z-zvec_d_ptr[kloc-1])/dz;
84  arr(i,j,k,0) = rho_d_ptr[idx1] + fac*(rho_d_ptr[idx2]-rho_d_ptr[idx1]);
85  arr(i,j,k,1) = uvel_d_ptr[idx1] + fac*(uvel_d_ptr[idx2]-uvel_d_ptr[idx1]);
86  arr(i,j,k,2) = vvel_d_ptr[idx1] + fac*(vvel_d_ptr[idx2]-vvel_d_ptr[idx1]);
87  arr(i,j,k,3) = wvel_d_ptr[idx1] + fac*(wvel_d_ptr[idx2]-wvel_d_ptr[idx1]);
88  arr(i,j,k,4) = theta_d_ptr[idx1] + fac*(theta_d_ptr[idx2]-theta_d_ptr[idx1]);
89  arr(i,j,k,5) = qv_d_ptr[idx1] + fac*(qv_d_ptr[idx2]-qv_d_ptr[idx1]);
90  arr(i,j,k,6) = qc_d_ptr[idx1] + fac*(qc_d_ptr[idx2]-qc_d_ptr[idx1]);
91  arr(i,j,k,7) = qr_d_ptr[idx1] + fac*(qr_d_ptr[idx2]-qr_d_ptr[idx1]);
92  idx1 = get_single_index(i,j,0,nx_d,ny_d);
93  arr(i,j,k,8) = latvec_d_ptr[idx1];
94  arr(i,j,k,9) = lonvec_d_ptr[idx1];
95  });
96  }
97 }
AMREX_FORCE_INLINE AMREX_GPU_HOST_DEVICE int get_single_index(int i, int j, int k, int nx, int ny)
Definition: ERF_Interpolation_Bilinear.H:13

Referenced by ERF::init_coarse_weather_data().

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

◆ find_bound_idx()

IntVect find_bound_idx ( const Real &  x,
const Real &  y,
const Real &  z,
const BoxList &  bl_weather,
const Geometry &  geom_weather,
BoundType  bound_type 
)
236 {
237  const auto prob_lo_weather = geom_weather.ProbLoArray();
238  const auto dx_weather = geom_weather.CellSizeArray();
239 
240  int i, j, k;
241 
242  if (bound_type == BoundType::Lo) {
243  i = static_cast<int>(std::floor((x - prob_lo_weather[0]) / dx_weather[0]));
244  j = static_cast<int>(std::floor((y - prob_lo_weather[1]) / dx_weather[1]));
245  k = static_cast<int>(std::floor((z - prob_lo_weather[2]) / dx_weather[2]));
246  } else { // BoundType::Hi
247  i = static_cast<int>(std::ceil((x - prob_lo_weather[0]) / dx_weather[0]));
248  j = static_cast<int>(std::ceil((y - prob_lo_weather[1]) / dx_weather[1]));
249  k = static_cast<int>(std::ceil((z - prob_lo_weather[2]) / dx_weather[2]));
250  }
251 
252  IntVect idx(i, j, k);
253 
254  for (const auto& b : bl_weather) {
255  if (b.contains(idx)) {
256  return idx;
257  }
258  }
259 
260  amrex::Abort("Bound index not found in any box in BoxList!");
261  return IntVect::TheZeroVector(); // unreachable if Abort
262 }

Referenced by ERF::interp_weather_data_onto_mesh().

Here is the caller graph for this function:

◆ PlotMultiFab()

void PlotMultiFab ( const MultiFab &  mf,
const Geometry &  geom_mf,
const std::string  plotfilename,
MultiFabType  mftype 
)
107 {
108 
109  Vector<std::string> varnames = {
110  "rho", "uvel", "vvel", "wvel", "theta", "qv", "qc", "qr", "latitude", "longitude"
111  }; // Customize variable names
112 
113  const Real time = 0.0;
114 
115 
116  // Assume weather_mf is nodal in all directions
117  if(mftype == MultiFabType::NC) {
118  BoxArray cba = mf.boxArray();
119  cba = amrex::convert(mf.boxArray(), IntVect::TheCellVector());
120 
121  MultiFab cc_mf(cba, mf.DistributionMap(),
122  mf.nComp(), 0);
123 
124  amrex::average_node_to_cellcenter(cc_mf, 0, mf, 0, mf.nComp());
125 
126  WriteSingleLevelPlotfile(
127  plotfilename,
128  cc_mf,
129  varnames,
130  geom_mf,
131  time,
132  0 // level
133  );
134  } else {
135  WriteSingleLevelPlotfile(
136  plotfilename,
137  mf,
138  varnames,
139  geom_mf,
140  time,
141  0 // level
142  );
143  }
144 }

Referenced by ERF::init_coarse_weather_data(), and ERF::interp_weather_data_onto_mesh().

Here is the caller graph for this function: