From 77cb4a3a950f46e37ad2b085a1140d78d275ec7d Mon Sep 17 00:00:00 2001
From: sbiri <sbiri@noc.ac.uk>
Date: Wed, 30 Jun 2021 12:26:46 +0100
Subject: [PATCH] removed strict constraint in get_L renamed L options "S80" ->
"tsrv" and "ecmwf" -> "Rb" changed default L option relative to method used
in get_init
---
AirSeaFluxCode.py | 7 ++++---
flux_subs.py | 6 +++---
get_init.py | 10 ++++++----
3 files changed, 13 insertions(+), 10 deletions(-)
diff --git a/AirSeaFluxCode.py b/AirSeaFluxCode.py
index de19904..9997d78 100644
--- a/AirSeaFluxCode.py
+++ b/AirSeaFluxCode.py
@@ -80,9 +80,10 @@ def AirSeaFluxCode(spd, T, SST, lat=None, hum=None, P=None, hin=18, hout=10,
set 1 to keep points
L : str
Monin-Obukhov length definition options
- "S80" : default for S80, S88, LP82, YT96 and LY04
- "ecmwf" : following ecmwf (IFS Documentation cy46r1), default for
- ecmwf
+ "tsrv" : default for "S80", "S88", "LP82", "YT96", "UA", "LY04",
+ "C30", "C35"
+ "Rb" : following ecmwf (IFS Documentation cy46r1), default for
+ "ecmwf", "Beljaars"
Returns
-------
res : array that contains
diff --git a/flux_subs.py b/flux_subs.py
index 771cf07..ecda6f9 100755
--- a/flux_subs.py
+++ b/flux_subs.py
@@ -935,6 +935,7 @@ def get_L(L, lat, usr, tsr, qsr, hin, Ta, sst, qair, qsea, wind, monob, zo,
Rb = np.empty(sst.shape)
# as in aerobulk One_on_L in mod_phymbl.f90
tsrv = tsr*(1+0.6077*qair)+0.6077*Ta*qsr
+ # tsrv = tsr+0.6077*Ta*qsr
# from eq. 3.24 ifs Cy46r1 pp. 37
thvs = sst*(1+0.6077*qsea) # virtual SST
dthv = (Ta-sst)*(1+0.6077*qair)+0.6077*Ta*(qair-qsea)
@@ -943,13 +944,12 @@ def get_L(L, lat, usr, tsr, qsr, hin, Ta, sst, qair, qsea, wind, monob, zo,
uz = (wind-usr/kappa*(np.log(hin[0]/hin[1])-psim +
psim_calc(hin[1]/monob, meth)))
Rb = g*dthv*hin[1]/(tv*uz*uz)
- if (L == "S80"):
+ if (L == "tsrv"):
temp = (g*kappa*tsrv /
np.maximum(np.power(usr, 2)*Ta*(1+0.6077*qair), 1e-9))
temp = np.minimum(np.abs(temp), 200)*np.sign(temp)
- temp = np.minimum(np.abs(temp), 10/hin[0])*np.sign(temp)
monob = 1/np.copy(temp)
- elif (L == "ecmwf"):
+ elif (L == "Rb"):
zol = (Rb*(np.power(np.log((hin[1]+zo)/zo)-psim_calc((hin[1]+zo) /
monob, meth) +
psim_calc(zo/monob, meth), 2) /
diff --git a/get_init.py b/get_init.py
index 43038e9..401e668 100644
--- a/get_init.py
+++ b/get_init.py
@@ -146,10 +146,12 @@ def get_init(spd, T, SST, lat, hum, P, Rl, Rs, cskin, skin, wl, gust, L, tol,
gust = [0, 0, 0]
elif (np.size(gust) < 3):
sys.exit("gust input must be a 3x1 array")
- if (L not in [None, "S80", "ecmwf"]):
- sys.exit("L input must be either None, S80 or ecmwf")
- if (L == None):
- L = "ecmwf"
+ if (L not in [None, "tsrv", "Rb"]):
+ sys.exit("L input must be either None, tsrv or Rb")
+ if ((L == None) and ((meth != "ecmwf") and (meth != "Beljaars"))):
+ L = "tsrv"
+ elif ((L == None) and ((meth == "ecmwf") or (meth == "Beljaars"))):
+ L = "Rb"
if (tol == None):
tol = ['all', 0.01, 0.01, 1e-05, 1e-3, 0.1, 0.1]
elif (tol[0] not in ['flux', 'ref', 'all']):
--
GitLab