81 auto dz_ptr = stretched_dz_d.data();
83 for (
int n(0); n<num_comp; ++n) {
84 const int qty_index = start_comp + n;
88 ParallelFor(
xbx, [=] AMREX_GPU_DEVICE (
int i,
int j,
int k) noexcept
90 const int prim_index = qty_index - 1;
93 Real rhoFace = 0.5 * ( cell_data(i, j, k,
Rho_comp) + cell_data(i-1, j, k,
Rho_comp) );
94 Real rhoAlpha = rhoFace *
d_alpha_eff[prim_scal_index];
102 bool SurfLayer_on_zlo = ( use_SurfLayer && rotate && k ==
dom_lo.z);
104 Real GradCx =
dx_inv * ( cell_prim(i, j, k , prim_index) - cell_prim(i-1, j, k , prim_index) );
107 xflux(i,j,k) = hfx_x(i,j,0);
108 }
else if (SurfLayer_on_zlo && (qty_index ==
RhoQ1_comp)) {
109 xflux(i,j,k) = qfx1_x(i,j,0);
111 xflux(i,j,k) = -rhoAlpha * mf_ux(i,j,0) * GradCx;
114 ParallelFor(
ybx, [=] AMREX_GPU_DEVICE (
int i,
int j,
int k) noexcept
116 const int prim_index = qty_index - 1;
119 Real rhoFace = 0.5 * ( cell_data(i, j, k,
Rho_comp) + cell_data(i, j-1, k,
Rho_comp) );
120 Real rhoAlpha = rhoFace *
d_alpha_eff[prim_scal_index];
121 rhoAlpha += 0.5 * ( mu_turb(i, j , k,
d_eddy_diff_idy[prim_scal_index])
127 bool SurfLayer_on_zlo = ( use_SurfLayer && rotate && k ==
dom_lo.z);
129 Real GradCy =
dy_inv * ( cell_prim(i, j, k , prim_index) - cell_prim(i, j-1, k , prim_index) );
132 yflux(i,j,k) = hfx_y(i,j,0);
133 }
else if (SurfLayer_on_zlo && (qty_index ==
RhoQ1_comp)) {
134 yflux(i,j,k) = qfx1_y(i,j,0);
136 yflux(i,j,k) = -rhoAlpha * mf_vy(i,j,0) * GradCy;
139 ParallelFor(
zbx, [=] AMREX_GPU_DEVICE (
int i,
int j,
int k) noexcept
141 const int prim_index = qty_index - 1;
144 Real rhoFace = 0.5 * ( cell_data(i, j, k,
Rho_comp) + cell_data(i, j, k-1,
Rho_comp) );
145 Real rhoAlpha = rhoFace *
d_alpha_eff[prim_scal_index];
146 rhoAlpha += 0.5 * ( mu_turb(i, j, k ,
d_eddy_diff_idz[prim_scal_index])
159 bool SurfLayer_on_zlo = ( use_SurfLayer && k ==
dom_lo.z);
161 if (ext_dir_on_zlo) {
163 Real dz0 = dz_ptr[k];
164 Real dz1 = dz_ptr[k+1];
165 Real idz0 = 1.0 / dz0;
166 Real f = (dz1 / dz0) + 2.0;
168 Real c3 = 2.0 / (f - f2);
170 Real
c1 = -(1.0-f2)*c3;
171 GradCz = idz0 * (
c1 * cell_prim(i, j, k-1, prim_index)
172 +
c2 * cell_prim(i, j, k , prim_index)
173 + c3 * cell_prim(i, j, k+1, prim_index) );
174 }
else if (ext_dir_on_zhi) {
176 Real dz0 = dz_ptr[k-1];
177 Real dz1 = dz_ptr[k-2];
178 Real idz0 = 1.0 / dz0;
179 Real f = (dz1 / dz0) + 2.0;
181 Real c3 = 2.0 / (f - f2);
183 Real
c1 = -(1.0-f2)*c3;
184 GradCz = idz0 * ( -(
c1 * cell_prim(i, j, k , prim_index)
185 +
c2 * cell_prim(i, j, k-1, prim_index)
186 + c3 * cell_prim(i, j, k-2, prim_index) ) );
188 Real dzk_inv = 2.0 / (dz_ptr[k] + dz_ptr[k-1]);
189 GradCz = dzk_inv * ( cell_prim(i, j, k, prim_index) - cell_prim(i, j, k-1, prim_index) );
193 zflux(i,j,k) = hfx_z(i,j,0);
194 }
else if (SurfLayer_on_zlo && (qty_index ==
RhoQ1_comp)) {
195 zflux(i,j,k) = qfx1_z(i,j,0);
197 zflux(i,j,k) = -rhoAlpha * GradCz;
201 if (!SurfLayer_on_zlo) {
202 hfx_z(i,j,k) = zflux(i,j,k);
205 if (!SurfLayer_on_zlo) {
206 qfx1_z(i,j,k) = zflux(i,j,k);
209 qfx2_z(i,j,k) = zflux(i,j,k);
214 ParallelFor(
xbx, [=] AMREX_GPU_DEVICE (
int i,
int j,
int k) noexcept
216 const int prim_index = qty_index - 1;
225 bool SurfLayer_on_zlo = ( use_SurfLayer && rotate && k ==
dom_lo.z);
227 Real GradCx =
dx_inv * ( cell_prim(i, j, k , prim_index) - cell_prim(i-1, j, k , prim_index) );
230 xflux(i,j,k) = hfx_x(i,j,0);
231 }
else if (SurfLayer_on_zlo && (qty_index ==
RhoQ1_comp)) {
232 xflux(i,j,k) = qfx1_x(i,j,0);
234 xflux(i,j,k) = -rhoAlpha * mf_ux(i,j,0) * GradCx;
237 ParallelFor(
ybx, [=] AMREX_GPU_DEVICE (
int i,
int j,
int k) noexcept
239 const int prim_index = qty_index - 1;
248 bool SurfLayer_on_zlo = ( use_SurfLayer && rotate && k ==
dom_lo.z);
250 Real GradCy =
dy_inv * ( cell_prim(i, j, k , prim_index) - cell_prim(i, j-1, k , prim_index) );
253 yflux(i,j,k) = hfx_y(i,j,0);
254 }
else if (SurfLayer_on_zlo && (qty_index ==
RhoQ1_comp)) {
255 yflux(i,j,k) = qfx1_y(i,j,0);
257 yflux(i,j,k) = -rhoAlpha * mf_vy(i,j,0) * GradCy;
260 ParallelFor(
zbx, [=] AMREX_GPU_DEVICE (
int i,
int j,
int k) noexcept
262 const int prim_index = qty_index - 1;
279 bool SurfLayer_on_zlo = ( use_SurfLayer && k ==
dom_lo.z);
281 if (ext_dir_on_zlo) {
283 Real dz0 = dz_ptr[k];
284 Real dz1 = dz_ptr[k+1];
285 Real idz0 = 1.0 / dz0;
286 Real f = (dz1 / dz0) + 2.0;
288 Real c3 = 2.0 / (f - f2);
290 Real
c1 = -(1.0-f2)*c3;
291 GradCz = idz0 * (
c1 * cell_prim(i, j, k-1, prim_index)
292 +
c2 * cell_prim(i, j, k , prim_index)
293 + c3 * cell_prim(i, j, k+1, prim_index) );
294 }
else if (ext_dir_on_zhi) {
296 Real dz0 = dz_ptr[k-1];
297 Real dz1 = dz_ptr[k-2];
298 Real idz0 = 1.0 / dz0;
299 Real f = (dz1 / dz0) + 2.0;
301 Real c3 = 2.0 / (f - f2);
303 Real
c1 = -(1.0-f2)*c3;
304 GradCz = idz0 * ( -(
c1 * cell_prim(i, j, k , prim_index)
305 +
c2 * cell_prim(i, j, k-1, prim_index)
306 + c3 * cell_prim(i, j, k-2, prim_index) ) );
308 Real dzk_inv = 2.0 / (dz_ptr[k] + dz_ptr[k-1]);
309 GradCz = dzk_inv * ( cell_prim(i, j, k, prim_index) - cell_prim(i, j, k-1, prim_index) );
313 zflux(i,j,k) = hfx_z(i,j,0);
314 }
else if (SurfLayer_on_zlo && (qty_index ==
RhoQ1_comp)) {
315 zflux(i,j,k) = qfx1_z(i,j,0);
317 zflux(i,j,k) = -rhoAlpha * GradCz;
321 if (!SurfLayer_on_zlo) {
322 hfx_z(i,j,k) = zflux(i,j,k);
325 if (!SurfLayer_on_zlo) {
326 qfx1_z(i,j,k) = zflux(i,j,k);
329 qfx2_z(i,j,k) = zflux(i,j,k);
334 ParallelFor(
xbx, [=] AMREX_GPU_DEVICE (
int i,
int j,
int k) noexcept
336 const int prim_index = qty_index - 1;
338 Real rhoFace = 0.5 * ( cell_data(i, j, k,
Rho_comp) + cell_data(i-1, j, k,
Rho_comp) );
344 bool SurfLayer_on_zlo = ( use_SurfLayer && rotate && k ==
dom_lo.z);
346 Real GradCx =
dx_inv * ( cell_prim(i, j, k , prim_index) - cell_prim(i-1, j, k , prim_index) );
349 xflux(i,j,k) = hfx_x(i,j,0);
350 }
else if (SurfLayer_on_zlo && (qty_index ==
RhoQ1_comp)) {
351 xflux(i,j,k) = qfx1_x(i,j,0);
353 xflux(i,j,k) = -rhoAlpha * mf_ux(i,j,0) * GradCx;
356 ParallelFor(
ybx, [=] AMREX_GPU_DEVICE (
int i,
int j,
int k) noexcept
358 const int prim_index = qty_index - 1;
360 Real rhoFace = 0.5 * ( cell_data(i, j, k,
Rho_comp) + cell_data(i, j-1, k,
Rho_comp) );
366 bool SurfLayer_on_zlo = ( use_SurfLayer && rotate && k ==
dom_lo.z);
368 Real GradCy =
dy_inv * ( cell_prim(i, j, k , prim_index) - cell_prim(i, j-1, k , prim_index) );
371 yflux(i,j,k) = hfx_y(i,j,0);
372 }
else if (SurfLayer_on_zlo && (qty_index ==
RhoQ1_comp)) {
373 yflux(i,j,k) = qfx1_y(i,j,0);
375 yflux(i,j,k) = -rhoAlpha * mf_vy(i,j,0) * GradCy;
378 ParallelFor(
zbx, [=] AMREX_GPU_DEVICE (
int i,
int j,
int k) noexcept
380 const int prim_index = qty_index - 1;
382 Real rhoFace = 0.5 * ( cell_data(i, j, k,
Rho_comp) + cell_data(i, j, k-1,
Rho_comp) );
395 bool SurfLayer_on_zlo = ( use_SurfLayer && k ==
dom_lo.z);
397 if (ext_dir_on_zlo) {
399 Real dz0 = dz_ptr[k];
400 Real dz1 = dz_ptr[k+1];
401 Real idz0 = 1.0 / dz0;
402 Real f = (dz1 / dz0) + 2.0;
404 Real c3 = 2.0 / (f - f2);
406 Real
c1 = -(1.0-f2)*c3;
407 GradCz = idz0 * (
c1 * cell_prim(i, j, k-1, prim_index)
408 +
c2 * cell_prim(i, j, k , prim_index)
409 + c3 * cell_prim(i, j, k+1, prim_index) );
410 }
else if (ext_dir_on_zhi) {
412 Real dz0 = dz_ptr[k-1];
413 Real dz1 = dz_ptr[k-2];
414 Real idz0 = 1.0 / dz0;
415 Real f = (dz1 / dz0) + 2.0;
417 Real c3 = 2.0 / (f - f2);
419 Real
c1 = -(1.0-f2)*c3;
420 GradCz = idz0 * ( -(
c1 * cell_prim(i, j, k , prim_index)
421 +
c2 * cell_prim(i, j, k-1, prim_index)
422 + c3 * cell_prim(i, j, k-2, prim_index) ) );
424 Real dzk_inv = 2.0 / (dz_ptr[k] + dz_ptr[k-1]);
425 GradCz = dzk_inv * ( cell_prim(i, j, k, prim_index) - cell_prim(i, j, k-1, prim_index) );
429 zflux(i,j,k) = hfx_z(i,j,0);
430 }
else if (SurfLayer_on_zlo && (qty_index ==
RhoQ1_comp)) {
431 zflux(i,j,k) = qfx1_z(i,j,0);
433 zflux(i,j,k) = -rhoAlpha * GradCz;
437 if (!SurfLayer_on_zlo) {
438 hfx_z(i,j,k) = zflux(i,j,k);
441 if (!SurfLayer_on_zlo) {
442 qfx1_z(i,j,k) = zflux(i,j,k);
445 qfx2_z(i,j,k) = zflux(i,j,k);
450 ParallelFor(
xbx, [=] AMREX_GPU_DEVICE (
int i,
int j,
int k) noexcept
452 const int prim_index = qty_index - 1;
459 bool SurfLayer_on_zlo = ( use_SurfLayer && rotate && k ==
dom_lo.z);
461 Real GradCx =
dx_inv * ( cell_prim(i, j, k , prim_index) - cell_prim(i-1, j, k , prim_index) );
464 xflux(i,j,k) = hfx_x(i,j,0);
465 }
else if (SurfLayer_on_zlo && (qty_index ==
RhoQ1_comp)) {
466 xflux(i,j,k) = qfx1_x(i,j,0);
468 xflux(i,j,k) = -rhoAlpha * mf_ux(i,j,0) * GradCx;
471 ParallelFor(
ybx, [=] AMREX_GPU_DEVICE (
int i,
int j,
int k) noexcept
473 const int prim_index = qty_index - 1;
480 bool SurfLayer_on_zlo = ( use_SurfLayer && rotate && k ==
dom_lo.z);
482 Real GradCy =
dy_inv * ( cell_prim(i, j, k , prim_index) - cell_prim(i, j-1, k , prim_index) );
485 yflux(i,j,k) = hfx_y(i,j,0);
486 }
else if (SurfLayer_on_zlo && (qty_index ==
RhoQ1_comp)) {
487 yflux(i,j,k) = qfx1_y(i,j,0);
489 yflux(i,j,k) = -rhoAlpha * mf_vy(i,j,0) * GradCy;
492 ParallelFor(
zbx, [=] AMREX_GPU_DEVICE (
int i,
int j,
int k) noexcept
494 const int prim_index = qty_index - 1;
508 bool SurfLayer_on_zlo = ( use_SurfLayer && k ==
dom_lo.z);
510 if (ext_dir_on_zlo) {
512 Real dz0 = dz_ptr[k];
513 Real dz1 = dz_ptr[k+1];
514 Real idz0 = 1.0 / dz0;
515 Real f = (dz1 / dz0) + 2.0;
517 Real c3 = 2.0 / (f - f2);
519 Real
c1 = -(1.0-f2)*c3;
520 GradCz = idz0 * (
c1 * cell_prim(i, j, k-1, prim_index)
521 +
c2 * cell_prim(i, j, k , prim_index)
522 + c3 * cell_prim(i, j, k+1, prim_index) );
523 }
else if (ext_dir_on_zhi) {
525 Real dz0 = dz_ptr[k-1];
526 Real dz1 = dz_ptr[k-2];
527 Real idz0 = 1.0 / dz0;
528 Real f = (dz1 / dz0) + 2.0;
530 Real c3 = 2.0 / (f - f2);
532 Real
c1 = -(1.0-f2)*c3;
533 GradCz = idz0 * ( -(
c1 * cell_prim(i, j, k , prim_index)
534 +
c2 * cell_prim(i, j, k-1, prim_index)
535 + c3 * cell_prim(i, j, k-2, prim_index) ) );
537 Real dzk_inv = 2.0 / (dz_ptr[k] + dz_ptr[k-1]);
538 GradCz = dzk_inv * ( cell_prim(i, j, k, prim_index) - cell_prim(i, j, k-1, prim_index) );
542 zflux(i,j,k) = hfx_z(i,j,0);
543 }
else if (SurfLayer_on_zlo && (qty_index ==
RhoQ1_comp)) {
544 zflux(i,j,k) = qfx1_z(i,j,0);
546 zflux(i,j,k) = -rhoAlpha * GradCz;
550 if (!SurfLayer_on_zlo) {
551 hfx_z(i,j,k) = zflux(i,j,k);
554 if (!SurfLayer_on_zlo) {
555 qfx1_z(i,j,k) = zflux(i,j,k);
558 qfx2_z(i,j,k) = zflux(i,j,k);
564 ParallelFor(
xbx, [=] AMREX_GPU_DEVICE (
int i,
int j,
int k) noexcept
566 xflux(i,j,k) /= mf_uy(i,j,0);
568 ParallelFor(
ybx, [=] AMREX_GPU_DEVICE (
int i,
int j,
int k) noexcept
570 yflux(i,j,k) /= mf_vx(i,j,0);
576 ParallelFor(bx,[=] AMREX_GPU_DEVICE (
int i,
int j,
int k) noexcept
578 Real mfsq = mf_mx(i,j,0) * mf_my(i,j,0);
579 Real dzk_inv = 1.0 / dz_ptr[k];
580 Real stateContrib = (xflux(i+1,j ,k ) - xflux(i, j, k)) *
dx_inv * mfsq
581 +(yflux(i ,j+1,k ) - yflux(i, j, k)) *
dy_inv * mfsq
582 +(zflux(i ,j ,k+1) - zflux(i, j, k)) * dzk_inv;
584 cell_rhs(i,j,k,qty_index) -= stateContrib;
const auto & dom_hi
Definition: ERF_DiffSetup.H:10
bool l_turb
Definition: ERF_DiffSetup.H:19
const auto & dom_lo
Definition: ERF_DiffSetup.H:9
bool l_consA
Definition: ERF_DiffSetup.H:18
int * d_eddy_diff_idx
Definition: ERF_DiffSetup.H:120
const Box zbx
Definition: ERF_DiffSetup.H:23
int * d_eddy_diff_idz
Definition: ERF_DiffSetup.H:122
const Real dx_inv
Definition: ERF_DiffSetup.H:6
const Real dy_inv
Definition: ERF_DiffSetup.H:7
int * d_eddy_diff_idy
Definition: ERF_DiffSetup.H:121
const Box xbx
Definition: ERF_DiffSetup.H:21
Real * d_alpha_eff
Definition: ERF_DiffSetup.H:119
const Box ybx
Definition: ERF_DiffSetup.H:22
void DiffusionSrcForState_S(const Box &bx, const Box &domain, int start_comp, int num_comp, const bool &rotate, const Array4< const Real > &u, const Array4< const Real > &v, const Array4< const Real > &cell_data, const Array4< const Real > &cell_prim, const Array4< Real > &cell_rhs, const Array4< Real > &xflux, const Array4< Real > &yflux, const Array4< Real > &zflux, const Gpu::DeviceVector< Real > &stretched_dz_d, const GpuArray< Real, AMREX_SPACEDIM > &cellSizeInv, const Array4< const Real > &SmnSmn_a, const Array4< const Real > &mf_mx, const Array4< const Real > &mf_ux, const Array4< const Real > &mf_vx, const Array4< const Real > &mf_my, const Array4< const Real > &mf_uy, const Array4< const Real > &mf_vy, Array4< Real > &hfx_x, Array4< Real > &hfx_y, Array4< Real > &hfx_z, Array4< Real > &qfx1_x, Array4< Real > &qfx1_y, Array4< Real > &qfx1_z, Array4< Real > &qfx2_z, Array4< Real > &diss, const Array4< const Real > &mu_turb, const SolverChoice &solverChoice, const int level, const Array4< const Real > &tm_arr, const GpuArray< Real, AMREX_SPACEDIM > grav_gpu, const BCRec *bc_ptr, const bool use_SurfLayer)
Definition: ERF_DiffusionSrcForState_S.cpp:41
#define RhoScalar_comp
Definition: ERF_IndexDefines.H:40
#define Rho_comp
Definition: ERF_IndexDefines.H:36
#define RhoTheta_comp
Definition: ERF_IndexDefines.H:37
#define RhoQ2_comp
Definition: ERF_IndexDefines.H:43
#define NSCALARS
Definition: ERF_IndexDefines.H:16
#define RhoQ1_comp
Definition: ERF_IndexDefines.H:42
#define PrimScalar_comp
Definition: ERF_IndexDefines.H:52
@ RhoScalar_bc_comp
Definition: ERF_IndexDefines.H:80
@ ext_dir
Definition: ERF_IndexDefines.H:209
@ ext_dir_prim
Definition: ERF_IndexDefines.H:211
real(c_double), parameter c2
Definition: ERF_module_model_constants.F90:35
real(c_double), private c1
Definition: ERF_module_mp_morr_two_moment.F90:211