*IF DEF,A03_6A KMKHZ6A.2
C *****************************COPYRIGHT****************************** KMKHZ6A.3
C (c) CROWN COPYRIGHT 1997, METEOROLOGICAL OFFICE, All Rights Reserved. KMKHZ6A.4
C KMKHZ6A.5
C Use, duplication or disclosure of this code is subject to the KMKHZ6A.6
C restrictions as set forth in the contract. KMKHZ6A.7
C KMKHZ6A.8
C Meteorological Office KMKHZ6A.9
C London Road KMKHZ6A.10
C BRACKNELL KMKHZ6A.11
C Berkshire UK KMKHZ6A.12
C RG12 2SZ KMKHZ6A.13
C KMKHZ6A.14
C If no contract has been raised with this copy of the code, the use, KMKHZ6A.15
C duplication or disclosure of it is strictly prohibited. Permission KMKHZ6A.16
C to do so must first be obtained in writing from the Head of Numerical KMKHZ6A.17
C Modelling at the above address. KMKHZ6A.18
C ******************************COPYRIGHT****************************** KMKHZ6A.19
! SUBROUTINE KMKHZ -------------------------------------------------- KMKHZ6A.20
!!! KMKHZ6A.21
!!! Purpose: To calculate the turbulent mixing coefficients KM and KH KMKHZ6A.22
!!! KMKHZ6A.23
!!! KMKHZ6A.25
!!! Model Modification history KMKHZ6A.26
!!! version date KMKHZ6A.27
!!! KMKHZ6A.28
!!! 4.4 05/09/97 New deck R.N.B.Smith KMKHZ6A.29
!!! 4.5 Jul. 98 Kill the IBM specific lines (JCThil) AJC1F405.131
!!! KMKHZ6A.30
!!! Programming standard: KMKHZ6A.31
!!! KMKHZ6A.32
!!! System component covered: Part of P243. KMKHZ6A.33
!!! KMKHZ6A.34
!!! Documentation: UMDP No.24 KMKHZ6A.35
!!! KMKHZ6A.36
!!!--------------------------------------------------------------------- KMKHZ6A.37
KMKHZ6A.38
! Arguments :- KMKHZ6A.39
KMKHZ6A.40
SUBROUTINE KMKHZ ( 1,6KMKHZ6A.41
& P_FIELD,P1,P_POINTS,BL_LEVELS KMKHZ6A.42
&,P,P_HALF,T,Q,QCL,QCF,BT,BQ,CF,DZL KMKHZ6A.43
&,RDZ,DELTAP,FTL,FQW KMKHZ6A.44
&,Z0M,Z_FULL,Z_HALF,Z_UV,Z_TQ,V_S,FB_SURF KMKHZ6A.45
&,QW,RHOKM,DB_SVL,RHOKH,TL,ZH,TV1_SD,T1_SD,Q1_SD,NTML ARN0F405.875
&,GRAD_T_ADJ,GRAD_Q_ADJ KMKHZ6A.47
&,BTM,BQM,DQSDT,BTM_CLD,BQM_CLD,A_QSM,A_DQSDTM,RHO_TQ,RHO_UV KMKHZ6A.48
&,RAD_HR,RADHR_DIM1,CUMULUS,Z_LCL,RHOKM_TOP,RHOKH_TOP,ZHT ARN0F405.876
&,BL_TYPE_1,BL_TYPE_2,BL_TYPE_3,BL_TYPE_4,BL_TYPE_5,BL_TYPE_6 ARN0F405.877
&,UNSTABLE,NTDSC,DSC ARN0F405.878
&,LTIMER KMKHZ6A.50
& ) KMKHZ6A.51
KMKHZ6A.52
IMPLICIT NONE KMKHZ6A.53
KMKHZ6A.54
LOGICAL LTIMER ! IN Flag for TIMER diagnostics KMKHZ6A.55
KMKHZ6A.56
INTEGER KMKHZ6A.57
& P_FIELD ! IN No. of P-grid points in whole field KMKHZ6A.58
&,P1 ! IN First P-grid point to be processed. KMKHZ6A.59
&,P_POINTS ! IN No. of P-grid points to be KMKHZ6A.60
! processed. KMKHZ6A.61
! ( = P_POINTS = 1 for SCM.) KMKHZ6A.63
&,BL_LEVELS ! IN No. of atmospheric levels for KMKHZ6A.65
! which boundary layer fluxes are KMKHZ6A.66
! calculated. KMKHZ6A.67
&,RADHR_DIM1 ! IN Dimension of Radiative heating KMKHZ6A.68
! ! rate (P_FIELD but used for KMKHZ6A.69
! ! dynamic allocation) KMKHZ6A.70
KMKHZ6A.71
REAL KMKHZ6A.72
& BQ(P_FIELD,BL_LEVELS) ! IN A buoyancy parameter for clear air KMKHZ6A.73
! ! on p,T,q-levels (full levels). KMKHZ6A.74
&,BT(P_FIELD,BL_LEVELS) ! IN A buoyancy parameter for clear air KMKHZ6A.75
! ! on p,T,q-levels (full levels). KMKHZ6A.76
&,BQM(P_FIELD,BL_LEVELS) ! IN A buoyancy parameter for clear air KMKHZ6A.77
! ! on intermediate levels (half levels): KMKHZ6A.78
! ! (*,K) elements are k+1/2 values. KMKHZ6A.79
&,BTM(P_FIELD,BL_LEVELS) ! IN A buoyancy parameter for clear air KMKHZ6A.80
! ! on intermediate levels (half levels): KMKHZ6A.81
! ! (*,K) elements are k+1/2 values. KMKHZ6A.82
&,BQM_CLD(P_FIELD,BL_LEVELS) KMKHZ6A.83
! ! IN A buoyancy parameter for cloudy air KMKHZ6A.84
! ! on intermediate levels (half levels): KMKHZ6A.85
! ! (*,K) elements are k+1/2 values. KMKHZ6A.86
&,BTM_CLD(P_FIELD,BL_LEVELS) KMKHZ6A.87
! ! IN A buoyancy parameter for cloudy air KMKHZ6A.88
! ! on intermediate levels (half levels): KMKHZ6A.89
! ! (*,K) elements are k+1/2 values. KMKHZ6A.90
&,A_QSM(P_FIELD,BL_LEVELS) KMKHZ6A.91
! ! IN Saturated lapse rate factor KMKHZ6A.92
! ! on intermediate levels (half levels): KMKHZ6A.93
! ! (*,K) elements are k+1/2 values. KMKHZ6A.94
&,A_DQSDTM(P_FIELD,BL_LEVELS) KMKHZ6A.95
! ! IN Saturated lapse rate factor KMKHZ6A.96
! ! on intermediate levels (half levels): KMKHZ6A.97
! ! (*,K) elements are k+1/2 values. KMKHZ6A.98
&,P(P_FIELD,BL_LEVELS) ! IN P(*,K) is pressure at full level k. KMKHZ6A.99
&,P_HALF(P_FIELD,BL_LEVELS) KMKHZ6A.100
! ! IN P_HALF(*,K) is pressure at half KMKHZ6A.101
! ! level k-1/2. KMKHZ6A.102
&,T(P_FIELD,BL_LEVELS) ! IN Temperature (K). KMKHZ6A.103
&,Q(P_FIELD,BL_LEVELS) ! IN Sp humidity (kg water per kg air). KMKHZ6A.104
&,QCF(P_FIELD,BL_LEVELS) ! IN Cloud ice (kg per kg air). KMKHZ6A.105
&,QCL(P_FIELD,BL_LEVELS) ! IN Cloud liquid water (kg per kg air). KMKHZ6A.106
&,CF(P_FIELD,BL_LEVELS) ! IN Cloud fractions for boundary levs. KMKHZ6A.107
&,DZL(P_FIELD,BL_LEVELS) ! IN Layer depths (m). DZL(,K) is the KMKHZ6A.108
! distance from layer boundary K-1/2 KMKHZ6A.109
! to layer boundary K+1/2. For K=1 KMKHZ6A.110
! the lower boundary is the surface. KMKHZ6A.111
KMKHZ6A.112
REAL KMKHZ6A.113
& RDZ(P_FIELD,BL_LEVELS) ! IN Reciprocal of distance between KMKHZ6A.114
! full levels (m-1). 1/RDZ(,K) is KMKHZ6A.115
! the vertical distance from level KMKHZ6A.116
! K-1 to level K, except that for KMKHZ6A.117
! K=1 it is the height of the KMKHZ6A.118
! lowest atmospheric full level. KMKHZ6A.119
&,Z0M(P_FIELD) ! IN Roughness length for momentum (m). KMKHZ6A.120
&,Z_FULL(P_FIELD,BL_LEVELS) ! IN Z_FULL(*,K) is the height of the KMKHZ6A.121
! ! k-th full level above the surface. KMKHZ6A.122
&,Z_HALF(P_FIELD,BL_LEVELS) ! IN Z_HALF(*,K) is the height of level KMKHZ6A.123
! k-1/2 above the surface (m). KMKHZ6A.124
&,DELTAP(P_FIELD,BL_LEVELS) ! IN Pressure difference between KMKHZ6A.125
! ! half-levels (Pa). KMKHZ6A.126
&,RAD_HR(RADHR_DIM1,BL_LEVELS) KMKHZ6A.127
! ! IN Radiative heating rate KMKHZ6A.128
! ! (Kelvins/s) KMKHZ6A.129
&,V_S(P_FIELD) ! IN Surface friction velocity (m/s) KMKHZ6A.130
&,FB_SURF(P_FIELD) ! IN Surface buoyancy flux over density KMKHZ6A.131
! (m^2/s^3). KMKHZ6A.132
&,FQW(P_FIELD,BL_LEVELS) ! IN "Explicit" flux of QW (i.e. KMKHZ6A.133
! evaporation) from layer below KMKHZ6A.134
! on P-grid (kg per sq m per s). KMKHZ6A.135
&,FTL(P_FIELD,BL_LEVELS) ! IN "Explicit" flux of TL = H/CP KMKHZ6A.136
! (sensible heat/CP) from layer KMKHZ6A.137
! below, on P-grid. KMKHZ6A.138
&,TV1_SD(P_FIELD) ! IN Standard Deviation of level 1 KMKHZ6A.139
! ! virtual temperature (K). KMKHZ6A.140
&,T1_SD(P_FIELD) ! IN Standard Deviation of level 1 KMKHZ6A.141
! ! temperature (K). KMKHZ6A.142
&,Q1_SD(P_FIELD) ! IN Standard Deviation of level 1 KMKHZ6A.143
! ! specific humidity (kg/kg). KMKHZ6A.144
&,DQSDT(P_FIELD,BL_LEVELS) ! IN Partial derivative of QSAT w.r.t. KMKHZ6A.145
! ! temperature. KMKHZ6A.146
&,Z_UV(P_FIELD,BL_LEVELS) ! IN For a vertically staggered grid KMKHZ6A.147
! ! with a u,v-level first above the KMKHZ6A.148
! ! surface, Z_UV(*,K) is the height KMKHZ6A.149
! ! of the k-th u,v-level (half level KMKHZ6A.150
! ! k-1/2) above the surface; KMKHZ6A.151
! ! for an unstaggered grid the KMKHZ6A.152
! ! heights of the half-levels KMKHZ6A.153
! ! 0.5 to BL_LEVELS-0.5 should be KMKHZ6A.154
! ! input to elements 1 to BL_LEVELS. KMKHZ6A.155
! ! (1st value not used in either KMKHZ6A.156
! ! case.) KMKHZ6A.157
&,Z_TQ(P_FIELD,BL_LEVELS) ! IN For a vertically staggered grid KMKHZ6A.158
! ! with a u,v-level first above the KMKHZ6A.159
! ! surface, Z_TQ(*,K) is the height KMKHZ6A.160
! ! of the k-th T,q-level (full level KMKHZ6A.161
! ! k) above the surface; KMKHZ6A.162
! ! for an unstaggered grid the KMKHZ6A.163
! ! heights of the half levels KMKHZ6A.164
! ! 1.5 to BL_LEVELS+0.5 should be KMKHZ6A.165
! ! input to elements 1 to BL_LEVELS. KMKHZ6A.166
! ! (value for BL_LEVELS not used KMKHZ6A.167
! ! in either case.) KMKHZ6A.168
! KMKHZ6A.169
&,RHO_UV(P_FIELD,BL_LEVELS) ! IN For a vertically staggered grid KMKHZ6A.170
! ! with a u,v-level first above the KMKHZ6A.171
! ! surface, RHO_UV(*,K) is the KMKHZ6A.172
! ! density at the k-th u,v-level KMKHZ6A.173
! ! above the surface; KMKHZ6A.174
! ! for an unstaggered grid the KMKHZ6A.175
! ! densities at the layer interfaces KMKHZ6A.176
! ! (half-levels) 0.5 to BL_LEVELS-0.5 KMKHZ6A.177
! ! should be input to elements 1 to KMKHZ6A.178
! ! BL_LEVELS. KMKHZ6A.179
! ! (1st value not used in either KMKHZ6A.180
! ! case.) KMKHZ6A.181
&,RHO_TQ(P_FIELD,BL_LEVELS) ! IN For a vertically staggered grid KMKHZ6A.182
! ! with a u,v-level first above the KMKHZ6A.183
! ! surface, RHO_TQ(*,K) is the KMKHZ6A.184
! ! density of the k-th T,q-level KMKHZ6A.185
! ! above the surface; KMKHZ6A.186
! ! for an unstaggered grid the KMKHZ6A.187
! ! densities at the layer interfaces KMKHZ6A.188
! ! (half-levels) 1.5 to BL_LEVELS+0.5 KMKHZ6A.189
! ! should be input to elements 1 to KMKHZ6A.190
! ! BL_LEVELS. KMKHZ6A.191
! ! (value for BL_LEVELS not used KMKHZ6A.192
! ! in either case.) KMKHZ6A.193
! KMKHZ6A.194
&,QW(P_FIELD,BL_LEVELS) ! IN Total water content (kg per kg KMKHZ6A.195
! air). KMKHZ6A.196
&,TL(P_FIELD,BL_LEVELS) ! IN Liquid/frozen water temperature KMKHZ6A.197
! (K). KMKHZ6A.198
REAL KMKHZ6A.199
& ZH(P_FIELD) ! INOUT Height of the top of the surface ARN0F405.879
! ! based turbulently mixed layer (m). ARN0F405.880
KMKHZ6A.208
REAL KMKHZ6A.209
& DB_SVL(P_FIELD,2:BL_LEVELS) ARN0F405.881
! ! OUT Buoyancy jump across layer ARN0F405.882
! interface based on SVL. ARN0F405.883
&,GRAD_T_ADJ(P_FIELD) ! OUT Temperature gradient adjustment KMKHZ6A.213
! for non-local mixing in unstable KMKHZ6A.214
! turbulent boundary layer. KMKHZ6A.215
&,GRAD_Q_ADJ(P_FIELD) ! OUT Humidity gradient adjustment KMKHZ6A.216
! for non-local mixing in unstable KMKHZ6A.217
! turbulent boundary layer. KMKHZ6A.218
&,Z_LCL(P_FIELD) ! OUT Height of lifting condensation ARN0F405.884
! level. ARN0F405.885
&,RHOKM(P_FIELD,2:BL_LEVELS) ARN0F405.886
! ! OUT Non-local turbulent mixing ARN0F405.887
! ! coefficient for momentum. ARN0F405.888
&,RHOKH(P_FIELD,2:BL_LEVELS) ARN0F405.889
! ! OUT Non-local turbulent mixing ARN0F405.890
! ! coefficient for scalars. ARN0F405.891
&,RHOKM_TOP(P_FIELD,2:BL_LEVELS) ARN0F405.892
! ! OUT Top-down turbulent mixing ARN0F405.893
! ! coefficient for momentum. ARN0F405.894
&,RHOKH_TOP(P_FIELD,2:BL_LEVELS) ARN0F405.895
! ! OUT Top-down turbulent mixing ARN0F405.896
! ! coefficient for scalars. ARN0F405.897
&,ZHT(P_FIELD) ! OUT Height below which there may be ARN0F405.898
! ! turbulent mixing (m). ARN0F405.899
&,BL_TYPE_1(P_FIELD) ! OUT Indicator set to 1.0 if stable ARN0F405.900
! ! b.l. diagnosed, 0.0 otherwise. ARN0F405.901
&,BL_TYPE_2(P_FIELD) ! OUT Indicator set to 1.0 if Sc over ARN0F405.902
! ! stable surface layer diagnosed, ARN0F405.903
! ! 0.0 otherwise. ARN0F405.904
&,BL_TYPE_3(P_FIELD) ! OUT Indicator set to 1.0 if well mixed ARN0F405.905
! ! b.l. diagnosed, 0.0 otherwise. ARN0F405.906
&,BL_TYPE_4(P_FIELD) ! OUT Indicator set to 1.0 if decoupled ARN0F405.907
! ! Sc layer (not over cumulus) ARN0F405.908
! ! diagnosed, 0.0 otherwise. ARN0F405.909
&,BL_TYPE_5(P_FIELD) ! OUT Indicator set to 1.0 if decoupled ARN0F405.910
! ! Sc layer over cumulus diagnosed, ARN0F405.911
! ! 0.0 otherwise. ARN0F405.912
&,BL_TYPE_6(P_FIELD) ! OUT Indicator set to 1.0 if a cumulus ARN0F405.913
! ! capped b.l. diagnosed, ARN0F405.914
! ! 0.0 otherwise. ARN0F405.915
ARN0F405.916
INTEGER KMKHZ6A.219
& NTML(P_FIELD) ! OUT Number of model levels in the KMKHZ6A.220
! turbulently mixed layer. KMKHZ6A.221
&,NTDSC(P_FIELD) ! OUT Top level for turbulent mixing in ARN0F405.917
! ! cloud layer. ARN0F405.918
KMKHZ6A.222
LOGICAL ARN0F405.919
& CUMULUS(P_FIELD) ! OUT Flag for cumulus in the b.l. ARN0F405.920
&,UNSTABLE(P_FIELD) ! OUT Flag to indicate an unstable ARN0F405.921
! boundary layer forced from ARN0F405.922
! the surface. ARN0F405.923
&,DSC(P_FIELD) ! OUT Flag set if decoupled stratocumulus ARN0F405.924
! ! layer found. ARN0F405.925
KMKHZ6A.223
ARN0F405.926
!!---------------------------------------------------------------------- KMKHZ6A.224
! External references :- KMKHZ6A.225
EXTERNAL TIMER, QSAT, EXCF_NL KMKHZ6A.226
KMKHZ6A.227
!!---------------------------------------------------------------------- KMKHZ6A.228
! Local and other symbolic constants :- KMKHZ6A.229
*CALL C_LHEAT 
KMKHZ6A.230
*CALL C_R_CP KMKHZ6A.231
*CALL C_G KMKHZ6A.232
*CALL C_0_DG_C KMKHZ6A.233
*CALL C_EPSLON KMKHZ6A.234
*CALL C_GAMMA KMKHZ6A.235
KMKHZ6A.236
REAL ETAR,GRCP,LCRCP,LFRCP,LS,LSRCP KMKHZ6A.237
PARAMETER ( KMKHZ6A.238
& ETAR=1.0/(1.0-EPSILON) ! Used in buoyancy parameter BETAC. KMKHZ6A.239
&,GRCP=G/CP ! Adiabatic lapse rate. KMKHZ6A.240
&,LCRCP=LC/CP ! Latent heat of condensation / CP. KMKHZ6A.241
&,LFRCP=LF/CP ! Latent heat of fusion / CP. KMKHZ6A.242
&,LS=LC+LF ! Latent heat of sublimation. KMKHZ6A.243
&,LSRCP=LS/CP ! Latent heat of sublimation / CP. KMKHZ6A.244
&) KMKHZ6A.245
REAL A_PLUME,B_PLUME,A_GRAD_ADJ,MAX_T_GRAD KMKHZ6A.246
& ,MIN_SVL_GRAD,SC_CFTOL,CT_RESID,LGF,FGF ARN0F405.927
PARAMETER ( KMKHZ6A.247
& A_PLUME=0.4 ARN0F405.928
&,B_PLUME=3.26 KMKHZ6A.249
&,A_GRAD_ADJ=3.26 KMKHZ6A.250
&,MAX_T_GRAD=1.0E-3 ARN0F405.929
&,MIN_SVL_GRAD=1.0E-3 ! Minimum SVL gradient in a mixed layer. ARN0F405.930
&,SC_CFTOL=0.1 ! Cloud fraction required for a ARN0F405.931
! stratocumulus layer to be diagnosed. ARN0F405.932
&,CT_RESID=500.0 ! Approximate parcel cloud top residence ARN0F405.933
! time (s). ARN0F405.934
&,LGF=1.0 ! Adiabatic gradient factor for cloud ARN0F405.935
! ! liquid water. ARN0F405.936
&,FGF=0.0) ! Adiabatic gradient factor for cloud ARN0F405.937
! ! frozen water. ARN0F405.938
KMKHZ6A.252
! KMKHZ6A.253
! Define local storage. KMKHZ6A.254
! KMKHZ6A.255
! (a) Workspace. 6*BL_LEVELS-1 blocks of real workspace are required KMKHZ6A.256
! plus 1 block of logical. KMKHZ6A.257
! KMKHZ6A.258
! KMKHZ6A.259
LOGICAL KMKHZ6A.260
& TOPBL(P_FIELD) ! Flag set when top of boundary layer KMKHZ6A.261
! ! is reached. KMKHZ6A.262
&,CLOUD_BASE(P_FIELD) ! Flag set when cloud base is reached. KMKHZ6A.263
! KMKHZ6A.264
&,ABOVE_LCL(P_FIELD) ! Flag set when parcel above LCL. ARN0F405.939
&,COUPLED(P_FIELD) ! Flag to indicate stratocumulus layer ARN0F405.940
! ! weakly coupled to surface (i.e. ARN0F405.941
! ! weakly decoupled). ARN0F405.942
ARN0F405.943
REAL KMKHZ6A.265
& QS(P_FIELD) ! Saturated sp humidity at pressure and KMKHZ6A.266
! ! temperature of sucessive levels. KMKHZ6A.267
&,QCL_IC_BOT(P_FIELD,BL_LEVELS) ! In-cloud liquid water content KMKHZ6A.268
! ! at the bottom of the model layer KMKHZ6A.269
&,QCF_IC_BOT(P_FIELD,BL_LEVELS) ! In-cloud frozen water content KMKHZ6A.270
! ! at the bottom of the model layer KMKHZ6A.271
&,QCL_IC_TOP(P_FIELD,BL_LEVELS) ! In-cloud liquid water content KMKHZ6A.272
! ! at the top of the model layer. KMKHZ6A.273
&,QCF_IC_TOP(P_FIELD,BL_LEVELS) ! In-cloud frozen water content KMKHZ6A.274
! ! at the top of the model layer. KMKHZ6A.275
&,CFL(P_FIELD,BL_LEVELS) ! Liquid cloud fraction. KMKHZ6A.276
&,CFF(P_FIELD,BL_LEVELS) ! Frozen cloud fraction. KMKHZ6A.277
&,DQCLDZ(P_FIELD,BL_LEVELS) ! Vertical gradient of in-cloud KMKHZ6A.278
! ! liquid cloud water in a KMKHZ6A.279
! ! well mixed layer. ARN0F405.944
&,DQCFDZ(P_FIELD,BL_LEVELS) ! Vertical gradient of in-cloud KMKHZ6A.281
! ! frozen cloud water in a KMKHZ6A.282
! ! well-mixed layer. KMKHZ6A.283
&,CHI_S(P_FIELD,2:BL_LEVELS) ! Mixing fraction of just saturate KMKHZ6A.284
! ! mixture. KMKHZ6A.285
&,DB(P_FIELD,2:BL_LEVELS) ! Buoyancy jump across layer ARN0F405.945
! ! interface. ARN0F405.946
&,BFLUX_SURF(P_FIELD) ! Surface buoyancy flux using ARN0F405.947
! ! clear air buoyancy parameters. ARN0F405.948
&,BFLUX_SAT_SURF(P_FIELD) ! Surface buoyancy flux using ARN0F405.949
! ! cloudy buoyancy parameters. ARN0F405.950
&,DB_TOP(P_FIELD) ! Buoyancy jump at the top of the KMKHZ6A.287
! ! boundary layer. KMKHZ6A.288
&,DB_CLD(P_FIELD,2:BL_LEVELS) ! In-cloud buoyancy jump across KMKHZ6A.289
! ! layer interface. KMKHZ6A.290
&,DF_OVER_CP(P_FIELD,BL_LEVELS) ! Radiative flux change over layer KMKHZ6A.291
! ! divided by c_P. KMKHZ6A.292
&,DF_TOP_OVER_CP(P_FIELD) ! Radiative flux change at cloud t KMKHZ6A.293
! ! divided by c_P. KMKHZ6A.294
&,SVL_PLUME(P_FIELD) ! Liquid/frozen water static ARN0F405.951
! ! energy over CP for a plume ARN0F405.952
! ! parcel rising without dilution ARN0F405.953
! ! from level 1. KMKHZ6A.301
&,ENV_SVL_KM1(P_FIELD) ! Density (virtual) static energy ARN0F405.954
! ! over CP for last layer. ARN0F405.955
! ! considered. ARN0F405.956
&,PAR_SVL_KM1(P_FIELD) ! Density (virtual) static energy ARN0F405.957
! over CP of parcel at last level ARN0F405.958
! ! considered. ARN0F405.959
&,Z_CLD(P_FIELD) ! Cloud fraction weighted KMKHZ6A.305
! ! thickness of b.l. cloud. KMKHZ6A.306
&,BT_TOP(P_FIELD) ! Buoyancy parameter at the top of KMKHZ6A.307
! ! the turbulently mixed layer. ARN0F405.960
&,BTT_TOP(P_FIELD) ! In-cloud buoyancy parameter at KMKHZ6A.309
! ! the top of the turbulently ARN0F405.961
! ! mixed layer. ARN0F405.962
&,BTC_TOP(P_FIELD) ! Cloud fraction weighted buoyancy KMKHZ6A.311
! ! parameter at the top of the b.l. KMKHZ6A.312
&,DB_TOP_CLD(P_FIELD) ! In-cloud buoyancy jump at the KMKHZ6A.313
! ! top of the b.l. KMKHZ6A.314
&,CLD_FACTOR(P_FIELD) ! Fraction of grid box potentially KMKHZ6A.315
! ! giving evaporative entrainment. KMKHZ6A.316
&,CHI_S_TOP(P_FIELD) ! Mixing fraction of just saturate KMKHZ6A.317
! ! mixture at top of the b.l. KMKHZ6A.318
&,ZETA_S(P_FIELD) ! Non-cloudy fraction of mixing KMKHZ6A.319
! ! layer for surface forced KMKHZ6A.320
! ! entrainment term. KMKHZ6A.321
&,ZETA_R(P_FIELD) ! Non-cloudy fraction of mixing KMKHZ6A.322
! ! layer for cloud top radiative KMKHZ6A.323
! ! cooling entrainment term. KMKHZ6A.324
&,ALPHA_R(P_FIELD) ! Fraction of the cloud top KMKHZ6A.325
! ! radiative cooling assumed to KMKHZ6A.326
! ! act above the minimum turbulent KMKHZ6A.327
! ! flux level. KMKHZ6A.328
&,ZC(P_FIELD) ! Cloud depth (not cloud fraction KMKHZ6A.329
! ! weighted). KMKHZ6A.330
&,T_LCL(P_FIELD) ! Temperature at liftng condensati KMKHZ6A.331
! ! level. KMKHZ6A.332
&,P_LCL(P_FIELD) ! Pressure at lifting condensation KMKHZ6A.333
! ! level. KMKHZ6A.334
&,ZHSC(P_FIELD) ! Height of cloud layer top above ARN0F405.963
! ! surface (m). ARN0F405.964
&,DSVL_TOP(P_FIELD) ! s_VL jump at cloud layer top(K) ARN0F405.965
&,DSCDEPTH(P_FIELD) ! Depth of cloud layer (m). ARN0F405.966
&,DB_DSCT(P_FIELD) ! Buoyancy jump at the top of the ARN0F405.967
! ! surface mixed layer. ARN0F405.968
&,SVL(P_FIELD,BL_LEVELS) ! Liquid/frozen water virtual ARN0F405.969
! ! static energy over CP. ARN0F405.970
&,DSVL_KP2(P_FIELD) ! Liquid/frozen water virtual ARN0F405.971
! ! static energy gradient over CP. ARN0F405.972
&,DF_DSCT_OVER_CP(P_FIELD) ! Radiative flux change at ARN0F405.973
! ! decoupled stratocumulus top ARN0F405.974
! ! divided by c_P. ARN0F405.975
&,BT_DSCT(P_FIELD) ! Buoyancy parameter at the top of ARN0F405.976
! ! the decoupled Sc. ARN0F405.977
&,BTT_DSCT(P_FIELD) ! In-cloud buoyancy parameter at ARN0F405.978
! ! the top of the decoupled Sc. ARN0F405.979
&,BTC_DSCT(P_FIELD) ! Cloud fraction weighted buoyancy ARN0F405.980
! ! parameter at the top of the ARN0F405.981
! ! decoupled Sc. ARN0F405.982
&,DB_DSCT_CLD(P_FIELD) ! In-cloud buoyancy jump at the ARN0F405.983
! ! top of the decoupled Sc. ARN0F405.984
&,CHI_S_DSCT(P_FIELD) ! Mixing fraction of just ARN0F405.985
! ! saturated mixture at top of ARN0F405.986
! ! the decoupled Sc. ARN0F405.987
&,ZETA_R_DSC(P_FIELD) ! Non-cloudy fraction of decoupled ARN0F405.988
! ! Sc layer for cloud top radiative ARN0F405.989
! ! cooling entrainment term. ARN0F405.990
&,ALPHA_R_DSC(P_FIELD) ! Fraction of the cloud top ARN0F405.991
! ! radiative cooling assumed to ARN0F405.992
! ! act above the minimum turbulent ARN0F405.993
! ! flux level for the decoupled Sc. ARN0F405.994
&,ZC_DSC(P_FIELD) ! Cloud depth (not cloud fraction ARN0F405.995
! ! weighted). ARN0F405.996
&,Z_CLD_DSC(P_FIELD) ! Cloud fraction weighted ARN0F405.997
! ! thickness of the decoupled Sc ARN0F405.998
! ! layer. ARN0F405.999
&,CLD_FACTOR_DSC(P_FIELD) ! Fraction of grid box potentially ARN0F405.1000
! ! giving evaporative entrainment. ARN0F405.1001
&,ZHPAR(P_FIELD) ! Temporary store of ZH before LCL ARN0F405.1002
! ! test. ARN0F405.1003
! KMKHZ6A.337
INTEGER KMKHZ6A.338
& K_CLOUD_TOP(P_FIELD) ! Level number of top of b.l. cloud. KMKHZ6A.339
&,NLCL(P_FIELD) ! No. of model layers below the lifting KMKHZ6A.340
! ! condensation level. KMKHZ6A.341
&,NTPAR(P_FIELD) ! Temporary location to save initial ARN0F405.1004
! ! NTML before LCL test. ARN0F405.1005
KMKHZ6A.342
! KMKHZ6A.343
! (b) Scalars. KMKHZ6A.344
! KMKHZ6A.345
REAL KMKHZ6A.346
& VIRT_FACTOR ! Temporary in calculation of buoyancy KMKHZ6A.347
! ! parameters. KMKHZ6A.348
&,DTLDZ ! Vertical gradient of TL in a well-mixed KMKHZ6A.349
! ! layer. KMKHZ6A.350
&,DQW ! Total water content change across layer KMKHZ6A.351
! ! interface. KMKHZ6A.352
&,DTL ! Liquid/ice static energy change across KMKHZ6A.353
! ! layer interface. KMKHZ6A.354
&,DQCL ! Cloud liquid water change across layer KMKHZ6A.355
! ! interface. KMKHZ6A.356
&,DQCF ! Cloud frozen water change across layer KMKHZ6A.357
! ! interface. KMKHZ6A.358
&,VAP_PRESS ! Vapour pressure. KMKHZ6A.359
&,GRAD_CLD ! Humidity gradient in layer above LCL. KMKHZ6A.360
&,GRAD_SUB ! Humidity gradient in layer below LCL. KMKHZ6A.361
&,Q_VAP_PARC ! Vapour content of parcel. ARN0F405.1006
&,Q_LIQ_PARC ! Condensed water content of parcel. ARN0F405.1007
&,T_PARC ! Temperature of parcel. ARN0F405.1008
&,T_DENS_PARC ! Density potential temperature of parcel. ARN0F405.1009
&,T_DENS_ENV ! Density potential temperature of environment. ARN0F405.1010
&,DENV_BYDZ ! Gradient of density potential ARN0F405.1011
! ! temperature in the environment. ARN0F405.1012
&,DPAR_BYDZ ! Gradient of density potential ARN0F405.1013
! ! temperature of the parcel. ARN0F405.1014
&,D_SIEMS ! Siems et al. (1990) cloud-top entrainment ARN0F405.1015
! ! instability parameter. ARN0F405.1016
&,DSVL_KP1 ! s_VL gradient between previous layers. ARN0F405.1017
! KMKHZ6A.364
INTEGER KMKHZ6A.365
& I ! Loop counter (horizontal field index). KMKHZ6A.366
&,J ! Offset counter in certain I loops. KMKHZ6A.367
&,K ! Loop counter (vertical level index). KMKHZ6A.368
&,KM1 ! K-1 KMKHZ6A.369
&,MBL ! Maximum number of model layers allowed in the rapidly KMKHZ6A.370
! ! mixing layer; set to BL_LEVELS-1. KMKHZ6A.371
&,INV_DROP ! Set to 1 if inversion spread over 2 levels. ARN0F405.1018
&,K_RAD_SMLT ! Highest SML level for radiative divergence calc. ARN0F405.1019
&,K_RAD_LIM ARN0F405.1020
! KMKHZ6A.372
!----------------------------------------------------------------------- KMKHZ6A.373
!! 0. Check that the scalars input to define the grid are consistent. KMKHZ6A.374
!! See comments to routine SF_EXCH for details. KMKHZ6A.375
!----------------------------------------------------------------------- KMKHZ6A.376
KMKHZ6A.377
IF (LTIMER) THEN KMKHZ6A.378
CALL TIMER('KMKHZ ',3) KMKHZ6A.379
ENDIF KMKHZ6A.380
KMKHZ6A.381
! Set MBL, "maximum number of boundary levels" for the purposes of KMKHZ6A.382
! boundary layer height calculation. KMKHZ6A.383
KMKHZ6A.384
MBL = BL_LEVELS - 1 KMKHZ6A.385
KMKHZ6A.386
! KMKHZ6A.387
!----------------------------------------------------------------------- KMKHZ6A.388
!! 0.1 Set TOPBL to .FALSE. and calculate boundary layer top using KMKHZ6A.389
!! a non-local, plume method. ARN0F405.1021
! ---------------------------------------------------------------------- KMKHZ6A.391
! KMKHZ6A.392
DO I = P1,P1-1+P_POINTS KMKHZ6A.393
SVL_PLUME(I) = TL(I,1) + GRCP * 0.5 * DZL(I,1) + A_PLUME + ARN0F405.1022
& MIN( MAX_T_GRAD * ZH(I) , B_PLUME * TV1_SD(I) ) KMKHZ6A.396
ARN0F405.1023
! KMKHZ6A.397
!----------------------------------------------------------------------- KMKHZ6A.398
!! 0.2 Calculate temperature and pressure of lifting condensation level KMKHZ6A.399
!! using approximations from Bolton (1980) KMKHZ6A.400
!----------------------------------------------------------------------- KMKHZ6A.401
! KMKHZ6A.402
VAP_PRESS = Q(I,1) * P(I,1) / (100.0*EPSILON) KMKHZ6A.403
IF (VAP_PRESS .GT .0.0) THEN KMKHZ6A.404
T_LCL(I) = 2840.0/(ALOG(T(I,1))/KAPPA-ALOG(VAP_PRESS)-4.805) KMKHZ6A.405
& + 55.0 KMKHZ6A.406
P_LCL(I) = P(I,1) * (T_LCL(I)/T(I,1))**(1.0/KAPPA) KMKHZ6A.407
ELSE KMKHZ6A.408
P_LCL(I) = P_HALF(I,1) KMKHZ6A.409
ENDIF KMKHZ6A.410
PAR_SVL_KM1(I) = SVL_PLUME(I) * ( 1.0 + C_VIRTUAL*Q(I,1) - ARN0F405.1024
& QCL(I,1) - QCF(I,1) ) ARN0F405.1025
ENV_SVL_KM1(I) = TL(I,1) + GRCP * 0.5 * DZL(I,1) ARN0F405.1026
CUMULUS(I)= .FALSE. ARN0F405.1027
TOPBL(I) = .FALSE. KMKHZ6A.412
NTML(I) = 1 KMKHZ6A.413
NLCL(I) = 1 KMKHZ6A.414
Z_LCL(I) = Z_HALF(I,1) ARN0F405.1028
ENDDO KMKHZ6A.415
! KMKHZ6A.416
!----------------------------------------------------------------------- KMKHZ6A.417
!! 0.3 Now compare plume s_VL with each model layer s_VL in turn to KMKHZ6A.418
!! find the first time that plume has negative buoyancy. KMKHZ6A.419
!----------------------------------------------------------------------- KMKHZ6A.420
! KMKHZ6A.421
DO K = 2,BL_LEVELS KMKHZ6A.422
CALL QSAT(QS(P1),T(P1,K),P(P1,K),P_POINTS) ARN0F405.1029
DO I = P1,P1-1+P_POINTS KMKHZ6A.423
IF ( P_LCL(I) .LT. P_HALF(I,K) ) THEN KMKHZ6A.424
! !----------------------------------------------------------- ARN0F405.1030
! ! Below the lifting condensation level ARN0F405.1031
! ! potential temperature is constant. ARN0F405.1032
! !----------------------------------------------------------- ARN0F405.1033
T_PARC = SVL_PLUME(I) - GRCP * Z_FULL(I,K) ARN0F405.1034
Q_VAP_PARC = Q(I,1) ARN0F405.1035
Q_LIQ_PARC = 0.0 ARN0F405.1036
Z_LCL(I) = Z_HALF(I,K) KMKHZ6A.425
NLCL(I) = K-1 KMKHZ6A.426
ABOVE_LCL(I) = .FALSE. ARN0F405.1037
ELSE ARN0F405.1038
! !----------------------------------------------------------- ARN0F405.1039
! ! Above the lifting condenstion level ARN0F405.1040
! ! calculate parcel potential temperature by linearising ARN0F405.1041
! ! q_sat about the environmental temperature. ARN0F405.1042
! !----------------------------------------------------------- ARN0F405.1043
T_PARC = ( SVL_PLUME(I) - GRCP * Z_FULL(I,K) ARN0F405.1044
& + LCRCP * (QW(I,1) - QS(I) + DQSDT(I,K)*T(I,K)) ARN0F405.1045
& ) / (1.0 + LCRCP*DQSDT(I,K)) ARN0F405.1046
IF (T_PARC .LT. TM) THEN ARN0F405.1047
T_PARC = ( SVL_PLUME(I) - GRCP * Z_FULL(I,K) ARN0F405.1048
& + LSRCP * (QW(I,1) - QS(I) + DQSDT(I,K)*T(I,K)) ARN0F405.1049
& ) / (1.0 + LSRCP*DQSDT(I,K)) ARN0F405.1050
ENDIF KMKHZ6A.427
Q_VAP_PARC = QS(I) + DQSDT(I,K) * (T_PARC - T(I,K)) ARN0F405.1051
Q_LIQ_PARC = QW(I,1) - Q_VAP_PARC ARN0F405.1052
ABOVE_LCL(I) = .TRUE. ARN0F405.1053
ENDIF ARN0F405.1054
! KMKHZ6A.428
T_DENS_PARC = T_PARC * ARN0F405.1055
& (1.0 + C_VIRTUAL*Q_VAP_PARC - Q_LIQ_PARC) ARN0F405.1056
T_DENS_ENV = T(I,K) * ARN0F405.1057
& (1.0 + C_VIRTUAL*Q(I,K) - QCL(I,K) - QCF(I,K)) ARN0F405.1058
! !------------------------------------------------------------- ARN0F405.1059
! ! Find vertical gradients in parcel and environment SVL ARN0F405.1060
! ! (using values from level below (i.e. K-1)). ARN0F405.1061
! !------------------------------------------------------------- ARN0F405.1062
DPAR_BYDZ = (T_DENS_PARC + GRCP*Z_FULL(I,K) - ARN0F405.1063
& PAR_SVL_KM1(I)) / (Z_FULL(I,K) - Z_FULL(I,K-1)) ARN0F405.1064
DENV_BYDZ = (T_DENS_ENV + GRCP*Z_FULL(I,K) - ARN0F405.1065
& ENV_SVL_KM1(I)) / (Z_FULL(I,K) - Z_FULL(I,K-1)) ARN0F405.1066
IF ( .NOT.TOPBL(I) .AND. KMKHZ6A.432
& ( ( T_DENS_PARC-T_DENS_ENV .LE. 0.0 ) .OR. ARN0F405.1067
! ARN0F405.1068
! plume non buoyant ARN0F405.1069
! ARN0F405.1070
& (ABOVE_LCL(I) .AND. (DENV_BYDZ .GT. 1.25*DPAR_BYDZ)) .OR. ARN0F405.1071
! ARN0F405.1072
! or environmental virtual temperature gradient ARN0F405.1073
! significantly larger than parcel gradient ARN0F405.1074
! above lifting condensation level ARN0F405.1075
! ARN0F405.1076
& (K .GT. MBL) KMKHZ6A.438
! or max no. of model layers in b.l. reached ARN0F405.1077
& ) KMKHZ6A.440
& ) THEN KMKHZ6A.441
! KMKHZ6A.442
TOPBL(I) = .TRUE. KMKHZ6A.443
ZH(I) = Z_HALF(I,K) KMKHZ6A.444
NTML(I) = K-1 KMKHZ6A.445
ENDIF KMKHZ6A.446
ENV_SVL_KM1(I) = T_DENS_ENV + GRCP*Z_FULL(I,K) ARN0F405.1078
PAR_SVL_KM1(I) = T_DENS_PARC + GRCP*Z_FULL(I,K) ARN0F405.1079
ENDDO KMKHZ6A.448
ENDDO KMKHZ6A.449
!----------------------------------------------------------------------- KMKHZ6A.450
!! 0.3a Look for decoupled cloudy mixed layer above b.l. top at ZH: ARN0F405.1080
!! find cloud base above surface mixed layer inversion, ARN0F405.1081
!! i.e. above NTML+1, ARN0F405.1082
!! then cloud top (i.e. CF < SC_CFTOL) ARN0F405.1083
!! and finally check that cloud is well mixed. ARN0F405.1084
!----------------------------------------------------------------------- ARN0F405.1085
! Initialise variables and calculate SVL: conserved variable used ARN0F405.1086
! to test for well mixed layers. ARN0F405.1087
! ARN0F405.1088
DO I = P1,P1-1+P_POINTS ARN0F405.1089
CLOUD_BASE(I) = .FALSE. ARN0F405.1090
DSC(I) = .FALSE. ARN0F405.1091
COUPLED(I) = .FALSE. ARN0F405.1092
ZHSC(I) = 0.0 ARN0F405.1093
NTDSC(I) = 0 ARN0F405.1094
ENDDO ARN0F405.1095
ARN0F405.1096
DO K = 1,BL_LEVELS ARN0F405.1097
DO I = P1,P1-1+P_POINTS ARN0F405.1098
SVL(I,K) = TL(I,K) * ( 1.0 + C_VIRTUAL*QW(I,K) ) + ARN0F405.1099
& GRCP * Z_FULL(I,K) ARN0F405.1100
ENDDO ARN0F405.1101
ENDDO ARN0F405.1102
! ARN0F405.1103
DO K=3,MBL ARN0F405.1104
DO I = P1,P1-1+P_POINTS ARN0F405.1105
! !-------------------------------------------------------------- ARN0F405.1106
! ! Find cloud base (where cloud here means CF > SC_CFTOL). ARN0F405.1107
! !-------------------------------------------------------------- ARN0F405.1108
IF ( K .GT. NTML(I)+1 .AND. CF(I,K) .GT. SC_CFTOL ARN0F405.1109
& .AND. .NOT.CLOUD_BASE(I) ARN0F405.1110
! not yet found cloud base ARN0F405.1111
& .AND. .NOT.DSC(I) ) THEN ARN0F405.1112
! not yet found a Sc layer ARN0F405.1113
CLOUD_BASE(I) = .TRUE. ARN0F405.1114
ENDIF ARN0F405.1115
IF ( CLOUD_BASE(I) .AND. .NOT.DSC(I) .AND. ( ARN0F405.1116
! found cloud base but not yet reached cloud top ARN0F405.1117
& (CF(I,K+1).LT.SC_CFTOL) ) ARN0F405.1118
! got to cloud top ARN0F405.1119
& ) THEN ARN0F405.1120
! !----------------------------------------------------------- ARN0F405.1121
! ! Look to see if at least top of cloud is well mixed: ARN0F405.1122
! ! test SVL-gradient for top 2 pairs of levels, in case ARN0F405.1123
! ! cloud top extends into the inversion. ARN0F405.1124
! ! Parcel descent in Section 4.0 below will determine depth ARN0F405.1125
! ! of mixed layer. ARN0F405.1126
! !---------------------------------------------------------- ARN0F405.1127
IF ( 2.0*(SVL(I,K)-SVL(I,K-1)) / ARN0F405.1128
& (DZL(I,K)+DZL(I,K-1)) .LT. MIN_SVL_GRAD ) THEN ARN0F405.1129
DSC(I) = .TRUE. ARN0F405.1130
NTDSC(I) = K ARN0F405.1131
ZHSC(I) = Z_HALF(I,K+1) ARN0F405.1132
ELSEIF ( 2.0*(SVL(I,K-1)-SVL(I,K-2)) / ARN0F405.1133
& (DZL(I,K-1)+DZL(I,K-2)) .LT. MIN_SVL_GRAD ) THEN ARN0F405.1134
DSC(I) = .TRUE. ARN0F405.1135
NTDSC(I) = K-1 ARN0F405.1136
ZHSC(I) = Z_HALF(I,K) ARN0F405.1137
ELSE ARN0F405.1138
CLOUD_BASE(I) = .FALSE. ARN0F405.1139
ENDIF ARN0F405.1140
ENDIF ARN0F405.1141
ENDDO ARN0F405.1142
ENDDO ARN0F405.1143
!----------------------------------------------------------------------- ARN0F405.1144
!! 0.4 Save parcel ascent top: this will be used to allow mixing and ARN0F405.1145
!! entrainment into decoupled Sc of single layer thickness when it ARN0F405.1146
!! occurs above Cu. ARN0F405.1147
!----------------------------------------------------------------------- ARN0F405.1148
DO I = P1,P1-1+P_POINTS ARN0F405.1149
ZHPAR(I) = ZH(I) ARN0F405.1150
NTPAR(I) = NTML(I) ARN0F405.1151
ENDDO ARN0F405.1152
! ARN0F405.1153
!----------------------------------------------------------------------- ARN0F405.1154
! Test height derived above against lifting condensation level KMKHZ6A.451
!----------------------------------------------------------------------- KMKHZ6A.452
DO I = P1,P1-1+P_POINTS KMKHZ6A.453
!----------------------------------------------------------------------- KMKHZ6A.454
! Check lifting condensation levels against height of parcel ascent, KMKHZ6A.455
! if lifting condensation level lower than parcel ascent then decide KMKHZ6A.456
! on type of cloudy layer. If lifting condensation level at or below KMKHZ6A.457
! low grid point, assume fog layer and turbulent mixing. For KMKHZ6A.458
! gradient tests assume any if LCL and top of parcel ascent is less KMKHZ6A.459
! than two levels then stratocumulus. KMKHZ6A.460
!----------------------------------------------------------------------- KMKHZ6A.461
IF ( (NTML(I)-NLCL(I) .GE. 2) .AND. (NLCL(I) .GT. 1) ) THEN KMKHZ6A.462
!----------------------------------------------------------------------- KMKHZ6A.463
! Cloudy boundary layer, diagnose whether stratocumulus or cumulus. KMKHZ6A.464
! For stratocumulus top of mixed layer = ZH KMKHZ6A.465
! For cumulus top of mixed layer = ZLCL KMKHZ6A.466
! Diagnosis is done by comparing gradients KMKHZ6A.467
!----------------------------------------------------------------------- KMKHZ6A.468
IF (NTML(I) .EQ. MBL) THEN ARN0F405.1155
CUMULUS(I) = .TRUE. ARN0F405.1156
ZH(I) = Z_LCL(I) KMKHZ6A.470
NTML(I) = NLCL(I) KMKHZ6A.471
ELSE KMKHZ6A.472
GRAD_CLD = ABS( QW(I,NTML(I)) - QW(I,NLCL(I)) ) / KMKHZ6A.473
& ( Z_FULL(I,NTML(I)) - Z_FULL(I,NLCL(I)) ) KMKHZ6A.474
GRAD_SUB = ABS( QW(I,NLCL(I)) - QW(I,1) ) / KMKHZ6A.475
& ( Z_FULL(I,NLCL(I)) - Z_FULL(I,1) ) KMKHZ6A.476
IF (GRAD_CLD .GT. 1.10*GRAD_SUB) THEN KMKHZ6A.477
!----------------------------------------------------------------------- KMKHZ6A.478
! Not well mixed, however it is possible that the depth of a well KMKHZ6A.479
! mixed boundary layer has increased but not yet been mixed yet so KMKHZ6A.480
! test gradient from next level down. KMKHZ6A.481
!----------------------------------------------------------------------- KMKHZ6A.482
KMKHZ6A.483
GRAD_CLD = ABS( QW(I,NTML(I)-1) - QW(I,NLCL(I)) ) / KMKHZ6A.485
& ( Z_FULL(I,NTML(I)-1) - Z_FULL(I,NLCL(I)) ) KMKHZ6A.486
ARN0F405.1157
IF ( GRAD_CLD .GT. 1.10*GRAD_SUB) THEN KMKHZ6A.487
!----------------------------------------------------------------------- KMKHZ6A.488
! Diagnoss a cumulus layer and set mixed layer depth to Z_LCL KMKHZ6A.489
!----------------------------------------------------------------------- KMKHZ6A.490
IF (P_LCL(I) .LT. (P(I,NLCL(I)))) THEN ARN0F405.1158
!----------------------------------------------------------------------- ARN0F405.1159
! If LCL is diagnosed in the upper half of the layer set Z_LCL to ARN0F405.1160
! the height of the upper layer interface ARN0F405.1161
! (in code above LCL is always set to the lower interface). ARN0F405.1162
!----------------------------------------------------------------------- ARN0F405.1163
NLCL(I) = NLCL(I)+1 ARN0F405.1164
Z_LCL(I) = Z_HALF(I,NLCL(I)+1) ARN0F405.1165
ENDIF KMKHZ6A.493
CUMULUS(I) = .TRUE. ARN0F405.1166
ZH(I) = Z_LCL(I) ARN0F405.1167
NTML(I) = NLCL(I) ARN0F405.1168
ENDIF KMKHZ6A.494
ENDIF KMKHZ6A.495
ENDIF KMKHZ6A.496
ENDIF ARN0F405.1169
ENDDO KMKHZ6A.497
! KMKHZ6A.498
CALL QSAT(QS(P1),T(P1,1),P(P1,1),P_POINTS) KMKHZ6A.499
! KMKHZ6A.500
!----------------------------------------------------------------------- ARN0F405.1170
!! 0.4a Check whether the inversion at NTML is spread over 2 levels. ARN0F405.1171
!! If it is, we want the lower as it may indicate a deepening ARN0F405.1172
!! layer which will deepen too fast if we jump up a level too ARN0F405.1173
!! soon. Note that NTML does not change if it is at NLCL (ie. Cu). ARN0F405.1174
!----------------------------------------------------------------------- ARN0F405.1175
! ARN0F405.1176
DO I = P1,P1-1+P_POINTS ARN0F405.1177
IF ( NTML(I).GT.1 .AND. .NOT.CUMULUS(I) ) THEN ARN0F405.1178
K=NTML(I) ARN0F405.1179
IF ( 2.0*(SVL(I,K)-SVL(I,K-1)) / ARN0F405.1180
& (DZL(I,K)+DZL(I,K-1)) .GT. MIN_SVL_GRAD ) THEN ARN0F405.1181
NTML(I) = K-1 ARN0F405.1182
ZH(I) = Z_HALF(I,K) ARN0F405.1183
ENDIF ARN0F405.1184
ENDIF ARN0F405.1185
ENDDO ARN0F405.1186
! ARN0F405.1187
!----------------------------------------------------------------------- ARN0F405.1188
!! 0.4aa If no decoupled cloud layer above ZHPAR was found but the layer ARN0F405.1189
!! to ZHPAR is cloudy, declare this layer a decoupled cloud layer ARN0F405.1190
!! and set ZHSC and NTDSC accordingly. Note that if the ARN0F405.1191
!! subsequent tests indicate a well mixed cloud layer to the ARN0F405.1192
!! surface, the decoupled layer will be switched off. ARN0F405.1193
!----------------------------------------------------------------------- ARN0F405.1194
! ARN0F405.1195
DO I = P1,P1-1+P_POINTS ARN0F405.1196
IF ( .NOT.DSC(I) .AND. CF(I,NTPAR(I)) .GT. SC_CFTOL ) THEN ARN0F405.1197
DSC(I) = .TRUE. ARN0F405.1198
ARN0F405.1199
IF ( ZHPAR(I) .GT. ZH(I) .AND. NTPAR(I) .GT. 1 ) THEN ARN0F405.1200
! ! May still need to lower this inversion, as at section 0.4a ARN0F405.1201
K=NTPAR(I) ARN0F405.1202
IF ( 2.0*(SVL(I,K)-SVL(I,K-1)) / ARN0F405.1203
& (DZL(I,K)+DZL(I,K-1)) .GT. MIN_SVL_GRAD ) THEN ARN0F405.1204
NTPAR(I) = K - 1 ARN0F405.1205
ZHPAR(I) = Z_HALF(I,K) ARN0F405.1206
ENDIF ARN0F405.1207
ENDIF ARN0F405.1208
ZHSC(I) = ZHPAR(I) ARN0F405.1209
NTDSC(I)= NTPAR(I) ARN0F405.1210
ENDIF ARN0F405.1211
ENDDO ARN0F405.1212
!----------------------------------------------------------------------- ARN0F405.1213
!! 0.4b Calculate the radiative flux changes across cloud top for the ARN0F405.1214
!! stratocumulus layer and thence a first guess for the top-down ARN0F405.1215
!! mixing depth of this layer, DSCDEPTH. ARN0F405.1216
!----------------------------------------------------------------------- ARN0F405.1217
! Initialise variables ARN0F405.1218
!------------------------------ ARN0F405.1219
DO I = P1,P1-1+P_POINTS ARN0F405.1220
K_CLOUD_TOP(I) = 0 ARN0F405.1221
DF_DSCT_OVER_CP(I) = 0.0 ARN0F405.1222
ENDDO ARN0F405.1223
ARN0F405.1224
DO K = 1,BL_LEVELS-1 ARN0F405.1225
DO I = P1,P1-1+P_POINTS ARN0F405.1226
ARN0F405.1227
DF_OVER_CP(I,K) = DELTAP(I,K)/G * RAD_HR(I,K) ARN0F405.1228
! !------------------------------------------------------------- ARN0F405.1229
! ! Find the layer with the greatest radiative flux jump and ARN0F405.1230
! ! assume that this is the top decoupled Sc layer. If this is ARN0F405.1231
! ! above the surface mixed layer (SML) (which implies strong ARN0F405.1232
! ! decoupling), limit the search to levels above the SML. ARN0F405.1233
! !------------------------------------------------------------- ARN0F405.1234
K_RAD_LIM = 0 ARN0F405.1235
IF ( NTDSC(I) .GT. NTML(I) ) K_RAD_LIM = NTML(I)+1 ARN0F405.1236
ARN0F405.1237
IF ( DSC(I) .AND. K .GT. K_RAD_LIM .AND. ARN0F405.1238
& K .LE. NTDSC(I)+2 .AND. ARN0F405.1239
& DF_OVER_CP(I,K) .GT. DF_DSCT_OVER_CP(I) ) THEN ARN0F405.1240
K_CLOUD_TOP(I) = K ARN0F405.1241
DF_DSCT_OVER_CP(I) = DF_OVER_CP(I,K) ARN0F405.1242
ENDIF ARN0F405.1243
ARN0F405.1244
ENDDO ARN0F405.1245
ENDDO ARN0F405.1246
!----------------------------------------------------------------------- ARN0F405.1247
! If cloud top extends into the level above NTDSC (so that the ARN0F405.1248
! radiative divergence maximum is at NTDSC+1) but there is ARN0F405.1249
! additional `cloud top' cooling at NTDSC, add this on to ARN0F405.1250
! DF_DSCT_OVER_CP. ARN0F405.1251
!----------------------------------------------------------------------- ARN0F405.1252
DO I = P1,P1-1+P_POINTS ARN0F405.1253
IF ( DSC(I) .AND. K_CLOUD_TOP(I) .EQ. NTDSC(I)+1 ) THEN ARN0F405.1254
K=NTDSC(I) ARN0F405.1255
DF_OVER_CP(I,K) = DELTAP(I,K)/G * RAD_HR(I,K) ARN0F405.1256
IF (DF_OVER_CP(I,K) .GT. 0.0 ) ARN0F405.1257
& DF_DSCT_OVER_CP(I) = DF_DSCT_OVER_CP(I) + DF_OVER_CP(I,K) ARN0F405.1258
ENDIF ARN0F405.1259
ENDDO ARN0F405.1260
!----------------------------------------------------------------------- ARN0F405.1261
! Set DSCDEPTH, the depth of the DSC layer. Note that this estimate ARN0F405.1262
! will be the length-scale used to calculate the entrainment rate ARN0F405.1263
! (although the dependence is only weak) but that a more accurate ARN0F405.1264
! plume descent (dependent on w_e) is subsequently used to determine ARN0F405.1265
! the depth over which the top-down mixing profiles will be applied. ARN0F405.1266
! If DSC is FALSE, DSCDEPTH = 0. The plume descent here uses ARN0F405.1267
! a radiative perturbation to the cloud-top SVL based roughly on a ARN0F405.1268
! typical cloud-top residence time. If the plume does not sink and ARN0F405.1269
! the cloud is decoupled from the surface, then it is assumed to be ARN0F405.1270
! stable, i.e. St rather than Sc, and no mixing or entrainment is ARN0F405.1271
! applied to it. ARN0F405.1272
!----------------------------------------------------------------------- ARN0F405.1273
DO I = P1,P1-1+P_POINTS ARN0F405.1274
DSCDEPTH(I) = 0.0 ARN0F405.1275
IF ( DSC(I) ) THEN ARN0F405.1276
IF (K_CLOUD_TOP(I) .GT. 0) THEN ARN0F405.1277
SVL_PLUME(I) = SVL(I,NTDSC(I)) ARN0F405.1278
& + CT_RESID * MIN( RAD_HR(I,K_CLOUD_TOP(I)), 0.0) ARN0F405.1279
ELSE ARN0F405.1280
SVL_PLUME(I) = SVL(I,NTDSC(I)) ARN0F405.1281
ENDIF ARN0F405.1282
ARN0F405.1283
IF ( K_CLOUD_TOP(I) .EQ. NTDSC(I)+1 ) ARN0F405.1284
& SVL_PLUME(I)=SVL_PLUME(I) ARN0F405.1285
& + CT_RESID * MIN( RAD_HR(I,NTDSC(I)), 0.0) ARN0F405.1286
ARN0F405.1287
ELSE ARN0F405.1288
SVL_PLUME(I)=0.0 ARN0F405.1289
ENDIF ARN0F405.1290
ENDDO ARN0F405.1291
ARN0F405.1292
DO K=BL_LEVELS-1,1,-1 ARN0F405.1293
DO I = P1,P1-1+P_POINTS ARN0F405.1294
IF ( K.LT.NTDSC(I) ! layer must at least 2 levels deep ARN0F405.1295
& .AND. SVL_PLUME(I)-SVL(I,K) .LT. 0.01 ) THEN ARN0F405.1296
DSCDEPTH(I) = ZHSC(I) - Z_HALF(I,K) ARN0F405.1297
ENDIF ARN0F405.1298
ENDDO ARN0F405.1299
ENDDO ARN0F405.1300
!----------------------------------------------------------------------- ARN0F405.1301
!! 0.4c If the layer to ZH is cloud capped and there is no Cu or ARN0F405.1302
!! decoupled Sc above (indicated by NTML=NTDSC), set spare logical ARN0F405.1303
!! CLOUD_BASE to true and try looking for a weak, ARN0F405.1304
!! possibly diurnally decoupling, inversion below ZH using local ARN0F405.1305
!! gradients method. If one is found, lower ZH to this decoupling ARN0F405.1306
!! inversion leaving ZHSC at cloud-top. ARN0F405.1307
!! If DSCDEPTH was found to be zero above, set DSCDEPTH=ZHSC-ZH. ARN0F405.1308
!----------------------------------------------------------------------- ARN0F405.1309
K=1 ARN0F405.1310
DO I = P1,P1-1+P_POINTS ARN0F405.1311
CLOUD_BASE(I) = .FALSE. ARN0F405.1312
IF ( NTML(I) .EQ. NTDSC(I) ) THEN ARN0F405.1313
!----------------------------------------------------------------------- ARN0F405.1314
! Well mixed cloudy layer mixing up to ZH: reset decoupled layer ARN0F405.1315
! parameters ready to look for decoupling inversion below ZH. ARN0F405.1316
!----------------------------------------------------------------------- ARN0F405.1317
DSC(I) = .FALSE. ARN0F405.1318
ZHSC(I) = 0.0 ARN0F405.1319
NTDSC(I)= 0 ARN0F405.1320
CLOUD_BASE(I) = .TRUE. ARN0F405.1321
ENDIF ARN0F405.1322
DSVL_KP2(I)= 2.0*(SVL(I,K+2) - SVL(I,K+1)) / ARN0F405.1323
& (DZL(I,K+2) + DZL(I,K+1)) ARN0F405.1324
ENDDO ARN0F405.1325
ARN0F405.1326
DO K=2,MBL-1 ARN0F405.1327
DO I = P1,P1-1+P_POINTS ARN0F405.1328
!---------------------------------------------------------------------- ARN0F405.1329
!! 0.4d Local gradients method: find where s_VL gradient decreases from ARN0F405.1330
!! a non-trivial positive value (ie. greater than MAX_T_GRAD). ARN0F405.1331
!! The inversion is then at the max s_VL-gradient height. ARN0F405.1332
!---------------------------------------------------------------------- ARN0F405.1333
DSVL_KP1 = DSVL_KP2(I) ARN0F405.1334
DSVL_KP2(I) = 2.0*(SVL(I,K+2) - SVL(I,K+1)) / ARN0F405.1335
& (DZL(I,K+2) + DZL(I,K+1)) ARN0F405.1336
IF ( K.LT.NTML(I) .AND. .NOT.DSC(I) .AND. CLOUD_BASE(I) .AND. ARN0F405.1337
! Below apparently well-mixed cloud-layer top ARN0F405.1338
& (DSVL_KP1 .GT. MAX_T_GRAD .AND. DSVL_KP2(I) .LT. DSVL_KP1) ARN0F405.1339
! and gradient conditions satisfied ARN0F405.1340
& ) THEN ARN0F405.1341
DSC(I) = .TRUE. ARN0F405.1342
NTDSC(I) = NTML(I) ARN0F405.1343
ZHSC(I) = Z_HALF(I,NTDSC(I)+1) ARN0F405.1344
NTML(I) = K ARN0F405.1345
ZH(I) = Z_HALF(I,NTML(I)+1) ARN0F405.1346
! !----------------------------------------------------------- ARN0F405.1347
! ! As in section 0.4c, check whether this ZH is at the base ARN0F405.1348
! ! of the inversion and if it is not, lower it. ARN0F405.1349
! !----------------------------------------------------------- ARN0F405.1350
IF ( 2.0*(SVL(I,K) - SVL(I,K-1)) / ARN0F405.1351
& (DZL(I,K) + DZL(I,K-1)) .GT. MIN_SVL_GRAD ) THEN ARN0F405.1352
NTML(I) = K-1 ARN0F405.1353
ZH(I) = Z_HALF(I,K) ARN0F405.1354
ENDIF ARN0F405.1355
! !----------------------------------------------------------- ARN0F405.1356
! ! Reset DSCDEPTH, if top-down mixing is only weak at present ARN0F405.1357
! !----------------------------------------------------------- ARN0F405.1358
IF ( DSCDEPTH(I) .EQ. 0.0 ) ARN0F405.1359
& DSCDEPTH(I) = ZHSC(I) - Z_HALF(I,NTDSC(I)-1) ARN0F405.1360
ENDIF ARN0F405.1361
ENDDO ARN0F405.1362
ENDDO ARN0F405.1363
!---------------------------------------------------------------------- ARN0F405.1364
!! 0.4e Tidy up variables associated with decoupled layer. ARN0F405.1365
!---------------------------------------------------------------------- ARN0F405.1366
DO I = P1,P1-1+P_POINTS ARN0F405.1367
IF (CUMULUS(I) .AND. DSC(I)) ARN0F405.1368
& DSCDEPTH(I)=MIN( DSCDEPTH(I),ZHSC(I)-Z_HALF(I,NTML(I)+2) ) ARN0F405.1369
IF ( DSCDEPTH(I) .EQ. 0.0 ) THEN ARN0F405.1370
IF (NTDSC(I) .EQ. NTPAR(I)) THEN ARN0F405.1371
! !---------------------------------------------------------- ARN0F405.1372
! ! Indicates a Sc layer over Cu: force mixing over single ARN0F405.1373
! ! layer. ARN0F405.1374
! !---------------------------------------------------------- ARN0F405.1375
DSCDEPTH(I) = DZL(I,NTDSC(I)) ARN0F405.1376
ELSE ARN0F405.1377
DSC(I)=.FALSE. ARN0F405.1378
NTDSC(I)=0 ARN0F405.1379
ZHSC(I)=0.0 ARN0F405.1380
DF_DSCT_OVER_CP(I) = 0.0 ARN0F405.1381
ENDIF ARN0F405.1382
ENDIF ARN0F405.1383
! !--------------------------------------------------------------- ARN0F405.1384
! ! These next two tests should not be necessary but put in to be ARN0F405.1385
! ! fail-safe. ARN0F405.1386
! !--------------------------------------------------------------- ARN0F405.1387
IF (.NOT.DSC(I)) DSCDEPTH(I) = 0.0 ARN0F405.1388
IF (NTDSC(I) .EQ. NTML(I)) THEN ARN0F405.1389
! ! Clearly modelling the same layer ARN0F405.1390
DSC(I) = .FALSE. ARN0F405.1391
NTDSC(I) = 0 ARN0F405.1392
ZHSC(I) = 0.0 ARN0F405.1393
DF_DSCT_OVER_CP(I) = 0.0 ARN0F405.1394
DSCDEPTH(I) = 0.0 ARN0F405.1395
ENDIF ARN0F405.1396
ENDDO ARN0F405.1397
! ARN0F405.1398
!---------------------------------------------------------------------- ARN0F405.1399
! If decoupled cloud layer found test to see if it is, in fact, ARN0F405.1400
! coupled to the surface mixed-layer: if SVL difference between ARN0F405.1401
! NTML and NTDSC is less than 0.2K then assume coupled. This will ARN0F405.1402
! mean that the surface-driven entrainment term will be applied at ARN0F405.1403
! ZHSC, no entrainment will be applied at ZH and ZHSC will be the ARN0F405.1404
! lengthscale used in the entrainment inputs. ARN0F405.1405
!---------------------------------------------------------------------- ARN0F405.1406
DO I = P1,P1-1+P_POINTS ARN0F405.1407
COUPLED(I) = .FALSE. ARN0F405.1408
IF ( DSC(I) ) THEN ARN0F405.1409
! !------------------------------------------------------------ ARN0F405.1410
! ! Note this IF test structure is required because if DSC is ARN0F405.1411
! ! false then NTDSC = 0 and cannot be used to index SVL. ARN0F405.1412
! !------------------------------------------------------------ ARN0F405.1413
IF ( SVL(I,NTDSC(I)) - SVL(I,NTML(I)) .LT. 0.2 ) ARN0F405.1414
& COUPLED(I) = .TRUE. ARN0F405.1415
ENDIF ARN0F405.1416
ENDDO ARN0F405.1417
! ARN0F405.1418
DO I = P1,P1+P_POINTS-1 KMKHZ6A.501
! KMKHZ6A.502
!----------------------------------------------------------------------- KMKHZ6A.503
!! 0.5 Calculate the within-layer vertical gradients of cloud liquid KMKHZ6A.504
!! and frozen water for the layer 1 KMKHZ6A.505
!----------------------------------------------------------------------- KMKHZ6A.506
! KMKHZ6A.507
VIRT_FACTOR = 1.0 + C_VIRTUAL*Q(I,1) - QCL(I,1) - QCF(I,1) KMKHZ6A.508
! KMKHZ6A.509
GRAD_T_ADJ(I) = MIN( MAX_T_GRAD , KMKHZ6A.510
& A_GRAD_ADJ * T1_SD(I) / ZH(I) ) KMKHZ6A.511
! IF (T1_SD(I) .GT. 0.0) THEN KMKHZ6A.512
! GRAD_Q_ADJ(I) = (Q1_SD(I) / T1_SD(I)) * GRAD_T_ADJ(I) KMKHZ6A.513
! ELSE KMKHZ6A.514
GRAD_Q_ADJ(I) = 0.0 KMKHZ6A.515
! ENDIF KMKHZ6A.516
DTLDZ = -GRCP + GRAD_T_ADJ(I) KMKHZ6A.517
DQCLDZ(I,1) = -LGF * ( DTLDZ*DQSDT(I,1) + ARN0F405.1419
& G*QS(I)/(R*T(I,1)*VIRT_FACTOR) ) KMKHZ6A.519
& / (1.0 + LCRCP*DQSDT(I,1)) KMKHZ6A.520
DQCFDZ(I,1) = -FGF * ( DTLDZ*DQSDT(I,1) + ARN0F405.1420
& G*QS(I)/(R*T(I,1)*VIRT_FACTOR) ) KMKHZ6A.522
& / (1.0 + LSRCP*DQSDT(I,1)) KMKHZ6A.523
! KMKHZ6A.524
!----------------------------------------------------------------------- KMKHZ6A.525
!! 0.6 Calculate the cloud liquid and frozen water contents at the KMKHZ6A.526
!! top and bottom of layer 1 KMKHZ6A.527
!----------------------------------------------------------------------- KMKHZ6A.528
! KMKHZ6A.529
IF ( QCL(I,1) + QCF(I,1) .GT. 0.0 ) THEN KMKHZ6A.530
CFL(I,1) = CF(I,1) * QCL(I,1)/(QCL(I,1)+QCF(I,1)) KMKHZ6A.531
CFF(I,1) = CF(I,1) * QCF(I,1)/(QCL(I,1)+QCF(I,1)) KMKHZ6A.532
ELSE KMKHZ6A.533
CFL(I,1) = 0.0 KMKHZ6A.534
CFF(I,1) = 0.0 KMKHZ6A.535
ENDIF KMKHZ6A.536
! KMKHZ6A.537
IF (CFL(I,1) .GT. 0.0) THEN KMKHZ6A.538
QCL_IC_TOP(I,1) = QCL(I,1) / CFL(I,1) + KMKHZ6A.539
& 0.5*DZL(I,1)*DQCLDZ(I,1) KMKHZ6A.540
ELSE KMKHZ6A.541
QCL_IC_TOP(I,1) = 0.0 KMKHZ6A.542
ENDIF KMKHZ6A.543
! KMKHZ6A.544
IF (CFF(I,1) .GT. 0.0) THEN KMKHZ6A.545
QCF_IC_TOP(I,1) = QCF(I,1) / CFF(I,1) + KMKHZ6A.546
& 0.5*DZL(I,1)*DQCFDZ(I,1) KMKHZ6A.547
ELSE KMKHZ6A.548
QCF_IC_TOP(I,1) = 0.0 KMKHZ6A.549
ENDIF KMKHZ6A.550
! KMKHZ6A.551
QCL_IC_BOT(I,1) = 0.0 KMKHZ6A.552
QCF_IC_BOT(I,1) = 0.0 KMKHZ6A.553
! KMKHZ6A.554
ENDDO KMKHZ6A.566
! KMKHZ6A.567
!----------------------------------------------------------------------- KMKHZ6A.568
!! 1. First loop round boundary layer levels. KMKHZ6A.569
!----------------------------------------------------------------------- KMKHZ6A.570
! KMKHZ6A.571
DO K=2,BL_LEVELS KMKHZ6A.572
! KMKHZ6A.573
CALL QSAT(QS(P1),T(P1,K),P(P1,K),P_POINTS) KMKHZ6A.574
! KMKHZ6A.575
DO I=P1,P1+P_POINTS-1 KMKHZ6A.576
! KMKHZ6A.577
!----------------------------------------------------------------------- KMKHZ6A.578
!! 1.4 Calculate the within-layer vertical gradients of cloud liquid KMKHZ6A.579
!! and frozen water for the current layer KMKHZ6A.580
!----------------------------------------------------------------------- KMKHZ6A.581
! KMKHZ6A.582
VIRT_FACTOR = 1.0 + C_VIRTUAL*Q(I,K) - QCL(I,K) - QCF(I,K) KMKHZ6A.583
! KMKHZ6A.584
IF (K. LE. NTML(I)) THEN KMKHZ6A.585
DTLDZ = -GRCP + GRAD_T_ADJ(I) KMKHZ6A.586
ELSE KMKHZ6A.587
DTLDZ = -GRCP KMKHZ6A.588
ENDIF KMKHZ6A.589
! KMKHZ6A.590
DQCLDZ(I,K) = -LGF * ( DTLDZ*DQSDT(I,K) ARN0F405.1421
& + G*QS(I)/(R*T(I,K)*VIRT_FACTOR) ) KMKHZ6A.592
& / ( 1.0 + LCRCP*DQSDT(I,K) ) KMKHZ6A.593
DQCFDZ(I,K) = -FGF * ( DTLDZ*DQSDT(I,K) ARN0F405.1422
& + G*QS(I)/(R*T(I,K)*VIRT_FACTOR) ) KMKHZ6A.595
& / ( 1.0 + LSRCP*DQSDT(I,K) ) KMKHZ6A.596
! KMKHZ6A.597
!----------------------------------------------------------------------- KMKHZ6A.598
!! 1.5 Calculate the cloud liquid and frozen water contents at the KMKHZ6A.599
!! top and bottom of the current layer KMKHZ6A.600
!----------------------------------------------------------------------- KMKHZ6A.601
! KMKHZ6A.602
IF ( QCL(I,K) + QCF(I,K) .GT. 0.0 ) THEN KMKHZ6A.603
CFL(I,K) = CF(I,K) * QCL(I,K)/( QCL(I,K) + QCF(I,K) ) KMKHZ6A.604
CFF(I,K) = CF(I,K) * QCF(I,K)/( QCL(I,K) + QCF(I,K) ) KMKHZ6A.605
ELSE KMKHZ6A.606
CFL(I,K) = 0.0 KMKHZ6A.607
CFF(I,K) = 0.0 KMKHZ6A.608
ENDIF KMKHZ6A.609
! KMKHZ6A.610
IF (CFL(I,K) .GT. 0.0) THEN KMKHZ6A.611
QCL_IC_TOP(I,K) = QCL(I,K) / CFL(I,K) + KMKHZ6A.612
& 0.5*DZL(I,K)*DQCLDZ(I,K) KMKHZ6A.613
QCL_IC_BOT(I,K) = MAX( 0.0 , QCL(I,K) / CFL(I,K) - KMKHZ6A.614
& 0.5*DZL(I,K)*DQCLDZ(I,K) ) KMKHZ6A.615
ELSE KMKHZ6A.616
QCL_IC_TOP(I,K) = 0.0 KMKHZ6A.617
QCL_IC_BOT(I,K) = 0.0 KMKHZ6A.618
ENDIF KMKHZ6A.619
! KMKHZ6A.620
IF (CFF(I,K) .GT. 0.0) THEN KMKHZ6A.621
QCF_IC_TOP(I,K) = QCF(I,K) / CFF(I,K) + KMKHZ6A.622
& 0.5*DZL(I,K)*DQCFDZ(I,K) KMKHZ6A.623
QCF_IC_BOT(I,K) = MAX( 0.0 , QCF(I,K) / CFF(I,K) - KMKHZ6A.624
& 0.5*DZL(I,K)*DQCFDZ(I,K) ) KMKHZ6A.625
ELSE KMKHZ6A.626
QCF_IC_TOP(I,K) = 0.0 KMKHZ6A.627
QCF_IC_BOT(I,K) = 0.0 KMKHZ6A.628
ENDIF KMKHZ6A.629
ENDDO KMKHZ6A.630
ENDDO KMKHZ6A.631
! KMKHZ6A.632
!----------------------------------------------------------------------- KMKHZ6A.633
!! 2. Second loop round boundary layer levels. KMKHZ6A.634
!----------------------------------------------------------------------- KMKHZ6A.635
! KMKHZ6A.636
DO K=2,BL_LEVELS KMKHZ6A.637
KM1 = K-1 KMKHZ6A.638
DO I=P1,P1+P_POINTS-1 KMKHZ6A.639
! KMKHZ6A.640
!----------------------------------------------------------------------- KMKHZ6A.641
!! 2.1 Calculate the jumps of QW and TL across the layer interface KMKHZ6A.642
!! at level k-1/2. KMKHZ6A.643
!----------------------------------------------------------------------- KMKHZ6A.644
! KMKHZ6A.645
DQW = QW(I,K) - QW(I,KM1) ! Used in P243.C2 KMKHZ6A.646
!----------------------------------------------------------------------- KMKHZ6A.647
! ! TL_K_BOT = TL(I,K) - 0.5*DZL(I,K) *(-GRCP) KMKHZ6A.648
! ! TL_KM1_TOP = TL(I,KM1) + 0.5*DZL(I,KM1)*(-GRCP) KMKHZ6A.649
! ! DTL = TL_K_BOT - TL_KM1_TOP so therefore KMKHZ6A.650
!----------------------------------------------------------------------- KMKHZ6A.651
DTL = TL(I,K) - TL(I,KM1) + GRCP/RDZ(I,K) ! Used in P243.C2 KMKHZ6A.652
IF (K .LE. NTML(I)) THEN KMKHZ6A.653
DTL = DTL - GRAD_T_ADJ(I)/RDZ(I,K) KMKHZ6A.654
DQW = DQW - GRAD_Q_ADJ(I)/RDZ(I,K) KMKHZ6A.655
ENDIF KMKHZ6A.656
! KMKHZ6A.657
DQCL = CFL(I,K)*QCL_IC_BOT(I,K) - KMKHZ6A.658
& CFL(I,KM1)*QCL_IC_TOP(I,KM1) KMKHZ6A.659
DQCF = CFF(I,K)*QCF_IC_BOT(I,K) - KMKHZ6A.660
& CFF(I,KM1)*QCF_IC_TOP(I,KM1) KMKHZ6A.661
! KMKHZ6A.662
!----------------------------------------------------------------------- KMKHZ6A.663
!! 2.3 Calculate the buoyancy jumps across the interface between layers KMKHZ6A.664
!! k and k-1 KMKHZ6A.665
!----------------------------------------------------------------------- KMKHZ6A.666
! KMKHZ6A.667
DB(I,K) = G * ( BTM(I,KM1)*DTL + BQM(I,KM1)*DQW + KMKHZ6A.668
& (LCRCP*BTM(I,KM1) - ETAR*BQM(I,KM1)) * DQCL + KMKHZ6A.669
& (LSRCP*BTM(I,KM1) - ETAR*BQM(I,KM1)) * DQCF ) KMKHZ6A.670
! KMKHZ6A.671
DB_SVL(I,K) = G * ( BTM(I,KM1)*DTL + BQM(I,KM1)*DQW ) ARN0F405.1423
! ARN0F405.1424
DB_CLD(I,K) = G * ( BTM_CLD(I,KM1)*DTL + BQM_CLD(I,KM1)*DQW ) KMKHZ6A.672
! KMKHZ6A.673
CHI_S(I,K) = -QCL_IC_TOP(I,KM1) / KMKHZ6A.674
& (A_QSM(I,KM1)*DQW + A_DQSDTM(I,KM1)*DTL) KMKHZ6A.675
ENDDO KMKHZ6A.676
ENDDO KMKHZ6A.677
! KMKHZ6A.678
!----------------------------------------------------------------------- KMKHZ6A.679
!! 3.0a Calculate surface buoyancy flux ARN0F405.1425
!----------------------------------------------------------------------- KMKHZ6A.681
! KMKHZ6A.682
DO I = P1,P1-1+P_POINTS KMKHZ6A.684
ARN0F405.1426
BFLUX_SURF(I) = G * ( BTM(I,NTML(I))*FTL(I,1) + ARN0F405.1427
& BQM(I,NTML(I))*FQW(I,1) ) ARN0F405.1428
ARN0F405.1429
IF ( BFLUX_SURF(I) . GT. 0.0 ) THEN ARN0F405.1430
BFLUX_SAT_SURF(I) = G * ( BTM_CLD(I,NTML(I))*FTL(I,1) + ARN0F405.1431
& BQM_CLD(I,NTML(I))*FQW(I,1) ) ARN0F405.1432
IF ( COUPLED(I) ) THEN ARN0F405.1433
BFLUX_SAT_SURF(I) = G * ( BTM_CLD(I,NTDSC(I))*FTL(I,1) + ARN0F405.1434
& BQM_CLD(I,NTDSC(I))*FQW(I,1) ) ARN0F405.1435
ENDIF KMKHZ6A.696
ELSE ARN0F405.1436
BFLUX_SAT_SURF(I) = 0.0 ARN0F405.1437
ENDIF KMKHZ6A.704
ARN0F405.1438
UNSTABLE(I) = (BFLUX_SURF(I) .GT. 0.0) .AND. (NTML(I) .GT. 1) ARN0F405.1439
ARN0F405.1440
ENDDO KMKHZ6A.709
!----------------------------------------------------------------------- KMKHZ6A.713
!! 3.0aa Calculate Sc layer cloud depth (not cloud fraction weighted). ARN0F405.1441
!! (If DSC(I)=.FALSE. then NTDSC=0 and ZC_DSC remains equal to 0.) ARN0F405.1442
!----------------------------------------------------------------------- KMKHZ6A.715
! KMKHZ6A.716
! Initialise variables ARN0F405.1443
! KMKHZ6A.718
DO I = P1,P1-1+P_POINTS ARN0F405.1444
ZC_DSC(I) = 0.0 ARN0F405.1445
CLOUD_BASE(I) = .TRUE. ARN0F405.1446
ENDDO ARN0F405.1447
! KMKHZ6A.726
DO K=BL_LEVELS,2,-1 ARN0F405.1448
KM1=K-1 ARN0F405.1449
DO I = P1,P1-1+P_POINTS ARN0F405.1450
IF ( (K .EQ. NTDSC(I)) .AND. (CF(I,K) .GT. SC_CFTOL) ) ARN0F405.1451
& CLOUD_BASE(I) = .FALSE. ARN0F405.1452
IF ( (K .LE. NTDSC(I)) .AND. (CF(I,K) .GT. SC_CFTOL) .AND. ARN0F405.1453
& .NOT.CLOUD_BASE(I) ) THEN ARN0F405.1454
ZC_DSC(I) = ZC_DSC(I) + 0.5*DZL(I,K) ARN0F405.1455
IF (CF(I,KM1) .GT. SC_CFTOL) THEN ARN0F405.1456
ZC_DSC(I) = ZC_DSC(I) + 0.5*DZL(I,K) ARN0F405.1457
ELSE KMKHZ6A.736
ZC_DSC(I) = ZC_DSC(I) + ARN0F405.1458
& MIN( 0.5*(DZL(I,K) + DZL(I,KM1)) * CFL(I,K) , ARN0F405.1459
& QCL(I,K) / DQCLDZ(I,K) ) / CF(I,K) ARN0F405.1460
IF (DQCFDZ(I,K) .GT. 0.0) THEN ARN0F405.1461
ZC_DSC(I) = ZC_DSC(I) + ARN0F405.1462
& MIN( 0.5*(DZL(I,K) + DZL(I,KM1)) * CFF(I,K) , ARN0F405.1463
& QCF(I,K) / DQCFDZ(I,K) ) / CF(I,K) ARN0F405.1464
ELSE ARN0F405.1465
ZC_DSC(I) = ZC_DSC(I) + ARN0F405.1466
& 0.5*(DZL(I,K) + DZL(I,KM1)) * CFF(I,K) / CF(I,K) ARN0F405.1467
ENDIF ARN0F405.1468
CLOUD_BASE(I) = .TRUE. ARN0F405.1469
ENDIF KMKHZ6A.744
ENDIF ARN0F405.1470
ENDDO KMKHZ6A.745
ENDDO ARN0F405.1471
ARN0F405.1472
DO I = P1,P1-1+P_POINTS ARN0F405.1473
IF ( DSC(I) .AND. .NOT.CLOUD_BASE(I) ARN0F405.1474
& .AND. (CF(I,1) .GT. SC_CFTOL) ) THEN ARN0F405.1475
ZC_DSC(I) = ZC_DSC(I) + 0.5*DZL(I,1) + ARN0F405.1476
& MIN( 0.5*DZL(I,1)*CFL(I,1) , ARN0F405.1477
& QCL(I,1) / DQCLDZ(I,1) ) / CF(I,1) ARN0F405.1478
IF(DQCFDZ(I,1) .GT. 0.0) THEN ARN0F405.1479
ZC_DSC(I) = ZC_DSC(I) + ARN0F405.1480
& MIN( 0.5*DZL(I,1)*CFF(I,1) , ARN0F405.1481
& QCF(I,1) / DQCFDZ(I,1) ) / CF(I,1) ARN0F405.1482
ELSE ARN0F405.1483
ZC_DSC(I) = ZC_DSC(I) + 0.5*DZL(I,1)*CFF(I,1)/CF(I,1) ARN0F405.1484
ENDIF ARN0F405.1485
CLOUD_BASE(I) = .TRUE. ARN0F405.1486
ENDIF ARN0F405.1487
ENDDO ARN0F405.1488
! !----------------------------------------------------------------- ARN0F405.1489
! ! Check for cloud within inversion. ARN0F405.1490
! !----------------------------------------------------------------- ARN0F405.1491
DO I = P1,P1-1+P_POINTS ARN0F405.1492
IF ( DSC(I) .AND. (CF(I,NTDSC(I)+1) .GT. SC_CFTOL) ) THEN ARN0F405.1493
K=NTDSC(I)+1 ARN0F405.1494
ZC_DSC(I) = ZC_DSC(I) + 0.5*DZL(I,K) ARN0F405.1495
IF (CF(I,K-1) .GT. SC_CFTOL) THEN ARN0F405.1496
ZC_DSC(I) = ZC_DSC(I) + 0.5*DZL(I,K) ARN0F405.1497
ELSE ARN0F405.1498
ZC_DSC(I) = ZC_DSC(I) + ARN0F405.1499
& MIN( 0.5*(DZL(I,K) + DZL(I,K-1)) * CFL(I,K) , ARN0F405.1500
& QCL(I,K) / DQCLDZ(I,K) ) / CF(I,K) ARN0F405.1501
IF (DQCFDZ(I,K) .GT. 0.0) THEN ARN0F405.1502
ZC_DSC(I) = ZC_DSC(I) + ARN0F405.1503
& MIN( 0.5*(DZL(I,K) + DZL(I,K-1)) * CFF(I,K) , ARN0F405.1504
& QCF(I,K) / DQCFDZ(I,K) ) / CF(I,K) ARN0F405.1505
ELSE ARN0F405.1506
ZC_DSC(I) = ZC_DSC(I) + ARN0F405.1507
& 0.5*(DZL(I,K) + DZL(I,K-1)) * CFF(I,K) / CF(I,K) ARN0F405.1508
ENDIF ARN0F405.1509
ENDIF ARN0F405.1510
ENDIF ARN0F405.1511
ENDDO ARN0F405.1512
! !----------------------------------------------------------------- ARN0F405.1513
! ! Layer cloud depth cannot be > the layer depth itself. ARN0F405.1514
! !----------------------------------------------------------------- ARN0F405.1515
DO I = P1,P1-1+P_POINTS ARN0F405.1516
ZC_DSC(I) = MIN( ZC_DSC(I), DSCDEPTH(I) ) ARN0F405.1517
ENDDO ARN0F405.1518
!----------------------------------------------------------------------- KMKHZ6A.747
!! 3.0b Calculate surface mixed layer cloud depth (not cloud fraction ARN0F405.1519
!! weighted) ARN0F405.1520
!----------------------------------------------------------------------- KMKHZ6A.749
! KMKHZ6A.750
DO I = P1,P1-1+P_POINTS KMKHZ6A.753
ZC(I) = 0.0 KMKHZ6A.754
CLOUD_BASE(I) = .TRUE. KMKHZ6A.755
ENDDO KMKHZ6A.756
! KMKHZ6A.757
DO K=BL_LEVELS-1,2,-1 ARN0F405.1521
KM1=K-1 KMKHZ6A.759
DO I = P1,P1-1+P_POINTS KMKHZ6A.760
IF ( (K .EQ. NTML(I)) .AND. (CF(I,K) .GT. 0.0) ) KMKHZ6A.761
& CLOUD_BASE(I) = .FALSE. KMKHZ6A.762
IF ( (K .LE. NTML(I)) .AND. (CF(I,K) .GT. 0.0) .AND. KMKHZ6A.763
& .NOT.CLOUD_BASE(I) ) THEN KMKHZ6A.764
ZC(I) = ZC(I) + 0.5*DZL(I,K) KMKHZ6A.765
IF (CF(I,KM1) .GT. 0.0) THEN KMKHZ6A.766
ZC(I) = ZC(I) + 0.5*DZL(I,K) KMKHZ6A.767
ELSE KMKHZ6A.768
ZC(I) = ZC(I) + KMKHZ6A.769
& MIN( 0.5*(DZL(I,K)+DZL(I,KM1))*CFL(I,K) , ARN0F405.1522
& QCL(I,K) / DQCLDZ(I,K) ) / CF(I,K) ARN0F405.1523
IF (DQCFDZ(I,K) .GT. 0.0) THEN ARN0F405.1524
ZC(I) = ZC(I) + ARN0F405.1525
& MIN( 0.5*(DZL(I,K)+DZL(I,KM1))*CFF(I,K) , ARN0F405.1526
& QCF(I,K) / DQCFDZ(I,K) ) / CF(I,K) ARN0F405.1527
ELSE ARN0F405.1528
ZC(I) = ZC(I) + ARN0F405.1529
& 0.5*(DZL(I,K)+DZL(I,KM1))*CFF(I,K)/CF(I,K) ARN0F405.1530
ENDIF ARN0F405.1531
CLOUD_BASE(I) = .TRUE. KMKHZ6A.775
ENDIF KMKHZ6A.776
ENDIF KMKHZ6A.777
ENDDO KMKHZ6A.778
ENDDO KMKHZ6A.779
DO I = P1,P1-1+P_POINTS KMKHZ6A.781
IF ( (CF(I,1) .GT. 0.0) .AND. .NOT.CLOUD_BASE(I) ) THEN KMKHZ6A.782
ZC(I) = ZC(I) + 0.5*DZL(I,1) + KMKHZ6A.783
& MIN( 0.5*DZL(I,1)*CFL(I,1) , ARN0F405.1532
& QCL(I,1) / DQCLDZ(I,1) ) / CF(I,1) ARN0F405.1533
IF (DQCFDZ(I,1) .GT. 0.0) THEN ARN0F405.1534
ZC(I) = ZC(I) + ARN0F405.1535
& MIN( 0.5*DZL(I,1)*CFF(I,1) , ARN0F405.1536
& QCF(I,1) / DQCFDZ(I,1) ) / CF(I,1) ARN0F405.1537
ELSE ARN0F405.1538
ZC(I) = ZC(I) + 0.5*DZL(I,1)*CFF(I,1) / CF(I,1) ARN0F405.1539
ENDIF ARN0F405.1540
CLOUD_BASE(I) = .TRUE. KMKHZ6A.789
ENDIF KMKHZ6A.790
ENDDO KMKHZ6A.791
! !----------------------------------------------------------------- ARN0F405.1541
! ! Check for cloud within inversion. ARN0F405.1542
! !----------------------------------------------------------------- ARN0F405.1543
DO I = P1,P1-1+P_POINTS ARN0F405.1544
IF ( CF(I,NTML(I)+1) .GT. SC_CFTOL ) THEN ARN0F405.1545
K=NTML(I)+1 ARN0F405.1546
ZC(I) = ZC(I) + 0.5*DZL(I,K) ARN0F405.1547
IF (CF(I,K-1) .GT. SC_CFTOL) THEN ARN0F405.1548
ZC(I) = ZC(I) + 0.5*DZL(I,K) ARN0F405.1549
ELSE ARN0F405.1550
ZC(I) = ZC(I) + ARN0F405.1551
& MIN( 0.5*(DZL(I,K) + DZL(I,K-1)) * CFL(I,K) , ARN0F405.1552
& QCL(I,K) / DQCLDZ(I,K) ) / CF(I,K) ARN0F405.1553
IF (DQCFDZ(I,K) .GT. 0.0) THEN ARN0F405.1554
ZC(I) = ZC(I) + ARN0F405.1555
& MIN( 0.5*(DZL(I,K) + DZL(I,K-1)) * CFF(I,K) , ARN0F405.1556
& QCF(I,K) / DQCFDZ(I,K) ) / CF(I,K) ARN0F405.1557
ELSE ARN0F405.1558
ZC(I) = ZC(I) + ARN0F405.1559
& 0.5*(DZL(I,K) + DZL(I,K-1)) * CFF(I,K) / CF(I,K) ARN0F405.1560
ENDIF ARN0F405.1561
ENDIF ARN0F405.1562
ENDIF ARN0F405.1563
ENDDO ARN0F405.1564
! KMKHZ6A.792
!----------------------------------------------------------------------- KMKHZ6A.793
!! 3.1 Calculate inputs for the top of b.l. entrainment parametrization. ARN0F405.1565
!----------------------------------------------------------------------- KMKHZ6A.796
! KMKHZ6A.797
DO I = P1,P1-1+P_POINTS ARN0F405.1566
ZETA_R_DSC(I) = 0.0 ARN0F405.1567
ALPHA_R_DSC(I) = 0.0 ARN0F405.1568
CHI_S_DSCT(I) = 0.0 ARN0F405.1569
CLD_FACTOR_DSC(I) = 0.0 ARN0F405.1570
BT_DSCT(I) = 0.0 ARN0F405.1571
BTT_DSCT(I) = 0.0 ARN0F405.1572
BTC_DSCT(I) = 0.0 ARN0F405.1573
DB_DSCT(I) = 0.0 ARN0F405.1574
DB_DSCT_CLD(I) = 0.0 ARN0F405.1575
CHI_S_TOP(I) = 0.0 ARN0F405.1576
CLD_FACTOR(I) = 0.0 ARN0F405.1577
BT_TOP(I) = 0.0 ARN0F405.1578
BTT_TOP(I) = 0.0 ARN0F405.1579
BTC_TOP(I) = 0.0 ARN0F405.1580
DB_TOP(I) = 0.0 ARN0F405.1581
DB_TOP_CLD(I) = 0.0 ! default required if COUPLED ARN0F405.1582
Z_CLD(I) = 0.0 ARN0F405.1583
Z_CLD_DSC(I) = 0.0 ARN0F405.1584
ENDDO ARN0F405.1585
! KMKHZ6A.799
DO K = 1,BL_LEVELS ARN0F405.1586
DO I = P1,P1-1+P_POINTS KMKHZ6A.800
IF ( K .LE. NTML(I)+1 ) THEN ARN0F405.1587
! !--------------------------------------------------- ARN0F405.1588
! ! Calculation of cloud fraction weighted ARN0F405.1589
! ! thickness of cloud in the surface mixed layer ARN0F405.1590
! !--------------------------------------------------- ARN0F405.1591
Z_CLD(I) = Z_CLD(I) + ARN0F405.1592
& CF(I,K) * 0.5 * DZL(I,K) + ARN0F405.1593
& MIN( CFL(I,K) * 0.5 * DZL(I,K) , ARN0F405.1594
& QCL(I,K) / DQCLDZ(I,K) ) ARN0F405.1595
IF (DQCFDZ(I,K) .GT. 0.0) THEN ARN0F405.1596
Z_CLD(I) = Z_CLD(I) + ARN0F405.1597
& MIN( CFF(I,K) * 0.5 * DZL(I,K) , ARN0F405.1598
& QCF(I,K) / DQCFDZ(I,K) ) ARN0F405.1599
ELSE ARN0F405.1600
Z_CLD(I) = Z_CLD(I) + ARN0F405.1601
& CFF(I,K) * 0.5 * DZL(I,K) ARN0F405.1602
ENDIF ARN0F405.1603
ENDIF ARN0F405.1604
! ARN0F405.1605
IF ( DSC(I) .AND. K.LE.NTDSC(I)+1 .AND. ARN0F405.1606
& ( COUPLED(I) .OR. ARN0F405.1607
& Z_HALF(I,K+1).GE.ZHSC(I)-ZC_DSC(I) ) ) THEN ARN0F405.1608
! !---------------------------------------------------- ARN0F405.1609
! ! Calculation of cloud fraction weighted thickness of ARN0F405.1610
! ! cloud in the DSC layer (or to the surface if COUPLED) ARN0F405.1611
! !---------------------------------------------------- ARN0F405.1612
Z_CLD_DSC(I) = Z_CLD_DSC(I) + ARN0F405.1613
& CF(I,K) * 0.5 * DZL(I,K) + ARN0F405.1614
& MIN( CFL(I,K) * 0.5 * DZL(I,K) , ARN0F405.1615
& QCL(I,K) / DQCLDZ(I,K) ) ARN0F405.1616
IF (DQCFDZ(I,K) .GT. 0.0) THEN ARN0F405.1617
Z_CLD_DSC(I) = Z_CLD_DSC(I) + ARN0F405.1618
& MIN( CFF(I,K) * 0.5 * DZL(I,K) , ARN0F405.1619
& QCF(I,K) / DQCFDZ(I,K) ) ARN0F405.1620
ELSE ARN0F405.1621
Z_CLD_DSC(I) = Z_CLD_DSC(I) + ARN0F405.1622
& CFF(I,K) * 0.5 * DZL(I,K) ARN0F405.1623
ENDIF ARN0F405.1624
ENDIF ARN0F405.1625
! ARN0F405.1626
IF (K .EQ. NTML(I) .AND. .NOT.COUPLED(I) ) THEN ARN0F405.1627
! !------------------------------------------------------ ARN0F405.1628
! ! Calculation of SML inputs. If COUPLED then these are ARN0F405.1629
! ! not used (as no entrainment is then applied at ZH) ARN0F405.1630
! !------------------------------------------------------ ARN0F405.1631
CHI_S_TOP(I) = MAX( 0.0, MIN( CHI_S(I,K+1), 1.) ) ARN0F405.1632
CLD_FACTOR(I) = MAX( 0.0 , CF(I,K)-CF(I,K+1) ) ARN0F405.1633
BT_TOP(I) = G * BTM(I,K) ARN0F405.1634
BTT_TOP(I) = G * BTM_CLD(I,K) ARN0F405.1635
BTC_TOP(I) = BTT_TOP(I) ARN0F405.1636
DB_TOP(I) = DB(I,K+1) ARN0F405.1637
DB_TOP_CLD(I) = DB_CLD(I,K+1) ARN0F405.1638
ENDIF ARN0F405.1639
IF (K .EQ. NTDSC(I)) THEN ARN0F405.1640
! !--------------------------------------------------- ARN0F405.1641
! ! Calculation of DSC inputs ARN0F405.1642
! ! (if DSC=.FALSE. then K never equals NTDSC(=0)) ARN0F405.1643
! !--------------------------------------------------- ARN0F405.1644
CHI_S_DSCT(I) = MAX( 0.0, MIN( CHI_S(I,K+1), 1.) ) ARN0F405.1645
CLD_FACTOR_DSC(I) = MAX( 0.0 , CF(I,K)-CF(I,K+1) ) ARN0F405.1646
BT_DSCT(I) = G * BTM(I,K) ARN0F405.1647
BTT_DSCT(I) = G * BTM_CLD(I,K) ARN0F405.1648
BTC_DSCT(I) = BTT_DSCT(I) ARN0F405.1649
DB_DSCT(I) = DB(I,K+1) ARN0F405.1650
DB_DSCT_CLD(I) = DB_CLD(I,K+1) ARN0F405.1651
ENDIF ARN0F405.1652
ENDDO ARN0F405.1653
ENDDO ARN0F405.1654
! !----------------------------------------------------------------- ARN0F405.1655
! ! Next those terms dependent on the presence of buoyancy reversal. ARN0F405.1656
! !"---------------------------------------------------------------- ARN0F405.1657
DO I = P1,P1-1+P_POINTS ARN0F405.1658
Z_CLD(I) = MIN( Z_CLD(I), ZH(I) ) ARN0F405.1659
Z_CLD_DSC(I) = MIN( Z_CLD_DSC(I), ZHSC(I) ) ARN0F405.1660
! !--------------------------------------------------------------- ARN0F405.1661
! ! First the surface mixed layer. ARN0F405.1662
! !--------------------------------------------------------------- ARN0F405.1663
IF ( COUPLED(I) ) THEN ARN0F405.1664
ZETA_S(I) = 1.0 - Z_CLD_DSC(I) / ZHSC(I) ARN0F405.1665
ZETA_R(I) = 1.0 - ZC_DSC(I) / ZHSC(I) ARN0F405.1666
ELSE ARN0F405.1667
ZETA_S(I) = 1.0 - Z_CLD(I) / ZH(I) ARN0F405.1668
ZETA_R(I) = 1.0 - ZC(I) / ZH(I) ARN0F405.1669
ENDIF ARN0F405.1670
! ARN0F405.1671
IF (DB_TOP_CLD(I) .GE. 0.0) THEN ARN0F405.1672
! !-------------------------------------------------- ARN0F405.1673
! ! i.e. no buoyancy reversal (or default if COUPLED) ARN0F405.1674
! !-------------------------------------------------- ARN0F405.1675
ALPHA_R(I) = 0.2 ARN0F405.1676
DB_TOP_CLD(I) = 0.0 ARN0F405.1677
ELSE ARN0F405.1678
! !---------------------------- ARN0F405.1679
! ! IF (DB_TOP_CLD(I) .LT. 0.0) ARN0F405.1680
! ! i.e. buoyancy reversal ARN0F405.1681
! !---------------------------- ARN0F405.1682
DB_TOP_CLD(I) = -DB_TOP_CLD(I) * CLD_FACTOR(I) ARN0F405.1683
D_SIEMS = MAX( 0.0, ARN0F405.1684
& CHI_S_TOP(I) * DB_TOP_CLD(I) / (DB_TOP(I)+1.E-14) ) ARN0F405.1685
ZETA_R(I) = MIN( ZETA_R(I)+10.0*(1.0-ZETA_R(I)) ARN0F405.1686
& *D_SIEMS, 1.0 ) ARN0F405.1687
ALPHA_R(I) = MIN( 0.2+10.0*(1.0-0.2)*D_SIEMS, 1.0 ) ARN0F405.1688
ENDIF ARN0F405.1689
! !--------------------------------------------------------------- ARN0F405.1690
! ! Now the decoupled Sc layer (DSC). ARN0F405.1691
! !--------------------------------------------------------------- ARN0F405.1692
IF (DSC(I)) THEN ARN0F405.1693
IF ( COUPLED(I) ) THEN ARN0F405.1694
ZETA_R_DSC(I) = 1.0 - ZC_DSC(I) / ZHSC(I) ARN0F405.1695
ELSE ARN0F405.1696
ZETA_R_DSC(I) = 1.0 - ZC_DSC(I) / DSCDEPTH(I) ARN0F405.1697
ENDIF ARN0F405.1698
ARN0F405.1699
IF (DB_DSCT_CLD(I) .GE. 0.0) THEN ARN0F405.1700
! !---------------------------- ARN0F405.1701
! ! i.e. no buoyancy reversal ARN0F405.1702
! !---------------------------- ARN0F405.1703
ALPHA_R_DSC(I) = 0.2 ARN0F405.1704
DB_DSCT_CLD(I) = 0.0 ARN0F405.1705
ELSE ARN0F405.1706
! !---------------------------- ARN0F405.1707
! ! IF (DB_DSCT_CLD(I) .LT. 0.0) ARN0F405.1708
! ! i.e. buoyancy reversal ARN0F405.1709
! !---------------------------- ARN0F405.1710
DB_DSCT_CLD(I) = -DB_DSCT_CLD(I) * CLD_FACTOR_DSC(I) ARN0F405.1711
D_SIEMS = MAX( 0.0, ARN0F405.1712
& CHI_S_DSCT(I) * DB_DSCT_CLD(I) / (DB_DSCT(I)+1.E-14) ) ARN0F405.1713
ZETA_R_DSC(I) = MIN( ZETA_R_DSC(I)+10.0*(1.0-ZETA_R_DSC(I)) ARN0F405.1714
& *D_SIEMS, 1.0 ) ARN0F405.1715
ALPHA_R_DSC(I) = MIN( 0.2+10.0*(1.0-0.2)*D_SIEMS, 1.0 ) ARN0F405.1716
ENDIF ARN0F405.1717
ENDIF ARN0F405.1718
ENDDO ARN0F405.1719
! ARN0F405.1720
!----------------------------------------------------------------------- ARN0F405.1721
!! 4. Calculate the radiative flux change across cloud top for mixed ARN0F405.1722
!! layer to ZH. If there is a decoupled layer above, restrict ARN0F405.1723
!! search for maximum divergence to below NTML+1. This may ARN0F405.1724
!! introduce errors if NTML changes during radiative timestep but ARN0F405.1725
!! can't be helped. ARN0F405.1726
!----------------------------------------------------------------------- ARN0F405.1727
! Initialise variables ARN0F405.1728
!------------------------------ ARN0F405.1729
DO I = P1,P1-1+P_POINTS ARN0F405.1730
K_CLOUD_TOP(I) = 0 KMKHZ6A.801
DF_TOP_OVER_CP(I) = 0.0 KMKHZ6A.802
ENDDO KMKHZ6A.804
! KMKHZ6A.805
DO K = 1,BL_LEVELS-1 KMKHZ6A.806
DO I = P1,P1-1+P_POINTS KMKHZ6A.807
ARN0F405.1731
DF_OVER_CP(I,K) = DELTAP(I,K)/G * RAD_HR(I,K) ARN0F405.1732
ARN0F405.1733
K_RAD_SMLT = NTML(I)+1 ARN0F405.1734
IF ( .NOT.DSC(I) ) K_RAD_SMLT = NTML(I)+2 ARN0F405.1735
! !------------------------------------------------------------- KMKHZ6A.808
! ! Find the layer with the greatest radiative flux jump below ARN0F405.1736
! ! K_RAD_SMLT and assume that this is the top of the SML. ARN0F405.1737
! !------------------------------------------------------------- KMKHZ6A.811
IF (DF_OVER_CP(I,K) .GT. DF_TOP_OVER_CP(I) ARN0F405.1738
& .AND. K .LE. K_RAD_SMLT ) THEN ARN0F405.1739
K_CLOUD_TOP(I) = K KMKHZ6A.814
DF_TOP_OVER_CP(I) = DF_OVER_CP(I,K) KMKHZ6A.815
ENDIF KMKHZ6A.816
ARN0F405.1740
ENDDO KMKHZ6A.817
ENDDO KMKHZ6A.818
! !----------------------------------------------------------------- ARN0F405.1741
! ! If cloud top extends into the level above NTML (so that the ARN0F405.1742
! ! radiative divergence maximum is at NTML+1) add on RAD_HR from ARN0F405.1743
! ! NTML if there is additional cloud top cooling there. ARN0F405.1744
! !----------------------------------------------------------------- ARN0F405.1745
DO I = P1,P1-1+P_POINTS KMKHZ6A.820
IF ( K_CLOUD_TOP(I) .EQ. NTML(I)+1 ) THEN ARN0F405.1746
K=NTML(I) ARN0F405.1747
DF_OVER_CP(I,K) = DELTAP(I,K)/G * RAD_HR(I,K) ARN0F405.1748
IF (DF_OVER_CP(I,K) .GT. 0.0 ) ARN0F405.1749
& DF_TOP_OVER_CP(I) = DF_TOP_OVER_CP(I) + DF_OVER_CP(I,K) ARN0F405.1750
ENDIF KMKHZ6A.823
ENDDO KMKHZ6A.824
! KMKHZ6A.826
!----------------------------------------------------------------------- KMKHZ6A.827
!! 5.1 Subroutine EXCF_NL. KMKHZ6A.828
!----------------------------------------------------------------------- KMKHZ6A.829
! KMKHZ6A.830
CALL EXCF_NL ( KMKHZ6A.831
& P_FIELD,P1,P_POINTS,BL_LEVELS, KMKHZ6A.832
& NTML,CF, ARN0F405.1751
& RDZ,ZH,Z_UV,Z_TQ,RHO_UV,RHO_TQ, KMKHZ6A.834
& Z0M,V_S,FB_SURF,DB_TOP, KMKHZ6A.835
& BFLUX_SURF,BFLUX_SAT_SURF,ZETA_S,BT_TOP,BTT_TOP, ARN0F405.1752
& DF_TOP_OVER_CP,ZETA_R,ALPHA_R,BTC_TOP, ARN0F405.1753
& DB_TOP_CLD,CHI_S_TOP,ZC, KMKHZ6A.838
& RHOKM(1,2),RHOKH(1,2),RHOKM_TOP(1,2),RHOKH_TOP(1,2), ARN0F405.1754
& ZHSC,DSCDEPTH,NTDSC,DB_DSCT,SVL, ARN0F405.1755
& BT_DSCT,BTT_DSCT, ARN0F405.1756
& DF_DSCT_OVER_CP,ZETA_R_DSC,ALPHA_R_DSC,BTC_DSCT, ARN0F405.1757
& DB_DSCT_CLD,CHI_S_DSCT,ZC_DSC,COUPLED, ARN0F405.1758
& LTIMER ARN0F405.1759
&) KMKHZ6A.840
! ARN0F405.1760
!----------------------------------------------------------------------- ARN0F405.1761
!! 5.2 Diagnose boundary layer type. ARN0F405.1762
! Six different types are considered: ARN0F405.1763
! 1 - Stable b.l. ARN0F405.1764
! 2 - Stratocumulus over a stable surface layer. ARN0F405.1765
! 3 - Well mixed b.l. (possibly with stratocumulus). ARN0F405.1766
! 4 - Decoupled stratocumulus (not over cumulus). ARN0F405.1767
! 5 - Decoupled stratocumulus over cumulus. ARN0F405.1768
! 6 - Cumulus capped b.l. ARN0F405.1769
!----------------------------------------------------------------------- ARN0F405.1770
! First initialise the type variables and set the diagnostic ZHT. ARN0F405.1771
! ARN0F405.1772
DO I = P1,P1-1+P_POINTS ARN0F405.1773
BL_TYPE_1(I) = 0.0 ARN0F405.1774
BL_TYPE_2(I) = 0.0 ARN0F405.1775
BL_TYPE_3(I) = 0.0 ARN0F405.1776
BL_TYPE_4(I) = 0.0 ARN0F405.1777
BL_TYPE_5(I) = 0.0 ARN0F405.1778
BL_TYPE_6(I) = 0.0 ARN0F405.1779
ZHT(I) = MAX( ZH(I) , ZHSC(I) ) ARN0F405.1780
ENDDO ARN0F405.1781
DO I = P1,P1-1+P_POINTS ARN0F405.1782
IF (.NOT.UNSTABLE(I) .AND. .NOT.DSC(I) .AND. ARN0F405.1783
& .NOT.CUMULUS(I)) THEN ARN0F405.1784
! Stable b.l. ARN0F405.1785
BL_TYPE_1(I) = 1.0 ARN0F405.1786
ELSEIF (.NOT.UNSTABLE(I) .AND. DSC(I) .AND. ARN0F405.1787
& .NOT.CUMULUS(I)) THEN ARN0F405.1788
! Stratocumulus over a stable surface layer ARN0F405.1789
BL_TYPE_2(I) = 1.0 ARN0F405.1790
ELSEIF (UNSTABLE(I) .AND. .NOT.CUMULUS(I) .AND. ARN0F405.1791
& (.NOT.DSC(I) .OR. COUPLED(I))) THEN ARN0F405.1792
! Well mixed b.l. (possibly with stratocumulus) ARN0F405.1793
BL_TYPE_3(I) = 1.0 ARN0F405.1794
ELSEIF (UNSTABLE(I) .AND. DSC(I) .AND. .NOT.CUMULUS(I)) THEN ARN0F405.1795
! Decoupled stratocumulus (not over cumulus) ARN0F405.1796
BL_TYPE_4(I) = 1.0 ARN0F405.1797
ELSEIF (DSC(I) .AND. CUMULUS(I)) THEN ARN0F405.1798
! Decoupled stratocumulus over cumulus ARN0F405.1799
BL_TYPE_5(I) = 1.0 ARN0F405.1800
ELSEIF (.NOT.DSC(I) .AND. CUMULUS(I)) THEN ARN0F405.1801
! Cumulus capped b.l. ARN0F405.1802
BL_TYPE_6(I) = 1.0 ARN0F405.1803
ENDIF ARN0F405.1804
ENDDO ARN0F405.1805
KMKHZ6A.841
IF (LTIMER) THEN KMKHZ6A.842
CALL TIMER('KMKHZ ',4) KMKHZ6A.843
ENDIF KMKHZ6A.844
KMKHZ6A.845
RETURN KMKHZ6A.846
END KMKHZ6A.847
*ENDIF KMKHZ6A.848