ERF
Energy Research and Forecasting: An Atmospheric Modeling Code
Functions
ERF_InitCustomTerrain.cpp File Reference
#include "AMReX_ParmParse.H"
#include "ERF_Constants.H"
#include "ERF_TerrainMetrics.H"
Include dependency graph for ERF_InitCustomTerrain.cpp:

Functions

void init_my_custom_terrain (const Geometry &geom, FArrayBox &terrain_fab, const Real &time)
 

Function Documentation

◆ init_my_custom_terrain()

void init_my_custom_terrain ( const Geometry &  geom,
FArrayBox &  terrain_fab,
const Real &  time 
)
15 {
16  //
17  // We put this here as a convenience for testing the map factor implementation
18  // Note that these factors must match those in Source/ERF_MakeNewArrays.cpp
19  //
20  ParmParse pp("erf");
21  bool test_mapfactor = false;
22  pp.query("test_mapfactor",test_mapfactor);
23 
24  Real mf_m;
25  if (test_mapfactor) {
26  mf_m = 0.5;
27  } else {
28  mf_m = 1.;
29  }
30 
31  // Domain cell size and real bounds
32  auto dx = geom.CellSizeArray();
33  auto ProbLoArr = geom.ProbLoArray();
34  auto ProbHiArr = geom.ProbHiArray();
35 
36  const amrex::Box& domain = geom.Domain();
37  int domlo_x = domain.smallEnd(0); int domhi_x = domain.bigEnd(0) + 1;
38  int domlo_y = domain.smallEnd(1); int domhi_y = domain.bigEnd(1) + 1;
39  int domlo_z = domain.smallEnd(2);
40 
41  // User function parameters
42  Real a = 0.5;
43  Real num = 8. * a * a * a;
44  Real xcen = 0.5 * (ProbLoArr[0] + ProbHiArr[0]) / mf_m;
45  Real ycen = 0.5 * (ProbLoArr[1] + ProbHiArr[1]) / mf_m;
46 
47  // Populate bottom plane
48  int k0 = domlo_z;
49 
50  std::string custom_terrain_type = "None";
51  ParmParse pp_prob("prob"); pp_prob.query("custom_terrain_type", custom_terrain_type);
52  amrex::Print() << "IN CUSTOM TERRAIN WITH TYPE = " << custom_terrain_type << std::endl;
53 
54  amrex::Box zbx = terrain_fab.box();
55  if (zbx.smallEnd(2) <= k0)
56  {
57  amrex::Array4<Real> const& z_arr = terrain_fab.array();
58 
59  if (custom_terrain_type == "WoA") {
60 
61  // Default to x-direction
62  int dir = 0; pp_prob.query("dir", dir);
63 
64  Real L = 100.0; pp_prob.query("L" , L);
65  Real z_offset = 0.0; pp_prob.query("z_offset" , z_offset);
66 
67  // If hm is nonzero, then use alternate hill definition
68  Real hm = 0.0; pp_prob.query("hmax" , hm);
69 
70  // This is a 2D hill with variation in only the x-direction
71  if (dir == 0) {
72  ParallelFor(zbx, [=] AMREX_GPU_DEVICE (int i, int j, int)
73  {
74  // Clip indices for ghost-cells
75  int ii = amrex::min(amrex::max(i,domlo_x),domhi_x);
76 
77  // Location of nodes
78  Real x = (ProbLoArr[0] + ii * dx[0] - xcen) * mf_m;
79 
80  // WoA Hill in x-direction
81  if (hm==0) {
82  z_arr(i,j,k0) = num / (x*x + 4 * a * a);
83  } else {
84  Real x_L = x / L;
85  z_arr(i,j,k0) = hm / (1 + x_L*x_L) + z_offset;
86  }
87  });
88  } else if (dir == 1) {
89  ParallelFor(zbx, [=] AMREX_GPU_DEVICE (int i, int j, int)
90  {
91  // Clip indices for ghost-cells
92  int jj = amrex::min(amrex::max(j,domlo_y),domhi_y);
93 
94  // Location of nodes
95  Real y = (ProbLoArr[1] + jj * dx[1] - ycen) * mf_m;
96 
97  // WoA Hill in y-direction
98  if (hm==0) {
99  z_arr(i,j,k0) = num / (y*y + 4.0 * a * a);
100  } else {
101  Real y_L = y / L;
102  z_arr(i,j,k0) = hm / (1.0 + y_L*y_L) + z_offset;
103  }
104  });
105  } else if (dir == 2) {
106  ParallelFor(zbx, [=] AMREX_GPU_DEVICE (int i, int j, int)
107  {
108  // Clip indices for ghost-cells
109  int ii = amrex::min(amrex::max(i,domlo_x),domhi_x);
110  int jj = amrex::min(amrex::max(j,domlo_y),domhi_y);
111 
112  // Location of nodes
113  Real x = (ProbLoArr[0] + ii * dx[0] - xcen) * mf_m;
114  Real y = (ProbLoArr[1] + jj * dx[1] - ycen) * mf_m;
115  Real r = std::sqrt(x*x + y*y);
116 
117  // WoA Hill in radial direction
118  if (hm==0) {
119  z_arr(i,j,k0) = num / (r*r + 4.0 * a * a);
120  } else {
121  Real r_L = r / L;
122  z_arr(i,j,k0) = hm / (1.0 + r_L*r_L) + z_offset;
123  }
124  });
125  } else {
126  amrex::Abort("Unknown dir in ERF_Prob.cpp");
127  }
128 
129  } else if (custom_terrain_type == "ScharMountain") {
130 
131  Real asq = 5000.0 * 5000.0;
132  Real Hm = 250.0;
133  Real lambda = 4000.0;
134 
135  ParallelFor(zbx, [=] AMREX_GPU_DEVICE (int i, int j, int)
136  {
137  // Clip indices for ghost-cells
138  int ii = amrex::min(amrex::max(i,domlo_x),domhi_x);
139 
140  // Location of nodes
141  Real x = (ProbLoArr[0] + ii * dx[0] - xcen);
142 
143  Real cosx = cos(PI * x / lambda);
144 
145  z_arr(i,j,k0) = Hm * std::exp(-x*x/asq) * cosx * cosx;
146  });
147 
148  } else if (custom_terrain_type == "HalfCylinder") {
149 
150  Real asq = 0.5 * 0.5;
151 
152  ParallelFor(zbx, [=] AMREX_GPU_DEVICE (int i, int j, int)
153  {
154  // Clip indices for ghost-cells
155  int ii = amrex::min(amrex::max(i,domlo_x),domhi_x);
156 
157  // Location of nodes
158  Real x = (ProbLoArr[0] + ii * dx[0] - xcen);
159 
160  Real rsq = x*x;
161 
162  if (rsq < asq) {
163  z_arr(i,j,k0) = pow(asq - rsq, 0.5);
164  } else {
165  z_arr(i,j,k0) = 0.0;
166  }
167  });
168 
169  } else if (custom_terrain_type == "Hemisphere") {
170 
171  Real asq = 0.5 * 0.5;
172 
173  ParallelFor(zbx, [=] AMREX_GPU_DEVICE (int i, int j, int)
174  {
175  // Clip indices for ghost-cells
176  int ii = amrex::min(amrex::max(i,domlo_x),domhi_x);
177  int jj = amrex::min(amrex::max(j,domlo_y),domhi_y);
178 
179  // Location of nodes
180  Real x = (ProbLoArr[0] + ii * dx[0] - xcen);
181  Real y = (ProbLoArr[1] + jj * dx[1] - ycen);
182 
183  Real rsq = x*x + y*y;
184 
185  if (rsq < asq) {
186  z_arr(i,j,k0) = std::pow(asq-rsq, 0.5);
187  } else {
188  z_arr(i,j,k0) = 0.0;
189  }
190  });
191 
192  } else if (custom_terrain_type == "MovingSineWave") {
193 
194  Real Ampl = 0.0; pp_prob.query("Ampl", Ampl);
195  Real wavelength = 100.; pp_prob.query("wavelength", wavelength);
196 
197  Real kp = 2.0 * PI / wavelength;
198  Real g = CONST_GRAV;
199  Real omega = std::sqrt(g * kp);
200 
201  ParallelFor(zbx, [=] AMREX_GPU_DEVICE (int i, int j, int)
202  {
203  // Clip indices for ghost-cells
204  int ii = amrex::min(amrex::max(i,domlo_x),domhi_x);
205 
206  // Location of nodes
207  Real x = ii * dx[0];
208 
209  // Wave height
210  Real height = Ampl * std::sin(kp * x - omega * time);
211 
212  // Populate terrain height
213  z_arr(i,j,0) = height;
214  });
215 
216  } else if (custom_terrain_type == "WindFarmTest") {
217 
218  ParallelFor(zbx, [=] AMREX_GPU_DEVICE (int i, int j, int)
219  {
220  // Clip indices for ghost-cells
221  int ii = amrex::min(amrex::max(i,domlo_x),domhi_x);
222  int jj = amrex::min(amrex::max(j,domlo_y),domhi_y);
223 
224  // Location of nodes
225  Real x = (ProbLoArr[0] + ii * dx[0] - xcen);
226  Real y = (ProbLoArr[1] + jj * dx[1] - ycen);
227 
228  Real x_L = x/100.0;
229  Real y_L = y/100.0;
230 
231  z_arr(i,j,k0) = 100.0 / (1.0 + x_L*x_L + y_L*y_L);
232  });
233 
234  } else if (custom_terrain_type == "RaisedFlat") {
235  Real z_offset = 0.0; pp_prob.query("z_offset" , z_offset);
236  ParallelFor(zbx, [=] AMREX_GPU_DEVICE (int i, int j, int)
237  {
238  z_arr(i,j,k0) = z_offset;
239  });
240 
241  } else if (custom_terrain_type == "None") {
242  ParallelFor(zbx, [=] AMREX_GPU_DEVICE (int i, int j, int)
243  {
244  z_arr(i,j,k0) = 0.0;
245  });
246  } else {
247  Abort("Don't know this custom_terrain_type");
248  }
249  }
250 }
PI
constexpr amrex::Real PI
Definition: ERF_Constants.H:6
CONST_GRAV
constexpr amrex::Real CONST_GRAV
Definition: ERF_Constants.H:21
num
@ num
Definition: ERF_DataStruct.H:24
Coord::y
@ y
Coord::x
@ x
pp
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE amrex::Real pp(amrex::Real y)
Definition: ERF_MicrophysicsUtils.H:233
zbx
const Box zbx
Definition: ERF_SetupDiff.H:9
Real
amrex::Real Real
Definition: ERF_ShocInterface.H:19
MicVar_Morr::omega
@ omega
Definition: ERF_Morrison.H:53
module_model_constants::g
real(c_double), parameter g
Definition: ERF_module_model_constants.F90:19
Here is the call graph for this function: