ERF
Energy Research and Forecasting: An Atmospheric Modeling Code
Classes | Functions
ERF_EddyViscosity.H File Reference
#include "AMReX_BCRec.H"
#include "ERF_SurfaceLayer.H"
#include "ERF_DataStruct.H"
#include "ERF_IndexDefines.H"
#include "ERF_NumericalConstants.H"
#include "ERF_EOS.H"
#include "ERF_EB.H"
Include dependency graph for ERF_EddyViscosity.H:
This graph shows which files directly or indirectly include this file:

Go to the source code of this file.

Classes

struct  StrainArrays
 

Functions

void ComputeTurbulentViscosity (double dt, const amrex::MultiFab &xvel, const amrex::MultiFab &yvel, amrex::Vector< std::unique_ptr< amrex::MultiFab >> &Tau_lev, amrex::MultiFab &cons_in, const amrex::MultiFab &wdist, amrex::MultiFab &eddyViscosity, amrex::MultiFab &Hfx1, amrex::MultiFab &Hfx2, amrex::MultiFab &Hfx3, amrex::MultiFab &Diss, const amrex::Geometry &geom, amrex::Vector< std::unique_ptr< amrex::MultiFab >> &mapfac, const std::unique_ptr< amrex::MultiFab > &z_phys_nd, const std::unique_ptr< amrex::MultiFab > &z_phys_cc, const SolverChoice &solverChoice, std::unique_ptr< SurfaceLayer > &SurfLayer, const amrex::MultiFab *z_0, const bool &use_terrain_fitted_coords, const bool &use_moisture, int level, const amrex::BCRec *bc_ptr, const eb_ &ebfact, bool vert_only=false, const amrex::MultiFab *qheating_rates=nullptr, const amrex::MultiFab *terrain_blank=nullptr, const amrex::MultiFab *forest_drag=nullptr)
 
void ComputeTurbulentViscosityLES (amrex::Vector< std::unique_ptr< amrex::MultiFab >> &Tau_lev, const amrex::MultiFab &cons_in, amrex::MultiFab &eddyViscosity, amrex::MultiFab &Diss, const amrex::Geometry &geom, bool use_terrain_fitted_coords, amrex::Vector< std::unique_ptr< amrex::MultiFab >> &mapfac, const std::unique_ptr< amrex::MultiFab > &z_phys_nd, const TurbChoice &turbChoice, const amrex::Real const_grav, std::unique_ptr< SurfaceLayer > &SurfLayer, const MoistureComponentIndices &moisture_indices, const amrex::MultiFab *xvel=nullptr, const amrex::MultiFab *yvel=nullptr)
 
void ComputeTurbulentViscosityRANS (int level, amrex::Vector< std::unique_ptr< amrex::MultiFab >> &Tau_lev, const amrex::MultiFab &cons_in, const amrex::MultiFab &wdist, amrex::MultiFab &eddyViscosity, amrex::MultiFab &Diss, const amrex::Geometry &geom, bool use_terrain_fitted_coords, const std::unique_ptr< amrex::MultiFab > &z_phys_nd, const TurbChoice &turbChoice, const amrex::Real const_grav, std::unique_ptr< SurfaceLayer > &SurfLayer, const amrex::MultiFab *z_0, const CanopyTurbulence &canopy, const amrex::MultiFab *forest_drag, const amrex::MultiFab *ib_blank=nullptr, bool flat_wall_diss=false, amrex::Real ib_small_volfrac=amrex::Real(0.005), amrex::Real ib_z0=amrex::Real(0.0))
 
void ComputeTurbulentViscosityLES_EB (amrex::Vector< std::unique_ptr< amrex::MultiFab >> &Tau_lev, const amrex::MultiFab &cons_in, amrex::MultiFab &eddyViscosity, amrex::MultiFab &Hfx1, amrex::MultiFab &Hfx2, amrex::MultiFab &Hfx3, const amrex::Geometry &geom, amrex::Vector< std::unique_ptr< amrex::MultiFab >> &mapfac, const TurbChoice &turbChoice, const amrex::Real const_grav, const SolverChoice &solverChoice, std::unique_ptr< SurfaceLayer > &SurfLayer, const MoistureComponentIndices &moisture_indices, const eb_ &ebfact, const amrex::MultiFab *xvel=nullptr, const amrex::MultiFab *yvel=nullptr)
 
AMREX_GPU_DEVICE AMREX_FORCE_INLINE amrex::Real ComputeSmnSmn (int &i, int &j, int &k, const amrex::Array4< amrex::Real const > &tau11, const amrex::Array4< amrex::Real const > &tau22, const amrex::Array4< amrex::Real const > &tau33, const amrex::Array4< amrex::Real const > &tau12, const amrex::Array4< amrex::Real const > &tau13, const amrex::Array4< amrex::Real const > &tau23)
 
StrainArrays strain_arrays (amrex::Vector< std::unique_ptr< amrex::MultiFab >> &Tau_lev, const amrex::MFIter &mfi)
 
AMREX_GPU_DEVICE AMREX_FORCE_INLINE amrex::Real ComputeSmnSmn_EB (int &i, int &j, int &k, const amrex::Array4< amrex::Real const > &tau11, const amrex::Array4< amrex::Real const > &tau22, const amrex::Array4< amrex::Real const > &tau33, const amrex::Array4< amrex::Real const > &tau12, const amrex::Array4< amrex::Real const > &tau13, const amrex::Array4< amrex::Real const > &tau23, const amrex::Array4< const amrex::EBCellFlag > &c_cflag, const amrex::Array4< const amrex::EBCellFlag > &u_cflag, const amrex::Array4< const amrex::EBCellFlag > &v_cflag, const amrex::Array4< const amrex::EBCellFlag > &w_cflag)
 
AMREX_GPU_DEVICE AMREX_FORCE_INLINE amrex::Real ComputeSmnSmn2D (int &i, int &j, int &k, const amrex::Array4< amrex::Real const > &tau11, const amrex::Array4< amrex::Real const > &tau22, const amrex::Array4< amrex::Real const > &tau12)
 

Function Documentation

◆ ComputeSmnSmn()

AMREX_GPU_DEVICE AMREX_FORCE_INLINE amrex::Real ComputeSmnSmn ( int &  i,
int &  j,
int &  k,
const amrex::Array4< amrex::Real const > &  tau11,
const amrex::Array4< amrex::Real const > &  tau22,
const amrex::Array4< amrex::Real const > &  tau33,
const amrex::Array4< amrex::Real const > &  tau12,
const amrex::Array4< amrex::Real const > &  tau13,
const amrex::Array4< amrex::Real const > &  tau23 
)

Compute twice-contracted strain-rate magnitude from stress-like strain components.

Parameters
[in]icell-centered i-index
[in]jcell-centered j-index
[in]kcell-centered k-index
[in]tau1111 strain component
[in]tau2222 strain component
[in]tau3333 strain component
[in]tau1212 strain component
[in]tau1313 strain component
[in]tau2323 strain component
Returns
S_mn S_mn at cell centers
140 {
141  amrex::Real s11bar = tau11(i,j,k);
142  amrex::Real s22bar = tau22(i,j,k);
143  amrex::Real s33bar = tau33(i,j,k);
144  amrex::Real s12bar = fourth * ( tau12(i , j , k ) + tau12(i , j+1, k )
145  + tau12(i+1, j , k ) + tau12(i+1, j+1, k ) );
146  amrex::Real s13bar = fourth * ( tau13(i , j , k ) + tau13(i , j , k+1)
147  + tau13(i+1, j , k ) + tau13(i+1, j , k+1) );
148  amrex::Real s23bar = fourth * ( tau23(i , j , k ) + tau23(i , j , k+1)
149  + tau23(i , j+1, k ) + tau23(i , j+1, k+1) );
150 
151  amrex::Real SmnSmn = s11bar*s11bar + s22bar*s22bar + s33bar*s33bar
152  + two*s12bar*s12bar + two*s13bar*s13bar + two*s23bar*s23bar;
153 
154  return SmnSmn;
155 }
tau12
@ tau12
Definition: ERF_DataStruct.H:49
tau23
@ tau23
Definition: ERF_DataStruct.H:49
tau33
@ tau33
Definition: ERF_DataStruct.H:49
tau22
@ tau22
Definition: ERF_DataStruct.H:49
tau11
@ tau11
Definition: ERF_DataStruct.H:49
tau13
@ tau13
Definition: ERF_DataStruct.H:49
two
constexpr amrex::Real two
Definition: ERF_NumericalConstants.H:38
fourth
constexpr amrex::Real fourth
Definition: ERF_NumericalConstants.H:42
Real
amrex::Real Real
Definition: ERF_ShocInterface.H:19

Referenced by ComputeTurbulentViscosityLES(), ComputeTurbulentViscosityLES_EB(), ComputeTurbulentViscosityRANS(), and erf_make_tau_terms().

Here is the caller graph for this function:

◆ ComputeSmnSmn2D()

AMREX_GPU_DEVICE AMREX_FORCE_INLINE amrex::Real ComputeSmnSmn2D ( int &  i,
int &  j,
int &  k,
const amrex::Array4< amrex::Real const > &  tau11,
const amrex::Array4< amrex::Real const > &  tau22,
const amrex::Array4< amrex::Real const > &  tau12 
)

Compute horizontal two-dimensional strain-rate magnitude.

Parameters
[in]icell-centered i-index
[in]jcell-centered j-index
[in]kcell-centered k-index
[in]tau1111 strain component
[in]tau2222 strain component
[in]tau1212 strain component
Returns
horizontal S_mn S_mn contribution
301 {
302  amrex::Real sdiff = tau11(i,j,k) - tau22(i,j,k);
303  amrex::Real s12bar = fourth * ( tau12(i , j , k ) + tau12(i , j+1, k )
304  + tau12(i+1, j , k ) + tau12(i+1, j+1, k ) );
305  return myhalf * sdiff*sdiff + two*s12bar*s12bar;
306 }
myhalf
constexpr amrex::Real myhalf
Definition: ERF_NumericalConstants.H:41

Referenced by ComputeTurbulentViscosityLES().

Here is the caller graph for this function:

◆ ComputeSmnSmn_EB()

AMREX_GPU_DEVICE AMREX_FORCE_INLINE amrex::Real ComputeSmnSmn_EB ( int &  i,
int &  j,
int &  k,
const amrex::Array4< amrex::Real const > &  tau11,
const amrex::Array4< amrex::Real const > &  tau22,
const amrex::Array4< amrex::Real const > &  tau33,
const amrex::Array4< amrex::Real const > &  tau12,
const amrex::Array4< amrex::Real const > &  tau13,
const amrex::Array4< amrex::Real const > &  tau23,
const amrex::Array4< const amrex::EBCellFlag > &  c_cflag,
const amrex::Array4< const amrex::EBCellFlag > &  u_cflag,
const amrex::Array4< const amrex::EBCellFlag > &  v_cflag,
const amrex::Array4< const amrex::EBCellFlag > &  w_cflag 
)

Compute twice-contracted strain-rate magnitude for embedded-boundary cells.

Parameters
[in]icell-centered i-index
[in]jcell-centered j-index
[in]kcell-centered k-index
[in]tau1111 strain component
[in]tau2222 strain component
[in]tau3333 strain component
[in]tau1212 strain component
[in]tau1313 strain component
[in]tau2323 strain component
[in]c_cflagcell-centered EB flags
[in]u_cflagx-face EB flags
[in]v_cflagy-face EB flags
[in]w_cflagz-face EB flags
Returns
S_mn S_mn at uncovered cell centers
203 {
204  amrex::Real SmnSmn = zero;
205  if (!c_cflag(i,j,k).isCovered()) {
206  amrex::Real s11bar = tau11(i,j,k);
207  amrex::Real s22bar = tau22(i,j,k);
208  amrex::Real s33bar = tau33(i,j,k);
209  amrex::Real s12bar = zero;
210  amrex::Real s13bar = zero;
211  amrex::Real s23bar = zero;
212 
213  amrex::Real count_s12 = zero;
214  if (!u_cflag(i,j,k).isCovered() || !v_cflag(i,j,k).isCovered()) {
215  s12bar += tau12(i,j,k);
216  count_s12 += one;
217  }
218  if (!u_cflag(i+1,j,k).isCovered() || !v_cflag(i,j,k).isCovered()) {
219  s12bar += tau12(i+1,j,k);
220  count_s12 += one;
221  }
222  if (!u_cflag(i,j,k).isCovered() || !v_cflag(i,j+1,k).isCovered()) {
223  s12bar += tau12(i,j+1,k);
224  count_s12 += one;
225  }
226  if (!u_cflag(i+1,j,k).isCovered() || !v_cflag(i,j+1,k).isCovered()) {
227  s12bar += tau12(i+1,j+1,k);
228  count_s12 += one;
229  }
230  if (count_s12 > zero) {
231  s12bar /= count_s12;
232  }
233 
234  amrex::Real count_s13 = zero;
235  if (!u_cflag(i,j,k).isCovered() || !w_cflag(i,j,k).isCovered()) {
236  s13bar += tau13(i, j, k);
237  count_s13 += one;
238  }
239  if (!u_cflag(i,j,k).isCovered() || !w_cflag(i,j,k+1).isCovered()) {
240  s13bar += tau13(i, j, k+1);
241  count_s13 += one;
242  }
243  if (!u_cflag(i+1,j,k).isCovered() || !w_cflag(i,j,k).isCovered()) {
244  s13bar += tau13(i+1, j, k);
245  count_s13 += one;
246  }
247  if (!u_cflag(i+1,j,k).isCovered() || !w_cflag(i,j,k+1).isCovered()) {
248  s13bar += tau13(i+1, j, k+1);
249  count_s13 += one;
250  }
251  if (count_s13 > zero) {
252  s13bar /= count_s13;
253  }
254 
255  amrex::Real count_s23 = zero;
256  if (!v_cflag(i,j,k).isCovered() || !w_cflag(i,j,k).isCovered()) {
257  s23bar += tau23(i, j, k);
258  count_s23 += one;
259  }
260  if (!v_cflag(i,j,k).isCovered() || !w_cflag(i,j,k+1).isCovered()) {
261  s23bar += tau23(i, j, k+1);
262  count_s23 += one;
263  }
264  if (!v_cflag(i,j+1,k).isCovered() || !w_cflag(i,j,k).isCovered()) {
265  s23bar += tau23(i, j+1, k);
266  count_s23 += one;
267  }
268  if (!v_cflag(i,j+1,k).isCovered() || !w_cflag(i,j,k+1).isCovered()) {
269  s23bar += tau23(i, j+1, k+1);
270  count_s23 += one;
271  }
272  if (count_s23 > zero) {
273  s23bar /= count_s23;
274  }
275 
276  SmnSmn = s11bar*s11bar + s22bar*s22bar + s33bar*s33bar
277  + two*s12bar*s12bar + two*s13bar*s13bar + two*s23bar*s23bar;
278  }
279 
280  return SmnSmn;
281 }
one
constexpr amrex::Real one
Definition: ERF_NumericalConstants.H:37
zero
constexpr amrex::Real zero
Definition: ERF_NumericalConstants.H:36

Referenced by ComputeTurbulentViscosityLES_EB().

Here is the caller graph for this function:

◆ ComputeTurbulentViscosity()

void ComputeTurbulentViscosity ( double  dt,
const amrex::MultiFab &  xvel,
const amrex::MultiFab &  yvel,
amrex::Vector< std::unique_ptr< amrex::MultiFab >> &  Tau_lev,
amrex::MultiFab &  cons_in,
const amrex::MultiFab &  wdist,
amrex::MultiFab &  eddyViscosity,
amrex::MultiFab &  Hfx1,
amrex::MultiFab &  Hfx2,
amrex::MultiFab &  Hfx3,
amrex::MultiFab &  Diss,
const amrex::Geometry &  geom,
amrex::Vector< std::unique_ptr< amrex::MultiFab >> &  mapfac,
const std::unique_ptr< amrex::MultiFab > &  z_phys_nd,
const std::unique_ptr< amrex::MultiFab > &  z_phys_cc,
const SolverChoice &  solverChoice,
std::unique_ptr< SurfaceLayer > &  SurfLayer,
const amrex::MultiFab *  z_0,
const bool &  use_terrain_fitted_coords,
const bool &  use_moisture,
int  level,
const amrex::BCRec *  bc_ptr,
const eb_ &  ebfact,
bool  vert_only = false,
const amrex::MultiFab *  qheating_rates = nullptr,
const amrex::MultiFab *  terrain_blank = nullptr,
const amrex::MultiFab *  forest_drag = nullptr 
)

Compute turbulent viscosity and related flux fields for the configured turbulence model.

◆ ComputeTurbulentViscosityLES()

void ComputeTurbulentViscosityLES ( amrex::Vector< std::unique_ptr< amrex::MultiFab >> &  Tau_lev,
const amrex::MultiFab &  cons_in,
amrex::MultiFab &  eddyViscosity,
amrex::MultiFab &  Diss,
const amrex::Geometry &  geom,
bool  use_terrain_fitted_coords,
amrex::Vector< std::unique_ptr< amrex::MultiFab >> &  mapfac,
const std::unique_ptr< amrex::MultiFab > &  z_phys_nd,
const TurbChoice &  turbChoice,
const amrex::Real  const_grav,
std::unique_ptr< SurfaceLayer > &  SurfLayer,
const MoistureComponentIndices &  moisture_indices,
const amrex::MultiFab *  xvel = nullptr,
const amrex::MultiFab *  yvel = nullptr 
)

Compute LES turbulent viscosity and dissipation.

No subgrid heat flux is stored: Hfx3 is z-nodal and the theta diffusion of each RK stage writes its face fluxes there (ERF_ComputeTurbulentViscosity.cpp, ERF_TurbKESources.H).

◆ ComputeTurbulentViscosityLES_EB()

void ComputeTurbulentViscosityLES_EB ( amrex::Vector< std::unique_ptr< amrex::MultiFab >> &  Tau_lev,
const amrex::MultiFab &  cons_in,
amrex::MultiFab &  eddyViscosity,
amrex::MultiFab &  Hfx1,
amrex::MultiFab &  Hfx2,
amrex::MultiFab &  Hfx3,
const amrex::Geometry &  geom,
amrex::Vector< std::unique_ptr< amrex::MultiFab >> &  mapfac,
const TurbChoice &  turbChoice,
const amrex::Real  const_grav,
const SolverChoice &  solverChoice,
std::unique_ptr< SurfaceLayer > &  SurfLayer,
const MoistureComponentIndices &  moisture_indices,
const eb_ &  ebfact,
const amrex::MultiFab *  xvel = nullptr,
const amrex::MultiFab *  yvel = nullptr 
)

Compute LES turbulent viscosity and heat fluxes for embedded-boundary cells.

◆ ComputeTurbulentViscosityRANS()

void ComputeTurbulentViscosityRANS ( int  level,
amrex::Vector< std::unique_ptr< amrex::MultiFab >> &  Tau_lev,
const amrex::MultiFab &  cons_in,
const amrex::MultiFab &  wdist,
amrex::MultiFab &  eddyViscosity,
amrex::MultiFab &  Diss,
const amrex::Geometry &  geom,
bool  use_terrain_fitted_coords,
const std::unique_ptr< amrex::MultiFab > &  z_phys_nd,
const TurbChoice &  turbChoice,
const amrex::Real  const_grav,
std::unique_ptr< SurfaceLayer > &  SurfLayer,
const amrex::MultiFab *  z_0,
const CanopyTurbulence &  canopy,
const amrex::MultiFab *  forest_drag,
const amrex::MultiFab *  ib_blank = nullptr,
bool  flat_wall_diss = false,
amrex::Real  ib_small_volfrac = amrex::Real(0.005),
amrex::Real  ib_z0 = amrex::Real(0.0) 
)

Compute the RANS eddy viscosity, length scale and dissipation (k-eqn, Axell & Liungman 2001). ib_blank is the solid fraction of the fraction-stress immersed wall law (else nullptr), and flat_wall_diss (erf.if_flat_wall_dissipation) makes its partial wall cells dissipate with the length of flat ground's first cell. Declared here for the unit tests; ComputeTurbulentViscosity is the caller in the code.

◆ strain_arrays()

StrainArrays strain_arrays ( amrex::Vector< std::unique_ptr< amrex::MultiFab >> &  Tau_lev,
const amrex::MFIter &  mfi 
)
inline

Strain arrays of the box of mfi from the Tau_lev container (tau11..tau23 of TauType).

165 {
166  return { Tau_lev[TauType::tau11]->const_array(mfi), Tau_lev[TauType::tau22]->const_array(mfi),
167  Tau_lev[TauType::tau33]->const_array(mfi), Tau_lev[TauType::tau12]->const_array(mfi),
168  Tau_lev[TauType::tau13]->const_array(mfi), Tau_lev[TauType::tau23]->const_array(mfi) };
169 }

Referenced by ComputeTurbulentViscosityLES(), ComputeTurbulentViscosityLES_EB(), and ComputeTurbulentViscosityRANS().

Here is the caller graph for this function: