REMORA
Regional Modeling of Oceans Refined Adaptively
Loading...
Searching...
No Matches
REMORA_advance_2d.cpp
Go to the documentation of this file.
1#include <REMORA.H>
2
3using namespace amrex;
4/** Nonlinear shallow-water rpimitive equations predictor (Leap-frog) and
5 * corrector (Adams-Moulton) time-stepping engine. Corresponds to Nonlinear/step2d_LF_AM3.h
6 * in ROMS.
7 *
8 * @param[in ] lev level of refinement (coarsest level is 0)
9 * @param[in ] mf_rhoS density perturbation
10 * @param[in ] mf_rhoA vertically-averaged density
11 * @param[inout] mf_ru2d RHS contributions to 2D u-momentum
12 * @param[inout] mf_rv2d RHS contribtuions to 2D v-momentum
13 * @param[inout] mf_rufrc before first predictor, vertical integral of 3D RHS for uvel, converted to forcing terms
14 * @param[inout] mf_rvfrc before first predictor, vertical integral of 3D RHS for vvel, converted to forcing term
15 * @param[inout] mf_Zt_avg1 average of sea surface height over all fast steps
16 * @param[inout] mf_DU_avg1 time-averaged u-flux for 2D equations
17 * @param[inout] mf_DU_avg2 time-averaged u-flux for 3D equation coupling
18 * @param[inout] mf_DV_avg1 time-averaged v-flux for 2D equations
19 * @param[inout] mf_DV_avg2 time-averaged v-flux for 3D equation coupling
20 * @param[inout] mf_rubar RHS of vertically integrated u-momentum
21 * @param[inout] mf_rvbar RHS of vertically integrated v-momentum
22 * @param[inout] mf_rzeta RHS of sea surface height
23 * @param[inout] mf_ubar vertically integrated u-momentum
24 * @param[inout] mf_vbar vertically integrated v-momentum
25 * @param[inout] mf_zeta Sea-surface height
26 * @param[in ] mf_h Bathymetry
27 * @param[in ] mf_pm 1 / dx
28 * @param[in ] mf_pn 1 / dy
29 * @param[in ] mf_fcor Coriolis factor
30 * @param[inout] mf_visc2_p Harmonic viscosity at psi points
31 * @param[inout] mf_visc2_r Harmoic viscosity at rho points
32 * @param[in ] mf_mskr Land-sea mask at rho-points
33 * @param[in ] mf_msku Land-sea mask at u-points
34 * @param[in ] mf_mskv Land-sea mask at v-points
35 * @param[in ] mf_mskp Land-sea mask at psi-points
36 * @param[in ] dtfast_lev Length of current barotropic step
37 * @param[in ] predictor_2d_step Is this a predictor step?
38 * @param[in ] first_2d_step Is this the first barotropic step?
39 * @param[in ] my_iif Which barotropic predictor-corrector pair?
40 * @param[inout] next_indx1 Cached index for
41 */
42
43void
44REMORA::advance_2d (int lev,
45 MultiFab const* mf_rhoS,
46 MultiFab const* mf_rhoA,
47 MultiFab * mf_ru2d,
48 MultiFab * mf_rv2d,
49 MultiFab * mf_rufrc,
50 MultiFab * mf_rvfrc,
51 MultiFab * mf_Zt_avg1,
52 std::unique_ptr<MultiFab>& mf_DU_avg1,
53 std::unique_ptr<MultiFab>& mf_DU_avg2,
54 std::unique_ptr<MultiFab>& mf_DV_avg1,
55 std::unique_ptr<MultiFab>& mf_DV_avg2,
56 std::unique_ptr<MultiFab>& mf_rubar,
57 std::unique_ptr<MultiFab>& mf_rvbar,
58 std::unique_ptr<MultiFab>& mf_rzeta,
59 std::unique_ptr<MultiFab>& mf_ubar,
60 std::unique_ptr<MultiFab>& mf_vbar,
61 MultiFab * mf_zeta,
62 MultiFab const* mf_h,
63 MultiFab const* mf_pm,
64 MultiFab const* mf_pn,
65 MultiFab const* mf_fcor,
66 MultiFab const* mf_visc2_p,
67 MultiFab const* mf_visc2_r,
68 MultiFab const* mf_mskr,
69 MultiFab const* mf_msku,
70 MultiFab const* mf_mskv,
71 MultiFab const* mf_mskp,
72 Real dtfast_lev,
73 bool predictor_2d_step,
74 bool first_2d_step, int my_iif,
75 int & next_indx1)
76{
77 int iic = istep[lev];
78 const int nnew = 0;
79 const int nstp = 0;
80 int ntfirst = 0;
81
82 int knew = 3;
83 int krhs = (my_iif + iic) % 2 + 1;
84 int kstp = my_iif <=1 ? iic % 2 + 1 : (iic % 2 + my_iif % 2 + 1) % 2 + 1;
85 int indx1 = krhs;
86 if (predictor_2d_step) {
87 next_indx1 = 3 - indx1;
88 } else {
89 knew = next_indx1;
90 kstp = 3 - knew;
91 krhs = 3;
92 //If it's not the auxiliary time step, set indx1 to next_indx1
93 // NOTE: should this ever not execute?
94 // Include indx1 updates for diagnostic purposes?
95 // if (my_iif<nfast+1)
96 // indx1=next_indx1;
97 }
98 int ptsk = 3-kstp;
99 knew-=1;
100 krhs-=1;
101 kstp-=1;
102 // Include indx1 updates for diagnostic purposes?
103 //indx1-=1;
104 ptsk-=1;
105 auto ba = mf_h->boxArray();
106 auto dm = mf_h->DistributionMap();
107
108 MultiFab mf_DUon(convert(ba,IntVect(1,0,0)),dm,1,IntVect(NGROW,NGROW,0));
109 MultiFab mf_DVom(convert(ba,IntVect(0,1,0)),dm,1,IntVect(NGROW,NGROW,0));
110
111 const auto dlo = amrex::lbound(Geom(lev).Domain());
112 const auto dhi = amrex::ubound(Geom(lev).Domain());
113
114 for ( MFIter mfi(*mf_rhoS, TilingIfNotGPU()); mfi.isValid(); ++mfi )
115 {
116 Array4<Real > const& ubar = mf_ubar->array(mfi);
117 Array4<Real > const& vbar = mf_vbar->array(mfi);
118 Array4<Real > const& zeta = mf_zeta->array(mfi);
119 Array4<Real const> const& h = mf_h->const_array(mfi);
120
121 Array4<Real const> const& pm = mf_pm->const_array(mfi);
122 Array4<Real const> const& pn = mf_pn->const_array(mfi);
123
124 Box bx = mfi.tilebox();
125 Box gbx = mfi.growntilebox();
126 Box gbx1 = mfi.growntilebox(IntVect(NGROW-1,NGROW-1,0));
127 Box gbx2 = mfi.growntilebox(IntVect(NGROW,NGROW,0));
128 Box xgbx2 = mfi.grownnodaltilebox(0, IntVect(NGROW,NGROW,0));
129 Box ygbx2 = mfi.grownnodaltilebox(1, IntVect(NGROW,NGROW,0));
130
131 Box tbxp1 = bx;
132 Box tbxp11 = bx;
133 Box tbxp2 = bx;
134 Box tbxp3 = bx;
135 tbxp1.grow(IntVect(NGROW-1,NGROW-1,0));
136 tbxp2.grow(IntVect(NGROW,NGROW,0));
137 tbxp11.grow(IntVect(NGROW-1,NGROW-1,NGROW-1));
138 tbxp3.grow(IntVect(NGROW+1,NGROW+1,0));
139
140 Box bxD = bx ; bxD.makeSlab(2,0);
141 Box gbxD = gbx ; gbxD.makeSlab(2,0);
142 Box gbx1D = gbx1; gbx1D.makeSlab(2,0);
143 Box gbx2D = gbx2; gbx2D.makeSlab(2,0);
144
145 Box tbxp2D = tbxp2;
146 tbxp2D.makeSlab(2,0);
147
148 // step2d work arrays
149 FArrayBox fab_Drhs(makeSlab(tbxp3,2,0),1,The_Async_Arena());
150 auto Drhs=fab_Drhs.array();
151
152 auto DUon = mf_DUon.array(mfi);
153 auto DVom = mf_DVom.array(mfi);
154
155 ParallelFor(makeSlab(tbxp3,2,0), [=] AMREX_GPU_DEVICE (int i, int j, int)
156 {
157 Drhs(i,j,0)=zeta(i,j,0,krhs)+h(i,j,0);
158 });
159
160 ParallelFor(makeSlab(xgbx2,2,0), [=] AMREX_GPU_DEVICE (int i, int j, int)
161 {
162 Real on_u = 2.0_rt / (pn(i,j,0)+pn(i-1,j,0));
163 Real cff1= 0.5_rt * on_u *(Drhs(i,j,0)+Drhs(i-1,j,0));
164 DUon(i,j,0)=ubar(i,j,0,krhs)*cff1;
165 });
166
167 ParallelFor(makeSlab(ygbx2,2,0), [=] AMREX_GPU_DEVICE (int i, int j, int)
168 {
169 Real om_v = 2.0_rt / (pm(i,j,0)+pm(i,j-1,0));
170 Real cff1= 0.5_rt * om_v * (Drhs(i,j,0)+Drhs(i,j-1,0));
171 DVom(i,j,0)=vbar(i,j,0,krhs)*cff1;
172 });
173 }
174
175 // These are needed to pass the tests with bathymetry but I don't quite see why
176 mf_DUon.FillBoundary(geom[lev].periodicity());
177 mf_DVom.FillBoundary(geom[lev].periodicity());
178
179#ifdef REMORA_USE_NETCDF
180 if (solverChoice.do_m2_clim_nudg) {
181 ubar_clim_data_from_file->update_interpolated_to_time(t_new[lev]);
182 vbar_clim_data_from_file->update_interpolated_to_time(t_new[lev]);
183 }
184#endif
185
186 for ( MFIter mfi(*mf_rhoS, TilingIfNotGPU()); mfi.isValid(); ++mfi )
187 {
188 Array4<Real const> const& rhoS = mf_rhoS->const_array(mfi);
189 Array4<Real const> const& rhoA = mf_rhoA->const_array(mfi);
190 Array4<Real const> const& h = mf_h->const_array(mfi);
191
192 Array4<Real > const& rufrc = mf_rufrc->array(mfi);
193 Array4<Real > const& rvfrc = mf_rvfrc->array(mfi);
194 Array4<Real > const& Zt_avg1 = mf_Zt_avg1->array(mfi);
195 Array4<Real > const& ubar = mf_ubar->array(mfi);
196 Array4<Real > const& vbar = mf_vbar->array(mfi);
197 Array4<Real > const& ubar_krhs = mf_ubar->array(mfi, krhs);
198 Array4<Real > const& vbar_krhs = mf_vbar->array(mfi, krhs);
199 Array4<Real > const& zeta = mf_zeta->array(mfi);
200 Array4<Real > const& DU_avg1 = (mf_DU_avg1)->array(mfi);
201 Array4<Real > const& DU_avg2 = (mf_DU_avg2)->array(mfi);
202 Array4<Real > const& DV_avg1 = (mf_DV_avg1)->array(mfi);
203 Array4<Real > const& DV_avg2 = (mf_DV_avg2)->array(mfi);
204 Array4<Real > const& ru2d = (mf_ru2d)->array(mfi);
205 Array4<Real > const& rv2d = (mf_rv2d)->array(mfi);
206 Array4<Real > const& rubar = (mf_rubar)->array(mfi);
207 Array4<Real > const& rvbar = (mf_rvbar)->array(mfi);
208 Array4<Real > const& rzeta = (mf_rzeta)->array(mfi);
209 Array4<Real const> const& visc2_p = mf_visc2_p->const_array(mfi);
210 Array4<Real const> const& visc2_r = mf_visc2_r->const_array(mfi);
211
212 Array4<Real const> const& pm = mf_pm->const_array(mfi);
213 Array4<Real const> const& pn = mf_pn->const_array(mfi);
214 Array4<Real const> const& fcor = mf_fcor->const_array(mfi);
215
216 Array4<Real const> const& mskr = mf_mskr->const_array(mfi);
217 Array4<Real const> const& msku = mf_msku->const_array(mfi);
218 Array4<Real const> const& mskv = mf_mskv->const_array(mfi);
219 Array4<Real const> const& mskp = mf_mskp->const_array(mfi);
220
221 Box bx = mfi.tilebox();
222 Box gbx = mfi.growntilebox();
223 Box gbx1 = mfi.growntilebox(IntVect(NGROW-1,NGROW-1,0));
224 Box gbx2 = mfi.growntilebox(IntVect(NGROW,NGROW,0));
225 Box gbx3 = mfi.growntilebox(IntVect(NGROW+1,NGROW+1,0));
226 Box xgbx2 = mfi.grownnodaltilebox(0, IntVect(NGROW,NGROW,0));
227 Box ygbx2 = mfi.grownnodaltilebox(1, IntVect(NGROW,NGROW,0));
228
229 Box xbxD = mfi.nodaltilebox(0);
230 xbxD.makeSlab(2,0);
231
232 Box ybxD = mfi.nodaltilebox(1);
233 ybxD.makeSlab(2,0);
234
235 Box xbxD_adj = mfi.nodaltilebox(0);
236 xbxD_adj.makeSlab(2,0);
237 Box ybxD_adj = mfi.nodaltilebox(1);
238 ybxD_adj.makeSlab(2,0);
239
240 auto xbxD_lo = lbound(xbxD_adj);
241 auto xbxD_hi = ubound(xbxD_adj);
242 auto ybxD_lo = lbound(ybxD_adj);
243 auto ybxD_hi = ubound(ybxD_adj);
244
245 if (xbxD_lo.x == dlo.x) {
246 xbxD_adj.growLo(0,-1);
247 } else if (xbxD_hi.x == dhi.x) {
248 xbxD_adj.growHi(0,-1);
249 }
250
251 if (ybxD_lo.y == dlo.y) {
252 ybxD_adj.growLo(1,-1);
253 } else if (ybxD_hi.y == dhi.y) {
254 ybxD_adj.growHi(1,-1);
255 }
256
257 Box tbxp1 = bx; tbxp1.grow(IntVect(NGROW-1,NGROW-1,0));
258 Box tbxp11 = bx; tbxp11.grow(IntVect(NGROW-1,NGROW-1,NGROW-1));
259 Box tbxp2 = bx; tbxp2.grow(IntVect(NGROW,NGROW,0));
260 Box tbxp3 = bx; tbxp3.grow(IntVect(NGROW+1,NGROW+1,0));
261
262 Box bxD = bx; bxD.makeSlab(2,0);
263 Box gbxD = gbx; gbxD.makeSlab(2,0);
264 Box gbx1D = gbx1; gbx1D.makeSlab(2,0);
265 Box gbx2D = gbx2; gbx2D.makeSlab(2,0);
266
267 Box tbxp2D = tbxp2;
268 tbxp2D.makeSlab(2,0);
269
270 FArrayBox fab_fomn(tbxp2,1,The_Async_Arena());
271
272 //step2d work arrays
273 FArrayBox fab_Drhs(tbxp3,1,The_Async_Arena());
274 FArrayBox fab_Dnew(tbxp2,1,The_Async_Arena());
275 FArrayBox fab_zwrk(tbxp1,1,The_Async_Arena());
276 FArrayBox fab_gzeta(tbxp1,1,The_Async_Arena());
277 FArrayBox fab_gzeta2(tbxp1,1,The_Async_Arena());
278 FArrayBox fab_gzetaSA(tbxp1,1,The_Async_Arena());
279 FArrayBox fab_Dstp(tbxp3,1,The_Async_Arena());
280 FArrayBox & fab_DUon=mf_DUon[mfi];
281 FArrayBox & fab_DVom=mf_DVom[mfi];
282 FArrayBox fab_rhs_ubar(xbxD,1,The_Async_Arena());
283 FArrayBox fab_rhs_vbar(ybxD,1,The_Async_Arena());
284 FArrayBox fab_rhs_zeta(tbxp1,1,The_Async_Arena());
285 FArrayBox fab_zeta_new(tbxp1,1,The_Async_Arena());
286
287 auto fomn=fab_fomn.array();
288
289 auto Drhs = fab_Drhs.array();
290 auto Drhs_const = fab_Drhs.const_array();
291 auto Dnew=fab_Dnew.array();
292 auto zwrk=fab_zwrk.array();
293 auto gzeta=fab_gzeta.array();
294 auto gzeta2=fab_gzeta2.array();
295 auto gzetaSA=fab_gzetaSA.array();
296 auto Dstp=fab_Dstp.array();
297 auto DUon=fab_DUon.array();
298 auto DVom=fab_DVom.array();
299 auto rhs_ubar=fab_rhs_ubar.array();
300 auto rhs_vbar=fab_rhs_vbar.array();
301 auto rhs_zeta=fab_rhs_zeta.array();
302 auto zeta_new=fab_zeta_new.array();
303
304 auto weight1 = vec_weight1.dataPtr();
305 auto weight2 = vec_weight2.dataPtr();
306
307 //From ana_grid.h and metrics.F
308 ParallelFor(xbxD, [=] AMREX_GPU_DEVICE (int i, int j, int)
309 {
310 rhs_ubar(i,j,0)=0.0_rt;
311 });
312
313 ParallelFor(ybxD, [=] AMREX_GPU_DEVICE (int i, int j, int)
314 {
315 rhs_vbar(i,j,0)=0.0_rt;
316 });
317
318 if (solverChoice.use_coriolis) {
319 ParallelFor(tbxp2D, [=] AMREX_GPU_DEVICE (int i, int j, int )
320 {
321 fomn(i,j,0) = fcor(i,j,0)*(1.0_rt/(pm(i,j,0)*pn(i,j,0)));
322 });
323 }
324
325 ParallelFor(makeSlab(tbxp3,2,0), [=] AMREX_GPU_DEVICE (int i, int j, int)
326 {
327 Drhs(i,j,0)=zeta(i,j,0,krhs)+h(i,j,0);
328 });
329
330 if(predictor_2d_step)
331 {
332 if(first_2d_step) {
333 Real cff2=(-1.0_rt/12.0_rt)*weight2[my_iif+1];
334 ParallelFor(makeSlab(gbx3,2,0),
335 [=] AMREX_GPU_DEVICE (int i, int j, int)
336 {
337 Zt_avg1(i,j,0)=0.0_rt;
338 });
339 ParallelFor(makeSlab(xgbx2,2,0),
340 [=] AMREX_GPU_DEVICE (int i, int j, int)
341 {
342 DU_avg1(i,j,0)=0.0_rt;
343 DU_avg2(i,j,0)=cff2*DUon(i,j,0);
344 });
345 ParallelFor(makeSlab(ygbx2,2,0),
346 [=] AMREX_GPU_DEVICE (int i, int j, int)
347 {
348 DV_avg1(i,j,0)=0.0_rt;
349 DV_avg2(i,j,0)=cff2*DVom(i,j,0);
350 });
351 }
352 else {
353 Real cff1_wt1 = weight1[my_iif-1];
354 Real cff2_wt1 = (8.0_rt/12.0_rt)*weight2[my_iif]-
355 (1.0_rt/12.0_rt)*weight2[my_iif+1];
356
357 ParallelFor(makeSlab(gbx3,2,0), [=] AMREX_GPU_DEVICE (int i, int j, int)
358 {
359 Zt_avg1(i,j,0) += cff1_wt1*zeta(i,j,0,krhs);
360 });
361
362 ParallelFor(makeSlab(xgbx2,2,0), [=] AMREX_GPU_DEVICE (int i, int j, int)
363 {
364 DU_avg1(i,j,0) += cff1_wt1*DUon(i,j,0);
365 DU_avg2(i,j,0) += cff2_wt1*DUon(i,j,0);
366 });
367
368 ParallelFor(makeSlab(ygbx2,2,0), [=] AMREX_GPU_DEVICE (int i, int j, int)
369 {
370 DV_avg1(i,j,0) += cff1_wt1*DVom(i,j,0);
371 DV_avg2(i,j,0) += cff2_wt1*DVom(i,j,0);
372 });
373 }
374 }
375 else {
376 Real cff2_wt2;
377
378 if (first_2d_step) {
379 cff2_wt2=weight2[my_iif];
380 } else {
381 cff2_wt2=5.0_rt/12.0_rt*weight2[my_iif];
382 }
383
384 ParallelFor(makeSlab(xgbx2,2,0), [=] AMREX_GPU_DEVICE (int i, int j, int)
385 {
386 DU_avg2(i,j,0)=DU_avg2(i,j,0)+cff2_wt2*DUon(i,j,0);
387 });
388
389 ParallelFor(makeSlab(ygbx2,2,0), [=] AMREX_GPU_DEVICE (int i, int j, int)
390 {
391 DV_avg2(i,j,0)=DV_avg2(i,j,0)+cff2_wt2*DVom(i,j,0);
392 });
393 }
394 //
395 // Do not perform the actual time stepping during the auxiliary
396 // (nfast(ng)+1) time step. Jump to next box
397 //
398
399 if (my_iif>=nfast) {
400 continue; }
401 //Load new free-surface values into shared array at both predictor
402 //and corrector steps
403 //
404 //=======================================================================
405 // Time step free-surface equation.
406 //=======================================================================
407 //
408 // During the first time-step, the predictor step is Forward-Euler
409 // and the corrector step is Backward-Euler. Otherwise, the predictor
410 // step is Leap-frog and the corrector step is Adams-Moulton.
411 //
412
413 // todo: gzeta
414
415 // todo: HACKHACKHACK Should use rho0 from prob.H
416 Real fac=1000.0_rt/1025.0_rt;
417
418 if (my_iif==0) {
419 Real cff1=dtfast_lev;
420
421 ParallelFor(makeSlab(tbxp1,2,0), [=] AMREX_GPU_DEVICE (int i, int j, int )
422 {
423 rhs_zeta(i,j,0) = (DUon(i,j,0)-DUon(i+1,j,0))+
424 (DVom(i,j,0)-DVom(i,j+1,0));
425 zeta_new(i,j,0) = (zeta(i,j,0,kstp)+ pm(i,j,0)*pn(i,j,0)*cff1*rhs_zeta(i,j,0)) * mskr(i,j,0);
426 Dnew(i,j,0) = zeta_new(i,j,0)+h(i,j,0);
427
428 //Pressure gradient terms:
429 zwrk(i,j,0)=0.5_rt*(zeta(i,j,0,kstp)+zeta_new(i,j,0));
430 gzeta(i,j,0)=(fac+rhoS(i,j,0))*zwrk(i,j,0);
431 gzeta2(i,j,0)=gzeta(i,j,0)*zwrk(i,j,0);
432 gzetaSA(i,j,0)=zwrk(i,j,0)*(rhoS(i,j,0)-rhoA(i,j,0));
433 });
434
435 } else if (predictor_2d_step) {
436
437 Real cff1=2.0_rt * dtfast_lev;
438 Real cff4=4.0_rt / 25.0_rt;
439 Real cff5=1.0_rt - 2.0_rt*cff4;
440
441 ParallelFor(makeSlab(tbxp1,2,0), [=] AMREX_GPU_DEVICE (int i, int j, int )
442 {
443 rhs_zeta(i,j,0)=(DUon(i,j,0)-DUon(i+1,j,0))+
444 (DVom(i,j,0)-DVom(i,j+1,0));
445 zeta_new(i,j,0)=(zeta(i,j,0,kstp)+
446 pm(i,j,0)*pn(i,j,0)*cff1*rhs_zeta(i,j,0)) * mskr(i,j,0);
447 Dnew(i,j,0)=zeta_new(i,j,0)+h(i,j,0);
448 //Pressure gradient terms
449 zwrk(i,j,0)=cff5*zeta(i,j,0,krhs)+
450 cff4*(zeta(i,j,0,kstp)+zeta_new(i,j,0));
451 gzeta(i,j,0)=(fac+rhoS(i,j,0))*zwrk(i,j,0);
452 gzeta2(i,j,0)=gzeta(i,j,0)*zwrk(i,j,0);
453 gzetaSA(i,j,0)=zwrk(i,j,0)*(rhoS(i,j,0)-rhoA(i,j,0));
454 });
455
456 } else if (!predictor_2d_step) { //AKA if(corrector_2d_step)
457
458 Real cff1=dtfast_lev * 5.0_rt/12.0_rt;
459 Real cff2=dtfast_lev * 8.0_rt/12.0_rt;
460 Real cff3=dtfast_lev * 1.0_rt/12.0_rt;
461 Real cff4=2.0_rt/5.0_rt;
462 Real cff5=1.0_rt-cff4;
463
464 ParallelFor(makeSlab(tbxp1,2,0), [=] AMREX_GPU_DEVICE (int i, int j, int )
465 {
466 Real cff=cff1*((DUon(i,j,0)-DUon(i+1,j,0))+
467 (DVom(i,j,0)-DVom(i,j+1,0)));
468 zeta_new(i,j,0)=zeta(i,j,0,kstp)+
469 pm(i,j,0)*pn(i,j,0)*(cff+
470 cff2*rzeta(i,j,0,kstp)-
471 cff3*rzeta(i,j,0,ptsk));
472 zeta_new(i,j,0) *= mskr(i,j,0);
473 Dnew(i,j,0)=zeta_new(i,j,0)+h(i,j,0);
474 //Pressure gradient terms
475 zwrk(i,j,0)=cff5*zeta_new(i,j,0)+cff4*zeta(i,j,0,krhs);
476 gzeta(i,j,0)=(fac+rhoS(i,j,0))*zwrk(i,j,0);
477 gzeta2(i,j,0)=gzeta(i,j,0)*zwrk(i,j,0);
478 gzetaSA(i,j,0)=zwrk(i,j,0)*(rhoS(i,j,0)-rhoA(i,j,0));
479 });
480 }
481
482 //
483 // Load new free-surface values into shared array at both predictor
484 // and corrector steps.
485 //
486 //// zeta(knew) only valid at zeta_new, i.e. tbxp1
487 ParallelFor(makeSlab(gbx1,2,0),
488 [=] AMREX_GPU_DEVICE (int i, int j, int )
489 {
490 zeta(i,j,0,knew) = zeta_new(i,j,0);
491 });
492
493 //
494 // If predictor step, load right-side-term into shared array.
495 //
496 if (predictor_2d_step) {
497 ParallelFor(makeSlab(gbx1,2,0), [=] AMREX_GPU_DEVICE (int i, int j, int )
498 {
499 rzeta(i,j,0,krhs)=rhs_zeta(i,j,0);
500 });
501 }
502
503 //
504 //=======================================================================
505 // Compute right-hand-side for the 2D momentum equations.
506 //=======================================================================
507 //
508/*
509!
510!-----------------------------------------------------------------------
511! Compute pressure gradient terms.
512!-----------------------------------------------------------------------
513!
514*/
515
516 Real cff1 = 0.5_rt * g;
517 Real cff2 = 1.0_rt / 3.0_rt;
518 ParallelFor(xbxD,
519 [=] AMREX_GPU_DEVICE (int i, int j, int )
520 {
521 Real on_u = 2.0_rt / (pn(i,j,0)+pn(i-1,j,0));
522 rhs_ubar(i,j,0)=cff1 * on_u *
523 (( h(i-1,j,0) + h(i,j,0))*
524 (gzeta(i-1,j,0) - gzeta(i,j,0))+
525 ( h(i-1,j,0) - h(i,j,0))*
526 ( gzetaSA(i-1,j,0) + gzetaSA(i,j,0)+
527 cff2*( rhoA(i-1,j,0) - rhoA(i,j,0))*
528 ( zwrk(i-1,j,0) - zwrk(i,j,0)))+
529 (gzeta2(i-1,j,0)- gzeta2(i ,j,0)));
530 });
531
532 ParallelFor(ybxD,
533 [=] AMREX_GPU_DEVICE (int i, int j, int )
534 {
535 Real om_v = 2.0_rt / (pm(i,j,0)+pm(i,j-1,0));
536 rhs_vbar(i,j,0) = cff1*om_v *
537 (( h(i,j-1,0) + h(i,j,0))*
538 (gzeta(i,j-1,0) - gzeta(i,j,0))+
539 ( h(i,j-1,0) - h(i,j,0))*
540 (gzetaSA(i,j-1,0)+ gzetaSA(i,j ,0)+
541 cff2*(rhoA(i,j-1,0)- rhoA(i,j ,0))*
542 (zwrk(i,j-1,0)- zwrk(i,j ,0)))+
543 (gzeta2(i,j-1,0)- gzeta2(i,j ,0)));
544 });
545
546 // Advection terms for 2d ubar, vbar added to rhs_ubar and rhs_vbar
547 //
548 //-----------------------------------------------------------------------
549 // rhs_uv_2d
550 //-----------------------------------------------------------------------
551 //
552 Array4<Real const> const& ubar_const = mf_ubar->const_array(mfi);
553 Array4<Real const> const& vbar_const = mf_vbar->const_array(mfi);
554
555 rhs_uv_2d(xbxD, ybxD, ubar_const, vbar_const, rhs_ubar, rhs_vbar, DUon, DVom, krhs);
556
557 //-----------------------------------------------------------------------
558 // Add Coriolis forcing
559 //-----------------------------------------------------------------------
560 if (solverChoice.use_coriolis) {
561 // Coriolis terms for 2d ubar, vbar added to rhs_ubar and rhs_vbar
562 //
563 //-----------------------------------------------------------------------
564 // coriolis
565 //-----------------------------------------------------------------------
566 //
567 coriolis(xbxD, ybxD, ubar_const, vbar_const, rhs_ubar, rhs_vbar, Drhs, fomn, krhs, 0);
568 }
569 //-----------------------------------------------------------------------
570 //Add in horizontal harmonic viscosity.
571 // Consider generalizing or copying uv3dmix, where Drhs is used instead of Hz and u=>ubar v=>vbar, drop dt terms
572 //-----------------------------------------------------------------------
573 uv3dmix(xbxD, ybxD, ubar, vbar, ubar, vbar, rhs_ubar, rhs_vbar,
574 visc2_p, visc2_r, Drhs_const,
575 pm, pn, mskp, krhs, nnew, 0.0_rt);
576
577#ifdef REMORA_USE_NETCDF
578 if (solverChoice.do_m2_clim_nudg) {
579 Array4<const Real> const& ubar_clim = ubar_clim_data_from_file->mf_interpolated->const_array(mfi);
580 Array4<const Real> const& vbar_clim = vbar_clim_data_from_file->mf_interpolated->const_array(mfi);
581 Array4<const Real> const& ubar_nudg_coeff = vec_nudg_coeff[BdyVars::ubar][lev]->const_array(mfi);
582 Array4<const Real> const& vbar_nudg_coeff = vec_nudg_coeff[BdyVars::vbar][lev]->const_array(mfi);
583 // Boxes are like this to match ROMS
584 apply_clim_nudg(xbxD_adj, 1, 0, rhs_ubar, ubar_krhs, ubar_clim, ubar_nudg_coeff, Drhs_const, pm, pn);
585 apply_clim_nudg(ybxD_adj, 0, 1, rhs_vbar, vbar_krhs, vbar_clim, vbar_nudg_coeff, Drhs_const, pm, pn);
586 }
587#endif
588
589 //-----------------------------------------------------------------------
590 // Coupling from 3d to 2d
591 //-----------------------------------------------------------------------
592 if (first_2d_step&&predictor_2d_step)
593 {
594 if (iic==ntfirst) {
595 ParallelFor(xbxD, [=] AMREX_GPU_DEVICE (int i, int j, int )
596 {
597 rufrc(i,j,0) -= rhs_ubar(i,j,0);
598 rhs_ubar(i,j,0) += rufrc(i,j,0);
599 ru2d(i,j,0,nstp) = rufrc(i,j,0);
600 });
601
602 ParallelFor(ybxD, [=] AMREX_GPU_DEVICE (int i, int j, int )
603 {
604 rvfrc(i,j,0) -= rhs_vbar(i,j,0);
605 rhs_vbar(i,j,0) += rvfrc(i,j,0);
606 rv2d(i,j,0,nstp) = rvfrc(i,j,0);
607 });
608
609 } else if (iic==(ntfirst+1)) {
610
611 ParallelFor(xbxD, [=] AMREX_GPU_DEVICE (int i, int j, int )
612 {
613 rufrc(i,j,0)=rufrc(i,j,0)-rhs_ubar(i,j,0);
614 rhs_ubar(i,j,0)=rhs_ubar(i,j,0)+1.5_rt*rufrc(i,j,0)-0.5_rt*ru2d(i,j,0,0);
615 ru2d(i,j,0,1)=rufrc(i,j,0);
616 Real r_swap= ru2d(i,j,0,1);
617 ru2d(i,j,0,1) = ru2d(i,j,0,0);
618 ru2d(i,j,0,0) = r_swap;
619 });
620
621 ParallelFor(ybxD, [=] AMREX_GPU_DEVICE (int i, int j, int )
622 {
623 rvfrc(i,j,0)=rvfrc(i,j,0)-rhs_vbar(i,j,0);
624 rhs_vbar(i,j,0)=rhs_vbar(i,j,0)+1.5_rt*rvfrc(i,j,0)-0.5_rt*rv2d(i,j,0,0);
625 rv2d(i,j,0,1)=rvfrc(i,j,0);
626 Real r_swap= rv2d(i,j,0,1);
627 rv2d(i,j,0,1) = rv2d(i,j,0,0);
628 rv2d(i,j,0,0) = r_swap;
629 });
630
631 } else {
632 cff1=23.0_rt/12.0_rt;
633 cff2=16.0_rt/12.0_rt;
634 Real cff3= 5.0_rt/12.0_rt;
635
636 ParallelFor(xbxD, [=] AMREX_GPU_DEVICE (int i, int j, int )
637 {
638 rufrc(i,j,0)=rufrc(i,j,0)-rhs_ubar(i,j,0);
639 rhs_ubar(i,j,0)=rhs_ubar(i,j,0)+
640 cff1*rufrc(i,j,0)-
641 cff2*ru2d(i,j,0,0)+
642 cff3*ru2d(i,j,0,1);
643 ru2d(i,j,0,1)=rufrc(i,j,0);
644 Real r_swap= ru2d(i,j,0,1);
645 ru2d(i,j,0,1) = ru2d(i,j,0,0);
646 ru2d(i,j,0,0) = r_swap;
647 });
648
649 ParallelFor(ybxD, [=] AMREX_GPU_DEVICE (int i, int j, int )
650 {
651 rvfrc(i,j,0)=rvfrc(i,j,0)-rhs_vbar(i,j,0);
652 rhs_vbar(i,j,0)=rhs_vbar(i,j,0)+
653 cff1*rvfrc(i,j,0)-
654 cff2*rv2d(i,j,0,0)+
655 cff3*rv2d(i,j,0,1);
656 rv2d(i,j,0,1)=rvfrc(i,j,0);
657
658 Real r_swap= rv2d(i,j,0,1);
659 rv2d(i,j,0,1) = rv2d(i,j,0,0);
660 rv2d(i,j,0,0) = r_swap;
661 });
662 }
663 } else {
664 ParallelFor(xbxD, [=] AMREX_GPU_DEVICE (int i, int j, int )
665 {
666 rhs_ubar(i,j,0) += rufrc(i,j,0);
667 });
668
669 ParallelFor(ybxD, [=] AMREX_GPU_DEVICE (int i, int j, int )
670 {
671 rhs_vbar(i,j,0) += rvfrc(i,j,0);
672 });
673 }
674
675 //
676 //=======================================================================
677 // Time step 2D momentum equations.
678 //=======================================================================
679 //
680 // Compute total water column depth.
681 //
682 ParallelFor(makeSlab(tbxp3,2,0), [=] AMREX_GPU_DEVICE (int i, int j, int )
683 {
684 Dstp(i,j,0)=zeta(i,j,0,kstp)+h(i,j,0);
685 });
686
687 //
688 // During the first time-step, the predictor step is Forward-Euler
689 // and the corrector step is Backward-Euler. Otherwise, the predictor
690 // step is Leap-frog and the corrector step is Adams-Moulton.
691 //
692 if (my_iif==0) {
693 cff1=0.5_rt*dtfast_lev;
694 ParallelFor(xbxD,
695 [=] AMREX_GPU_DEVICE (int i, int j, int )
696 {
697 Real cff=(pm(i,j,0)+pm(i-1,j,0))*(pn(i,j,0)+pn(i-1,j,0));
698 Real Dnew_avg =1.0_rt/(Dnew(i,j,0)+Dnew(i-1,j,0));
699 ubar(i,j,0,knew)=(ubar(i,j,0,kstp)*
700 (Dstp(i,j,0)+Dstp(i-1,j,0))+
701 cff*cff1*rhs_ubar(i,j,0))*Dnew_avg * msku(i,j,0);
702 });
703 ParallelFor(ybxD,
704 [=] AMREX_GPU_DEVICE (int i, int j, int )
705 {
706 Real cff=(pm(i,j,0)+pm(i,j-1,0))*(pn(i,j,0)+pn(i,j-1,0));
707 Real Dnew_avg=1.0_rt/(Dnew(i,j,0)+Dnew(i,j-1,0));
708 vbar(i,j,0,knew)=(vbar(i,j,0,kstp)*
709 (Dstp(i,j,0)+Dstp(i,j-1,0))+
710 cff*cff1*rhs_vbar(i,j,0))*Dnew_avg * mskv(i,j,0);
711 });
712
713 } else if (predictor_2d_step) {
714
715 cff1=dtfast_lev;
716 ParallelFor(xbxD,
717 [=] AMREX_GPU_DEVICE (int i, int j, int )
718 {
719 Real cff=(pm(i,j,0)+pm(i-1,j,0))*(pn(i,j,0)+pn(i-1,j,0));
720 Real Dnew_avg=1.0_rt/(Dnew(i,j,0)+Dnew(i-1,j,0));
721 ubar(i,j,0,knew)=(ubar(i,j,0,kstp)*
722 (Dstp(i,j,0)+Dstp(i-1,j,0))+
723 cff*cff1*rhs_ubar(i,j,0))*Dnew_avg * msku(i,j,0);
724 });
725 ParallelFor(ybxD,
726 [=] AMREX_GPU_DEVICE (int i, int j, int )
727 {
728 Real cff=(pm(i,j,0)+pm(i,j-1,0))*(pn(i,j,0)+pn(i,j-1,0));
729 Real Dnew_avg=1.0_rt/(Dnew(i,j,0)+Dnew(i,j-1,0));
730 vbar(i,j,0,knew)=(vbar(i,j,0,kstp)*
731 (Dstp(i,j,0)+Dstp(i,j-1,0))+
732 cff*cff1*rhs_vbar(i,j,0))*Dnew_avg * mskv(i,j,0);
733 });
734
735 } else if ((!predictor_2d_step)) {
736
737 cff1=0.5_rt*dtfast_lev*5.0_rt/12.0_rt;
738 cff2=0.5_rt*dtfast_lev*8.0_rt/12.0_rt;
739 Real cff3=0.5_rt*dtfast_lev*1.0_rt/12.0_rt;
740 ParallelFor(xbxD,
741 [=] AMREX_GPU_DEVICE (int i, int j, int )
742 {
743 Real cff=(pm(i,j,0)+pm(i-1,j,0))*(pn(i,j,0)+pn(i-1,j,0));
744 Real Dnew_avg=1.0_rt/(Dnew(i,j,0)+Dnew(i-1,j,0));
745 ubar(i,j,0,knew)=(ubar(i,j,0,kstp)*
746 (Dstp(i,j,0)+Dstp(i-1,j,0))+
747 cff*(cff1*rhs_ubar(i,j,0)+
748 cff2*rubar(i,j,0,kstp)-
749 cff3*rubar(i,j,0,ptsk)))*Dnew_avg * msku(i,j,0);
750 });
751 ParallelFor(ybxD,
752 [=] AMREX_GPU_DEVICE (int i, int j, int )
753 {
754 Real cff=(pm(i,j,0)+pm(i,j-1,0))*(pn(i,j,0)+pn(i,j-1,0));
755 Real Dnew_avg=1.0_rt/(Dnew(i,j,0)+Dnew(i,j-1,0));
756 vbar(i,j,0,knew)=(vbar(i,j,0,kstp)*
757 (Dstp(i,j,0)+Dstp(i,j-1,0))+
758 cff*(cff1*rhs_vbar(i,j,0)+
759 cff2*rvbar(i,j,0,kstp)-
760 cff3*rvbar(i,j,0,ptsk)))*Dnew_avg * mskv(i,j,0);
761 });
762 }
763
764 //store rhs_ubar and rhs_vbar to save later
765 //
766 // If predictor step, load right-side-term into shared arrays for
767 // future use during the subsequent corrector step.
768 //
769
770 if (predictor_2d_step) {
771 ParallelFor(xbxD, [=] AMREX_GPU_DEVICE (int i, int j, int )
772 {
773 rubar(i,j,0,krhs)=rhs_ubar(i,j,0);
774 });
775 ParallelFor(ybxD, [=] AMREX_GPU_DEVICE (int i, int j, int )
776 {
777 rvbar(i,j,0,krhs)=rhs_vbar(i,j,0);
778 });
779 }
780 }
781
782 // Don't do the FillPatch or rivers at the last truncated predictor step.
783 // We may need to move the zeta FillPatch further up
784 if (my_iif<nfast) {
785 int know;
786 Real dt2d;
787 if (my_iif==0) {
788 know = krhs;
789 dt2d = dtfast_lev;
790 } else if (predictor_2d_step) {
791 know = krhs;
792 dt2d = 2.0_rt * dtfast_lev;
793 } else {
794 know = kstp;
795 dt2d = dtfast_lev;
796 }
797
798 FillPatch(lev, t_old[lev], *vec_ubar[lev], GetVecOfPtrs(vec_ubar), BCVars::ubar_bc, BdyVars::ubar,
799 knew, false,true, 0,know, dt2d);
800 FillPatch(lev, t_old[lev], *vec_vbar[lev], GetVecOfPtrs(vec_vbar), BCVars::vbar_bc, BdyVars::vbar,
801 knew, false,true, 0,know, dt2d);
802 FillPatch(lev, t_old[lev], *vec_zeta[lev], GetVecOfPtrs(vec_zeta), BCVars::zeta_bc, BdyVars::zeta,
803 knew, false,false, 0,know, dt2d);
804
805#ifdef REMORA_USE_NETCDF
806 if (solverChoice.do_rivers) {
807 river_source_transportbar->update_interpolated_to_time(t_old[lev]);
808 int* river_direction_d = river_direction.data();
809 for ( MFIter mfi(*mf_rhoS, TilingIfNotGPU()); mfi.isValid(); ++mfi )
810 {
811 Array4<const int > const& river_pos = vec_river_position[lev]->const_array(mfi);
812 Array4<const Real> const& river_transportbar = river_source_transportbar->fab_interp->array();
813 Array4<Real > const& ubar = mf_ubar->array(mfi);
814 Array4<Real > const& vbar = mf_vbar->array(mfi);
815 Array4<Real const> const& zeta = mf_zeta->const_array(mfi);
816 Array4<Real const> const& h = mf_h->const_array(mfi);
817 Array4<Real const> const& pm = mf_pm->const_array(mfi);
818 Array4<Real const> const& pn = mf_pn->const_array(mfi);
819
820 Box gbx1D = mfi.growntilebox(IntVect(NGROW-1,NGROW-1,0));
821 gbx1D.makeSlab(2,0);
822
823 ParallelFor(gbx1D, [=] AMREX_GPU_DEVICE (int i, int j, int )
824 {
825 int iriver = river_pos(i,j,0);
826 if (iriver >= 0) {
827 if (river_direction_d[iriver] == 0) {
828 Real on_u = 2.0_rt / (pn(i,j,0)+pn(i-1,j,0));
829 Real cff = 1.0_rt / (on_u * 0.5_rt * (zeta(i-1,j,0,knew) + h(i-1,j,0) +
830 zeta(i,j,0,knew) + h(i,j,0)));
831 ubar(i,j,0,knew) = river_transportbar(iriver,0,0) * cff;
832 } else {
833 Real om_v = 2.0_rt / (pm(i,j,0)+pm(i,j-1,0));
834 Real cff = 1.0_rt / (om_v * 0.5_rt * (zeta(i,j-1,0,knew) + h(i,j-1,0) +
835 zeta(i,j,0,knew) + h(i,j,0)));
836 vbar(i,j,0,knew) = river_transportbar(iriver,0,0) * cff;
837 }
838 }
839 });
840 }
841 }
842 FillPatchNoBC(lev, t_old[lev], *vec_ubar[lev], GetVecOfPtrs(vec_ubar), BdyVars::ubar,
843 knew, false,false);
844 FillPatchNoBC(lev, t_old[lev], *vec_vbar[lev], GetVecOfPtrs(vec_vbar), BdyVars::vbar,
845 knew, false,false);
846#endif
847 }
848}
g
constexpr amrex::Real g
Definition REMORA_Constants.H:16
NGROW
#define NGROW
Definition REMORA_IndexDefines.H:13
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_transportbar
NCTimeSeriesRiver * river_source_transportbar
Data container for vertically integrated momentum transport in rivers.
Definition REMORA.H:1049
REMORA::nfast
int nfast
Number of fast steps to take.
Definition REMORA.H:1196
REMORA::rhs_uv_2d
void rhs_uv_2d(const amrex::Box &xbx, const amrex::Box &ybx, const amrex::Array4< amrex::Real const > &uold, const amrex::Array4< amrex::Real const > &vold, const amrex::Array4< amrex::Real > &ru, const amrex::Array4< amrex::Real > &rv, const amrex::Array4< amrex::Real const > &Duon, const amrex::Array4< amrex::Real const > &Dvom, const int nrhs)
RHS terms for 2D momentum.
Definition REMORA_rhs_uv_2d.cpp:17
REMORA::vbar_clim_data_from_file
NCTimeSeries * vbar_clim_data_from_file
Data container for vbar climatology data read from file.
Definition REMORA.H:1034
REMORA::FillPatchNoBC
void FillPatchNoBC(int lev, amrex::Real time, amrex::MultiFab &mf_to_be_filled, amrex::Vector< amrex::MultiFab * > const &mfs, const int bdy_var_type=BdyVars::null, const int icomp=0, const bool fill_all=true, const bool fill_set=true)
Fill a new MultiFab by copying in phi from valid region and filling ghost cells without applying boun...
Definition REMORA_FillPatch.cpp:195
REMORA::ubar_clim_data_from_file
NCTimeSeries * ubar_clim_data_from_file
Data container for ubar climatology data read from file.
Definition REMORA.H:1032
REMORA::vec_weight2
amrex::Vector< amrex::Real > vec_weight2
Weights for calculating avg2 in 2D advance.
Definition REMORA.H:402
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:1110
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::uv3dmix
void uv3dmix(const amrex::Box &xbx, const amrex::Box &ybx, const amrex::Array4< amrex::Real > &u, const amrex::Array4< amrex::Real > &v, const amrex::Array4< amrex::Real const > &uold, const amrex::Array4< amrex::Real const > &vold, const amrex::Array4< amrex::Real > &rufrc, const amrex::Array4< amrex::Real > &rvfrc, const amrex::Array4< amrex::Real const > &visc2_p, const amrex::Array4< amrex::Real const > &visc2_r, const amrex::Array4< amrex::Real const > &Hz, const amrex::Array4< amrex::Real const > &pm, const amrex::Array4< amrex::Real const > &pn, const amrex::Array4< amrex::Real const > &mskp, int nrhs, int nnew, const amrex::Real dt_lev)
Harmonic viscosity.
Definition REMORA_uv3dmix.cpp:25
REMORA::t_new
amrex::Vector< amrex::Real > t_new
new time at each level
Definition REMORA.H:1114
REMORA::solverChoice
static SolverChoice solverChoice
Container for algorithmic choices.
Definition REMORA.H:1237
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:1052
REMORA::vec_weight1
amrex::Vector< amrex::Real > vec_weight1
Weights for calculating avg1 in 2D advance.
Definition REMORA.H:400
REMORA::vec_zeta
amrex::Vector< std::unique_ptr< amrex::MultiFab > > vec_zeta
free surface height (2D)
Definition REMORA.H:346
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:1054
REMORA::vec_vbar
amrex::Vector< std::unique_ptr< amrex::MultiFab > > vec_vbar
barotropic y velocity (2D)
Definition REMORA.H:344
REMORA::advance_2d
void advance_2d(int lev, amrex::MultiFab const *mf_rhoS, amrex::MultiFab const *mf_rhoA, amrex::MultiFab *mf_ru2d, amrex::MultiFab *mf_rv2d, amrex::MultiFab *mf_rufrc, amrex::MultiFab *mf_rvfrc, amrex::MultiFab *mf_Zt_avg1, 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_rubar, std::unique_ptr< amrex::MultiFab > &mf_rvbar, std::unique_ptr< amrex::MultiFab > &mf_rzeta, std::unique_ptr< amrex::MultiFab > &mf_ubar, std::unique_ptr< amrex::MultiFab > &mf_vbar, amrex::MultiFab *mf_zeta, amrex::MultiFab const *mf_h, amrex::MultiFab const *mf_pm, amrex::MultiFab const *mf_pn, amrex::MultiFab const *mf_fcor, amrex::MultiFab const *mf_visc2_p, amrex::MultiFab const *mf_visc2_r, amrex::MultiFab const *mf_mskr, amrex::MultiFab const *mf_msku, amrex::MultiFab const *mf_mskv, amrex::MultiFab const *mf_mskp, amrex::Real dtfast_lev, bool predictor_2d_step, bool first_2d_step, int my_iif, int &next_indx1)
Perform a 2D predictor (predictor_2d_step=True) or corrector (predictor_2d_step=False) step.
Definition REMORA_advance_2d.cpp:44
REMORA::vec_ubar
amrex::Vector< std::unique_ptr< amrex::MultiFab > > vec_ubar
barotropic x velocity (2D)
Definition REMORA.H:342
REMORA::vec_nudg_coeff
amrex::Vector< amrex::Vector< std::unique_ptr< amrex::MultiFab > > > vec_nudg_coeff
Climatology nudging coefficients.
Definition REMORA.H:417
REMORA::t_old
amrex::Vector< amrex::Real > t_old
old time at each level
Definition REMORA.H:1116
REMORA::coriolis
void coriolis(const amrex::Box &xbx, const amrex::Box &ybx, const amrex::Array4< amrex::Real const > &uold, const amrex::Array4< amrex::Real const > &vold, const amrex::Array4< amrex::Real > &ru, const amrex::Array4< amrex::Real > &rv, const amrex::Array4< amrex::Real const > &Hz, const amrex::Array4< amrex::Real const > &fomn, int nrhs, int nr)
Calculate Coriolis terms.
Definition REMORA_coriolis.cpp:19