*IF DEF,A03_3A,OR,DEF,A03_5A AJS1F401.1483
C ******************************COPYRIGHT****************************** GTS2F400.2611
C (c) CROWN COPYRIGHT 1995, METEOROLOGICAL OFFICE, All Rights Reserved. GTS2F400.2612
C GTS2F400.2613
C Use, duplication or disclosure of this code is subject to the GTS2F400.2614
C restrictions as set forth in the contract. GTS2F400.2615
C GTS2F400.2616
C Meteorological Office GTS2F400.2617
C London Road GTS2F400.2618
C BRACKNELL GTS2F400.2619
C Berkshire UK GTS2F400.2620
C RG12 2SZ GTS2F400.2621
C GTS2F400.2622
C If no contract has been raised with this copy of the code, the use, GTS2F400.2623
C duplication or disclosure of it is strictly prohibited. Permission GTS2F400.2624
C to do so must first be obtained in writing from the Head of Numerical GTS2F400.2625
C Modelling at the above address. GTS2F400.2626
C ******************************COPYRIGHT****************************** GTS2F400.2627
C GTS2F400.2628
C*LL SUBROUTINE EX_COEF------------------------------------------------ EXCOEF3A.3
CLL EXCOEF3A.4
CLL Purpose: To calculate exchange coefficients for boundary layer EXCOEF3A.5
CLL subroutine KMKH. EXCOEF3A.6
CLL EXCOEF3A.7
CLL Suitable for single column use (via *IF definition IBM). EXCOEF3A.8
CLL EXCOEF3A.9
CLL Version 1 written by Fiona Hewer, May 1992. EXCOEF3A.10
CLL EXCOEF3A.11
CLL Model Modification history: EXCOEF3A.12
CLL version Date EXCOEF3A.13
CLL EXCOEF3A.14
CLL 3.4 18/10/94 *DECK inserted into UM version 3.4. S Jackson EXCOEF3A.15
CLL EXCOEF3A.16
CLL 4.3 04/02/97 Reduction in mixing length above b.l. top ARN1F403.51
CLL under control of logical switch L_MIXLEN & ARN1F403.52
CLL removal of enhanced stable momentum mixing ARN1F403.53
CLL under control of logical switch L_MOM. ARN1F403.54
CLL R.N.B.Smith ARN1F403.55
!LL 4.5 23/10/98 Prevent TIMER from performing barrier P.Burton GPB8F405.66
CLL ARN1F403.56
CLL Programming standard: Unified Model Documentation Paper No 4, EXCOEF3A.17
CLL Version 2, dated 18/1/90. EXCOEF3A.18
CLL EXCOEF3A.19
CLL System component covered: Part of P243. EXCOEF3A.20
CLL EXCOEF3A.21
CLL Project task: EXCOEF3A.22
CLL EXCOEF3A.23
CLL Documentation: UMDP No.24 EXCOEF3A.24
CLL EXCOEF3A.25
CLL--------------------------------------------------------------------- EXCOEF3A.26
C* EXCOEF3A.27
C*L Arguments :- EXCOEF3A.28
SUBROUTINE EX_COEF ( 5,6EXCOEF3A.29
+ P_FIELD,U_FIELD,P1,P_POINTS,BL_LEVELS, EXCOEF3A.30
+ CCB,CCT,NTML,L_MOM,L_MIXLEN, ARN1F403.57
+ AKH,BKH,CCA,DZL,PSTAR,RDZ,RI,TV,U_P,V_P,ZH,ZLB,Z0M,H_BLEND, EXCOEF3A.32
+ RHOKM,RHOKH,LTIMER EXCOEF3A.33
+) EXCOEF3A.34
IMPLICIT NONE EXCOEF3A.35
LOGICAL LTIMER EXCOEF3A.36
C ARN1F403.58
LOGICAL ARN1F403.59
& L_MOM ! IN Switch for convective momentum transport. ARN1F403.60
&,L_MIXLEN ! IN Switch for reducing the turbulent mixing ARN1F403.61
C ! length above the top of the boundary layer. ARN1F403.62
C ARN1F403.63
INTEGER EXCOEF3A.37
+ P_FIELD ! IN No. of P-grid points in whole field. EXCOEF3A.38
+,U_FIELD ! IN No. of U-grid points in whole field. EXCOEF3A.39
+,P1 ! IN First P-grid point to be processed. EXCOEF3A.40
+,P_POINTS ! IN No. of P-grid points to be processed. EXCOEF3A.41
+,BL_LEVELS ! IN maximum number of boundary layer levels EXCOEF3A.42
INTEGER EXCOEF3A.43
+ CCB(P_FIELD) ! IN Convective Cloud Base. EXCOEF3A.44
+,CCT(P_FIELD) ! IN Convective Cloud Top. EXCOEF3A.45
&,NTML(P_FIELD) ! IN Number of turbulently mixed ARN1F403.64
C ! layers. ARN1F403.65
REAL EXCOEF3A.46
+ AKH(BL_LEVELS) ! IN Hybrid "A" for layer interfaces. EXCOEF3A.47
C ! AKH(K) is value for lower boundary EXCOEF3A.48
C ! of layer K. EXCOEF3A.49
+,BKH(BL_LEVELS) ! IN Hybrid "B" for layer interfaces. EXCOEF3A.50
C ! BKH(K) is value for lower boundary EXCOEF3A.51
C ! of layer K. EXCOEF3A.52
+,CCA(P_FIELD) ! IN Convective Cloud Amount. EXCOEF3A.53
+,DZL(P_FIELD,BL_LEVELS) ! IN Layer depths (m). DZL(,K) is the EXCOEF3A.54
C ! distance from layer boundary K-1/2 EXCOEF3A.55
C ! to layer boundary K+1/2. For K=1 EXCOEF3A.56
C ! the lower boundary is the surface. EXCOEF3A.57
+,PSTAR(P_FIELD) ! IN Surface pressure (Pa). EXCOEF3A.58
+,RDZ(P_FIELD,BL_LEVELS) ! IN Reciprocal of distance between EXCOEF3A.59
C ! hybrid levels (m-1). 1/RDZ(,K) is EXCOEF3A.60
C ! the vertical distance from level EXCOEF3A.61
C ! K-1 to level K, except that for EXCOEF3A.62
C ! K=1 it is just the height of the EXCOEF3A.63
C ! lowest atmospheric level. EXCOEF3A.64
+,RI(P_FIELD,2:BL_LEVELS) ! IN Richardson number for lower EXCOEF3A.65
C ! interface of layer. EXCOEF3A.66
+,TV(P_FIELD,BL_LEVELS) ! IN Virtual temperature (K) - from EXCOEF3A.67
C ! SUBROUTINE Z. EXCOEF3A.68
+,U_P(P_FIELD,BL_LEVELS) ! IN Westerly wind component on P-grid EXCOEF3A.69
C ! (m per s). EXCOEF3A.70
+,V_P(P_FIELD,BL_LEVELS) ! IN Southerly wind component on P-grid EXCOEF3A.71
C ! (m per s). EXCOEF3A.72
+,ZH(P_FIELD) ! IN Boundary layer height (m). EXCOEF3A.73
+,ZLB(P_FIELD,BL_LEVELS) ! IN ZLB(,K) is height above surface of EXCOEF3A.74
C ! lower boundary of layer K (m). EXCOEF3A.75
+,Z0M(P_FIELD) ! IN Roughness length for momentum (m). EXCOEF3A.76
+,H_BLEND(P_FIELD) ! IN Blending height for effective EXCOEF3A.77
C ! roughness length scheme EXCOEF3A.78
REAL EXCOEF3A.79
+ RHOKM(U_FIELD,BL_LEVELS) ! INOUT Layer K-1 - to - layer K EXCOEF3A.80
C ! exchange coefficient for EXCOEF3A.81
C ! momentum, on UV-grid with first EXCOEF3A.82
C ! and last rows set to "missing EXCOEF3A.83
C ! data".Output as GAMMA*RHOKM* EXCOEF3A.84
C ! RDZUV for P244 (IMPL_CAL). EXCOEF3A.85
+,RHOKH(P_FIELD,BL_LEVELS) ! INOUT Layer K-1 - to - layer K EXCOEF3A.86
C ! exchange coefficient for FTL, EXCOEF3A.87
C ! on P-grid.Output as GAMMA* EXCOEF3A.88
C ! *RHOKH*RDZ for P244(IMPL_CAL) EXCOEF3A.89
C* EXCOEF3A.90
C*L--------------------------------------------------------------------- EXCOEF3A.91
EXTERNAL TIMER EXCOEF3A.92
C* EXCOEF3A.93
C*L--------------------------------------------------------------------- EXCOEF3A.94
C Local and other symbolic constants :- EXCOEF3A.95
*CALL C_LHEAT 
EXCOEF3A.96
*CALL C_R_CP EXCOEF3A.97
*CALL C_EPSLON EXCOEF3A.98
*CALL C_VKMAN EXCOEF3A.99
REAL EH,EM,G0,DH,DM,LAMBDA_MIN,A_LAMBDA EXCOEF3A.100
PARAMETER ( EXCOEF3A.101
+ EH=25.0 ! Used in calc of stability function FH. EXCOEF3A.102
+,EM=4.0 ! Used in calc of stability function FM. EXCOEF3A.103
+,G0=10.0 ! Used in stability function calcs. EXCOEF3A.104
+,DH=G0/EH ! Used in calc of stability function FH. EXCOEF3A.105
+,DM=G0/EM ! Used in calc of stability function FM. EXCOEF3A.106
+,LAMBDA_MIN=40.0 ! Minimum value of length scale LAMBDA. EXCOEF3A.107
+,A_LAMBDA=2.0 ! used in calc of LAMBDA_EFF EXCOEF3A.108
+) EXCOEF3A.109
C* EXCOEF3A.110
C EXCOEF3A.111
C Define local storage. EXCOEF3A.112
C EXCOEF3A.113
C Scalars. EXCOEF3A.114
C EXCOEF3A.115
REAL EXCOEF3A.116
+ DVDZM ! Modulus of wind shear gradient across lower level bdy. EXCOEF3A.117
+,DZU ! Westerly wind shear between adjacent levels. EXCOEF3A.118
+,DZV ! Southerly wind shear between adjacent levels. EXCOEF3A.119
+,DVMOD2 ! Square of modulus of wind shear between adjacent levels EXCOEF3A.120
+,ELH ! Mixing length for heat and moisture at lower layer bdy. EXCOEF3A.121
+,ELM ! Mixing length for momentum at lower layer boundary. EXCOEF3A.122
+,FH ! (Value of) stability function for heat & moisture. EXCOEF3A.123
+,FM ! (Value of) stability function for momentum transport. EXCOEF3A.124
+,RHO ! Air density at lower layer boundary. EXCOEF3A.125
+,RTMRI ! Temporary in stability function calculation. EXCOEF3A.126
+,VKZ ! Temporary in calculation of ELH. EXCOEF3A.127
+,WK ! Temporary in calculation of RHO. EXCOEF3A.128
+,WKM1 ! Temporary in calculation of RHO. EXCOEF3A.129
+,LAMBDAM ! Asymptotic mixing length for turbulent transport EXCOEF3A.130
C ! of momentum. EXCOEF3A.131
+,LAMBDAH ! Asymptotic mixing length for turbulent transport EXCOEF3A.132
C ! of heat/moisture. EXCOEF3A.133
+,LAMBDA_EFF ! Effective mixing length used with effective EXCOEF3A.134
C ! roughness length scheme. EXCOEF3A.135
INTEGER EXCOEF3A.136
+ I ! Loop counter (horizontal field index). EXCOEF3A.137
+,K ! Loop counter (vertical level index). EXCOEF3A.138
+,KM1 ! K-1. EXCOEF3A.139
C EXCOEF3A.140
C Layer interface K_LOG_LAYR-1/2 is the highest which requires log EXCOEF3A.141
C profile correction factors to the vertical finite differences. EXCOEF3A.142
C The value should be reassessed if the vertical resolution is changed. EXCOEF3A.143
C We could set K_LOG_LAYR = BL_LEVELS and thus apply the correction EXCOEF3A.144
C factors for all the interfaces treated by the boundary layer scheme; EXCOEF3A.145
C this would be desirable theoretically but expensive computationally EXCOEF3A.146
C because of the use of the log function. EXCOEF3A.147
C EXCOEF3A.148
INTEGER K_LOG_LAYR EXCOEF3A.149
PARAMETER (K_LOG_LAYR=2) EXCOEF3A.150
C* EXCOEF3A.151
IF (LTIMER) THEN EXCOEF3A.152
CALL TIMER('EX_COEF ',103) GPB8F405.67
ENDIF EXCOEF3A.154
C EXCOEF3A.155
C----------------------------------------------------------------------- EXCOEF3A.156
CL 1. Loop round "boundary" levels; calculate the stability- EXCOEF3A.157
CL dependent turbulent mixing coefficients. EXCOEF3A.158
C----------------------------------------------------------------------- EXCOEF3A.159
C EXCOEF3A.160
DO 2 K=2,BL_LEVELS EXCOEF3A.161
KM1 = K-1 EXCOEF3A.162
CMIC$ DO ALL VECTOR SHARED(BL_LEVELS, P_POINTS, P1, DM, DH, EXCOEF3A.163
CMIC$1 P_FIELD, PSTAR, TV, DZL, RDZ, ZLB, Z0M, LAMBDA_MIN, ZH, EXCOEF3A.164
CMIC$2 U_P, V_P, RI, K, KM1, U_FIELD, RHOKM, RHOKH, AKH, BKH, CCA, EXCOEF3A.165
CMIC$3 CCB, CCT, H_BLEND, NTML) ARN1F403.66
CMIC$4 PRIVATE(RHO, VKZ, ELM, ELH, DZU, DZV, DVMOD2, DVDZM, EXCOEF3A.167
CMIC$5 RTMRI, FM, FH, I, WKM1, WK, LAMBDAM, LAMBDAH, LAMBDA_EFF, EXCOEF3A.168
CMIC$6 A_LAMBDA) EXCOEF3A.169
CDIR$ IVDEP EXCOEF3A.170
! Fujitsu vectorization directive GRB0F405.227
!OCL NOVREC GRB0F405.228
DO 21 I=P1,P1+P_POINTS-1 EXCOEF3A.171
C----------------------------------------------------------------------- EXCOEF3A.172
CL 2.1 Calculate asymptotic mixing lengths LAMBDAM and LAMBDAH EXCOEF3A.173
CL (currently assumed equal). EXCOEF3A.174
C----------------------------------------------------------------------- EXCOEF3A.175
C EXCOEF3A.176
LAMBDAM = MAX ( LAMBDA_MIN , 0.15*ZH(I) ) EXCOEF3A.177
LAMBDAH = LAMBDAM EXCOEF3A.178
LAMBDA_EFF = MAX (LAMBDAM, A_LAMBDA*H_BLEND(I) ) EXCOEF3A.179
IF ( L_MIXLEN ) THEN ARN1F403.67
IF ( K .GE. NTML(I) + 2 ) THEN ARN1F403.68
LAMBDAM = LAMBDA_MIN ARN1F403.69
LAMBDAH = LAMBDA_MIN ARN1F403.70
IF (ZLB(I,K) .GT. A_LAMBDA*H_BLEND(I)) ARN1F403.71
& LAMBDA_EFF = LAMBDA_MIN ARN1F403.72
ENDIF ARN1F403.73
ENDIF ARN1F403.74
C EXCOEF3A.180
C----------------------------------------------------------------------- EXCOEF3A.181
CL 2.2 Calculate mixing lengths ELH, ELM at layer interface K-1/2. EXCOEF3A.182
C----------------------------------------------------------------------- EXCOEF3A.183
C EXCOEF3A.184
C Incorporate log profile corrections to the vertical finite EXCOEF3A.185
C differences into the definitions of ELM and ELH. EXCOEF3A.186
C To save computing logarithms for all K, the values of ELM and ELH EXCOEF3A.187
C are unchanged for K > K_LOG_LAYR. EXCOEF3A.188
C EXCOEF3A.189
IF (K .LE. K_LOG_LAYR) THEN EXCOEF3A.190
VKZ = VKMAN / RDZ(I,K) EXCOEF3A.191
ELM = VKZ / ( LOG( ( ZLB(I,K) + Z0M(I) + 0.5*DZL(I,K) ) / EXCOEF3A.192
& ( ZLB(I,K) + Z0M(I) - 0.5*DZL(I,KM1) ) ) EXCOEF3A.193
& + VKZ / LAMBDA_EFF ) EXCOEF3A.194
ELH = VKZ / ( LOG( ( ZLB(I,K) + Z0M(I) + 0.5*DZL(I,K) ) / EXCOEF3A.195
& ( ZLB(I,K) + Z0M(I) - 0.5*DZL(I,KM1) ) ) EXCOEF3A.196
& + VKZ / LAMBDAH ) EXCOEF3A.197
ELSE EXCOEF3A.198
VKZ = VKMAN * ( ZLB(I,K) + Z0M(I) ) EXCOEF3A.199
ELM = VKZ / (1.0 + VKZ/LAMBDA_EFF ) EXCOEF3A.200
ELH = VKZ / (1.0 + VKZ/LAMBDAH ) EXCOEF3A.201
ENDIF EXCOEF3A.202
C EXCOEF3A.203
C----------------------------------------------------------------------- EXCOEF3A.204
CL 2.2 Calculate air density RHO = P/(R*TV) for interface between EXCOEF3A.205
CL "current" and previous (lower) layers (i.e. at level K-1/2). EXCOEF3A.206
C----------------------------------------------------------------------- EXCOEF3A.207
C Repeat of KMKH calculation, could be passed in from KMKH. EXCOEF3A.208
WKM1 = 0.5 * DZL(I,KM1) * RDZ(I,K) EXCOEF3A.209
WK = 0.5 * DZL(I,K) * RDZ(I,K) EXCOEF3A.210
RHO = ! Calculate rho at K-1/2, from P243.111 :- EXCOEF3A.211
+ ( AKH(K) + BKH(K)*PSTAR(I) ) ! Pressure at K-1/2, P243.112 EXCOEF3A.212
+ / ! divided by ... EXCOEF3A.213
+ ( R * ! R times ... EXCOEF3A.214
+ ( TV(I,KM1)*WK + TV(I,K)*WKM1 ) ! TV at K-1/2, from P243.113 EXCOEF3A.215
+ ) EXCOEF3A.216
C EXCOEF3A.217
C----------------------------------------------------------------------- EXCOEF3A.218
CL 2.3 Calculate wind shear and magnitude of gradient thereof across EXCOEF3A.219
CL interface K-1/2. EXCOEF3A.220
C----------------------------------------------------------------------- EXCOEF3A.221
C Repeat of KMKH calculation, could be passed in from KMKH. EXCOEF3A.222
DZU = U_P(I,K) - U_P(I,KM1) EXCOEF3A.223
DZV = V_P(I,K) - V_P(I,KM1) EXCOEF3A.224
DVMOD2 = MAX ( 1.0E-12 , DZU*DZU + DZV*DZV ) EXCOEF3A.225
DVDZM = SQRT (DVMOD2) * RDZ(I,K) EXCOEF3A.226
C EXCOEF3A.227
C----------------------------------------------------------------------- EXCOEF3A.228
CL 2.4 Calculate (values of) stability functions FH, FM. EXCOEF3A.229
C----------------------------------------------------------------------- EXCOEF3A.230
C EXCOEF3A.231
IF (RI(I,K) .GE. 0.0) THEN EXCOEF3A.232
RTMRI = 0.0 EXCOEF3A.233
FM = 1.0 / ( 1.0 + G0*RI(I,K) ) EXCOEF3A.234
FH = FM EXCOEF3A.235
C !----------------------------------------------------------- EXCOEF3A.236
C ! If convective cloud exists in layer K allow neutral mixing EXCOEF3A.237
C ! of momentum between layers K-1 and K. This is to ensure EXCOEF3A.238
C ! that a reasonable amount of momentum is mixed in the EXCOEF3A.239
C ! presence of convection; it is not be required when ARN1F403.75
C ! momentum transport is included in the convection scheme. EXCOEF3A.241
C !----------------------------------------------------------- EXCOEF3A.242
IF ( .NOT.L_MOM .AND. (CCA(I) .GT. 0.0) .AND. ARN1F403.76
& (K .GE. CCB(I)) .AND. (K .LT. CCT(I)) ) ARN1F403.77
& FM = 1.0 EXCOEF3A.244
ELSE EXCOEF3A.245
RTMRI = (ELM/ELH) * SQRT ( -RI(I,K) ) EXCOEF3A.246
FM = 1.0 - ( G0*RI(I,K) / ( 1.0 + DM*RTMRI ) ) EXCOEF3A.247
FH = 1.0 - ( G0*RI(I,K) / ( 1.0 + DH*RTMRI ) ) EXCOEF3A.248
ENDIF EXCOEF3A.249
C EXCOEF3A.250
C----------------------------------------------------------------------- EXCOEF3A.251
CL 2.5 Calculate exchange coefficients RHO*KM, RHO*KH for interface EXCOEF3A.252
CL K-1/2. EXCOEF3A.253
C----------------------------------------------------------------------- EXCOEF3A.254
C EXCOEF3A.255
RHOKM(I,K) = RHO * ELM * ELM * DVDZM * FM EXCOEF3A.256
RHOKH(I,K) = RHO * ELH * ELM * DVDZM * FH EXCOEF3A.257
21 CONTINUE EXCOEF3A.258
2 CONTINUE EXCOEF3A.259
IF (LTIMER) THEN EXCOEF3A.260
CALL TIMER('EX_COEF ',104) GPB8F405.68
ENDIF EXCOEF3A.262
RETURN EXCOEF3A.263
END EXCOEF3A.264
*ENDIF EXCOEF3A.265