*IF DEF,A03_6A EXCOEF6A.2
C *****************************COPYRIGHT****************************** EXCOEF6A.3
C (c) CROWN COPYRIGHT 1997, METEOROLOGICAL OFFICE, All Rights Reserved. EXCOEF6A.4
C EXCOEF6A.5
C Use, duplication or disclosure of this code is subject to the EXCOEF6A.6
C restrictions as set forth in the contract. EXCOEF6A.7
C EXCOEF6A.8
C Meteorological Office EXCOEF6A.9
C London Road EXCOEF6A.10
C BRACKNELL EXCOEF6A.11
C Berkshire UK EXCOEF6A.12
C RG12 2SZ EXCOEF6A.13
C EXCOEF6A.14
C If no contract has been raised with this copy of the code, the use, EXCOEF6A.15
C duplication or disclosure of it is strictly prohibited. Permission EXCOEF6A.16
C to do so must first be obtained in writing from the Head of Numerical EXCOEF6A.17
C Modelling at the above address. EXCOEF6A.18
C ******************************COPYRIGHT****************************** EXCOEF6A.19
!!! SUBROUTINE EX_COEF------------------------------------------------ EXCOEF6A.20
!!! EXCOEF6A.21
!!! Purpose: To calculate exchange coefficients for boundary layer EXCOEF6A.22
!!! subroutine KMKH. EXCOEF6A.23
!!! EXCOEF6A.24
!!! Suitable for single column use (via *IF definition IBM). EXCOEF6A.25
!!! EXCOEF6A.26
!!! Model Modification history: EXCOEF6A.27
!!! version Date EXCOEF6A.28
!!! EXCOEF6A.29
!!! 4.4 10/9/97 New deck R.N.B.Smith EXCOEF6A.30
!LL 4.5 23/10/98 Prevent TIMER from performing barrier P.Burton GPB8F405.69
!!! EXCOEF6A.31
!!! Programming standard: EXCOEF6A.32
!!! EXCOEF6A.33
!!! System component covered: Part of P243. EXCOEF6A.34
!!! EXCOEF6A.35
!!! Project task: EXCOEF6A.36
!!! EXCOEF6A.37
!!! Documentation: UMDP No.24 EXCOEF6A.38
!!! EXCOEF6A.39
!!!--------------------------------------------------------------------- EXCOEF6A.40
EXCOEF6A.41
!! Arguments :- EXCOEF6A.42
SUBROUTINE EX_COEF ( 5,6EXCOEF6A.43
& P_FIELD,P1,P_POINTS,BL_LEVELS EXCOEF6A.44
&,CCB,CCT,NTML,L_MOM EXCOEF6A.45
&,CCA,DZL,RDZ,DB,U_P,V_P EXCOEF6A.46
&,RHO,ZH,ZLB,Z0M,H_BLEND EXCOEF6A.47
&,CUMULUS,Z_LCL ARN0F405.727
&,RHOKM,RHOKH,LTIMER EXCOEF6A.48
& ) EXCOEF6A.49
EXCOEF6A.50
IMPLICIT NONE EXCOEF6A.51
EXCOEF6A.52
LOGICAL LTIMER EXCOEF6A.53
LOGICAL EXCOEF6A.54
& L_MOM ! IN Switch for convective momentum transport. EXCOEF6A.55
EXCOEF6A.56
INTEGER EXCOEF6A.57
& P_FIELD ! IN No. of P-grid points in whole field. EXCOEF6A.58
&,P1 ! IN First P-grid point to be processed. EXCOEF6A.59
&,P_POINTS ! IN No. of P-grid points to be processed. EXCOEF6A.60
&,BL_LEVELS ! IN maximum number of boundary layer levels EXCOEF6A.61
EXCOEF6A.62
INTEGER EXCOEF6A.63
& CCB(P_FIELD) ! IN Convective Cloud Base. EXCOEF6A.64
&,CCT(P_FIELD) ! IN Convective Cloud Top. EXCOEF6A.65
&,NTML(P_FIELD) ! IN Number of turbulently mixed EXCOEF6A.66
! layers. EXCOEF6A.67
EXCOEF6A.68
REAL EXCOEF6A.69
& CCA(P_FIELD) ! IN Convective Cloud Amount. EXCOEF6A.70
&,DZL(P_FIELD,BL_LEVELS) ! IN Layer depths (m). DZL(,K) is the EXCOEF6A.71
! distance from layer boundary K-1/2 EXCOEF6A.72
! to layer boundary K+1/2. For K=1 EXCOEF6A.73
! the lower boundary is the surface. EXCOEF6A.74
&,RDZ(P_FIELD,BL_LEVELS) ! IN Reciprocal of distance between EXCOEF6A.75
! hybrid levels (m-1). 1/RDZ(,K) is EXCOEF6A.76
! the vertical distance from level EXCOEF6A.77
! K-1 to level K, except that for EXCOEF6A.78
! K=1 it is just the height of the EXCOEF6A.79
! lowest atmospheric level. EXCOEF6A.80
&,RHO(P_FIELD,BL_LEVELS) ! IN Density at levels where Richardson EXCOEF6A.81
! number is calculated. EXCOEF6A.82
&,DB(P_FIELD,2:BL_LEVELS) ! IN Buoyancy jump across lower EXCOEF6A.83
! interface of layer. EXCOEF6A.84
&,U_P(P_FIELD,BL_LEVELS) ! IN Westerly wind component EXCOEF6A.85
! horizontally interpolated to EXCOEF6A.86
! P-grid. (m/s) EXCOEF6A.87
&,V_P(P_FIELD,BL_LEVELS) ! IN Southerly wind component EXCOEF6A.88
! horizontally interpolated to EXCOEF6A.89
! P-grid. (m/s) EXCOEF6A.90
&,ZLB(P_FIELD,BL_LEVELS) ! IN ZLB(,K) is height above surface of EXCOEF6A.92
! lower boundary of layer K (m). EXCOEF6A.93
&,Z0M(P_FIELD) ! IN Roughness length for momentum (m). EXCOEF6A.94
&,H_BLEND(P_FIELD) ! IN Blending height for effective EXCOEF6A.95
! roughness length scheme EXCOEF6A.96
&,Z_LCL(P_FIELD) ! IN Height of lifting condensation ARN0F405.728
! level (m). ARN0F405.729
LOGICAL ARN0F405.730
& CUMULUS(P_FIELD) ! IN Flag for boundary layer cumulus. ARN0F405.731
EXCOEF6A.97
REAL EXCOEF6A.98
& RHOKM(P_FIELD,BL_LEVELS) ! INOUT Layer K-1 - to - layer K EXCOEF6A.99
! exchange coefficient for EXCOEF6A.100
! momentum, on UV-grid with first EXCOEF6A.101
! and last rows set to "missing EXCOEF6A.102
! data". EXCOEF6A.103
&,RHOKH(P_FIELD,BL_LEVELS) ! INOUT Layer K-1 - to - layer K EXCOEF6A.104
! exchange coefficient for FTL, EXCOEF6A.105
! on P-grid. EXCOEF6A.106
&,ZH(P_FIELD) ! INOUT Mixing layer height (m). ARN0F405.732
EXCOEF6A.107
!----------------------------------------------------------------------- EXCOEF6A.108
EXTERNAL TIMER EXCOEF6A.109
EXCOEF6A.110
!----------------------------------------------------------------------- EXCOEF6A.111
EXCOEF6A.112
! Local and other symbolic constants :- EXCOEF6A.113
*CALL C_LHEAT
EXCOEF6A.114
*CALL C_R_CP
EXCOEF6A.115
*CALL C_EPSLON
EXCOEF6A.116
*CALL C_VKMAN
EXCOEF6A.117
EXCOEF6A.118
REAL EH,EM,G0,DH,DM,LAMBDA_MIN,A_LAMBDA EXCOEF6A.119
PARAMETER ( EXCOEF6A.120
& EH=25.0 ! Used in calc of stability function FH. EXCOEF6A.121
&,EM=4.0 ! Used in calc of stability function FM. EXCOEF6A.122
&,G0=10.0 ! Used in stability function calcs. EXCOEF6A.123
&,DH=G0/EH ! Used in calc of stability function FH. EXCOEF6A.124
&,DM=G0/EM ! Used in calc of stability function FM. EXCOEF6A.125
&,LAMBDA_MIN=40.0 ! Minimum value of length scale LAMBDA. EXCOEF6A.126
&,A_LAMBDA=2.0 ! used in calc of LAMBDA_EFF EXCOEF6A.127
&) EXCOEF6A.128
EXCOEF6A.129
EXCOEF6A.130
! Define local storage. EXCOEF6A.131
EXCOEF6A.132
! Arrays. ARN0F405.733
ARN0F405.734
REAL ARN0F405.735
& DVDZM(P_FIELD,2:BL_LEVELS) ! Modulus of wind shear. ARN0F405.736
&,RI(P_FIELD,2:BL_LEVELS) ! Local Richardson number. ARN0F405.737
&,ZH_LOCAL(P_FIELD) ! Mixed layer depth determined ARN0F405.738
! from local Richardson number ARN0F405.739
! profile. ARN0F405.740
ARN0F405.741
INTEGER ARN0F405.742
& NTML_LOCAL(P_FIELD) ! Number of model layers in the ARN0F405.743
! turbulently mixed layer as ARN0F405.744
! determined from the local ARN0F405.745
! Richardson number profile. ARN0F405.746
ARN0F405.747
LOGICAL ARN0F405.748
& TOPBL(P_FIELD) ! Flag for having reached the top ARN0F405.749
! of the turbulently mixed layer. ARN0F405.750
ARN0F405.751
ARN0F405.752
! Scalars. EXCOEF6A.133
EXCOEF6A.134
REAL EXCOEF6A.135
& DZU ! Westerly wind shear between adjacent levels. ARN0F405.753
&,DZV ! Southerly wind shear between adjacent levels. EXCOEF6A.138
&,DVMOD2 ! Square of modulus of wind shear between adjacent EXCOEF6A.139
! levels EXCOEF6A.140
&,ELH ! Mixing length for heat & moisture at lower layer bdy. EXCOEF6A.141
&,ELM ! Mixing length for momentum at lower layer boundary. EXCOEF6A.142
&,FH ! (Value of) stability function for heat & moisture. EXCOEF6A.143
&,FM ! (Value of) stability function for momentum transport. EXCOEF6A.144
&,RTMRI ! Temporary in stability function calculation. EXCOEF6A.146
&,VKZ ! Temporary in calculation of ELH. EXCOEF6A.147
&,LAMBDAM ! Asymptotic mixing length for turbulent transport EXCOEF6A.148
! of momentum. EXCOEF6A.149
&,LAMBDAH ! Asymptotic mixing length for turbulent transport EXCOEF6A.150
! of heat/moisture. EXCOEF6A.151
&,LAMBDA_EFF! Effective mixing length used with effective EXCOEF6A.152
! roughness length scheme. EXCOEF6A.153
INTEGER EXCOEF6A.154
& I ! Loop counter (horizontal field index). EXCOEF6A.155
&,K ! Loop counter (vertical level index). EXCOEF6A.156
&,KM1 ! K-1. EXCOEF6A.157
&,MBL ! Maximum number of model layers below mixed layer top. ARN0F405.754
EXCOEF6A.158
! Layer interface K_LOG_LAYR-1/2 is the highest which requires log EXCOEF6A.159
! profile correction factors to the vertical finite differences. EXCOEF6A.160
! The value should be reassessed if the vertical resolution is changed. EXCOEF6A.161
! We could set K_LOG_LAYR = BL_LEVELS and thus apply the correction EXCOEF6A.162
! factors for all the interfaces treated by the boundary layer scheme; EXCOEF6A.163
! this would be desirable theoretically but expensive computationally EXCOEF6A.164
! because of the use of the log function. EXCOEF6A.165
EXCOEF6A.166
INTEGER K_LOG_LAYR EXCOEF6A.167
PARAMETER (K_LOG_LAYR=2) EXCOEF6A.168
EXCOEF6A.169
IF (LTIMER) THEN EXCOEF6A.170
CALL TIMER
('EX_COEF ',103) GPB8F405.70
ENDIF EXCOEF6A.172
EXCOEF6A.173
! Set MBL, "maximum number of boundary layer levels" for the purposes ARN0F405.755
! of mixed layer height calculation. ARN0F405.756
ARN0F405.757
MBL = BL_LEVELS - 1 ARN0F405.758
ARN0F405.759
!----------------------------------------------------------------------- EXCOEF6A.174
!! 0. Initialise flag for having reached top of turbulently mixed layer ARN0F405.760
!! and also the number of turbulently mixed layers. ARN0F405.761
!----------------------------------------------------------------------- ARN0F405.762
ARN0F405.763
DO I=P1,P1+P_POINTS-1 ARN0F405.764
TOPBL(I) = .FALSE. ARN0F405.765
NTML_LOCAL(I) = 1 ARN0F405.766
ENDDO ARN0F405.767
ARN0F405.768
! ARN0F405.769
!----------------------------------------------------------------------- ARN0F405.770
! 1.1 Loop over levels calculating local wind shears and Richardson ARN0F405.771
! numbers. ARN0F405.772
!----------------------------------------------------------------------- ARN0F405.773
! ARN0F405.774
DO K=2,BL_LEVELS ARN0F405.775
KM1 = K-1 ARN0F405.776
DO I=P1,P1+P_POINTS-1 ARN0F405.777
ARN0F405.778
DZU = U_P(I,K) - U_P(I,KM1) ARN0F405.779
DZV = V_P(I,K) - V_P(I,KM1) ARN0F405.780
DVMOD2 = MAX ( 1.0E-12 , DZU*DZU + DZV*DZV ) ARN0F405.781
DVDZM(I,K) = SQRT (DVMOD2) * RDZ(I,K) ARN0F405.782
ARN0F405.783
RI(I,K) = DB(I,K) / ( RDZ(I,K) * DVMOD2 ) ARN0F405.784
ARN0F405.785
!----------------------------------------------------------------------- ARN0F405.786
!! 1.2 If either a stable layer (Ri > 1) or the maximum boundary layer ARN0F405.787
!! height has been reached, set boundary layer height (ZH) to ARN0F405.788
!! the height of the lower boundary of the current layer ARN0F405.789
!----------------------------------------------------------------------- ARN0F405.790
ARN0F405.791
IF ( .NOT.TOPBL(I) .AND. ( RI(I,K).GT.1.0 .OR. K.GT.MBL ARN0F405.792
& .OR. (CUMULUS(I) .AND. ZLB(I,K).GE.Z_LCL(I)) ) ) THEN ARN0F405.793
TOPBL(I) = .TRUE. ARN0F405.794
ZH_LOCAL(I) = ZLB(I,K) ARN0F405.795
NTML_LOCAL(I) = K-1 ARN0F405.796
ENDIF ARN0F405.797
ENDDO ! Loop over points ARN0F405.798
ENDDO ! Loop over levels ARN0F405.799
! ARN0F405.800
!----------------------------------------------------------------------- ARN0F405.801
!! 2. Loop round "boundary" levels; calculate the stability ARN0F405.802
!! dependent turbulent mixing coefficients. EXCOEF6A.176
!----------------------------------------------------------------------- EXCOEF6A.177
EXCOEF6A.178
DO K=2,BL_LEVELS EXCOEF6A.179
KM1 = K-1 EXCOEF6A.180
DO I=P1,P1+P_POINTS-1 EXCOEF6A.181
!----------------------------------------------------------------------- EXCOEF6A.182
!! 2.1 Calculate asymptotic mixing lengths LAMBDAM and LAMBDAH EXCOEF6A.183
!! (currently assumed equal). EXCOEF6A.184
!----------------------------------------------------------------------- EXCOEF6A.185
EXCOEF6A.186
LAMBDAM = MAX ( LAMBDA_MIN , 0.15*ZH_LOCAL(I) ) ARN0F405.803
LAMBDAH = LAMBDAM EXCOEF6A.188
LAMBDA_EFF = MAX (LAMBDAM, A_LAMBDA*H_BLEND(I) ) EXCOEF6A.189
IF ( K .GE. NTML_LOCAL(I) + 1) THEN ARN0F405.804
LAMBDAM = LAMBDA_MIN EXCOEF6A.191
LAMBDAH = LAMBDA_MIN EXCOEF6A.192
IF (ZLB(I,K) .GT. A_LAMBDA*H_BLEND(I)) LAMBDA_EFF = LAMBDA_MIN EXCOEF6A.193
ENDIF EXCOEF6A.194
EXCOEF6A.195
!----------------------------------------------------------------------- EXCOEF6A.196
!! 2.2 Calculate mixing lengths ELH, ELM at layer interface K-1/2. EXCOEF6A.197
!----------------------------------------------------------------------- EXCOEF6A.198
EXCOEF6A.199
! Incorporate log profile corrections to the vertical finite EXCOEF6A.200
! differences into the definitions of ELM and ELH. EXCOEF6A.201
! To save computing logarithms for all K, the values of ELM and ELH EXCOEF6A.202
! are unchanged for K > K_LOG_LAYR. EXCOEF6A.203
EXCOEF6A.204
IF (K .LE. K_LOG_LAYR) THEN EXCOEF6A.205
VKZ = VKMAN / RDZ(I,K) EXCOEF6A.206
ELM = VKZ / ( LOG( ( ZLB(I,K) + Z0M(I) + 0.5*DZL(I,K) ) / EXCOEF6A.207
& ( ZLB(I,K) + Z0M(I) - 0.5*DZL(I,KM1) ) ) EXCOEF6A.208
& + VKZ / LAMBDA_EFF ) EXCOEF6A.209
ELH = VKZ / ( LOG( ( ZLB(I,K) + Z0M(I) + 0.5*DZL(I,K) ) / EXCOEF6A.210
& ( ZLB(I,K) + Z0M(I) - 0.5*DZL(I,KM1) ) ) EXCOEF6A.211
& + VKZ / LAMBDAH ) EXCOEF6A.212
ELSE EXCOEF6A.213
VKZ = VKMAN * ( ZLB(I,K) + Z0M(I) ) EXCOEF6A.214
ELM = VKZ / (1.0 + VKZ/LAMBDA_EFF ) EXCOEF6A.215
ELH = VKZ / (1.0 + VKZ/LAMBDAH ) EXCOEF6A.216
ENDIF EXCOEF6A.217
EXCOEF6A.218
!----------------------------------------------------------------------- EXCOEF6A.219
!! 2.4 Calculate stability functions FH, FM. ARN0F405.805
!----------------------------------------------------------------------- EXCOEF6A.222
EXCOEF6A.223
IF (RI(I,K) .GE. 0.0) THEN ARN0F405.806
RTMRI = 0.0 EXCOEF6A.236
FM = 1.0 / ( 1.0 + G0 * RI(I,K) ) ARN0F405.807
FH = 1.0 / ( 1.0 + G0 * RI(I,K) ) ARN0F405.808
! ARN0F405.809
! !----------------------------------------------------------- EXCOEF6A.240
! ! If convective cloud exists in layer K allow neutral mixing EXCOEF6A.241
! ! of momentum between layers K-1 and K. This is to ensure EXCOEF6A.242
! ! that a reasonable amount of momentum is mixed in the EXCOEF6A.243
! ! presence of convection; it is not be required when EXCOEF6A.244
! ! momentum transport is included in the convection scheme. EXCOEF6A.245
! !----------------------------------------------------------- EXCOEF6A.246
! ARN0F405.810
IF ( .NOT.L_MOM .AND. (CCA(I) .GT. 0.0) .AND. EXCOEF6A.248
& (K .GE. CCB(I)) .AND. (K .LT. CCT(I)) ) EXCOEF6A.249
& FM = 1.0 EXCOEF6A.250
ELSE EXCOEF6A.251
RTMRI = (ELM/ELH) * SQRT ( -RI(I,K) ) ARN0F405.811
FM = 1.0 - ( G0*RI(I,K) / ( 1.0 + DM*RTMRI ) ) ARN0F405.812
FH = 1.0 - ( G0*RI(I,K) / ( 1.0 + DH*RTMRI ) ) ARN0F405.813
ENDIF EXCOEF6A.255
EXCOEF6A.256
!----------------------------------------------------------------------- EXCOEF6A.257
!! 2.5 Calculate exchange coefficients RHO*KM, RHO*KH for interface EXCOEF6A.258
!! K-1/2. EXCOEF6A.259
!----------------------------------------------------------------------- EXCOEF6A.260
EXCOEF6A.261
RHOKM(I,K) = RHO(I,K) * ELM * ELM * DVDZM(I,K) * FM ARN0F405.814
RHOKH(I,K) = RHO(I,K) * ELH * ELM * DVDZM(I,K) * FH ARN0F405.815
EXCOEF6A.264
ENDDO ! p_points EXCOEF6A.265
ENDDO ! bl_levels EXCOEF6A.266
ARN0F405.816
DO I = P1,P1+P_POINTS-1 ARN0F405.817
IF ( NTML_LOCAL(I) .GT. NTML(I)+1 .OR. ARN0F405.818
& ( NTML(I) .EQ. 1 .AND. NTML_LOCAL(I) .GT. NTML(I) ) ) THEN ARN0F405.819
ZH(I) = ZH_LOCAL(I) ARN0F405.820
NTML(I) = NTML_LOCAL(I) ARN0F405.821
ENDIF ARN0F405.822
ENDDO ARN0F405.823
EXCOEF6A.267
IF (LTIMER) THEN EXCOEF6A.268
CALL TIMER
('EX_COEF ',104) GPB8F405.71
ENDIF EXCOEF6A.270
EXCOEF6A.271
RETURN EXCOEF6A.272
END EXCOEF6A.273
*ENDIF EXCOEF6A.274