Changeset 1804 in lmdz_wrf

Timestamp:

Mar 15, 2018, 5:54:49 PM (7 years ago)

Author:

lfita

Message:

Adding:

• `evspsblpot': potential evapotranspiration following Penman-Monteith formulation implemented in ORCHIDEE

Location:

trunk/tools

Files:

• diag_tools.py (modified) (3 diffs)
• diagnostics.inf (modified) (1 diff)
• diagnostics.py (modified) (1 diff)
• module_ForDiagnostics.f90 (modified) (1 diff)
• module_ForDiagnosticsVars.f90 (modified) (4 diffs)
• variables_values.dat (modified) (1 diff)

trunk/tools/diag_tools.py

@@ -29 +29 @@
 # compute_WRFuasvas: Fucntion to compute geographical rotated WRF 2-meter winds
 # derivate_centered: Function to compute the centered derivate of a given field
-# Forcompute_cllmh: Function to compute cllmh: low/medium/hight cloud fraction following newmicro.F90 from LMDZ via Fortran subroutine
-# Forcompute_clt: Function to compute the total cloud fraction following 'newmicro.F90' from LMDZ via a Fortran module
-# Forcompute_psl_ptarget: Function to compute the sea-level pressure following target_pressure value found in `p_interp.F'
-# Forcompute_zmla_gen: Function to compute the boundary layer height following a generic method with Fortran
-# Forcompute_zwind: Function to compute the wind at a given height following the power law method
-# Forcompute_zwind_log: Function to compute the wind at a given height following the logarithmic law method
-# Forcompute_zwindMO: Function to compute the wind at a given height following the Monin-Obukhov theory
+# Forcompute_cllmh: Function to compute cllmh: low/medium/hight cloud fraction 
+#   following newmicro.F90 from LMDZ via Fortran subroutine
+# Forcompute_clt: Function to compute the total cloud fraction following 
+#   'newmicro.F90' from LMDZ via a Fortran module
+# Forcompute_potevap_orPM: Function to compute potential evapotranspiration following
+#   Penman-Monteith formulation implemented in ORCHIDEE
+# Forcompute_psl_ptarget: Function to compute the sea-level pressure following 
+#   target_pressure value found in `p_interp.F'
+# Forcompute_zmla_gen: Function to compute the boundary layer height following a 
+#   generic method with Fortran
+# Forcompute_zwind: Function to compute the wind at a given height following the 
+#   power law method
+# Forcompute_zwind_log: Function to compute the wind at a given height following the 
+#   logarithmic law method
+# Forcompute_zwindMO: Function to compute the wind at a given height following the 
+#   Monin-Obukhov theory
 # W_diagnostic: Class to compute WRF diagnostics variables
 
 # Others just providing variable values
 # var_cllmh: Fcuntion to compute cllmh on a 1D column
-# var_clt: Function to compute the total cloud fraction following 'newmicro.F90' from LMDZ using 1D vertical column values
-# var_hur: Function to compute relative humidity following 'August - Roche - Magnus' formula
-# var_hur_Uhus: Function to compute relative humidity following 'August - Roche - Magnus' formula using hus
+# var_clt: Function to compute the total cloud fraction following 'newmicro.F90' from
+#   LMDZ using 1D vertical column values
+# var_hur: Function to compute relative humidity following 'August - Roche - Magnus' 
+#   formula
+# var_hur_Uhus: Function to compute relative humidity following 'August-Roche-Magnus'
+#   formula using hus
 # var_mslp: Fcuntion to compute mean sea-level pressure
-# var_td: Function to compute dew-point air temperature from temperature and pressure values
-# var_td_Uhus: Function to compute dew-point air temperature from temperature and pressure values using hus
+# var_td: Function to compute dew-point air temperature from temperature and pressure
+#   values
+# var_td_Uhus: Function to compute dew-point air temperature from temperature and 
+#   pressure values using hus
 # var_virtualTemp: This function returns virtual temperature in K, 
 # var_WRFtime: Function to copmute CFtimes from WRFtime variable
@@ -51 +65 @@
 # var_wd: Function to compute the wind speed
 # rotational_z: z-component of the rotatinoal of horizontal vectorial field
-# turbulence_var: Function to compute the Taylor's decomposition turbulence term from a a given variable
+# turbulence_var: Function to compute the Taylor's decomposition turbulence term from
+#   a a given variable
 
 import numpy as np
@@ -1079 +1094 @@
 
     return var1, var2, vardims, varvdims
+
+def Forcompute_potevap_orPM(rho1, ust, uas, vas, tas, ps, qv1, dimns, dimvns):
+    """ Function to compute potential evapotranspiration following Penman-Monteith 
+      formulation implemented in ORCHIDEE
+    Forcompute_potevap_orPM(rho1, uas, vas, tas, ps, qv2, qv1, dimns, dimvns)
+      [rho1]= air-density at the first layer (assuming [[t],y,m]) [kgm-3]
+      [ust]= u* in similarity theory (assuming [[t],y,x]) [ms-1]
+      [uas]= x-component of unstaggered 10 m wind (assuming [[t],y,x]) [ms-1]
+      [vas]= y-component of unstaggered 10 m wind (assuming [[t],y,x]) [ms-1]
+      [tas]= 2m air temperature [K]
+      [ps]= surface pressure [Pa]
+      [qv1]= mixing ratio at the first atmospheric layer [kgkg-1]
+      [dimns]= list of the name of the dimensions of [uas]
+      [dimvns]= list of the name of the variables with the values of the 
+        dimensions of [uas]
+    """
+    fname = 'Forcompute_potevap_orPM'
+
+    vardims = dimns[:]
+    varvdims = dimvns[:]
+
+    if len(uas.shape) == 3:
+        var1= np.zeros((uas.shape[0],uas.shape[1],uas.shape[2]), dtype=np.float)
+        var2= np.zeros((uas.shape[0],uas.shape[1],uas.shape[2]), dtype=np.float)
+
+        dx = uas.shape[2]
+        dy = uas.shape[1]
+        dt = uas.shape[0]
+
+        pvar = fdin.module_fordiagnostics.compute_potevap_orpm3d(                    \
+          rho1=rho1.transpose(), ust=ust.transpose(), uas=uas.transpose(),           \
+          vas=vas.transpose(), tas=tas.transpose(), ps=ps.transpose(),               \
+          qv1=qv1.transpose(), d1=dx, d2=dy, d3=dt)
+        var = pvar.transpose()
+    else:
+        print errormsg
+        print '  ' + fname + ': rank', len(uas.shape), 'not ready !!'
+        print '  it only computes 3D [t,y,x] rank values'
+        quit(-1)
+
+    return var, vardims, varvdims
+
+####### ###### ##### #### ### ## # END Fortran diagnostics
 
 def compute_OMEGAw(omega, p, t, dimns, dimvns):

trunk/tools/diagnostics.inf

@@ -18 +18 @@
 mrso, WRFmrso,  SMOIS@DZS
 p, WRFp, P@PB
+evspsblpot, WRFpotevap_orPM, WRFdens@UST@U10@V10@T2@PSFC@QVAPOR
 pr, RAINTOT, RAINC@RAINNC@WRFtime
 pracc, ACRAINTOT, RAINC@RAINNC@WRFtime

trunk/tools/diagnostics.py

@@ -871 +871 @@
         ncvar.insert_variable(ncobj, 'zg', diagout, dnames, diagoutvd, newnc)
 
+# WRFpotevap_orPM potential evapotranspiration following Penman-Monteith formulation
+#   implemented in ORCHIDEE as: WRFdens, UST, U10, V10, T2, PSFC, QVAPOR
+    elif diagn == 'WRFpotevap_orPM':
+        var0 = WRFdens[:,0,:,:]
+        var1 = ncobj.variables[depvars[1]][:]
+        var2 = ncobj.variables[depvars[2]][:]
+        var3 = ncobj.variables[depvars[3]][:]
+        var4 = ncobj.variables[depvars[4]][:]
+        var5 = ncobj.variables[depvars[5]][:]
+        var6 = ncobj.variables[depvars[6]][:,0,:,:]
+
+        dnamesvar = list(ncobj.variables[depvars[1]].dimensions)
+        dvnamesvar = ncvar.var_dim_dimv(dnamesvar,dnames,dvnames)
+
+        diagout = np.zeros(var1.shape, dtype=np.float)
+        diagout, diagoutd, diagoutvd = diag.Forcompute_potevap_orPM(var0, var1, var2,\
+          var3, var4, var5, var6, dnamesvar, dvnamesvar)
+
+        # Removing the nonChecking variable-dimensions from the initial list
+        varsadd = []
+        for nonvd in NONchkvardims:
+            if gen.searchInlist(dvnames,nonvd): diagoutvd.remove(nonvd)
+            varsadd.append(nonvd)
+
+        ncvar.insert_variable(ncobj, 'evspsblpot', diagout, diagoutd, diagoutvd, newnc)
+
 # WRFmslp_ptarget sea-level pressure following ECMWF method as PSFC, HGT, WRFt, WRFp, ZNU, ZNW
     elif diagn == 'WRFpsl_ecmwf':

trunk/tools/module_ForDiagnostics.f90

@@ -795 +795 @@
   END SUBROUTINE compute_zwindMO3D
 
+  SUBROUTINE compute_potevap_orPM3D(d1, d2, d3, rho1, ust, uas, vas, tas, ps, qv1, potevap)
+! Subroutine to compute potential evapotranspiration Penman-Monteith formulation implemented in 
+!   ORCHIDEE
+
+    IMPLICIT NONE
+
+    INTEGER, INTENT(in)                                  :: d1, d2, d3
+    REAL(r_k), DIMENSION(d1,d2,d3), INTENT(in)           :: rho1, ust, uas, vas, tas, ps, qv1
+    REAL(r_k), DIMENSION(d1,d2,d3), INTENT(out)          :: potevap
+ 
+! Local
+    INTEGER                                              :: i, j, it
+
+!!!!!!! Variables
+! rho1: atsmophere density at the first layer [kgm-3]
+! ust: u* in similarity theory [ms-1]
+! uas: x-component 10-m wind speed [ms-1]
+! vas: y-component 10-m wind speed [ms-1]
+! tas: 2-m atmosphere temperature [K]
+! ps: surface pressure [Pa]
+! qv1: 1st layer atmospheric mixing ratio [kgkg-1]
+! potevap: potential evapo transpiration [kgm-2s-1]
+
+    fname = 'compute_potevap_orPM3D'
+
+    DO i=1, d1
+      DO j=1, d2
+        DO it=1, d3
+          CALL var_potevap_orPM(rho1(i,j,it), ust(i,j,it), uas(i,j,it), vas(i,j,it), tas(i,j,it),     &
+            ps(i,j,it), qv1(i,j,it), potevap(i,j,it))
+        END DO
+      END DO
+    END DO
+
+    RETURN
+
+  END SUBROUTINE compute_potevap_orPM3D
+
 END MODULE module_ForDiagnostics

trunk/tools/module_ForDiagnosticsVars.f90

@@ -15 +15 @@
 
 !!!!!!! Variables
+! Cdrag_0: Fuction to compute a first order generic approximation of the drag coefficient
 ! compute_psl_ptarget4d2: Compute sea level pressure using a target pressure. Similar to the Benjamin 
 !   and Miller (1990). Method found in p_interp.F90
@@ -28 +29 @@
 ! var_cllmh: low, medium, high-cloud [0,1]
 ! var_clt: total cloudiness [0,1]
+! var_potevap_orPM: potential evapotranspiration following Penman-Monteith formulation implemented in ORCHIDEE
 ! var_psl_ecmwf: sea level pressure using ECMWF method following Mats Hamrud and Philippe Courtier [Pa]
 ! var_zmla_generic: Subroutine to compute pbl-height following a generic method
@@ -1279 +1281 @@
   !   Springer (p376-383)
   ! Only usefull for newz < 80. m
+  ! Ackonwledgement: M. A. Jimenez, UIB
 
     IMPLICIT NONE
@@ -1378 +1381 @@
   END SUBROUTINE stabfunc_businger
 
+  REAL(r_k) FUNCTION Cdrag_0(ust,uas,vas)
+! Fuction to compute a first order generic approximation of the drag coefficient as
+!   CD = (ust/wss)**2
+!  after, Garratt, J.R., 1992.: The Atmospheric Boundary Layer. Cambridge Univ. Press, 
+!    Cambridge, U.K., 316 pp
+! Ackonwledgement: M. A. Jimenez, UIB
+!
+    IMPLICIT NONE
+
+    REAL(r_k), INTENT(in)                                :: ust, uas, vas
+
+!!!!!!! Variables
+! ust: u* in similarity theory [ms-1]
+! uas, vas: x/y-components of wind at 10 m
+
+    fname = 'Cdrag_0'
+
+    Cdrag_0 = ust**2/(uas**2+vas**2)
+
+  END FUNCTION Cdrag_0
+
+  SUBROUTINE var_potevap_orPM(rho1, ust, uas, vas, tas, ps, qv1, potevap)
+! Subroutine to compute the potential evapotranspiration following Penman-Monteith formulation 
+!   implemented in ORCHIDEE
+!      potevap = dt*rho1*qc*(q2sat-qv1)
+
+    IMPLICIT NONE
+
+    REAL(r_k), INTENT(in)                                :: rho1, ust, uas, vas, tas, ps, qv1
+    REAL(r_k), INTENT(out)                               :: potevap
+
+! Local
+    REAL(r_k)                                            :: q2sat, Cd, qc
+
+!!!!!!! Variables
+! rho1: atsmophere density at the first layer [kgm-3]
+! ust: u* in similarity theory [ms-1]
+! uas, vas: x/y-components of 10-m wind [ms-1]
+! tas: 2-m atmosphere temperature [K]
+! ps: surface pressure [Pa]
+! qv1: 1st layer atmospheric mixing ratio [kgkg-1]
+! potevap: potential evapo transpiration [kgm-2s-1]
+  fname = 'var_potevap_orPM'
+
+  ! q2sat: Saturated air at 2m (can be assumed to be q2 == qsfc?)
+  q2sat = SaturationMixingRatio(tas, ps)
+
+  ! Cd: drag coeffiecient
+  Cd = Cdrag_0(ust, uas, vas)
+
+  ! qc: surface drag coefficient
+  qc = SQRT(uas**2 + vas**2)*Cd
+
+  potevap = MAX(zeroRK, rho1*qc*(q2sat - qv1))
+
+  END SUBROUTINE var_potevap_orPM
+
 END MODULE module_ForDiagnosticsVars

trunk/tools/variables_values.dat

@@ -159 +159 @@
 evspsbl, evspsblac, water_evaporation_flux_ac, 0., 1.5e-4, accumulated|water|evaporation|flux, kgm-2, Blues, $evspsbl^{ac}$, evspsblac
 SFCEVPde, evspsblac, water_evaporation_flux_ac, 0., 1.5e-4, accumulated|water|evaporation|flux, kgm-2, Blues, $evspsbl^{ac}$, evspsblac
+evspsblpot, evspsblpot, potential_water_evapotranspiration_flux, 0., 1.5e-4, potential|evapotranspiration|flux, kgm-2s-1, Blues, $evspsblpot$, evspsblpot
+EVSPSBLPOT, evspsblpot, potential_water_evapotranspiration_flux, 0., 1.5e-4, potential|evapotranspiration|flux, kgm-2s-1, Blues, $evspsblpot$, evspsblpot
 g, g, grauepl_mixing_ratio, 0., 0.0003, graupel|mixing|ratio, kgkg-1, Purples, $qg$, qg
 QGRAUP, g, grauepl_mixing_ratio, 0., 0.0003, graupel|mixing|ratio, kgkg-1, Purples, $qg$, qg