REMORA
Regional Modeling of Oceans Refined Adaptively
Loading...
Searching...
No Matches
REMORA_advance_3d.cpp
Go to the documentation of this file.
1#include <REMORA.H>
2
3using namespace amrex;
4
5/** Corresponds to step3d_uv.F and step3d_t.F in ROMS
6 *
7 * @param[in ] lev level of refinement (coarsest level is 0)
8 * @param[inout] mf_cons scalar variables (temp, salt, scalar, etc)
9 * @param[inout] mf_u u-velocity
10 * @param[inout] mf_v v-velocity
11 * @param[inout] mf_sstore intermediate-step scalar variables
12 * @param[inout] mf_ru RHS of total u-velocity
13 * @param[inout] mf_rv RHS of total v-velocity
14 * @param[in ] mf_DU_avg1 time-averaged u-flux for 2D equations
15 * @param[in ] mf_DU_avg2 time-averaged u-flux for 3D equation coupling
16 * @param[in ] mf_DV_avg1 time-averaged v-flux for 2D equations
17 * @param[in ] mf_DV_avg2 time-averaged v-flux for 3D equation coupling
18 * @param[inout] mf_ubar vertically integrated u-momentum
19 * @param[inout] mf_vbar vertically integrated v-momentum
20 * @param[inout] mf_Akv vertical viscosity coefficient
21 * @param[inout] mf_Akt vertical diffusivity coefficient
22 * @param[inout] mf_Hz vertical height of cells
23 * @param[inout] mf_Huon u-volume flux
24 * @param[inout] mf_Huon v-volume flux
25 * @param[inout] mf_z_w vertical coordinates on w points
26 * @param[in ] mf_h bathymetry
27 * @param[in ] mf_pm 1 / dx
28 * @param[in ] mf_pn 1 / dy
29 * @param[in ] mf_msku land-sea mask at u-points
30 * @param[in ] mf_mskv land-sea mask at v-points
31 * @param[in ] N number of vertical levels
32 * @param[in ] dt_lev time step at this refinement level
33 */
34void
35REMORA::advance_3d (int lev, MultiFab& mf_cons,
36 MultiFab& mf_u , MultiFab& mf_v ,
37 MultiFab* mf_sstore,
38 MultiFab* mf_ru , MultiFab* mf_rv,
39 std::unique_ptr<MultiFab>& mf_DU_avg1,
40 std::unique_ptr<MultiFab>& mf_DU_avg2,
41 std::unique_ptr<MultiFab>& mf_DV_avg1,
42 std::unique_ptr<MultiFab>& mf_DV_avg2,
43 std::unique_ptr<MultiFab>& mf_ubar,
44 std::unique_ptr<MultiFab>& mf_vbar,
45 std::unique_ptr<MultiFab>& mf_Akv,
46 std::unique_ptr<MultiFab>& mf_Akt,
47 std::unique_ptr<MultiFab>& mf_Hz,
48 std::unique_ptr<MultiFab>& mf_Huon,
49 std::unique_ptr<MultiFab>& mf_Hvom,
50 std::unique_ptr<MultiFab>& mf_z_w,
51 MultiFab const* mf_h,
52 MultiFab const* mf_pm,
53 MultiFab const* mf_pn,
54 MultiFab const* mf_mskr,
55 MultiFab const* mf_msku,
56 MultiFab const* mf_mskv,
57 const int N, Real dt_lev)
58{
59 const int nrhs = 0;
60 int nnew = 0;
61
62 int iic = istep[lev];
63 int ntfirst = 0;
64
65 // Because zeta may have changed
66 stretch_transform(lev);
67
68 // These temporaries used to be made in advance_3d_ml and passed in;
69 // now we make them here
70
71 const BoxArray& ba = mf_cons.boxArray();
72 const DistributionMapping& dm = mf_cons.DistributionMap();
73
74 //Only used locally, probably should be rearranged into FArrayBox declaration
75 MultiFab mf_AK (convert(ba,IntVect(0,0,1)),dm,1,IntVect(NGROW,NGROW,0)); //2d missing j coordinate
76 MultiFab mf_DC (ba,dm,1,IntVect(NGROW,NGROW,NGROW-1)); //2d missing j coordinate
77 MultiFab mf_Hzk(ba,dm,1,IntVect(NGROW,NGROW,NGROW-1)); //2d missing j coordinate
78
79 for ( MFIter mfi(mf_cons, TilingIfNotGPU()); mfi.isValid(); ++mfi )
80 {
81 Array4<Real > const& u = mf_u.array(mfi);
82 Array4<Real > const& v = mf_v.array(mfi);
83
84 Array4<Real > const& ru = mf_ru->array(mfi);
85 Array4<Real > const& rv = mf_rv->array(mfi);
86
87 Array4<Real > const& AK = mf_AK.array(mfi);
88 Array4<Real > const& DC = mf_DC.array(mfi);
89
90 Array4<Real > const& Hzk = mf_Hzk.array(mfi);
91 Array4<Real > const& Akv = mf_Akv->array(mfi);
92
93 Array4<Real const> const& Hz = mf_Hz->const_array(mfi);
94
95 Array4<Real const> const& DU_avg1 = mf_DU_avg1->const_array(mfi);
96 Array4<Real const> const& DV_avg1 = mf_DV_avg1->const_array(mfi);
97
98 Array4<Real const> const& pm = mf_pm->const_array(mfi);
99 Array4<Real const> const& pn = mf_pn->const_array(mfi);
100
101 Array4<Real const> const& msku = mf_msku->const_array(mfi);
102 Array4<Real const> const& mskv = mf_mskv->const_array(mfi);
103
104 Box bx = mfi.tilebox();
105 Box gbx2 = mfi.growntilebox(IntVect(NGROW,NGROW,0));
106 Box gbx21 = mfi.growntilebox(IntVect(NGROW,NGROW,NGROW-1));
107
108 Box xbx = mfi.nodaltilebox(0);
109 Box ybx = mfi.nodaltilebox(1);
110
111 Box gbx2D = gbx2;
112 gbx2D.makeSlab(2,0);
113
114 Box tbxp1 = bx;
115 Box tbxp11 = bx;
116 Box tbxp2 = bx;
117 tbxp1.grow(IntVect(NGROW-1,NGROW-1,0));
118 tbxp2.grow(IntVect(NGROW,NGROW,0));
119 tbxp11.grow(IntVect(NGROW-1,NGROW-1,NGROW-1));
120
121 FArrayBox fab_FC(convert(gbx2,IntVect(0,0,1)),1,amrex::The_Async_Arena());
122 FArrayBox fab_BC(convert(gbx2,IntVect(0,0,1)),1,amrex::The_Async_Arena());
123 FArrayBox fab_CF(gbx21,1,amrex::The_Async_Arena());
124 FArrayBox fab_W(tbxp2,1,amrex::The_Async_Arena());
125
126 auto FC = fab_FC.array();
127 auto BC = fab_BC.array();
128 auto CF = fab_CF.array();
129
130 Real cff;
131 if (iic==ntfirst) {
132 cff=0.25_rt*dt_lev;
133 } else if (iic==ntfirst+1) {
134 cff=0.25_rt*dt_lev*3.0_rt/2.0_rt;
135 } else {
136 cff=0.25_rt*dt_lev*23.0_rt/12.0_rt;
137 }
138
139 ParallelFor(xbx, [=] AMREX_GPU_DEVICE (int i, int j, int k)
140 {
141 u(i,j,k) += cff * (pm(i,j,0)+pm(i-1,j,0)) * (pn(i,j,0)+pn(i-1,j,0)) * ru(i,j,k,nrhs);
142 u(i,j,k) *= 2.0_rt / (Hz(i-1,j,k) + Hz(i,j,k));
143 });
144
145 ParallelFor(ybx, [=] AMREX_GPU_DEVICE (int i, int j, int k)
146 {
147 v(i,j,k) += cff * (pm(i,j,0)+pm(i,j-1,0)) * (pn(i,j,0)+pn(i,j-1,0)) * rv(i,j,k,nrhs);
148 v(i,j,k) *= 2.0_rt / (Hz(i,j-1,k) + Hz(i,j,k));
149 });
150
151 // NOTE: DC is only used as scratch in vert_visc_3d -- no need to pass or return a value
152 // NOTE: may not actually need to set these to zero
153
154 // Reset to zero on the box on which they'll be used
155 mf_DC[mfi].template setVal<RunOn::Device>(0.,xbx);
156 fab_CF.template setVal<RunOn::Device>(0.,xbx);
157
158 vert_visc_3d(xbx,1,0,u,Hz,Hzk,AK,Akv,BC,DC,FC,CF,nnew,N,dt_lev);
159
160 // Reset to zero on the box on which they'll be used
161 mf_DC[mfi].template setVal<RunOn::Device>(0.,ybx);
162 fab_CF.template setVal<RunOn::Device>(0.,ybx);
163
164 vert_visc_3d(ybx,0,1,v,Hz,Hzk,AK,Akv,BC,DC,FC,CF,nnew,N,dt_lev);
165
166 // Reset to zero on the box on which they'll be used
167 mf_DC[mfi].template setVal<RunOn::Device>(0.,xbx);
168 fab_CF.template setVal<RunOn::Device>(0.,xbx);
169
170 vert_mean_3d(xbx,1,0,u,Hz,DU_avg1,DC,CF,pn,msku,nnew,N);
171
172 // Reset to zero on the box on which they'll be used
173 mf_DC[mfi].template setVal<RunOn::Device>(0.,ybx);
174 fab_CF.template setVal<RunOn::Device>(0.,ybx);
175
176 vert_mean_3d(ybx,0,1,v,Hz,DV_avg1,DC,CF,pm,mskv,nnew,N);
177 }
178
179 // Apply physical boundary conditions to u and v
180 (*physbcs[lev])(mf_u,*mf_msku,0,1,mf_u.nGrowVect(),t_old[lev],BCVars::xvel_bc,0,*xvel_old[lev]);
181 (*physbcs[lev])(mf_v,*mf_mskv,0,1,mf_v.nGrowVect(),t_old[lev],BCVars::yvel_bc,0,*yvel_old[lev]);
182
183#ifdef REMORA_USE_NETCDF
184 // Fill the data which is stored in the boundary data read from netcdf files
185 if ( (solverChoice.ic_bc_type == IC_BC_Type::netcdf) && (lev==0) )
186 {
187 fill_from_bdyfiles(mf_u,*mf_msku,t_old[lev],BCVars::xvel_bc,BdyVars::u,0,0,*xvel_old[lev],dt_lev);
188 fill_from_bdyfiles(mf_v,*mf_mskv,t_old[lev],BCVars::yvel_bc,BdyVars::v,0,0,*yvel_old[lev],dt_lev);
189 }
190
191 if (solverChoice.do_rivers) {
192 river_source_transport->update_interpolated_to_time(t_old[lev]);
193 }
194#endif
195 for ( MFIter mfi(mf_cons, TilingIfNotGPU()); mfi.isValid(); ++mfi )
196 {
197 Array4<Real > const& u = mf_u.array(mfi);
198 Array4<Real > const& v = mf_v.array(mfi);
199
200 Array4<Real > const& DC = mf_DC.array(mfi);
201
202 Array4<Real const> const& Hz = mf_Hz->const_array(mfi);
203
204 Array4<Real const> const& DU_avg1 = mf_DU_avg1->const_array(mfi);
205 Array4<Real const> const& DV_avg1 = mf_DV_avg1->const_array(mfi);
206
207 Array4<Real const> const& DU_avg2 = mf_DU_avg2->const_array(mfi);
208 Array4<Real const> const& DV_avg2 = mf_DV_avg2->const_array(mfi);
209
210 Array4<Real> const& ubar = mf_ubar->array(mfi);
211 Array4<Real> const& vbar = mf_vbar->array(mfi);
212
213 Array4<Real> const& Huon = mf_Huon->array(mfi);
214 Array4<Real> const& Hvom = mf_Hvom->array(mfi);
215
216 Array4<Real const> const& pm = mf_pm->const_array(mfi);
217 Array4<Real const> const& pn = mf_pn->const_array(mfi);
218 Array4<Real const> const& msku = mf_msku->const_array(mfi);
219 Array4<Real const> const& mskv = mf_mskv->const_array(mfi);
220
221 Array4<Real const> const& z_w = mf_z_w->const_array(mfi);
222
223 Box bx = mfi.tilebox();
224 Box gbx1 = mfi.growntilebox(IntVect(NGROW-1,NGROW-1,0));
225 Box gbx2 = mfi.growntilebox(IntVect(NGROW,NGROW,0));
226
227 FArrayBox fab_FC(gbx2,1,amrex::The_Async_Arena());
228 auto FC = fab_FC.array();
229
230#ifdef REMORA_USE_NETCDF
231 if (solverChoice.do_rivers) {
232 Array4<int const> const& river_pos = vec_river_position[lev]->const_array(mfi);
233 Array4<Real const> const& river_transport = river_source_transport->fab_interp->array();
234 int* river_direction_d = river_direction.data();
235 ParallelFor(gbx1, [=] AMREX_GPU_DEVICE(int i, int j, int k)
236 {
237 int iriver = river_pos(i,j,0);
238 if (iriver >= 0) {
239 if (river_direction_d[iriver] == 0) {
240 Real on_u = 2.0_rt / (pn(i,j,0)+pn(i-1,j,0));
241 Real cff = 1.0_rt / (on_u * 0.5_rt * (z_w(i-1,j,k+1) - z_w(i-1,j,k) + z_w(i,j,k+1) - z_w(i,j,k)));
242 u(i,j,k) = cff * river_transport(iriver,0,k);
243 } else {
244 Real om_v = 2.0_rt / (pm(i,j,0)+pm(i,j-1,0));
245 Real cff = 1.0_rt / (om_v * 0.5_rt * (z_w(i,j-1,k+1) - z_w(i,j-1,k) + z_w(i,j,k+1) - z_w(i,j,k)));
246 v(i,j,k) = cff * river_transport(iriver,0,k);
247 }
248 }
249 });
250 }
251#endif
252
253#if 0
254 // Reset to zero on the box on which they'll be used
255 mf_DC[mfi].template setVal<RunOn::Device>(0.,grow(xbx,IntVect(0,0,1)));
256 fab_CF.template setVal<RunOn::Device>(0.,grow(xbx,IntVect(0,0,1)));
257
258 update_massflux_3d(xbx,1,0,u,ubar,Huon,Hz,pn,DU_avg1,DU_avg2,DC,FC,nnew);
259
260 // Reset to zero on the box on which they'll be used
261 mf_DC[mfi].template setVal<RunOn::Device>(0.,grow(ybx,IntVect(0,0,1)));
262 fab_CF.template setVal<RunOn::Device>(0.,grow(ybx,IntVect(0,0,1)));
263
264 update_massflux_3d(ybx,0,1,v,vbar,Hvom,Hz,pm,DV_avg1,DV_avg2,DC,FC,nnew);
265
266#else
267
268 // Reset to zero on the box on which they'll be used
269 fab_FC.template setVal<RunOn::Device>(0.,gbx2);
270 mf_DC[mfi].template setVal<RunOn::Device>(0.,grow(gbx2,IntVect(0,0,1)));
271 update_massflux_3d(gbx2,1,0,u,ubar,Huon,Hz,pn,DU_avg1,DU_avg2,DC,FC,msku,nnew);
272
273 // Reset to zero on the box on which they'll be used
274 fab_FC.template setVal<RunOn::Device>(0.,gbx2);
275 mf_DC[mfi].template setVal<RunOn::Device>(0.,grow(gbx2,IntVect(0,0,1)));
276 update_massflux_3d(gbx2,0,1,v,vbar,Hvom,Hz,pm,DV_avg1,DV_avg2,DC,FC,mskv,nnew);
277#endif
278 }
279
280 // WE BELIEVE THESE VALUES SHOULD ALREADY BE FILLED
281 // mf_Huon->FillBoundary(geom[lev].periodicity());
282 // mf_Hvom->FillBoundary(geom[lev].periodicity());
283
284 // ************************************************************************
285 // This should fill both temp and salt with temp/salt currently in cons_old
286 // ************************************************************************
287
288 MultiFab mf_W(convert(ba,IntVect(0,0,1)),dm,1,IntVect(NGROW+1,NGROW+1,0));
289 mf_W.setVal(0.0_rt);
290 for ( MFIter mfi(mf_cons, TilingIfNotGPU()); mfi.isValid(); ++mfi )
291 {
292 Array4<Real> const& Huon = mf_Huon->array(mfi);
293 Array4<Real> const& Hvom = mf_Hvom->array(mfi);
294
295 Array4<Real const> const& z_w = mf_z_w->const_array(mfi);
296 Array4<Real const> const& h = mf_h->const_array(mfi);
297
298 Box bx = mfi.tilebox();
299 Box gbx1 = mfi.growntilebox(IntVect(NGROW-1,NGROW-1,0));
300 Box gbx2 = mfi.growntilebox(IntVect(NGROW,NGROW,0));
301 Box gbx21 = mfi.growntilebox(IntVect(NGROW,NGROW,NGROW-1));
302
303 Box tbxp1 = bx;
304 Box tbxp11 = bx;
305 Box tbxp2 = bx;
306 tbxp1.grow(IntVect(NGROW-1,NGROW-1,0));
307 tbxp2.grow(IntVect(NGROW,NGROW,0));
308 tbxp11.grow(IntVect(NGROW-1,NGROW-1,NGROW-1));
309
310 FArrayBox fab_FC(surroundingNodes(gbx2,2),1,amrex::The_Async_Arena());
311 FArrayBox fab_BC(gbx2,1,amrex::The_Async_Arena());
312 FArrayBox fab_CF(gbx21,1,amrex::The_Async_Arena());
313
314 auto W = mf_W.array(mfi);
315
316 //
317 //------------------------------------------------------------------------
318 // Vertically integrate horizontal mass flux divergence.
319 //------------------------------------------------------------------------
320 //
321 //Should really use gbx3uneven
322 //TODO: go over these boxes and compare to other spots where we do the same thing
323 Box gbx1D = gbx1;
324 gbx1D.makeSlab(2,0);
325
326 // Starting with zero vertical velocity at the bottom, integrate
327 // from the bottom (k=0) to the free-surface (k=N). The w(:,:,N(ng))
328 // contains the vertical velocity at the free-surface, d(zeta)/d(t).
329 // Notice that barotropic mass flux divergence is not used directly.
330 ParallelFor(gbx1D, [=] AMREX_GPU_DEVICE (int i, int j, int )
331 {
332 W(i,j,0) = 0.0_rt;
333 for (int k=1; k<=N+1; k++) {
334 W(i,j,k) = W(i,j,k-1) - (Huon(i+1,j,k-1)-Huon(i,j,k-1)) - (Hvom(i,j+1,k-1)-Hvom(i,j,k-1));
335 }
336 });
337
338 // Starting with zero vertical velocity at the bottom, integrate
339 // from the bottom (k=0) to the free-surface (k=N). The w(:,:,N(ng))
340 // contains the vertical velocity at the free-surface, d(zeta)/d(t).
341 // Notice that barotropic mass flux divergence is not used directly.
342 //
343 ParallelFor(gbx1D, [=] AMREX_GPU_DEVICE (int i, int j, int )
344 {
345 Real wrk_ij = W(i,j,N+1) / (z_w(i,j,N+1)+h(i,j,0,0));
346
347 for (int k=1; k<=N; k++) {
348 W(i,j,k) -= wrk_ij * (z_w(i,j,k)+h(i,j,0,0));
349 }
350 W(i,j,N+1) = 0.0_rt;
351 });
352 }
353
354 const int nstp = (iic) % 2;
355 nnew = 1-nstp;
356 if (solverChoice.vert_mixing_type == VertMixingType::GLS) {
357 gls_corrector(lev, vec_gls[lev].get(), vec_tke[lev].get(), mf_W, vec_Akv[lev].get(),
358 vec_Akt[lev].get(),vec_Akk[lev].get(), vec_Akp[lev].get(), vec_mskr[lev].get(),
359 vec_msku[lev].get(), vec_mskv[lev].get(),
360 nstp, nnew, N, dt_lev);
361 }
362 nnew = 0;
363
364 for ( MFIter mfi(mf_cons, TilingIfNotGPU()); mfi.isValid(); ++mfi )
365 {
366 Array4<Real> const& Hz = mf_Hz->array(mfi);
367
368 Array4<Real> const& Huon = mf_Huon->array(mfi);
369 Array4<Real> const& Hvom = mf_Hvom->array(mfi);
370
371 Array4<Real const> const& pm = mf_pm->const_array(mfi);
372 Array4<Real const> const& pn = mf_pn->const_array(mfi);
373 Array4<Real const> const& mskr = mf_mskr->const_array(mfi);
374 Array4<Real const> const& msku = mf_msku->const_array(mfi);
375 Array4<Real const> const& mskv = mf_mskv->const_array(mfi);
376
377 Box bx = mfi.tilebox();
378 Box gbx2 = mfi.growntilebox(IntVect(NGROW,NGROW,0));
379 Box gbx21 = mfi.growntilebox(IntVect(NGROW,NGROW,NGROW-1));
380
381 Box tbxp1 = bx;
382 Box tbxp11 = bx;
383 Box tbxp2 = bx;
384 tbxp1.grow(IntVect(NGROW-1,NGROW-1,0));
385 tbxp2.grow(IntVect(NGROW,NGROW,0));
386 tbxp11.grow(IntVect(NGROW-1,NGROW-1,NGROW-1));
387
388 FArrayBox fab_FC(surroundingNodes(gbx2,2),1,amrex::The_Async_Arena());
389 FArrayBox fab_BC(gbx2,1,amrex::The_Async_Arena());
390 FArrayBox fab_CF(gbx21,1,amrex::The_Async_Arena());
391
392 auto FC = fab_FC.array();
393 auto W = mf_W.array(mfi);
394 //
395 //-----------------------------------------------------------------------
396 // rhs_t_3d
397 //-----------------------------------------------------------------------
398 //
399 for (int i_comp=0; i_comp < NCONS; i_comp++)
400 {
401#ifdef REMORA_USE_NETCDF
402 FArrayBox* fab_river_source;
403 if (solverChoice.do_rivers_cons[i_comp]) {
404 fab_river_source = river_source_cons[i_comp]->fab_interp;
405 }
406 const Array4<const int>& river_pos = (solverChoice.do_rivers) ? vec_river_position[lev]->const_array(mfi) : Array4<const int>();
407 const Array4<const Real>& river_source = (solverChoice.do_rivers_cons[i_comp]) ? fab_river_source->const_array() : Array4<const Real>();
408#else
409 const Array4<const int >& river_pos = Array4<const int>();
410 const Array4<const Real>& river_source = Array4<const Real>();
411#endif
412 Array4<Real> const& sstore = mf_sstore->array(mfi, i_comp);
413 rhs_t_3d(bx, mf_cons.array(mfi,i_comp), sstore, Huon, Hvom,
414 Hz, pn, pm, W, FC, mskr, msku, mskv, river_pos, river_source, nrhs, nnew, N,dt_lev);
415 }
416
417 } // mfi
418
419 FillPatch(lev, t_old[lev], mf_cons, cons_new, BCVars::cons_bc, BdyVars::t,0,true,false,0,0,dt_lev,*cons_old[lev]);
420
421 for ( MFIter mfi(mf_cons, TilingIfNotGPU()); mfi.isValid(); ++mfi )
422 {
423 Array4<Real> const& AK = mf_AK.array(mfi);
424 Array4<Real> const& DC = mf_DC.array(mfi);
425
426 Array4<Real> const& Hzk = mf_Hzk.array(mfi);
427 Array4<Real const> const& Hz = mf_Hz->const_array(mfi);
428
429 Box bx = mfi.tilebox();
430
431 // Copy the tilebox
432 Box tbxp1 = bx;
433 Box tbxp11 = bx;
434 Box tbxp2 = bx;
435 Box tbxp21 = bx;
436 //make only gbx be grown to match multifabs
437 tbxp21.grow(IntVect(NGROW,NGROW,NGROW-1));
438 tbxp2.grow(IntVect(NGROW,NGROW,0));
439 tbxp1.grow(IntVect(NGROW-1,NGROW-1,0));
440 tbxp11.grow(IntVect(NGROW-1,NGROW-1,NGROW-1));
441
442 FArrayBox fab_FC(tbxp2,1,amrex::The_Async_Arena());
443 FArrayBox fab_BC(tbxp2,1,amrex::The_Async_Arena());
444 FArrayBox fab_CF(tbxp21,1,amrex::The_Async_Arena());
445 FArrayBox fab_W(tbxp2,1,amrex::The_Async_Arena());
446
447 auto FC = fab_FC.array();
448 auto BC = fab_BC.array();
449 auto CF = fab_CF.array();
450
451 for (int i_comp=0; i_comp < NCONS; i_comp++) {
452 vert_visc_3d(bx,0,0,mf_cons.array(mfi,i_comp),Hz,Hzk,
453 AK,mf_Akt->array(mfi,i_comp),BC,DC,FC,CF,nnew,N,dt_lev);
454 }
455 } // MFiter
456 FillPatch(lev, t_old[lev], *cons_new[lev], cons_new, BCVars::cons_bc, BdyVars::t,0,true,false,0,0,dt_lev,*cons_old[lev]);
457
458
459#ifdef REMORA_USE_NETCDF
460 if (solverChoice.do_temp_clim_nudg) {
461 temp_clim_data_from_file->update_interpolated_to_time(t_old[lev]);
462
463 for ( MFIter mfi(mf_cons, TilingIfNotGPU()); mfi.isValid(); ++mfi )
464 {
465 Box bx = mfi.growntilebox(IntVect(1,1,0));
466 Array4<const Real> const& temp_nudg_coeff = vec_nudg_coeff[BdyVars::t][lev]->const_array(mfi);
467 Array4<const Real> const& temp_clim = temp_clim_data_from_file->mf_interpolated->const_array(mfi);
468 Array4< Real> const& temp = cons_new[lev]->array(mfi, Temp_comp);
469 Array4<const Real> const& Hz = vec_Hz[lev]->const_array(mfi);
470 Array4<const Real> const& pm = vec_pm[lev]->const_array(mfi);
471 Array4<const Real> const& pn = vec_pn[lev]->const_array(mfi);
472
473 apply_clim_nudg(bx, 0, 0, temp, temp, temp_clim, temp_nudg_coeff, Hz, pm, pn, dt_lev);
474 }
475 }
476 if (solverChoice.do_salt_clim_nudg) {
477 salt_clim_data_from_file->update_interpolated_to_time(t_old[lev]);
478
479 for ( MFIter mfi(mf_cons, TilingIfNotGPU()); mfi.isValid(); ++mfi )
480 {
481 Box bx = mfi.growntilebox(IntVect(1,1,0));
482 Array4<const Real> const& salt_nudg_coeff = vec_nudg_coeff[BdyVars::s][lev]->const_array(mfi);
483 Array4<const Real> const& salt_clim = salt_clim_data_from_file->mf_interpolated->const_array(mfi);
484 Array4< Real> const& salt = cons_new[lev]->array(mfi, Salt_comp);
485 Array4<const Real> const& Hz = vec_Hz[lev]->const_array(mfi);
486 Array4<const Real> const& pm = vec_pm[lev]->const_array(mfi);
487 Array4<const Real> const& pn = vec_pn[lev]->const_array(mfi);
488
489 apply_clim_nudg(bx, 0, 0, salt, salt, salt_clim, salt_nudg_coeff, Hz, pm, pn, dt_lev);
490 }
491 }
492#endif
493}
IC_BC_Type::netcdf
@ netcdf
VertMixingType::GLS
@ GLS
NGROW
#define NGROW
Definition REMORA_IndexDefines.H:13
Temp_comp
#define Temp_comp
Definition REMORA_IndexDefines.H:8
Salt_comp
#define Salt_comp
Definition REMORA_IndexDefines.H:9
NCONS
#define NCONS
Definition REMORA_IndexDefines.H:11
NCTimeSeriesRiver::update_interpolated_to_time
void update_interpolated_to_time(amrex::Real time)
Calculate.
Definition REMORA_NCTimeSeriesRiver.cpp:109
NCTimeSeriesRiver::fab_interp
amrex::FArrayBox * fab_interp
Container for interpolated data; Only used if save_interpolated == true.
Definition REMORA_NCTimeSeriesRiver.H:26
NCTimeSeries::mf_interpolated
amrex::MultiFab * mf_interpolated
Definition REMORA_NCTimeSeries.H:35
NCTimeSeries::update_interpolated_to_time
void update_interpolated_to_time(amrex::Real time)
Calculate interpolated values at time, reading in data as necessary.
Definition REMORA_NCTimeSeries.cpp:96
REMORA::river_source_cons
amrex::Vector< NCTimeSeriesRiver * > river_source_cons
Vector of data containers for scalar data in rivers.
Definition REMORA.H:1029
REMORA::vec_pm
amrex::Vector< std::unique_ptr< amrex::MultiFab > > vec_pm
horizontal scaling factor: 1 / dx (2D)
Definition REMORA.H:355
REMORA::cons_new
amrex::Vector< amrex::MultiFab * > cons_new
multilevel data container for current step's scalar data: temperature, salinity, passive scalar
Definition REMORA.H:217
REMORA::stretch_transform
void stretch_transform(int lev)
Calculate vertical stretched coordinates.
Definition REMORA_DepthStretchTransform.H:15
REMORA::vec_mskr
amrex::Vector< std::unique_ptr< amrex::MultiFab > > vec_mskr
land/sea mask at cell centers (2D)
Definition REMORA.H:346
REMORA::advance_3d
void advance_3d(int lev, amrex::MultiFab &mf_cons, amrex::MultiFab &mf_u, amrex::MultiFab &mf_v, amrex::MultiFab *mf_sstore, amrex::MultiFab *mf_ru, amrex::MultiFab *mf_rv, std::unique_ptr< amrex::MultiFab > &mf_DU_avg1, std::unique_ptr< amrex::MultiFab > &mf_DU_avg2, std::unique_ptr< amrex::MultiFab > &mf_DV_avg1, std::unique_ptr< amrex::MultiFab > &mf_DV_avg2, std::unique_ptr< amrex::MultiFab > &mf_ubar, std::unique_ptr< amrex::MultiFab > &mf_vbar, std::unique_ptr< amrex::MultiFab > &mf_Akv, std::unique_ptr< amrex::MultiFab > &mf_Akt, std::unique_ptr< amrex::MultiFab > &mf_Hz, std::unique_ptr< amrex::MultiFab > &mf_Huon, std::unique_ptr< amrex::MultiFab > &mf_Hvom, std::unique_ptr< amrex::MultiFab > &mf_z_w, amrex::MultiFab const *mf_h, amrex::MultiFab const *mf_pm, amrex::MultiFab const *mf_pn, amrex::MultiFab const *mf_mskr, amrex::MultiFab const *mf_msku, amrex::MultiFab const *mf_mskv, const int N, const amrex::Real dt_lev)
Advance the 3D variables.
Definition REMORA_advance_3d.cpp:35
REMORA::rhs_t_3d
void rhs_t_3d(const amrex::Box &bx, const amrex::Array4< amrex::Real > &t, const amrex::Array4< amrex::Real const > &tempstore, const amrex::Array4< amrex::Real const > &Huon, const amrex::Array4< amrex::Real const > &Hvom, const amrex::Array4< amrex::Real const > &Hz, const amrex::Array4< amrex::Real const > &pn, const amrex::Array4< amrex::Real const > &pm, const amrex::Array4< amrex::Real const > &W, const amrex::Array4< amrex::Real > &FC, const amrex::Array4< amrex::Real const > &mskr, const amrex::Array4< amrex::Real const > &msku, const amrex::Array4< amrex::Real const > &mskv, const amrex::Array4< int const > &river_pos, const amrex::Array4< amrex::Real const > &river_source, int nrhs, int nnew, int N, const amrex::Real dt_lev)
RHS terms for tracer.
Definition REMORA_rhs_t_3d.cpp:27
REMORA::vec_tke
amrex::Vector< std::unique_ptr< amrex::MultiFab > > vec_tke
Turbulent kinetic energy.
Definition REMORA.H:403
REMORA::vec_gls
amrex::Vector< std::unique_ptr< amrex::MultiFab > > vec_gls
Turbulent generic length scale.
Definition REMORA.H:405
REMORA::vert_mean_3d
void vert_mean_3d(const amrex::Box &bx, const int ioff, const int joff, const amrex::Array4< amrex::Real > &phi, const amrex::Array4< amrex::Real const > &Hz, const amrex::Array4< amrex::Real const > &Dphi_avg1, const amrex::Array4< amrex::Real > &DC, const amrex::Array4< amrex::Real > &CF, const amrex::Array4< amrex::Real const > &pm_or_pn, const amrex::Array4< amrex::Real const > &msk, const int nnew, const int N)
Adjust 3D momentum variables based on vertical mean momentum.
Definition REMORA_vert_mean_3d.cpp:19
REMORA::vec_Hz
amrex::Vector< std::unique_ptr< amrex::MultiFab > > vec_Hz
Width of cells in the vertical (z-) direction (3D, Hz in ROMS)
Definition REMORA.H:232
REMORA::vec_Akt
amrex::Vector< std::unique_ptr< amrex::MultiFab > > vec_Akt
Vertical diffusion coefficient (3D)
Definition REMORA.H:252
REMORA::vec_msku
amrex::Vector< std::unique_ptr< amrex::MultiFab > > vec_msku
land/sea mask at x-faces (2D)
Definition REMORA.H:348
REMORA::xvel_old
amrex::Vector< amrex::MultiFab * > xvel_old
multilevel data container for last step's x velocities (u in ROMS)
Definition REMORA.H:210
REMORA::vec_mskv
amrex::Vector< std::unique_ptr< amrex::MultiFab > > vec_mskv
land/sea mask at y-faces (2D)
Definition REMORA.H:350
REMORA::physbcs
amrex::Vector< std::unique_ptr< REMORAPhysBCFunct > > physbcs
Vector (over level) of functors to apply physical boundary conditions.
Definition REMORA.H:1108
REMORA::FillPatch
void FillPatch(int lev, amrex::Real time, amrex::MultiFab &mf_to_be_filled, amrex::Vector< amrex::MultiFab * > const &mfs, const int bccomp, const int bdy_var_type=BdyVars::null, const int icomp=0, const bool fill_all=true, const bool fill_set=true, const int n_not_fill=0, const int icomp_calc=0, const amrex::Real dt=amrex::Real(0.0), const amrex::MultiFab &mf_calc=amrex::MultiFab())
Fill a new MultiFab by copying in phi from valid region and filling ghost cells.
Definition REMORA_FillPatch.cpp:30
REMORA::istep
amrex::Vector< int > istep
which step?
Definition REMORA.H:1094
REMORA::apply_clim_nudg
void apply_clim_nudg(const amrex::Box &bx, int ioff, int joff, const amrex::Array4< amrex::Real > &var, const amrex::Array4< amrex::Real const > &var_old, const amrex::Array4< amrex::Real const > &var_clim, const amrex::Array4< amrex::Real const > &clim_coeff, const amrex::Array4< amrex::Real const > &Hz, const amrex::Array4< amrex::Real const > &pm, const amrex::Array4< amrex::Real const > &pn, const amrex::Real dt_lev=amrex::Real(0.0))
Apply climatology nudging.
Definition REMORA_apply_clim_nudg.cpp:19
REMORA::fill_from_bdyfiles
void fill_from_bdyfiles(amrex::MultiFab &mf_to_fill, const amrex::MultiFab &mf_mask, const amrex::Real time, const int bccomp, const int bdy_var_type, const int icomp_to_fill, const int icomp_calc=0, const amrex::MultiFab &mf_calc=amrex::MultiFab(), const amrex::Real=amrex::Real(0.0))
Fill boundary data from netcdf file.
Definition REMORA_BoundaryConditions_netcdf.cpp:19
REMORA::update_massflux_3d
void update_massflux_3d(const amrex::Box &bx, const int ioff, const int joff, const amrex::Array4< amrex::Real > &phi, const amrex::Array4< amrex::Real > &phibar, const amrex::Array4< amrex::Real > &Hphi, const amrex::Array4< amrex::Real const > &Hz, const amrex::Array4< amrex::Real const > &pm_or_pn, const amrex::Array4< amrex::Real const > &Dphi1, const amrex::Array4< amrex::Real const > &Dphi2, const amrex::Array4< amrex::Real > &DC, const amrex::Array4< amrex::Real > &FC, const amrex::Array4< amrex::Real const > &msk, const int nnew)
Correct mass flux.
Definition REMORA_update_massflux_3d.cpp:22
REMORA::yvel_old
amrex::Vector< amrex::MultiFab * > yvel_old
multilevel data container for last step's y velocities (v in ROMS)
Definition REMORA.H:212
REMORA::river_source_transport
NCTimeSeriesRiver * river_source_transport
Data container for momentum transport in rivers.
Definition REMORA.H:1031
REMORA::temp_clim_data_from_file
NCTimeSeries * temp_clim_data_from_file
Data container for temperature climatology data read from file.
Definition REMORA.H:1024
REMORA::solverChoice
static SolverChoice solverChoice
Container for algorithmic choices.
Definition REMORA.H:1221
REMORA::vec_river_position
amrex::Vector< std::unique_ptr< amrex::iMultiFab > > vec_river_position
iMultiFab for river positions; contents are indices of rivers
Definition REMORA.H:1036
REMORA::vec_Akk
amrex::Vector< std::unique_ptr< amrex::MultiFab > > vec_Akk
Turbulent kinetic energy vertical diffusion coefficient.
Definition REMORA.H:409
REMORA::salt_clim_data_from_file
NCTimeSeries * salt_clim_data_from_file
Data container for salinity climatology data read from file.
Definition REMORA.H:1026
REMORA::vert_visc_3d
void vert_visc_3d(const amrex::Box &bx, const int ioff, const int joff, const amrex::Array4< amrex::Real > &phi, const amrex::Array4< amrex::Real const > &Hz, const amrex::Array4< amrex::Real > &Hzk, const amrex::Array4< amrex::Real > &AK, const amrex::Array4< amrex::Real const > &Akv, const amrex::Array4< amrex::Real > &BC, const amrex::Array4< amrex::Real > &DC, const amrex::Array4< amrex::Real > &FC, const amrex::Array4< amrex::Real > &CF, const int nnew, const int N, const amrex::Real dt_lev)
Calculate effects of vertical viscosity or diffusivity.
Definition REMORA_vert_visc_3d.cpp:23
REMORA::river_direction
amrex::Gpu::DeviceVector< int > river_direction
Vector over rivers of river direction: 0: u-face; 1: v-face; 2: w-face.
Definition REMORA.H:1038
REMORA::cons_old
amrex::Vector< amrex::MultiFab * > cons_old
multilevel data container for last step's scalar data: temperature, salinity, passive scalar
Definition REMORA.H:208
REMORA::vec_nudg_coeff
amrex::Vector< amrex::Vector< std::unique_ptr< amrex::MultiFab > > > vec_nudg_coeff
Climatology nudging coefficients.
Definition REMORA.H:414
REMORA::vec_pn
amrex::Vector< std::unique_ptr< amrex::MultiFab > > vec_pn
horizontal scaling factor: 1 / dy (2D)
Definition REMORA.H:357
REMORA::vec_Akv
amrex::Vector< std::unique_ptr< amrex::MultiFab > > vec_Akv
Vertical viscosity coefficient (3D)
Definition REMORA.H:250
REMORA::t_old
amrex::Vector< amrex::Real > t_old
old time at each level
Definition REMORA.H:1100
REMORA::vec_Akp
amrex::Vector< std::unique_ptr< amrex::MultiFab > > vec_Akp
Turbulent length scale vertical diffusion coefficient.
Definition REMORA.H:411
REMORA::gls_corrector
void gls_corrector(int lev, amrex::MultiFab *mf_gls, amrex::MultiFab *mf_tke, amrex::MultiFab &mf_W, amrex::MultiFab *mf_Akv, amrex::MultiFab *mf_Akt, amrex::MultiFab *mf_Akk, amrex::MultiFab *mf_Akp, amrex::MultiFab *mf_mskr, amrex::MultiFab *mf_msku, amrex::MultiFab *mf_mskv, const int nstp, const int nnew, const int N, const amrex::Real dt_lev)
Corrector step for GLS calculation.
Definition REMORA_gls.cpp:246
SolverChoice::ic_bc_type
IC_BC_Type ic_bc_type
Definition REMORA_DataStruct.H:552
SolverChoice::do_rivers_cons
amrex::Vector< int > do_rivers_cons
Definition REMORA_DataStruct.H:534
SolverChoice::vert_mixing_type
VertMixingType vert_mixing_type
Definition REMORA_DataStruct.H:570