ERF
Energy Research and Forecasting: An Atmospheric Modeling Code
ERF_MetgridUtils.H
Go to the documentation of this file.
1 #ifndef ERF_METGRIDUTIL_H_
2 #define ERF_METGRIDUTIL_H_
3 
4 #include <ERF.H>
5 #include <ERF_EOS.H>
6 #include <ERF_Utils.H>
7 #include <ERF_ProbCommon.H>
8 #include <ERF_HSEUtils.H>
9 
10 void
11 read_from_metgrid (int lev,
12  const amrex::Box& domain,
13  const std::string& fname,
14  std::string& NC_dateTime,
15  amrex::Real& NC_epochTime,
16  int& flag_psfc,
17  int& flag_msf,
18  int& flag_sst,
19  int& flag_lmask,
20  int& NC_nx,
21  int& NC_ny,
22  amrex::Real& NC_dx,
23  amrex::Real& NC_dy,
24  amrex::FArrayBox& NC_xvel_fab,
25  amrex::FArrayBox& NC_yvel_fab,
26  amrex::FArrayBox& NC_temp_fab,
27  amrex::FArrayBox& NC_rhum_fab,
28  amrex::FArrayBox& NC_pres_fab,
29  amrex::FArrayBox& NC_ght_fab,
30  amrex::FArrayBox& NC_hgt_fab,
31  amrex::FArrayBox& NC_psfc_fab,
32  amrex::FArrayBox& NC_msfu_fab,
33  amrex::FArrayBox& NC_msfv_fab,
34  amrex::FArrayBox& NC_msfm_fab,
35  amrex::FArrayBox& NC_sst_fab,
36  amrex::FArrayBox& NC_LAT_fab,
37  amrex::FArrayBox& NC_LON_fab,
38  amrex::IArrayBox& NC_lmask_iab,
39  amrex::Real& Latitude,
40  amrex::Real& Longitude,
41  amrex::Geometry& geom);
42 
43 
44 
45 void
46 init_terrain_from_metgrid (amrex::FArrayBox& z_phys_nd_fab,
47  const amrex::Vector<amrex::FArrayBox>& NC_hgt_fab);
48 
49 void
50 init_state_from_metgrid (const bool use_moisture,
51  const bool interp_theta,
52  const bool metgrid_debug_quiescent,
53  const bool metgrid_debug_isothermal,
54  const bool metgrid_debug_dry,
55  const bool metgrid_basic_linear,
56  const bool metgrid_use_below_sfc,
57  const bool metgrid_use_sfc,
58  const bool metgrid_retain_sfc,
59  const amrex::Real metgrid_proximity,
60  const int metgrid_order,
61  const int metgrid_metgrid_force_sfc_k,
62  const amrex::Real l_rdOcp,
63  amrex::Box& tbxc,
64  amrex::Box& tbxu,
65  amrex::Box& tbxv,
66  amrex::FArrayBox& state_fab,
67  amrex::FArrayBox& x_vel_fab,
68  amrex::FArrayBox& y_vel_fab,
69  amrex::FArrayBox& z_vel_fab,
70  amrex::FArrayBox& z_phys_nd_fab,
71  const amrex::Vector<amrex::FArrayBox>& NC_hgt_fab,
72  const amrex::Vector<amrex::FArrayBox>& NC_ght_fab,
73  const amrex::Vector<amrex::FArrayBox>& NC_xvel_fab,
74  const amrex::Vector<amrex::FArrayBox>& NC_yvel_fab,
75  const amrex::Vector<amrex::FArrayBox>& NC_temp_fab,
76  const amrex::Vector<amrex::FArrayBox>& NC_rhum_fab,
77  const amrex::Vector<amrex::FArrayBox>& NC_pres_fab,
78  amrex::FArrayBox& p_interp_fab,
79  amrex::FArrayBox& t_interp_fab,
80  amrex::FArrayBox& theta_fab,
81  amrex::FArrayBox& mxrat_fab,
82  amrex::Vector<amrex::Vector<amrex::FArrayBox>>& fabs_for_bcs_xlo,
83  amrex::Vector<amrex::Vector<amrex::FArrayBox>>& fabs_for_bcs_xhi,
84  amrex::Vector<amrex::Vector<amrex::FArrayBox>>& fabs_for_bcs_ylo,
85  amrex::Vector<amrex::Vector<amrex::FArrayBox>>& fabs_for_bcs_yhi,
86  const amrex::Array4<const int>& mask_c_arr,
87  const amrex::Array4<const int>& mask_u_arr,
88  const amrex::Array4<const int>& mask_v_arr);
89 
90 void
91 init_msfs_from_metgrid (const bool metgrid_debug_msf,
92  amrex::FArrayBox& msfu_fab,
93  amrex::FArrayBox& msfv_fab,
94  amrex::FArrayBox& msfm_fab,
95  const int& flag_msf,
96  const amrex::Vector<amrex::FArrayBox>& NC_MSFU_fab,
97  const amrex::Vector<amrex::FArrayBox>& NC_MSFV_fab,
98  const amrex::Vector<amrex::FArrayBox>& NC_MSFM_fab);
99 
100 void
101 init_base_state_from_metgrid (const bool use_moisture,
102  const bool metgrid_debug_psfc,
103  const amrex::Real l_rdOcp,
104  const amrex::Box& valid_bx,
105  const amrex::Vector<int>& flag_psfc,
106  amrex::FArrayBox& state,
107  amrex::FArrayBox& r_hse_fab,
108  amrex::FArrayBox& p_hse_fab,
109  amrex::FArrayBox& pi_hse_fab,
110  amrex::FArrayBox& th_hse_fab,
111  amrex::FArrayBox& z_phys_cc_fab,
112  const amrex::Vector<amrex::FArrayBox>& NC_psfc_fab);
113 
114 AMREX_FORCE_INLINE
115 AMREX_GPU_DEVICE
116 void
117 lagrange_interp (const int& order,
118  amrex::Real* x,
119  amrex::Real* y,
120  amrex::Real& new_x,
121  amrex::Real& new_y)
122 {
123  // Interpolation using Lagrange polynomials.
124  // P(x) = f(x0)Ln0(x) + ... + f(xn)Lnn(x)
125  // where Lnk(x) = (x -x0)(x -x1)...(x -xk-1)(x -xk+1)...(x -xn)
126  // ---------------------------------------------
127  // (xk-x0)(xk-x1)...(xk-xk-1)(xk-xk+1)...(xk-xn)
128  amrex::Real Px = 0.0;
129  for (int i=0; i <= order; i++) {
130  amrex::Real n = 1.0;
131  amrex::Real d = 1.0;
132  for (int k=0; k <= order; k++) {
133  if (k == i) continue;
134  n *= new_x-x[k];
135  d *= x[i]-x[k];
136  }
137  if (d != 0.0) {
138  Px += y[i]*n/d;
139  }
140  }
141  new_y = Px;
142 }
143 
144 AMREX_FORCE_INLINE
145 AMREX_GPU_DEVICE
146 void
147 lagrange_setup (char var_type,
148  const bool& exp_interp,
149  const int& orig_n,
150  const int& new_n,
151  const int& order,
152  amrex::Real* orig_x_z,
153  amrex::Real* orig_x_p,
154  amrex::Real* orig_y,
155  amrex::Real* new_x_z,
156  amrex::Real* new_x_p,
157  amrex::Real* new_y)
158 {
159  if (order < 1) amrex::Abort("metgrid initialization, we cannot go lower order than linear");
160 
161  amrex::Real CRC_const1 = 11880.516; // m
162  amrex::Real CRC_const2 = 0.1902632;
163  amrex::Real CRC_const3 = 0.0065; // K km-1
164  int vboundb = 4;
165  int vboundt = 0;
166 
167 #ifndef AMREX_USE_GPU
168  bool debug = false;
169 #endif
170 
171  for (int new_k=0; new_k < new_n; new_k++) {
172 #ifndef AMREX_USE_GPU
173  if (debug) amrex::Print() << "new_k=" << new_k;
174 #endif
175  // Find bounding x values and store the indices.
176  bool extrapolating = true;
177  int kl, kr;
178  for (int ko=0; ko < orig_n-1; ko++) {
179  amrex::Real a = new_x_z[new_k]-orig_x_z[ko];
180  amrex::Real b = new_x_z[new_k]-orig_x_z[ko+1];
181  if (a*b <= 0.0) {
182  kl = ko;
183  kr = ko+1;
184  extrapolating = false;
185  break;
186  }
187  }
188 
189  if (extrapolating) {
190  if (var_type == 'T') {
191  // Assume a standard atmosphere -6.5 K km-1 lapse rate.
192  // Comparable to the WRF default, t_extrap_type=2.
193  amrex::Real depth_of_extrap_in_p = new_x_p[new_k]-orig_x_p[1];
194  amrex::Real avg_of_extrap_p = 0.5*(new_x_p[new_k]+orig_x_p[1]);
195  amrex::Real temp_extrap_starting_point = orig_y[1]*std::pow(orig_x_p[1]/100000.0, R_d/Cp_d);
196  amrex::Real dZdP = CRC_const1*CRC_const2*std::pow(avg_of_extrap_p/100.0, CRC_const2-1.0);
197  amrex::Real dZ = dZdP*(depth_of_extrap_in_p/100.0);
198  amrex::Real dT = dZ*CRC_const3;
199  new_y[new_k] = (temp_extrap_starting_point+dT)*std::pow(100000.0/new_x_p[new_k], R_d/Cp_d);
200  } else {
201  // Use a constant value below ground.
202  // Comparable to the WRF default, extrap_type=2.
203  new_y[new_k] = orig_y[1];
204  }
205  continue;
206  }
207 
208  if (order%2 != 0) {
209  if ((kl-((order+1)/2-1) >= 0) && (kr+((order+1)/2-1) <= orig_n-1)) {
210  int ksta = kl-((order/2)-1);
211  int kend = ksta+order;
212  amrex::Real new_x;
213  amrex::GpuArray<amrex::Real,2> orig_x_sub;
214  amrex::Real* orig_x_sub_p = orig_x_sub.data();
215  if (exp_interp) {
216  new_x = new_x_p[new_k];
217  orig_x_sub_p[0] = orig_x_p[ksta];
218  orig_x_sub_p[1] = orig_x_p[kend];
219  } else {
220  new_x = new_x_z[new_k];
221  orig_x_sub_p[0] = orig_x_z[ksta];
222  orig_x_sub_p[1] = orig_x_z[kend];
223  }
224  amrex::GpuArray<amrex::Real,2> orig_y_sub;
225  amrex::Real* orig_y_sub_p = orig_y_sub.data();
226  orig_y_sub_p[0] = orig_y[ksta];
227  orig_y_sub_p[1] = orig_y[kend];
228  lagrange_interp(order, orig_x_sub_p, orig_y_sub_p, new_x, new_y[new_k]);
229  } else {
230  amrex::Abort("metgrid initialization, lost in lagrange_setup (odd order)");
231  }
232  } else if ((order%2 == 0) && (new_k >= 1+vboundb) && (new_k < new_n-vboundt)) {
233  if ((kl-(order/2) >= 0) && (kr+order/2 <= orig_n-1)) {
234  amrex::Real new_y_l, new_y_r;
235  {
236  int ksta = kl-(order/2-1);
237  int kend = ksta+order;
238 #ifndef AMREX_USE_GPU
239  int ksize = kend-ksta;
240  if (debug) amrex::Print() << " (1a) new_x=" << new_x_z[new_k] << " new_x_p=" << new_x_p[new_k] << " kl=" << kl << " kr=" << kr <<" ksta=" << ksta << " kend=" << kend << std::endl;
241 #endif
242  amrex::Real new_x;
243  amrex::GpuArray<amrex::Real,256> orig_x_sub;
244  amrex::GpuArray<amrex::Real,256> orig_y_sub;
245  amrex::Real* orig_x_sub_p = orig_x_sub.data();
246  amrex::Real* orig_y_sub_p = orig_y_sub.data();
247  if (exp_interp) {
248  new_x = new_x_p[new_k];
249  for (int k=ksta; k <= kend; k++) { orig_x_sub_p[k-ksta] = orig_x_p[k]; }
250  } else {
251  new_x = new_x_z[new_k];
252  for (int k=ksta; k <= kend; k++) { orig_x_sub_p[k-ksta] = orig_x_z[k]; }
253  }
254  for (int k=ksta; k <= kend; k++) { orig_y_sub_p[k-ksta] = orig_y[k]; }
255 #ifndef AMREX_USE_GPU
256  if (debug) {
257  amrex::Print() << " orig_x_sub = [";
258  for (int k=0; k < ksize; k++) amrex::Print() << " " << orig_x_sub_p[k];
259  amrex::Print() << "]" << std::endl;
260  amrex::Print() << " orig_y_sub = [";
261  for (int k=0; k < ksize; k++) amrex::Print() << " " << orig_y_sub_p[k];
262  amrex::Print() << "]" << std::endl;
263  }
264 #endif
265  lagrange_interp(order, orig_x_sub_p, orig_y_sub_p, new_x, new_y_l);
266  }
267  {
268  int ksta = kl-order/2;
269  int kend = ksta+order;
270 #ifndef AMREX_USE_GPU
271  int ksize = kend-ksta;
272  if (debug) amrex::Print() << "new_k=" << new_k << " (1b) new_x=" << new_x_z[new_k] << " new_x_p=" << new_x_p[new_k] << " kl=" << kl << " kr=" << kr <<" ksta=" << ksta << " kend=" << kend << std::endl;
273 #endif
274  amrex::Real new_x;
275  amrex::GpuArray<amrex::Real,256> orig_x_sub;
276  amrex::GpuArray<amrex::Real,256> orig_y_sub;
277  amrex::Real* orig_x_sub_p = orig_x_sub.data();
278  amrex::Real* orig_y_sub_p = orig_y_sub.data();
279  if (exp_interp) {
280  new_x = new_x_p[new_k];
281  for (int k=ksta; k <= kend; k++) { orig_x_sub_p[k-ksta] = orig_x_p[k]; }
282  } else {
283  new_x = new_x_z[new_k];
284  for (int k=ksta; k <= kend; k++) { orig_x_sub_p[k-ksta] = orig_x_z[k]; }
285  }
286  for (int k=ksta; k <= kend; k++) { orig_y_sub_p[k-ksta] = orig_y[k]; }
287 #ifndef AMREX_USE_GPU
288  if (debug) {
289  amrex::Print() << " orig_x_sub = [";
290  for (int k=0; k < ksize; k++) amrex::Print() << " " << orig_x_sub_p[k];
291  amrex::Print() << "]" << std::endl;
292  amrex::Print() << " orig_y_sub = [";
293  for (int k=0; k < ksize; k++) amrex::Print() << " " << orig_y_sub_p[k];
294  amrex::Print() << "]" << std::endl;
295  }
296 #endif
297  lagrange_interp(order, orig_x_sub_p, orig_y_sub_p, new_x, new_y_r);
298  }
299  new_y[new_k] = 0.5*(new_y_l+new_y_r);
300 #ifndef AMREX_USE_GPU
301  if (debug) amrex::Print() << " new_y=" << new_y[new_k] << std::endl;
302 #endif
303  } else if ((kl-(order/2-1) >= 0) && (kr+order/2 <= orig_n-1)) {
304  int ksta = kl-(order/2-1);
305  int kend = ksta+order;
306 #ifndef AMREX_USE_GPU
307  int ksize = kend-ksta;
308  if (debug) amrex::Print() << " (2) new_x=" << new_x_z[new_k] << " new_x_p=" << new_x_p[new_k] << " kl=" << kl << " kr=" << kr <<" ksta=" << ksta << " kend=" << kend << std::endl;
309 #endif
310  amrex::Real new_x;
311  amrex::GpuArray<amrex::Real,256> orig_x_sub;
312  amrex::GpuArray<amrex::Real,256> orig_y_sub;
313  amrex::Real* orig_x_sub_p = orig_x_sub.data();
314  amrex::Real* orig_y_sub_p = orig_y_sub.data();
315  if (exp_interp) {
316  new_x = new_x_p[new_k];
317  for (int k=ksta; k <= kend; k++) { orig_x_sub_p[k-ksta] = orig_x_p[k]; }
318  } else {
319  new_x = new_x_z[new_k];
320  for (int k=ksta; k <= kend; k++) { orig_x_sub_p[k-ksta] = orig_x_z[k]; }
321  }
322  for (int k=ksta; k <= kend; k++) { orig_y_sub_p[k-ksta] = orig_y[k]; }
323 #ifndef AMREX_USE_GPU
324  if (debug) {
325  amrex::Print() << " orig_x_sub = [";
326  for (int k=0; k < ksize; k++) amrex::Print() << " " << orig_x_sub_p[k];
327  amrex::Print() << "]" << std::endl;
328  amrex::Print() << " orig_y_sub = [";
329  for (int k=0; k < ksize; k++) amrex::Print() << " " << orig_y_sub_p[k];
330  amrex::Print() << "]" << std::endl;
331  }
332 #endif
333  lagrange_interp(order, orig_x_sub_p, orig_y_sub_p, new_x, new_y[new_k]);
334 #ifndef AMREX_USE_GPU
335  if (debug) amrex::Print() << " new_y=" << new_y[new_k] << std::endl;
336 #endif
337  } else if ((kl-order/2 >= 0) && (kr+(order/2-1) <= orig_n-1)) {
338  int ksta = kl-order/2;
339  int kend = ksta+order;
340 #ifndef AMREX_USE_GPU
341  int ksize = kend-ksta;
342  if (debug) amrex::Print() << " (3) new_x=" << new_x_z[new_k] << " new_x_p=" << new_x_p[new_k] << " kl=" << kl << " kr=" << kr <<" ksta=" << ksta << " kend=" << kend << std::endl;
343 #endif
344  amrex::Real new_x;
345  amrex::GpuArray<amrex::Real,256> orig_x_sub;
346  amrex::GpuArray<amrex::Real,256> orig_y_sub;
347  amrex::Real* orig_x_sub_p = orig_x_sub.data();
348  amrex::Real* orig_y_sub_p = orig_y_sub.data();
349  if (exp_interp) {
350  new_x = new_x_p[new_k];
351  for (int k=ksta; k <= kend; k++) { orig_x_sub_p[k-ksta] = orig_x_p[k]; }
352  } else {
353  new_x = new_x_z[new_k];
354  for (int k=ksta; k <= kend; k++) { orig_x_sub_p[k-ksta] = orig_x_z[k]; }
355  }
356  for (int k=ksta; k <= kend; k++) { orig_y_sub_p[k-ksta] = orig_y[k]; }
357 #ifndef AMREX_USE_GPU
358  if (debug) {
359  amrex::Print() << " orig_x_sub = [";
360  for (int k=0; k < ksize; k++) amrex::Print() << " " << orig_x_sub_p[k];
361  amrex::Print() << "]" << std::endl;
362  amrex::Print() << " orig_y_sub = [";
363  for (int k=0; k < ksize; k++) amrex::Print() << " " << orig_y_sub_p[k];
364  amrex::Print() << "]" << std::endl;
365  }
366 #endif
367  lagrange_interp(order, orig_x_sub_p, orig_y_sub_p, new_x, new_y[new_k]);
368 #ifndef AMREX_USE_GPU
369  if (debug) amrex::Print() << " new_y=" << new_y[new_k] << std::endl;
370 #endif
371  } else {
372  amrex::Abort("metgrid initialization, lost in lagrange_setup (even order)");
373  }
374  } else {
375  // Linear interpolation.
376  int ksta = kl;
377  int kend = kr;
378 #ifndef AMREX_USE_GPU
379  int ksize = kend-ksta;
380  if (debug) amrex::Print() << " (4) new_x=" << new_x_z[new_k] << " new_x_p=" << new_x_p[new_k] << " kl=" << kl << " kr=" << kr <<" ksta=" << ksta << " kend=" << kend << std::endl;
381 #endif
382  amrex::Real new_x;
383  amrex::GpuArray<amrex::Real,256> orig_x_sub;
384  amrex::GpuArray<amrex::Real,256> orig_y_sub;
385  amrex::Real* orig_x_sub_p = orig_x_sub.data();
386  amrex::Real* orig_y_sub_p = orig_y_sub.data();
387  if (exp_interp) {
388  new_x = new_x_p[new_k];
389  for (int k=ksta; k <= kend; k++) { orig_x_sub_p[k-ksta] = orig_x_p[k]; }
390  } else {
391  new_x = new_x_z[new_k];
392  for (int k=ksta; k <= kend; k++) { orig_x_sub_p[k-ksta] = orig_x_z[k]; }
393  }
394  for (int k=ksta; k <= kend; k++) { orig_y_sub_p[k-ksta] = orig_y[k]; }
395 #ifndef AMREX_USE_GPU
396  if (debug) {
397  amrex::Print() << " orig_x_sub = [";
398  for (int k=0; k < ksize; k++) amrex::Print() << " " << orig_x_sub_p[k];
399  amrex::Print() << "]" << std::endl;
400  amrex::Print() << " orig_y_sub = [";
401  for (int k=0; k < ksize; k++) amrex::Print() << " " << orig_y_sub_p[k];
402  amrex::Print() << "]" << std::endl;
403  }
404 #endif
405  lagrange_interp(1, orig_x_sub_p, orig_y_sub_p, new_x, new_y[new_k]);
406 #ifndef AMREX_USE_GPU
407  if (debug) amrex::Print() << " new_y=" << new_y[new_k] << std::endl;
408 #endif
409  }
410  }
411 }
412 
413 AMREX_FORCE_INLINE
414 AMREX_GPU_DEVICE
415 void
416 calc_p_isothermal (const amrex::Real& z,
417  amrex::Real& p)
418 {
419  p = p_0*exp(-CONST_GRAV*z/(290.0*R_d));
420 }
421 
422 AMREX_FORCE_INLINE
423 AMREX_GPU_DEVICE
424 void
425 interpolate_column_metgrid (const bool& metgrid_use_below_sfc,
426  const bool& metgrid_use_sfc,
427  const bool& exp_interp,
428  const bool& metgrid_retain_sfc,
429  const amrex::Real& metgrid_proximity,
430  const int& metgrid_order,
431  const int& metgrid_force_sfc_k,
432  const int& i,
433  const int& j,
434  const int& src_comp,
435  const int& itime,
436  char var_type,
437  char stag,
438  const amrex::Array4<amrex::Real const>& orig_z_full,
439  const amrex::Array4<amrex::Real const>& orig_data,
440  const amrex::Array4<amrex::Real const>& new_z_full,
441  const amrex::Array4<amrex::Real>& new_data_full,
442  const bool& update_bc_data,
443  const amrex::Array4<amrex::Real>& bc_data_xlo,
444  const amrex::Array4<amrex::Real>& bc_data_xhi,
445  const amrex::Array4<amrex::Real>& bc_data_ylo,
446  const amrex::Array4<amrex::Real>& bc_data_yhi,
447  const amrex::Box& bx_xlo,
448  const amrex::Box& bx_xhi,
449  const amrex::Box& bx_ylo,
450  const amrex::Box& bx_yhi,
451  const amrex::Array4<const int>& mask)
452 {
453  // Here we closely follow WRF's vert_interp from
454  // dyn_em/module_initialize_real.F, although changes have been
455  // made to accommodate interpolation relative to height instead of
456  // pressure.
457  int imax_orig = amrex::ubound(amrex::Box(orig_data)).x;
458  int jmax_orig = amrex::ubound(amrex::Box(orig_data)).y;
459  int kmax_orig = amrex::ubound(amrex::Box(orig_data)).z;
460  int kmax_new = amrex::ubound(amrex::Box(new_z_full)).z;
461 
462  AMREX_ASSERT(kmax_orig < 256);
463  AMREX_ASSERT(kmax_new < 256);
464 
465  amrex::GpuArray<amrex::Real,256> new_z;
466  amrex::GpuArray<amrex::Real,256> new_p;
467  amrex::GpuArray<amrex::Real,256> new_data;
468  amrex::Real* new_z_p = new_z.data();
469  amrex::Real* new_p_p = new_p.data();
470  amrex::Real* new_data_p = new_data.data();
471  for (int k=0; k < kmax_new; k++) {
472  if (stag == 'X') {
473  new_z_p[k] = 0.25*(new_z_full(i,j,k)+new_z_full(i,j+1,k)+new_z_full(i,j,k+1)+new_z_full(i,j+1,k+1));
474  } else if (stag == 'Y') {
475  new_z_p[k] = 0.25*(new_z_full(i,j,k)+new_z_full(i+1,j,k)+new_z_full(i,j,k+1)+new_z_full(i+1,j,k+1));
476  } else if (stag == 'M') {
477  new_z_p[k] = 0.125*(new_z_full(i,j,k )+new_z_full(i,j+1,k )+new_z_full(i+1,j,k )+new_z_full(i+1,j+1,k )+
478  new_z_full(i,j,k+1)+new_z_full(i,j+1,k+1)+new_z_full(i+1,j,k+1)+new_z_full(i+1,j+1,k+1));
479  }
480  calc_p_isothermal(new_z_p[k], new_p_p[k]);
481  }
482 
483  amrex::GpuArray<amrex::Real,256> orig_z;
484  amrex::Real* orig_z_p = orig_z.data();
485  for (int k=0; k < kmax_orig; k++) {
486  if (stag == 'M') {
487  orig_z_p[k] = orig_z_full(i,j,k);
488  } else if (stag == 'X') {
489  if (i == 0) {
490  orig_z_p[k] = orig_z_full(i,j,k);
491  } else if (i == imax_orig) {
492  orig_z_p[k] = orig_z_full(imax_orig-1,j,k);
493  } else {
494  orig_z_p[k] = 0.5*(orig_z_full(i,j,k)+orig_z_full(i-1,j,k));
495  }
496  } else if (stag == 'Y') {
497  if (j == 0) {
498  orig_z_p[k] = orig_z_full(i,j,k);
499  } else if (j == jmax_orig) {
500  orig_z_p[k] = orig_z_full(i,jmax_orig-1,k);
501  } else {
502  orig_z_p[k] = 0.5*(orig_z_full(i,j,k)+orig_z_full(i,j-1,k));
503  }
504  }
505  }
506 
507  // Check if the data is top-down instead of bottom-up.
508  bool flip_data_required = false;
509  if (orig_z[1] > orig_z[kmax_orig-1]) flip_data_required = true;
510  if (flip_data_required) amrex::Abort("metgrid initialization flip_data_required. Not yet implemented.");
511 
512  // Search for the first level above the surface in the origin data.
513  // This is needed since the origin model topography will be
514  // different than the topography processed by WPS.
515  int k_above_sfc = 0;
516  for (int k=1; k < kmax_orig; k++) {
517  if (orig_z_p[k] > orig_z_p[0]) {
518  k_above_sfc = k;
519  break;
520  }
521  }
522 
523  int zap = 0;
524  int zap_below = 0;
525  int zap_above = 0;
526  int kend_order;
527  amrex::ignore_unused(zap,zap_above);
528  amrex::GpuArray<amrex::Real,256> ordered_z;
529  amrex::GpuArray<amrex::Real,256> ordered_data;
530  amrex::Real* ordered_z_p = ordered_z.data();
531  amrex::Real* ordered_data_p = ordered_data.data();
532  if (k_above_sfc > 1) {
533  // The levels are not monotonically increasing in height, so
534  // we sort and then make "artistic" quality control choices.
535  int count = 0;
536 
537  // Insert levels that are below the surface.
538  for (int k=1; k < k_above_sfc; k++) {
539  ordered_z_p[count] = orig_z_p[k];
540  ordered_data_p[count] = orig_data(i,j,k);
541  count++;
542  }
543 
544  // Check if the level that is nearest to and below the surface
545  // is "too close". If so, we'll ignore the upper level and keep
546  // the lower. Origin data is likely to be on pressure levels
547  // with higher spatial resolution near-surface, which supports
548  // the choice of eliminating levels that are "too close" in
549  // pressure-space. For simplicity, calculate delta P assuming a
550  // baroclinic atmosphere.
551  amrex::Real Pu, Pl;
552  calc_p_isothermal(orig_z_p[0], Pu);
553  calc_p_isothermal(ordered_z_p[count-1], Pl);
554  if (Pl-Pu < metgrid_proximity) {
555  count--;
556  zap = 1;
557  zap_below = 1;
558  }
559 
560  // Insert the surface level.
561  ordered_z_p[count] = orig_z_p[0];
562  ordered_data_p[count] = orig_data(i,j,0);
563  count++;
564 
565  // Quoting WRF's comments, the next level to use is at,
566  // "... ta da, the first level above the surface. I know, wow."
567  int knext = k_above_sfc;
568  // Conditionally more strongly use the surface data by removing
569  // levels between the surface and the height corresponding to a
570  // set number of ERF levels from the surface. This forces the
571  // interpolator to use the surface data up through a number of
572  // ERF levels from the surface.
573  if (metgrid_force_sfc_k > 0) {
574  for (int k=k_above_sfc; k < kmax_orig; k++) {
575  if (orig_z_p[k] > new_z_p[metgrid_force_sfc_k-1]) {
576  knext = k;
577  break;
578  } else {
579  zap++;
580  zap_above++;
581  }
582  }
583  }
584 
585  // Check if the level that is nearest to and above the surface
586  // is "too close". If so, we'll ignore that level.
587  calc_p_isothermal(orig_z_p[knext], Pu);
588  calc_p_isothermal(ordered_z_p[count-1], Pl);
589  if (Pl-Pu < metgrid_proximity) {
590  knext++;
591  zap++;
592  zap_above++;
593  }
594 
595  // Insert levels that are above the surface.
596  for (int k=knext; k < kmax_orig; k++) {
597  ordered_z_p[count] = orig_z_p[k];
598  ordered_data_p[count] = orig_data(i,j,k);
599  count++;
600  }
601 
602  kend_order = count;
603  } else {
604  // The surface is the lowest level in the origin data.
605 
606  // Insert the surface.
607  ordered_z_p[0] = orig_z[0];
608  ordered_data_p[0] = orig_data(i,j,0);
609 
610  // Similar to above, conditionally more strongly use the
611  // surface data.
612  int count = 1;
613  int knext = count;
614  if (metgrid_force_sfc_k > 0) {
615  for (int k=knext; k < kmax_orig; k++) {
616  if (orig_z_p[k] > new_z_p[metgrid_force_sfc_k]) {
617  knext = k;
618  break;
619  } else {
620  zap++;
621  zap_above++;
622  }
623  }
624  }
625 
626  // Insert the remaining levels, again ignoring levels that are
627  // "too close" to the prior valid level.
628  for (int k=knext; k < kmax_orig; k++) {
629  amrex::Real Pu, Pl;
630  calc_p_isothermal(orig_z_p[k], Pu);
631  calc_p_isothermal(ordered_z_p[count-1], Pl);
632  if (Pl-Pu < metgrid_proximity) {
633  zap++;
634  zap_above++;
635  continue;
636  }
637  ordered_z_p[count] = orig_z_p[k];
638  ordered_data_p[count] = orig_data(i,j,k);
639  count++;
640  }
641  kend_order = count;
642  }
643 
644  int ksta(0), kend(0);
645  if (metgrid_use_below_sfc && metgrid_use_sfc) {
646  // Use all levels.
647  ksta = 0;
648  kend = ksta+kend_order-1;
649  } else if (metgrid_use_below_sfc && !metgrid_use_sfc) {
650  // Use all levels except for the surface.
651  int ksfc = 0;
652  for (int k=0; k < kmax_orig; k++) {
653  if (ordered_z_p[k] == orig_z_p[0]) {
654  ksfc = k;
655  break;
656  }
657  }
658  for (int k=ksfc; k < kmax_orig-1; k++) {
659  ordered_z_p[k] = ordered_z_p[k+1];
660  ordered_data_p[k] = ordered_data_p[k+1];
661  }
662  ksta = 0;
663  kend = ksta+kend_order-2;
664  } else if (!metgrid_use_below_sfc && metgrid_use_sfc) {
665  // Use all levels above and including the surface.
666  int kcount = k_above_sfc-1-zap_below;
667  int count = 0;
668  for (int k=0; k < kmax_orig; k++) {
669  if (ordered_z_p[kcount] == orig_z_p[k]) {
670  kcount++;
671  count++;
672  }
673  }
674  ksta = k_above_sfc-1-zap_below;
675  kend = ksta+count-1;
676  } else {
677  // We shouldn't be in here!
678  amrex::Abort("metgrid initialization, !use_levels below_ground && !metgrid_use_sfc");
679  }
680 
681  // Insert the level of maximum winds.
682 // amrex::Real maxw_above_this_level = 30000.0;
683 // amrex::Real maxw_horiz_pres_diff = 5000.0;
684 // if ((flag_maxw == 1) && (use_maxw)) {
685 // amrex::Abort("metgrid initialization, use_maxw not yet implemented");
686 // }
687 
688  // Insert the level of the tropopause.
689 // amrex::Real trop_horiz_pres_diff = 5000.0;
690 // if ((flag_trop == 1) && (use_trop)) {
691 // amrex::Abort("metgrid initialization, use_trop not yet implemented");
692 // }
693 
694  amrex::GpuArray<amrex::Real,256> ordered_p;
695  amrex::Real* ordered_p_p = ordered_p.data();
696  for (int k=0; k < kend_order; k++) {
697  calc_p_isothermal(ordered_z_p[k], ordered_p_p[k]);
698  }
699 
700  int new_n = 0;
701  int zap_final = 0;
702  amrex::GpuArray<amrex::Real,256> final_z;
703  amrex::GpuArray<amrex::Real,256> final_p;
704  amrex::GpuArray<amrex::Real,256> final_data;
705  amrex::Real* final_z_p = final_z.data();
706  amrex::Real* final_p_p = final_p.data();
707  amrex::Real* final_data_p = final_data.data();
708  final_z_p[0] = ordered_z[ksta];
709  final_p_p[0] = ordered_p[ksta];
710  final_data_p[0] = ordered_data[ksta];
711  for (int k=ksta+1; k <= kend; k++) {
712  if ((final_p_p[new_n]-ordered_p_p[k]) < metgrid_proximity) {
713  zap_final++;
714  } else {
715  new_n++;
716  final_z_p[new_n] = ordered_z_p[k];
717  final_p_p[new_n] = ordered_p_p[k];
718  final_data_p[new_n] = ordered_data_p[k];
719  }
720  }
721  kend -= zap_final;
722 
723  // Call the interpolator.
724  lagrange_setup(var_type,
725  exp_interp,
726  kend-ksta,
727  kmax_new,
728  metgrid_order,
729  final_z_p,
730  final_p_p,
731  final_data_p,
732  new_z_p,
733  new_p_p,
734  new_data_p);
735 
736  // Optionally replace the lowest level of data with the surface
737  // field from the origin data.
738  if (metgrid_retain_sfc) new_data[0] = ordered_data[0];
739 
740  // Save the interpolated data.
741  for (int k=0; k < kmax_new; k++) {
742  if (mask(i,j,k) && update_bc_data && bx_xlo.contains(i,j,k)) bc_data_xlo(i,j,k,0) = new_data[k];
743  if (mask(i,j,k) && update_bc_data && bx_xhi.contains(i,j,k)) bc_data_xhi(i,j,k,0) = new_data[k];
744  if (mask(i,j,k) && update_bc_data && bx_ylo.contains(i,j,k)) bc_data_ylo(i,j,k,0) = new_data[k];
745  if (mask(i,j,k) && update_bc_data && bx_yhi.contains(i,j,k)) bc_data_yhi(i,j,k,0) = new_data[k];
746  if (itime == 0) new_data_full(i,j,k,src_comp) = new_data[k];
747  }
748 
749 }
750 
751 AMREX_FORCE_INLINE
752 AMREX_GPU_DEVICE
753 amrex::Real
754 interpolate_column_metgrid_linear (const int& i,
755  const int& j,
756  const int& k,
757  char stag,
758  int src_comp,
759  const amrex::Array4<amrex::Real const>& orig_z,
760  const amrex::Array4<amrex::Real const>& orig_data,
761  const amrex::Array4<amrex::Real const>& new_z)
762 {
763  // This subroutine is a bit ham-handed and can be cleaned up later.
764  int imax_orig = amrex::ubound(amrex::Box(orig_data)).x;
765  int jmax_orig = amrex::ubound(amrex::Box(orig_data)).y;
766  int kmax_orig = amrex::ubound(amrex::Box(orig_data)).z;
767 
768  amrex::Real z;
769  if (stag == 'X') {
770  z = 0.25*(new_z(i,j,k)+new_z(i,j+1,k)+new_z(i,j,k+1)+new_z(i,j+1,k+1));
771  }
772  else if (stag == 'Y') {
773  z = 0.25*(new_z(i,j,k)+new_z(i+1,j,k)+new_z(i,j,k+1)+new_z(i+1,j,k+1));
774  }
775  else if (stag == 'M') {
776  z = 0.125*(new_z(i,j,k )+new_z(i,j+1,k )+new_z(i+1,j,k )+new_z(i+1,j+1,k )+
777  new_z(i,j,k+1)+new_z(i,j+1,k+1)+new_z(i+1,j,k+1)+new_z(i+1,j+1,k+1));
778  }
779 
780  amrex::Real z0, z1;
781  int klow = -1;
782  int khi0 = -1;
783  amrex::Real dzlow = 1.0e12;
784  amrex::Real dzhi0 = -1.0e12;
785  for (int kk = 0; kk < kmax_orig; kk++) {
786  amrex::Real orig_z_stag = 0.0;
787  if (stag == 'M') {
788  orig_z_stag = orig_z(i,j,kk);
789  }
790  if (stag == 'X') {
791  if (i == 0) {
792  orig_z_stag = orig_z(i,j,kk);
793  }
794  else if (i == imax_orig) {
795  orig_z_stag = orig_z(imax_orig-1,j,kk);
796  }
797  else {
798  orig_z_stag = 0.5*(orig_z(i,j,kk)+orig_z(i-1,j,kk));
799  }
800  }
801  else if (stag == 'Y') {
802  if (j == 0) {
803  orig_z_stag = orig_z(i,j,kk);
804  }
805  else if (j == jmax_orig) {
806  orig_z_stag = orig_z(i,jmax_orig-1,kk);
807  }
808  else {
809  orig_z_stag = 0.5*(orig_z(i,j,kk)+orig_z(i,j-1,kk));
810  }
811  }
812 
813  amrex::Real dz = z - orig_z_stag;
814  if ((dz < 0.0) && (dz > dzhi0)) {
815  dzhi0 = dz;
816  khi0 = kk;
817  z1 = orig_z_stag;
818  }
819  if ((dz >= 0.0) && (dz < dzlow)) {
820  dzlow = dz;
821  klow = kk;
822  z0 = orig_z_stag;
823  }
824  } // kk
825 
826  // extrapolate below the bottom surface
827  if (klow == -1) {
828  int khi1 = -1;
829  amrex::Real dzhi1 = -1.0e12;
830  for (int kk = 0; kk < kmax_orig; kk++) {
831  amrex::Real orig_z_stag = 0.0;
832  if (stag == 'M') {
833  orig_z_stag = orig_z(i,j,kk);
834  }
835  else if (stag == 'X') {
836  if (i == 0) {
837  orig_z_stag = orig_z(i,j,kk);
838  }
839  else if (i == imax_orig) {
840  orig_z_stag = orig_z(imax_orig-1,j,kk);
841  }
842  else {
843  orig_z_stag = 0.5*(orig_z(i,j,kk)+orig_z(i-1,j,kk));
844  }
845  }
846  else if (stag == 'Y') {
847  if (j == 0) {
848  orig_z_stag = orig_z(i,j,kk);
849  }
850  else if (j == jmax_orig) {
851  orig_z_stag = orig_z(i,jmax_orig-1,kk);
852  }
853  else {
854  orig_z_stag = 0.5*(orig_z(i,j,kk)+orig_z(i,j-1,kk));
855  }
856  }
857  amrex::Real dz = z - orig_z_stag;
858  if ((dz < 0.0) && (dz > dzhi1) && (kk != khi0)) {
859  dzhi1 = dz;
860  khi1 = kk;
861  z1 = orig_z_stag;
862  }
863  }
864  amrex::Real y0 = orig_data(i,j,khi0,src_comp);
865  amrex::Real y1 = orig_data(i,j,khi1,src_comp);
866  return ( y0-(y1-y0)/(z1-z0)*(z0-z) );
867 
868  // Extrapolate above the top surface
869  } else if (khi0 == -1) {
870  khi0 = klow - 1;
871  int khi1 = klow;
872  if (stag == 'M') {
873  z0 = orig_z(i,j,khi0);
874  }
875  else if (stag == 'X') {
876  if (i == 0) {
877  z0 = orig_z(i,j,khi0);
878  }
879  else if (i == imax_orig) {
880  z0 = orig_z(imax_orig-1,j,khi0);
881  }
882  else {
883  z0 = 0.5*(orig_z(i,j,khi0)+orig_z(i-1,j,khi0));
884  }
885  }
886  else if (stag == 'Y') {
887  if (j == 0) {
888  z0 = orig_z(i,j,khi0);
889  }
890  else if (j == jmax_orig) {
891  z0 = orig_z(i,jmax_orig-1,khi0);
892  }
893  else {
894  z0 = 0.5*(orig_z(i,j,khi0)+orig_z(i,j-1,khi0));
895  }
896  }
897  amrex::Real y0 = orig_data(i,j,khi0,src_comp);
898  amrex::Real y1 = orig_data(i,j,khi1,src_comp);
899  return ( y0+(y1-y0)/(z1-z0)*(z-z0) );
900  } else {
901  // interpolate
902  amrex::Real y0 = orig_data(i,j,klow,src_comp);
903  amrex::Real y1 = orig_data(i,j,khi0,src_comp);
904  return ( y0+(y1-y0)/(z1-z0)*(z-z0) );
905 
906  }
907 }
908 
909 AMREX_FORCE_INLINE
910 AMREX_GPU_DEVICE
911 void
912 rh_to_mxrat (int i,
913  int j,
914  int k,
915  const amrex::Array4<amrex::Real const>& rhum,
916  const amrex::Array4<amrex::Real const>& temp,
917  const amrex::Array4<amrex::Real const>& pres,
918  const amrex::Array4<amrex::Real>& mxrat)
919 {
920  amrex::Real qv_max_p_safe = 10000.0; // WRF default value
921  amrex::Real qv_max_flag = 1.0e-5; // WRF default value
922  amrex::Real qv_max_value = 3.0e-6; // WRF default value
923  amrex::Real qv_min_p_safe = 110000.0; // WRF default value
924  amrex::Real qv_min_flag = 1.0e-6; // WRF default value
925  amrex::Real qv_min_value = 1.0e-6; // WRF default value
926  amrex::Real eps = 0.622;
927  amrex::Real svp1 = 0.6112;
928  amrex::Real svp2 = 17.67;
929  amrex::Real svp3 = 29.65;
930  amrex::Real svpt0 = 273.15;
931  // WRF's method when model_config_rec%rh2qv_wrt_liquid=.true. (default behavior)
932  if (temp(i,j,k) != 0.0) {
933  amrex::Real es=0.01*rhum(i,j,k)*svp1*10.0*exp(svp2*(temp(i,j,k)-svpt0)/(temp(i,j,k)-svp3));
934  if (es >= pres(i,j,k)/100.0) {
935  // vapor pressure exceeds total pressure
936  mxrat(i,j,k) = std::pow(10.0, -6);
937  }
938  else {
939  mxrat(i,j,k) = amrex::max(eps*es/(pres(i,j,k)/100.0-es), 1.0e-6);
940  }
941  }
942  else {
943  // I don't know why there's a fringe case handled in WRF where T is absolute zero...
944  // Let's just deal with it here in case we also end up needing it.
945  mxrat(i,j,k) = 1.0e-6;
946  }
947  // See the below comment from WRF dyn_em/module_initialize_real.F rh_to_mxrat1.
948  // For pressures above a defined level, reasonable Qv values should be
949  // a certain value or smaller. If they are larger than this, the input data
950  // probably had "missing" RH, and we filled in some values. This is an
951  // attempt to catch those. Also, set the minimum value for the entire
952  // domain that is above the selected pressure level.
953  if (pres(i,j,k) < qv_max_p_safe) {
954  if (mxrat(i,j,k) > qv_max_flag) {
955  mxrat(i,j,k) = qv_max_value;
956  }
957  }
958  if (pres(i,j,k) < qv_min_p_safe) {
959  if (mxrat(i,j,k) < qv_min_flag) {
960  mxrat(i,j,k) = qv_min_value;
961  }
962  }
963 }
964 #endif
Cp_d
constexpr amrex::Real Cp_d
Definition: ERF_Constants.H:12
p_0
constexpr amrex::Real p_0
Definition: ERF_Constants.H:18
CONST_GRAV
constexpr amrex::Real CONST_GRAV
Definition: ERF_Constants.H:21
R_d
constexpr amrex::Real R_d
Definition: ERF_Constants.H:10
Coord::y
@ y
Coord::x
@ x
Coord::z
@ z
ERF_EOS.H
ERF_HSEUtils.H
interpolate_column_metgrid
AMREX_FORCE_INLINE AMREX_GPU_DEVICE void interpolate_column_metgrid(const bool &metgrid_use_below_sfc, const bool &metgrid_use_sfc, const bool &exp_interp, const bool &metgrid_retain_sfc, const amrex::Real &metgrid_proximity, const int &metgrid_order, const int &metgrid_force_sfc_k, const int &i, const int &j, const int &src_comp, const int &itime, char var_type, char stag, const amrex::Array4< amrex::Real const > &orig_z_full, const amrex::Array4< amrex::Real const > &orig_data, const amrex::Array4< amrex::Real const > &new_z_full, const amrex::Array4< amrex::Real > &new_data_full, const bool &update_bc_data, const amrex::Array4< amrex::Real > &bc_data_xlo, const amrex::Array4< amrex::Real > &bc_data_xhi, const amrex::Array4< amrex::Real > &bc_data_ylo, const amrex::Array4< amrex::Real > &bc_data_yhi, const amrex::Box &bx_xlo, const amrex::Box &bx_xhi, const amrex::Box &bx_ylo, const amrex::Box &bx_yhi, const amrex::Array4< const int > &mask)
Definition: ERF_MetgridUtils.H:425
init_base_state_from_metgrid
void init_base_state_from_metgrid(const bool use_moisture, const bool metgrid_debug_psfc, const amrex::Real l_rdOcp, const amrex::Box &valid_bx, const amrex::Vector< int > &flag_psfc, amrex::FArrayBox &state, amrex::FArrayBox &r_hse_fab, amrex::FArrayBox &p_hse_fab, amrex::FArrayBox &pi_hse_fab, amrex::FArrayBox &th_hse_fab, amrex::FArrayBox &z_phys_cc_fab, const amrex::Vector< amrex::FArrayBox > &NC_psfc_fab)
lagrange_interp
AMREX_FORCE_INLINE AMREX_GPU_DEVICE void lagrange_interp(const int &order, amrex::Real *x, amrex::Real *y, amrex::Real &new_x, amrex::Real &new_y)
Definition: ERF_MetgridUtils.H:117
interpolate_column_metgrid_linear
AMREX_FORCE_INLINE AMREX_GPU_DEVICE amrex::Real interpolate_column_metgrid_linear(const int &i, const int &j, const int &k, char stag, int src_comp, const amrex::Array4< amrex::Real const > &orig_z, const amrex::Array4< amrex::Real const > &orig_data, const amrex::Array4< amrex::Real const > &new_z)
Definition: ERF_MetgridUtils.H:754
read_from_metgrid
void read_from_metgrid(int lev, const amrex::Box &domain, const std::string &fname, std::string &NC_dateTime, amrex::Real &NC_epochTime, int &flag_psfc, int &flag_msf, int &flag_sst, int &flag_lmask, int &NC_nx, int &NC_ny, amrex::Real &NC_dx, amrex::Real &NC_dy, amrex::FArrayBox &NC_xvel_fab, amrex::FArrayBox &NC_yvel_fab, amrex::FArrayBox &NC_temp_fab, amrex::FArrayBox &NC_rhum_fab, amrex::FArrayBox &NC_pres_fab, amrex::FArrayBox &NC_ght_fab, amrex::FArrayBox &NC_hgt_fab, amrex::FArrayBox &NC_psfc_fab, amrex::FArrayBox &NC_msfu_fab, amrex::FArrayBox &NC_msfv_fab, amrex::FArrayBox &NC_msfm_fab, amrex::FArrayBox &NC_sst_fab, amrex::FArrayBox &NC_LAT_fab, amrex::FArrayBox &NC_LON_fab, amrex::IArrayBox &NC_lmask_iab, amrex::Real &Latitude, amrex::Real &Longitude, amrex::Geometry &geom)
calc_p_isothermal
AMREX_FORCE_INLINE AMREX_GPU_DEVICE void calc_p_isothermal(const amrex::Real &z, amrex::Real &p)
Definition: ERF_MetgridUtils.H:416
init_state_from_metgrid
void init_state_from_metgrid(const bool use_moisture, const bool interp_theta, const bool metgrid_debug_quiescent, const bool metgrid_debug_isothermal, const bool metgrid_debug_dry, const bool metgrid_basic_linear, const bool metgrid_use_below_sfc, const bool metgrid_use_sfc, const bool metgrid_retain_sfc, const amrex::Real metgrid_proximity, const int metgrid_order, const int metgrid_metgrid_force_sfc_k, const amrex::Real l_rdOcp, amrex::Box &tbxc, amrex::Box &tbxu, amrex::Box &tbxv, amrex::FArrayBox &state_fab, amrex::FArrayBox &x_vel_fab, amrex::FArrayBox &y_vel_fab, amrex::FArrayBox &z_vel_fab, amrex::FArrayBox &z_phys_nd_fab, const amrex::Vector< amrex::FArrayBox > &NC_hgt_fab, const amrex::Vector< amrex::FArrayBox > &NC_ght_fab, const amrex::Vector< amrex::FArrayBox > &NC_xvel_fab, const amrex::Vector< amrex::FArrayBox > &NC_yvel_fab, const amrex::Vector< amrex::FArrayBox > &NC_temp_fab, const amrex::Vector< amrex::FArrayBox > &NC_rhum_fab, const amrex::Vector< amrex::FArrayBox > &NC_pres_fab, amrex::FArrayBox &p_interp_fab, amrex::FArrayBox &t_interp_fab, amrex::FArrayBox &theta_fab, amrex::FArrayBox &mxrat_fab, amrex::Vector< amrex::Vector< amrex::FArrayBox >> &fabs_for_bcs_xlo, amrex::Vector< amrex::Vector< amrex::FArrayBox >> &fabs_for_bcs_xhi, amrex::Vector< amrex::Vector< amrex::FArrayBox >> &fabs_for_bcs_ylo, amrex::Vector< amrex::Vector< amrex::FArrayBox >> &fabs_for_bcs_yhi, const amrex::Array4< const int > &mask_c_arr, const amrex::Array4< const int > &mask_u_arr, const amrex::Array4< const int > &mask_v_arr)
rh_to_mxrat
AMREX_FORCE_INLINE AMREX_GPU_DEVICE void rh_to_mxrat(int i, int j, int k, const amrex::Array4< amrex::Real const > &rhum, const amrex::Array4< amrex::Real const > &temp, const amrex::Array4< amrex::Real const > &pres, const amrex::Array4< amrex::Real > &mxrat)
Definition: ERF_MetgridUtils.H:912
init_terrain_from_metgrid
void init_terrain_from_metgrid(amrex::FArrayBox &z_phys_nd_fab, const amrex::Vector< amrex::FArrayBox > &NC_hgt_fab)
lagrange_setup
AMREX_FORCE_INLINE AMREX_GPU_DEVICE void lagrange_setup(char var_type, const bool &exp_interp, const int &orig_n, const int &new_n, const int &order, amrex::Real *orig_x_z, amrex::Real *orig_x_p, amrex::Real *orig_y, amrex::Real *new_x_z, amrex::Real *new_x_p, amrex::Real *new_y)
Definition: ERF_MetgridUtils.H:147
init_msfs_from_metgrid
void init_msfs_from_metgrid(const bool metgrid_debug_msf, amrex::FArrayBox &msfu_fab, amrex::FArrayBox &msfv_fab, amrex::FArrayBox &msfm_fab, const int &flag_msf, const amrex::Vector< amrex::FArrayBox > &NC_MSFU_fab, const amrex::Vector< amrex::FArrayBox > &NC_MSFV_fab, const amrex::Vector< amrex::FArrayBox > &NC_MSFM_fab)
ERF_ProbCommon.H
ERF_Utils.H
MicVar_Kess::pres
@ pres
Definition: ERF_Kessler.H:33