!copyright (C) 2001 MSC-RPN COMM %%%RPNPHY%%%
*** S/P DRAG1
*
#include "phy_macros_f.h"
SUBROUTINE DRAG1 ( TS, WG, WR, THETAA, VMOD, HU, 1,1
1 PS, RS, VEG, Z0H, Z0TOT, WFC,
1 PSNG, PSNV, PSNZ0, LAI, ZA, FCOR,
1 RESA, ILMO, HST, UE2, FTEMP, FVAP,
1 CH, CD, HRSURF, HUSURF, HV, DEL, ZQS,
1 CTU, N )
*
#include "impnone.cdk"
*
INTEGER N
REAL TS(N), WG(N), WR(N), THETAA(N), VMOD(N), HU(N)
REAL PS(N), RS(N), VEG(N), Z0TOT(N), WFC(N)
REAL Z0H(N)
REAL PSNG(N), PSNV(N), PSNZ0(N), LAI(N), ZA(N)
REAL FCOR(N)
REAL RESA(N), ILMO(N), HST(N), UE2(N), FTEMP(N), FVAP(N)
REAL CH(N), CD(N), HUSURF(N), HV(N), DEL(N), ZQS(N)
REAL HRSURF(N), CTU(N)
*
*Author
* S. Belair (January 1997)
*Revisions
* 001 S. Belair (November 1998)
* Use FLXSURF1 to calculate the surface transfer
* coefficients instead of the MOMCOEF and HEATCOEF
* subroutines (Mascart method).
*
* 002 S. Belair (November 1998)
* Remove Z0 from the arguments (because now we
* use FLXSURF1.
*
* 003 B. Bilodeau (January 2001)
* Automatic arrays
*
* 004 S. Belair (September 2001)
* Add protection to calculation of DEL
*
*Object
*
* Calculates the drag coefficients for heat and momentum transfers
* over ground (i.e., Ch and Cd).
*
*
*Method
*
*
* 1) computes hu, hv, and DEL
*
* 2) use this to find qsoil, the grid-averaged relative humidity
* of the soil
*
* 3) find the transfer and resistance coefficients Ch, Cd, and Ra
* Calculate the surface fluxes of heat, moisture,
* and momentum over water surfaces.
*
*Arguments
*
* - Input/Output -
* RESA aerodynamical surface resistance
* ILMO
* HST
* UE2
* FTEMP
* FVAP
*
* - Input -
* TS surface temperature
* WG superficial volumetric water content
* WR water content retained by the vegetation canopy
* THETAA potential temperature at the lowest level
* VMOD module of the surface winds
* HU specific humidity of air at the lowest level
* PS surface pressure
* RS surface or stomatal resistance
* VEG fraction of a model grid area covered by vegetation
* Z0H roughness length for heat transfers
* Z0TOT roughness length including the effect of snow
* WFC volumetric water content at the field capacity
* PSNG fraction of bare ground covered by snow
* PSNV fraction of vegetation covered by snow
* PSNZ0 snow fraction for roughness length calculations
* LAI leaf area index
* ZA reference height
* FCOR Coriolis factor
*
*
* - Output -
* CH drag coefficient for heat
* CD drag coefficient for momentum
* HRSURF relative humidity of the surface
* HUSURF specific humidity of the surface
* HV Halstead coefficient (i.e., relative humidity of the
* vegetation canopy)
* DEL fraction of canopy covered by intercepted water
* ZQS area-averaged relative humidity of a model tile
* CTU homogeneous boundary condition term in the
* diffusion equation for Theta and Q
*
*
#include "consphy.cdk"
*
INTEGER I
REAL UE
*
*
*
EXTERNAL FLXSURF2
*
*
************************************************************************
* AUTOMATIC ARRAYS
************************************************************************
*
AUTOMATIC ( TEMP , REAL , (N) )
AUTOMATIC ( WRMAX , REAL , (N) )
AUTOMATIC ( QSAT , REAL , (N) )
AUTOMATIC ( COEF , REAL , (N) )
AUTOMATIC ( CMU , REAL , (N) )
AUTOMATIC ( SCR1 , REAL , (N) )
AUTOMATIC ( SCR6 , REAL , (N) )
AUTOMATIC ( SCR7 , REAL , (N) )
AUTOMATIC ( SCR8 , REAL , (N) )
AUTOMATIC ( SCR9 , REAL , (N) )
AUTOMATIC ( SCR10 , REAL , (N) )
AUTOMATIC ( SCR11 , REAL , (N) )
*
************************************************************************
*
*
#include "dintern.cdk"
#include "fintern.cdk"
*
*
*
*------------------------------------------------------------------------
*
*
*
*
** 1. RELATIVE AND SPECIFIC HUMIDITY OF THE GROUND (HU)
* -------------------------------------------------
*
* this relative humidity is related to
* the superficial soil moisture and the
* field capacity of the ground
*
DO I=1,N
TEMP(I) = PI*WG(I)/WFC(I)
HRSURF(I) = 0.5 * ( 1.-COS(TEMP(I)) )
END DO
*
* there is a specific treatment for dew
* (see Mahfouf and Noilhan, jam, 1991)
*
* first calculate the saturation vapor
* pressure and specific humidity
*
DO I=1,N
QSAT(I) = FOQST( TS(I), PS(I) )
END DO
*
*
DO I=1,N
*
*
* when hu*qsat < qa, there are two
* possibilities
*
* low-level air is dry, i.e.,
* qa < qsat
*
*
IF ( HRSURF(I)*QSAT(I).LT.HU(I).AND.QSAT(I).GT.HU(I) )
1 HRSURF(I) = HU(I) / QSAT(I)
*
*
* b) low-level air is humid, i.e.,
* qa >= qsat
*
IF ( HRSURF(I)*QSAT(I).LT.HU(I).AND.QSAT(I).LE.HU(I) )
1 HRSURF(I) = 1.0
*
* for very humid soil (i.e., wg > wfc ),
* we take hu=1
*
IF ( WG(I).GT.WFC(I) )
1 HRSURF(I) = 1.0
*
*
END DO
*
*
DO I=1,N
HUSURF(I) = HRSURF(I) * QSAT(I)
END DO
*
*
*
*
** 2. FRACTION OF THE FOLIAGE COVERED BY INTERCEPTED WATER (DEL)
* ------------------------------------------------------------
*
* first calculate the maximum value of
* equivalent water content in the
* vegetation canopy
*
DO I=1,N
*
WRMAX(I) = 0.2 * VEG(I) * LAI(I)
*
* calculate DEL
*
COEF(I) = 1. + 2.*LAI(I)
*
IF ( VEG(I).GT.0.0.AND.WRMAX(I).GT.0.0 ) THEN
DEL(I) = MIN(WR(I),WRMAX(I)) /
1 ( (1.-COEF(I))*MIN(WR(I),WRMAX(I)) + COEF(I)*WRMAX(I) )
ELSE
DEL(I) = 0.0
END IF
*
END DO
*
*
** 3. HALSTEAD COEFFICIENT (RELATIVE HUMIDITY OF THE VEGETATION) (HV)
* ---------------------------------------------------------------
*
DO I=1,N
HV(I) = 1. - MAX(0.,SIGN(1.,QSAT(I)-HU(I)))
1 *RS(I)*(1.-DEL(I)) / (RESA(I)+RS(I))
END DO
*
*
*
*
** 4. GRID-AVERAGED HUMIDITY OF THE SOIL (ZQS)
* ----------------------------------------
*
DO I=1,N
ZQS(I) = ( (1.-VEG(I))*(1.-PSNG(I))*HRSURF(I)
1 + (1.-VEG(I))* PSNG(I)
1 + VEG(I) *(1.-PSNV(I))*HV(I)
1 + VEG(I) * PSNV(I) )*QSAT(I)
1 + VEG(I) *(1.-PSNV(I))*(1.-HV(I))*HU(I)
END DO
*
*
*** 5. SURFACE TRANSFER COEFFICIENTS FOR HEAT AND MOMENTUM (CH and CD)
** ---------------------------------------------------------------
*
*
* Arbitrary values for the height
* of the boundary layer (for FLXSURF2)
* No impact on the results.
*
*
CALL FLXSURF2( CMU, CTU, SCR1, FTEMP, FVAP, ILMO,
1 UE2, FCOR, THETAA, HU, ZA, VMOD, TS,
1 ZQS, HST, Z0TOT, Z0H, SCR6, SCR7,
1 SCR8, SCR9, SCR10, SCR11, N )
*
*
DO I=1,N
UE = SQRT( UE2(I) )
CMU(I) = CMU(I) / UE
*
CD(I) = CMU(I) * CMU(I)
CH(I) = CMU(I) * CTU(I)/UE
*
RESA(I) = 1. / CH(I) / (VMOD(I)+0.001)
END DO
*
*
*
*
** 7. HALSTEAD COEFFICIENT (WITH THE NEW VALUES OF RESA)
* --------------------------------------------------
*
*
DO I=1,N
HV(I) = 1. - MAX(0.,SIGN(1.,QSAT(I)-HU(I)))
1 *RS(I)*(1.-DEL(I)) / (RESA(I)+RS(I))
END DO
*
*
*
*
RETURN
END