Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Turbulent viscosity corner cells #1281

Merged
merged 2 commits into from
Nov 2, 2023
Merged
Show file tree
Hide file tree
Changes from all commits
Commits
File filter

Filter by extension

Filter by extension

Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
27 changes: 17 additions & 10 deletions Source/BoundaryConditions/MOSTStress.H
Original file line number Diff line number Diff line change
Expand Up @@ -149,7 +149,7 @@ struct adiabatic_charnock
private:
most_data mdata;
similarity_funs sfuns;
amrex::Real tol = 1.0e-5;
const amrex::Real tol = 1.0e-5;
};


Expand Down Expand Up @@ -202,7 +202,7 @@ struct adiabatic_mod_charnock
private:
most_data mdata;
similarity_funs sfuns;
amrex::Real tol = 1.0e-5;
const amrex::Real tol = 1.0e-5;
};


Expand Down Expand Up @@ -260,7 +260,7 @@ struct surface_flux
private:
most_data mdata;
similarity_funs sfuns;
amrex::Real tol = 1.0e-5;
const amrex::Real tol = 1.0e-5;
};


Expand Down Expand Up @@ -323,7 +323,7 @@ struct surface_flux_charnock
private:
most_data mdata;
similarity_funs sfuns;
amrex::Real tol = 1.0e-5;
const amrex::Real tol = 1.0e-5;
};


Expand Down Expand Up @@ -387,7 +387,7 @@ struct surface_flux_mod_charnock
private:
most_data mdata;
similarity_funs sfuns;
amrex::Real tol = 1.0e-5;
const amrex::Real tol = 1.0e-5;
};


Expand Down Expand Up @@ -447,7 +447,7 @@ struct surface_temp
private:
most_data mdata;
similarity_funs sfuns;
amrex::Real tol = 1.0e-5;
const amrex::Real tol = 1.0e-5;
};


Expand Down Expand Up @@ -512,7 +512,7 @@ struct surface_temp_charnock
private:
most_data mdata;
similarity_funs sfuns;
amrex::Real tol = 1.0e-5;
const amrex::Real tol = 1.0e-5;
};


Expand Down Expand Up @@ -578,7 +578,7 @@ struct surface_temp_mod_charnock
private:
most_data mdata;
similarity_funs sfuns;
amrex::Real tol = 1.0e-5;
const amrex::Real tol = 1.0e-5;
};


Expand Down Expand Up @@ -639,6 +639,7 @@ struct moeng_flux
rho = cons_arr(ic,jc,zlo,Rho_comp);
theta = cons_arr(ic,jc,zlo,RhoTheta_comp) / rho;
eta = eta_arr(ie,je,zlo,EddyDiff::Theta_v); // == rho * alpha [kg/m^3 * m^2/s]
eta = amrex::max(eta,eps);

amrex::Real theta_mean = tm_arr(ic,jc,zlo);
amrex::Real wsp_mean = umm_arr(ic,jc,zlo);
Expand All @@ -649,7 +650,7 @@ struct moeng_flux
amrex::Real wsp = sqrt(velx*velx+vely*vely);
amrex::Real num1 = (theta-theta_mean)*wsp_mean;
amrex::Real num2 = (theta_mean-theta_surf)*wsp;
amrex::Real moflux = (std::abs(tstar) > EPS) ?
amrex::Real moflux = (std::abs(tstar) > eps) ?
tstar*ustar*(num1+num2)/((theta_mean-theta_surf)*wsp_mean) : 0.0;
amrex::Real deltaz = dz * (zlo - k);

Expand Down Expand Up @@ -700,6 +701,7 @@ struct moeng_flux
+ cons_arr(ic ,jc,zlo,Rho_comp) );
eta = 0.5 *( eta_arr(ie-1,je,zlo,EddyDiff::Mom_v)
+ eta_arr(ie ,je,zlo,EddyDiff::Mom_v) );
eta = amrex::max(eta,eps);

amrex::Real umean = um_arr(i,j,zlo);
amrex::Real wsp_mean = 0.5 * ( umm_arr(ic-1,jc,zlo) + umm_arr(ic,jc,zlo) );
Expand Down Expand Up @@ -763,6 +765,7 @@ struct moeng_flux
+ cons_arr(ic,jc ,zlo,Rho_comp) );
eta = 0.5*( eta_arr(ie,je-1,zlo,EddyDiff::Mom_v)
+ eta_arr(ie,je ,zlo,EddyDiff::Mom_v) );
eta = amrex::max(eta,eps);

amrex::Real vmean = vm_arr(i,j,zlo);
amrex::Real wsp_mean = 0.5 * ( umm_arr(ic,jc-1,zlo) + umm_arr(ic,jc,zlo) );
Expand All @@ -784,7 +787,7 @@ struct moeng_flux

private:
int zlo;
static constexpr amrex::Real EPS = 1e-16;
const amrex::Real eps = 1e-16;
};


Expand Down Expand Up @@ -833,6 +836,7 @@ struct donelan_flux
rho = cons_arr(ic,jc,zlo,Rho_comp);
theta = cons_arr(ic,jc,zlo,RhoTheta_comp) / rho;
eta = eta_arr(ie,je,zlo,EddyDiff::Theta_v); // == rho * alpha [kg/m^3 * m^2/s]
eta = amrex::max(eta,eps);

amrex::Real Cd = 0.0012;
amrex::Real wsp_mean = umm_arr(ic,jc,zlo);
Expand Down Expand Up @@ -888,6 +892,7 @@ struct donelan_flux
+ cons_arr(ic ,jc,zlo,Rho_comp) );
eta = 0.5 *( eta_arr(ie-1,je,zlo,EddyDiff::Mom_v)
+ eta_arr(ie ,je,zlo,EddyDiff::Mom_v) );
eta = amrex::max(eta,eps);

amrex::Real Cd = 0.001;
const amrex::Real c = 7e-5;
Expand Down Expand Up @@ -955,6 +960,7 @@ struct donelan_flux
+ cons_arr(ic,jc ,zlo,Rho_comp) );
eta = 0.5*( eta_arr(ie,je-1,zlo,EddyDiff::Mom_v)
+ eta_arr(ie,je ,zlo,EddyDiff::Mom_v) );
eta = amrex::max(eta,eps);

amrex::Real Cd = 0.001;
const amrex::Real c = 7e-5;
Expand All @@ -980,5 +986,6 @@ struct donelan_flux

private:
int zlo;
const amrex::Real eps = 1e-16;
};
#endif
20 changes: 11 additions & 9 deletions Source/Diffusion/ComputeTurbulentViscosity.cpp
Original file line number Diff line number Diff line change
Expand Up @@ -193,7 +193,7 @@ void ComputeTurbulentViscosityLES (const amrex::MultiFab& Tau11, const amrex::Mu
for ( amrex::MFIter mfi(eddyViscosity,amrex::TilingIfNotGPU()); mfi.isValid(); ++mfi)
{
Box bxcc = mfi.tilebox();
Box planex = bxcc; planex.setSmall(0, 1); planex.setBig(0, ngc);
Box planex = bxcc; planex.setSmall(0, 1); planex.setBig(0, ngc); planex.grow(1,1);
Box planey = bxcc; planey.setSmall(1, 1); planey.setBig(1, ngc);
int i_lo = bxcc.smallEnd(0); int i_hi = bxcc.bigEnd(0);
int j_lo = bxcc.smallEnd(1); int j_hi = bxcc.bigEnd(1);
Expand All @@ -202,19 +202,21 @@ void ComputeTurbulentViscosityLES (const amrex::MultiFab& Tau11, const amrex::Mu

const Array4<Real>& mu_turb = eddyViscosity.array(mfi);

// Extrapolate outside the domain in lateral directions
// Extrapolate outside the domain in lateral directions (planex owns corner cells)
if (i_lo == domain.smallEnd(0)) {
amrex::ParallelFor(planex, [=] AMREX_GPU_DEVICE (int i, int j, int k) noexcept
{
mu_turb(i_lo-i, j, k, EddyDiff::Mom_h) = mu_turb(i_lo, j, k, EddyDiff::Mom_h);
mu_turb(i_lo-i, j, k, EddyDiff::Mom_v) = mu_turb(i_lo, j, k, EddyDiff::Mom_v);
int lj = amrex::min(amrex::max(j, domain.smallEnd(1)), domain.bigEnd(1));
mu_turb(i_lo-i, j, k, EddyDiff::Mom_h) = mu_turb(i_lo, lj, k, EddyDiff::Mom_h);
mu_turb(i_lo-i, j, k, EddyDiff::Mom_v) = mu_turb(i_lo, lj, k, EddyDiff::Mom_v);
});
}
if (i_hi == domain.bigEnd(0)) {
amrex::ParallelFor(planex, [=] AMREX_GPU_DEVICE (int i, int j, int k) noexcept
{
mu_turb(i_hi+i, j, k, EddyDiff::Mom_h) = mu_turb(i_hi, j, k, EddyDiff::Mom_h);
mu_turb(i_hi+i, j, k, EddyDiff::Mom_v) = mu_turb(i_hi, j, k, EddyDiff::Mom_v);
int lj = amrex::min(amrex::max(j, domain.smallEnd(1)), domain.bigEnd(1));
mu_turb(i_hi+i, j, k, EddyDiff::Mom_h) = mu_turb(i_hi, lj, k, EddyDiff::Mom_h);
mu_turb(i_hi+i, j, k, EddyDiff::Mom_v) = mu_turb(i_hi, lj, k, EddyDiff::Mom_v);
});
}
if (j_lo == domain.smallEnd(1)) {
Expand Down Expand Up @@ -336,9 +338,9 @@ void ComputeTurbulentViscosityLES (const amrex::MultiFab& Tau11, const amrex::Mu
mu_turb(i, j, k_hi+k, indx_v) = mu_turb(i, j, k_hi, indx_v);
});
break;
}
}
}
}
}
}
}

/**
Expand Down