Creating the MAC velocities
To create the normal velocities on faces, we first extrapolate from the cell centers on each side using the slopes as computed earlier, and upwind the face value to define \(U^{pred}\)
. To compute the x-velocity on the x-faces of regular (ie not cut) cells, we call
AMREX_CUDA_HOST_DEVICE_FOR_3D(ubx, i, j, k, { // X-faces Real upls = ccvel_fab(i ,j,k,0) - 0.5 * xslopes_fab(i ,j,k,0); Real umns = ccvel_fab(i-1,j,k,0) + 0.5 * xslopes_fab(i-1,j,k,0); if ( umns < 0.0 && upls > 0.0 ) { umac_fab(i,j,k) = 0.0; } else { Real avg = 0.5 * ( upls + umns ); if ( std::abs(avg) < small_vel) { umac_fab(i,j,k) = 0.0; } else if (avg >= 0) { umac_fab(i,j,k) = umns; } else { umac_fab(i,j,k) = upls; } } });
For cut cells we test on whether the area fraction is non-zero:
AMREX_CUDA_HOST_DEVICE_FOR_3D(ubx, i, j, k, { // X-faces if (ax_fab(i,j,k) > 0.0) { Real upls = ccvel_fab(i ,j,k,0) - 0.5 * xslopes_fab(i ,j,k,0); Real umns = ccvel_fab(i-1,j,k,0) + 0.5 * xslopes_fab(i-1,j,k,0); if ( umns < 0.0 && upls > 0.0 ) { umac_fab(i,j,k) = 0.0; } else { Real avg = 0.5 * ( upls + umns ); if ( std::abs(avg) < small_vel) { umac_fab(i,j,k) = 0.0; } else if (avg >= 0) { umac_fab(i,j,k) = umns; } else { umac_fab(i,j,k) = upls; } } } else { umac_fab(i,j,k) = huge_vel; } });
We then perform a MAC projection on the face-centered velocities to enforce that they satisfy
\[\nabla \cdot (\varepsilon_g U^{MAC}) = 0\]
We do this by solving
\[\nabla \cdot \frac{\varepsilon_g}{\rho_g} \nabla \phi^{MAC} = \nabla \cdot \left( \varepsilon_g U^{pred} \right)\]
then defining
\[U^{MAC} = U^{pred} - \frac{1}{\rho_g} \nabla \phi^{MAC}\]