IndexRundeck: P1SoM40Created: Thu May 9 03:30:17 EDT 2024

SOCPBL

File: PBL.f
Summary: module SOCPBL defines subroutines and variables associated with the boundary layer physics. It sets up npbl(=8) sublayers between the surface (sublayer 1) and the middle of the first GCM layer (sublayer npbl), and integrates, over these sublayers, the dynamic equations for the mean turbulent variables using turbulence models, to find the surface values of these variables and related fluxes. t_pbl_args is a derived type structure which contains all input/output arguments for PBL. SOCPBL contains the following subroutines: advanc,stars,getl,dflux,simil,griddr,tfix ccoeff0,getk,e_eqn,t_eqn,q_eqn,uv_eqn, t_eqn_sta,q_eqn_sta,uv_eqn_sta, inits,tcheck,ucheck,check1,output,rtsafe.
Author : Ye Cheng/G. Hartke (modifications by G. Schmidt)
Version:

Subroutines:
advanc
time steps the solutions for the boundary layer variables. It is called from within the subroutine pbl (in PBL_DRV.f). All its outputs are contained in the structure pbl_args (an instance of t_pbl_args).
alloc_pbl_args

ccoeff0
ccoeff0 sets/calculates model coefficients for the GISS 2002 turbulence model (Cheng et al., 2002).
check1
check1 checks for NaN'S and INF'S in real 1-D arrays.
dealloc_pbl_args

dflux
dflux computes the dimensionless surface fluxes of momentum, heat and moisture (drag coefficient Cm , Stanton number Ch, and Dalton number Cq), with explicit Schmidt number (Sc) and Prandtl number (Pr) dependence and flexibility for water isotopes. It also computes the roughness lengths for momentum, z0m (for itype=1 or 2, i.e., surface type ocean or seaice), for temperature, z0h, and for water vapor, z0q, all in meters (Hartke and Rind, 1997). It is called from within subroutine stars.
e_eqn
e_eqn integrates differential eqns for the turbulent kinetic energy, e, using tridiagonal method over npbl-1(=7) sublayer edges (i.e., the secondary grids). The boundary condition near the bottom is: e(1)=(1/2)*B1**(2/3)*ustar**2. At the top secondary grid, nearest to the middle of the first GCM layer, e is prescribed.
e_les
e_les finds e according to the parameterization of les data
find_dpsih

find_dpsim

find_phih

find_phim

find_phim0

get_tv
get_tv converts temperature T to virtual temperature Tv
getchq
calculates the Stanton number or Dalton number (ch or cq) for heat and latent heat fluxes (Hartke and Rind, 1997; Zeng et al., 1998; updated by Y. Cheng, 2014). It is called from within subroutine dflux.
getcm
calculates the drag coefficient (cm) for momentum flux (Hartke and Rind, 1997; Zeng et al., 1998; updated by Y. Cheng, 2014). It is called from within subroutine dflux.
getk
getk computes the turbulent diffusivities for momentum, Km, for heat, Kh, for moisture, Kq and for kinetic energy, Ke, at the secondary grids, using the GISS second order closure model (Cheng et al., 2000). u,v,t,q,ke are calculated at the primary grid z, while e,lscale,km,kh,gm,gh are calculated at the secondary grid zhat.
getl
getl computes the master length scale, lscale, of the turbulence model on the secondary grid, zhat, using the formulas by Nakanishi (2001) from the LES data.
getl1
getl1 estimates the master length scale, lscale, of the turbulence model on the secondary grid, zhat.
getzhq
getzhq calculates the roughness lengths for heat (z0h) and for humidity (z0q), modified from Eqs 5.24, 5.27 and 5.35 in Brutsaert (1982). It is called from within subroutine dflux.
griddr
griddr computes altitudes on vertical grid. The xi coordinates are uniformly spaced and are mapped in a log-linear fashion onto the z grid. (The z's are the physical coords.) Also computes the altitudes on the secondary grid, zhat, and the derivatives dxi/dz evaluated at both all z and zhat. z and zhat are staggered: mean quantitied are calculated at z, turbulent kinetic enery and fluxes are calculated at zhat.
inits
inits initializes the winds, virtual potential temperature, and humidity by solving their differential equations for the static solutions, using tridiagonal method over npbl(=8) sublayers between the surface (sublayer 1) and the middle of the first GCM layer (sublayer npbl). (Cheng et a., 2002). It is called by subroutine init_pbl (in PBL_DRV.f), and the latter (init_pbl) is called by subroutine INPUT (in MODELE.f).
level2
level2 computes the turbulent kinetic energy e (Cheng et al., 2002).
output
output produces output for diagnostic purposes.
q_eqn
q_eqn integrates differential eqns for the specific humidity, Q, using tridiagonal method over npbl(=8) sublayers between the surface (sublayer 1) and the middle of the first GCM layer (sublayer npbl). The boundary condition at the bottom is: kq * dq/dz = min ( cq * usurf * (q1 - qgrnd) + cq * gusti * qprime , fr_sat * ( cq * usurf * (q1 - qgrnd) + cq * gusti * qprime ) - ( 1 - fr_sat ) * flux_max ), which includes the effects on the surface flux due to the moist convection wind gustiness and the downdraft specific humidity perturbation (Redelsperger et al. 2000; Emanuel and Zivkovic 1999), where qprime=qdns-qtop, qtop is Q at the first layer and qdns is the downdraft humidity in kg/kg, (i,j), which is calculated in subroutines CONDSE (in CLOUDS2_DRV.f) and PBL (in PBL_DRV.f). At the top, i.e., the middle of the first GCM layer, Q is prescribed.
q_eqn_sta
q_eqn_sta computes the static solutions of the specific humidity, Q, between the surface and the first GCM layer. The boundary condition at the bottom is: kq * dq/dz = cq * usurf * (q - qg). At the top, Q is prescribed. It is called only at the initialization (from within subroutine inits).
simil
simil calculates the similarity solutions for wind speed, virtual potential temperature, and moisture mixing ratio at height z.
stars
computes the friction speed, ustar, the virtual potential temperature scale, tstar, and the specific humidity scale, qstar. Note that surface momentum flux = ustar*ustar surface heat flux = ustar*tstar surface moisture flux = ustar*qstar It also calculates and outputs the Monin-Obukov length, lmonin, the roughness lengths (z0m,z0h,z0q), the drag coefficient (cm), the Stanton number (ch) and the Dalton number (cq) by calling subroutine dflux.
t_eqn
t_eqn integrates differential eqns for the virtual potential temperature, T, using tridiagonal method over npbl(=8) sublayers between the surface (sublayer 1) and the middle of the first GCM layer (sublayer npbl). The boundary condition at the bottom is: kh * dt/dz = ch * ( usurf*(t1 - tgrnd) +(1+xdelt*qtop)*gusti*tprime ) which includes the effects on the surface flux due to the moist convection wind gustiness and the downdraft temperature perturbation (Redelsperger et al. 2000; Emanuel and Zivkovic 1999), where tprime=tdns-ttop/(1+xdelt*qtop), t1, q1 are the T and Q at the surface, and tdns is the downdraft temperature in K at (i,j), which is calculated in subroutines CONDSE (in CLOUDS2_DRV.f) and PBL (in PBL_DRV.f). At the top, i.e., the middle of the first GCM layer, T is prescribed. ### the following term was removed from BC at the bottom ### + xdelt * t1/(1+xdelt*q1) * kq * dqdz
t_eqn_sta
t_eqn_sta computes the static solutions of the virtual potential temperature, T, between the surface and the first GCM layer. The boundary condition at the bottom is: kh * dt/dz = ch * usurf * (t - tg). At the top, T is prescribed. It is called only at the initialization (from within subroutine inits).
tcheck
tcheck checks for reasonable temperatures
tfix
tfix linearly interpolates between the ground temperature tgrnd and the virtual potential temperature at the middle of the first GCM layer to reset the T(z) profile. It is called when the T(z) profile becomes irregular.
tr_eqn
tr_eqn integrates differential eqns for the tracers, TR, using tridiagonal method over npbl(=8) sublayers between the surface (sublayer 1) and the middle of the first GCM layer (sublayer npbl). The boundary condition at the bottom is: kq * dtr/dz = sfac * trs - constflx, i.e. for moisture, sfac=cq*usurf, constflx=cq*usurf*qg, to get: kq * dq/dz = cq * usurf * (qs - qg); for new moisture (including downdraft effects), sfac=cq*(usurf-dusurf), constflx=cq*(usurf*qg + dusurf*qdns), or sfac=cq*usurf0, constflx=cq*(usurf*(qg+qdns)-usurf0*qdns), to get: kq * dq/dz = cq*(usurf*(qs-qg) + dusurf*(qdns-qs)). This should be flexible enough to deal with most situations. At the top, i.e., the middle of the first GCM layer, TR is prescribed.
ucheck
ucheck makes sure that the winds remain within reasonable bounds during the initialization process. (Sometimes the computed wind speed iterated out in left field someplace, *way* outside any reasonable range.) Tests and corrects both direction and magnitude of the wind rotation with altitude. Tests the total wind speed via comparison to similarity theory. Note that it works from the top down so that it can assume that at level (i), level (i+1) displays reasonable behavior.
uv_eqn
uv_eqn integrates differential eqns for mean velocity u and v using tridiagonal method over npbl(=8) sublayers between the surface (sublayer 1) and the middle of the first GCM layer (sublayer npbl). The boundary condition at the bottom is: km * du/dz = cm * usurf * u and km * dv/dz = cm * usurf * v. At the top, i.e., the middle of the first GCM layer, u, v are prescribed.
uv_eqn_sta
uv_eqn_sta computes the static solutions of the wind components, u and v, between the surface and the first GCM layer. The boundary conditions at the bottom are: km * du/dz = cm * usurf * u, km * dv/dz = cm * usurf * v. At the top, u and v are prescribed. It is called only at the initialization (from within subroutine inits).

Global Variables:
grav : used from constant

pi : used from constant

radian : used from constant

bygrav : used from constant

teeny : used from constant

deltx : used from constant

tf : used from constant

by3 : used from constant

lhe : used from constant

rgas : used from constant

rhows : used from constant

mair : used from constant

byrhows : used from constant

sha : used from constant

shv : used from constant

shw : used from constant

stbo : used from constant

visc_air : used from constant

tfrez : used from seaice

snmin : used from landice

tridiag : used from tridiag_mod

n : integer, parameter
no of pbl. layers
Initial Value = 8
Used by: | :diurn_defs | :init_pbl | :surface | DIAG_COM | PBLCOM | PBL_DRV:pbl | SOIL_DRV |
npbl : integer, parameter
Initial Value = n
maxntm : integer, public

Used by: | :init_pbl |
qg_sat :
Saturation vapor mixing ratio (kg vapor / kg air in a given volume)
wsgcm :
magnitude of the GCM surface wind - ocean currents [m/s]
wspdf :
mean surface wind calculated from PDF of wind speed [m/s]
wsubtke :
turbulent kinetic energy velocity scale [m/s]
wsubwd :
dry convective velocity scale [m/s]
wsubwm :
moist convective velocity scale [m/s]
mcfrac :
fractional area with moist convection in grid cell
xcdpbl : real*8
factor for momentum drag coefficient
Initial Value = 1d0
Used by: | :init_pbl |
rimax : real*8

Used by: | :e_gcm |
ghmin : real*8

Used by: | :k_gcm |
ghmax : real*8

Used by: | :k_gcm |
gmmax0 : real*8

gm_at_rimax : real*8

d1 : real*8

Used by: | :k_gcm |
d2 : real*8

Used by: | :k_gcm |
d3 : real*8

Used by: | :k_gcm |
d4 : real*8

Used by: | :k_gcm |
d5 : real*8

Used by: | :k_gcm |
s0 : real*8

Used by: | :k_gcm |
s1 : real*8

Used by: | :k_gcm |
s2 : real*8

Used by: | :k_gcm |
s4 : real*8

Used by: | :k_gcm |
s5 : real*8

Used by: | :k_gcm |
s6 : real*8

Used by: | :k_gcm |
s7 : real*8

Used by: | :k_gcm |
s8 : real*8

Used by: | :k_gcm |
c1 : real*8

Used by: | :e_gcm |
c2 : real*8

Used by: | :e_gcm |
c3 : real*8

Used by: | :e_gcm |
c4 : real*8

Used by: | :e_gcm |
c5 : real*8

Used by: | :e_gcm |
b1 : real*8

Used by: | :atm_diffus | :e_gcm | :k_gcm |
b123 : real*8

Used by: | :atm_diffus |
b2 : real*8

prt : real*8

Used by: | :atm_diffus | :k_gcm |
kappa : real*8, parameter
Von Karman constant
Initial Value = 0.40d0
Used by: | :atm_diffus | :e_gcm | :k_gcm | :l_gcm |
g0 : real*8

d1_3 : real*8

d2_3 : real*8

d3_3 : real*8

d4_3 : real*8

d5_3 : real*8

s0_3 : real*8

s1_3 : real*8

s2_3 : real*8

s3_3 : real*8

s4_3 : real*8

s5_3 : real*8

s6_3 : real*8

g1 : real*8

g2 : real*8

g3 : real*8

g4 : real*8

g5 : real*8

Used by: | :k_gcm |
g6 : real*8

g7 : real*8

g8 : real*8

zgs : real*8, parameter
height of surface layer (m)
Initial Value = 10.
Used by: | :atm_diffus | :get_dbl | :init_pbl | :zze | PBL_DRV:pbl |
zgs : real*8, parameter
height of surface layer (m)
Initial Value = 10.
Used by: | :atm_diffus | :get_dbl | :init_pbl | :zze | PBL_DRV:pbl |
sigma : real*8, parameter
Initial Value = 0.95d0
sigma1 : real*8, parameter
Initial Value = 1.-sigma
sigma2 : real*8, parameter
Initial Value = 1.+sigma
gamamu : real*8, parameter
Initial Value = 19.0d0
gamahu : real*8, parameter
Initial Value = 11.6d0
gamams : real*8, parameter
Initial Value = 5.3d0
gamahs : real*8, parameter
Initial Value = 8.d0/sigma
zet1 : real*8, parameter
Initial Value = 0.5d0
slope1 : real*8, parameter
Initial Value = 0.1d0
zetm : real*8, parameter
Initial Value = -1.464d0
zeth : real*8, parameter
Initial Value = -1.072d0
cmin : real*8, parameter
limits on drag coeffs.
Initial Value = smin*smin
cmax : real*8, parameter
limits on drag coeffs.
Initial Value = smax*smax
smin : real*8, parameter
limits on drag coeffs.
Initial Value = 0.005d0
smax : real*8, parameter
limits on drag coeffs.
Initial Value = 0.25d0
emax : real*8, parameter
limit on turbulent kinetic energy
Initial Value = 1.d5
Used by: | :e_gcm |
ustar_min : real*8, parameter
limit on surface friction speed
Initial Value = 1d-2
Used by: | :atm_diffus |
se : real*8, parameter
Initial Value = 0.1d0
k_max : real*8, parameter
Initial Value = 500.d0
kmmin : real*8, parameter
Initial Value = 1.5d-5
khmin : real*8, parameter
Initial Value = 2.5d-5
kqmin : real*8, parameter
Initial Value = 2.5d-5
kemin : real*8, parameter
Initial Value = 1.5d-5
emin : real*8, parameter
Initial Value = 1d-6
xdelt : real*8, parameter
When used in place of deltx in expressions involving
Initial Value = 0d0
Used by: | :init_pbl | PBL_DRV |
lmonin_min : real*8, parameter
Initial Value = 1d-20
Used by: | :atm_diffus | :e_gcm |
lmonin_max : real*8, parameter
Initial Value = 1d20
Used by: | :atm_diffus | :e_gcm |
twoby3 : real*8, parameter
2/3 constant
Initial Value = 2d0/3d0
skin_effect : integer
sets whether skin effects are used or not
Initial Value = 1
Used by: | :init_pbl |
calc_wspdf : integer
if calc_wspdf==1, calculate mean surface
Initial Value = 0
Used by: | :init_pbl |

Simplex Website Curator: Igor Aleinov — NASA Official: Gavin A. Schmidt

Contact GISS NASA Privacy PolicyAccessibility