*IF DEF,A03_3A,OR,DEF,A03_5A AJS1F401.1480
C ******************************COPYRIGHT****************************** GTS2F400.8875
C (c) CROWN COPYRIGHT 1995, METEOROLOGICAL OFFICE, All Rights Reserved. GTS2F400.8876
C GTS2F400.8877
C Use, duplication or disclosure of this code is subject to the GTS2F400.8878
C restrictions as set forth in the contract. GTS2F400.8879
C GTS2F400.8880
C Meteorological Office GTS2F400.8881
C London Road GTS2F400.8882
C BRACKNELL GTS2F400.8883
C Berkshire UK GTS2F400.8884
C RG12 2SZ GTS2F400.8885
C GTS2F400.8886
C If no contract has been raised with this copy of the code, the use, GTS2F400.8887
C duplication or disclosure of it is strictly prohibited. Permission GTS2F400.8888
C to do so must first be obtained in writing from the Head of Numerical GTS2F400.8889
C Modelling at the above address. GTS2F400.8890
C ******************************COPYRIGHT****************************** GTS2F400.8891
C GTS2F400.8892
CLL SUBROUTINE SFL_INT------------------------------------------------ SFLINT3A.3
CLL SFLINT3A.4
CLL Purpose: To calculate interpolation coefficients for 10m winds SFLINT3A.5
CLL and 1.5m temperature/specific humidity diagnostics SFLINT3A.6
CLL using a generalisation of the method of Geleyn (1988). SFLINT3A.7
CLL SFLINT3A.8
CLL Suitable for single column use (via *IF definition IBM). SFLINT3A.9
CLL SFLINT3A.10
CLL Model Modification history: SFLINT3A.11
CLL version Date SFLINT3A.12
CLL SFLINT3A.13
CLL 3.4 18/10/94 *DECK inserted into UM version 3.4. S Jackson SFLINT3A.14
CLL SFLINT3A.15
CLL 4.0 30/12/94 Modified calculation of 10m wind interpolation ARN1F400.83
CLL factor when effective roughness length used; ARN1F400.84
CLL 10m wind assumed to lie on "local" profile at ARN1F400.85
CLL height z0m+10 metres above the surface. ARN1F400.86
CLL R.N.B.Smith ARN1F400.87
CLL 4.2 Oct. 96 T3E migration - *DEF CRAY removed GSS1F402.81
CLL S J Swarbrick GSS1F402.82
CLL ARN1F400.88
CLL Programming standard: Unified Model Documentation Paper No 4, SFLINT3A.16
CLL Version 2, dated 18/1/90. SFLINT3A.17
CLL SFLINT3A.18
CLL Logical component covered: Part of P243. SFLINT3A.19
CLL SFLINT3A.20
CLL System Task: SFLINT3A.21
CLL SFLINT3A.22
CLL External Documentation: UMDP No.24 SFLINT3A.23
CLL SFLINT3A.24
CLL--------------------------------------------------------------------- SFLINT3A.25
C*L Arguments :- SFLINT3A.26
SUBROUTINE SFL_INT ( 3,8SFLINT3A.27
& P_POINTS,U_POINTS,RIB,Z1,Z0M,Z0M_EFF,Z0H,Z0F,CD_STD,CD,CH ARN1F400.89
&,RESFT,WIND_PROFILE_FACTOR ARN1F400.90
&,CDR10M,CHR1P5M,CER1P5M ARN1F400.91
+,SU10,SV10,ST1P5,SQ1P5,LTIMER SFLINT3A.30
+) SFLINT3A.31
IMPLICIT NONE SFLINT3A.32
INTEGER SFLINT3A.33
+ P_POINTS ! IN No. of P-grid points to be processed. SFLINT3A.34
+,U_POINTS ! IN No. of UV-grid points to be processed. SFLINT3A.35
REAL SFLINT3A.36
+ RIB(P_POINTS) ! IN Bulk Richardson number for SFLINT3A.37
C ! lowest layer. SFLINT3A.38
+,Z1(P_POINTS) ! IN Height of lowest model level (m). SFLINT3A.39
+,Z0M(P_POINTS) ! IN Roughness length for momentum (m). SFLINT3A.40
+,Z0M_EFF(P_POINTS) ! IN Effective roughness length for SFLINT3A.41
C ! momentum (m). SFLINT3A.42
+,Z0H(P_POINTS) ! IN Roughness length for heat and SFLINT3A.43
C ! moisture (m). SFLINT3A.44
+,Z0F(P_POINTS) ! IN Roughness length in the free SFLINT3A.45
C ! convective limit (m). SFLINT3A.46
&,CD_STD(P_POINTS) ! IN Surface drag coefficient for shear stress ARN1F400.92
C ! only, i.e. without orographic part of drag ARN1F400.93
&,CD(P_POINTS) ! IN Effective surface drag coefficient, ARN1F400.94
C ! including orographic part of drag ARN1F400.95
&,CH(P_POINTS) ! IN Surface transfer coefficient for heat and ARN1F400.96
C ! moisture. SFLINT3A.49
+,RESFT(P_POINTS) ! IN Total resistance factor for moisture SFLINT3A.50
C ! transfer from the surface. SFLINT3A.51
&,WIND_PROFILE_FACTOR(P_POINTS) ARN1F400.97
C ! IN Factor for converting effective friction ARN1F400.98
C ! velocity to local one. ARN1F400.99
LOGICAL SFLINT3A.52
+ SU10 ! IN 10m U-wind diagnostic flag SFLINT3A.53
+,SV10 ! IN 10m V-wind diagnostic flag SFLINT3A.54
+,ST1P5 ! IN screen temp diagnostic flag SFLINT3A.55
+,SQ1P5 ! IN screen specific humidity SFLINT3A.56
C ! diagnostic flag SFLINT3A.57
+,LTIMER ! IN TIMER diagnostics flag SFLINT3A.58
C Output variables SFLINT3A.59
C SFLINT3A.60
REAL SFLINT3A.61
+ CDR10M(U_POINTS) ! OUT interpolation coefficicent for 10m wind SFLINT3A.62
+,CHR1P5M(P_POINTS) ! OUT Interpolation coefficient for 1.5m SFLINT3A.63
C ! temperature SFLINT3A.64
+,CER1P5M(P_POINTS) ! OUT Interpolation coefficient for 1.5m SFLINT3A.65
C ! specific humidity SFLINT3A.66
C* SFLINT3A.67
C*L--------------------------------------------------------------------- SFLINT3A.68
EXTERNAL TIMER SFLINT3A.69
C* SFLINT3A.70
C*L--------------------------------------------------------------------- SFLINT3A.71
C Local and other symbolic constants :- SFLINT3A.72
*CALL C_VKMAN
SFLINT3A.73
REAL Z1P5M,Z10M SFLINT3A.74
PARAMETER ( SFLINT3A.75
+ Z1P5M = 1.5 ! for diagnosis of screen values of temperature SFLINT3A.76
C ! and humidity (assumed to be at 1.5m). SFLINT3A.77
+,Z10M = 10.0 ! for diagnosis of 10m winds. SFLINT3A.78
+) SFLINT3A.79
C SFLINT3A.80
C Define local storage. SFLINT3A.81
C SFLINT3A.82
C (a) Local work arrays. SFLINT3A.83
C SFLINT3A.84
REAL SFLINT3A.86
+ Z1E(P_POINTS) ! Level 1 height + Z0M_EFF SFLINT3A.87
+,SQRTCD(P_POINTS) ! Square root of CD SFLINT3A.88
C SFLINT3A.95
C (b) Scalars. SFLINT3A.96
C SFLINT3A.97
REAL SFLINT3A.98
+ Z1S ! Level 1 height + Z0M_EFF - Z0H SFLINT3A.99
+,Z1P5ME ! Z1P5M + Z0H ARN1F400.100
+,Z10ME ! Z10M + Z0M ARN1F400.101
+,SQRTCD_K ! Temporary storage in calc of 1.5 amd 10m diags SFLINT3A.102
+,Z_OVER_Z1 ! Temporary storage in calc of 1.5 amd 10m diags SFLINT3A.103
+,CDNZ ! Neutral drag coef. for momentum @ 10m SFLINT3A.104
+,CDNZ1 ! Neutral drag coef. for momentum @ level1 SFLINT3A.105
+,CHNZ ! Neutral drag coef. for heat/moisture @ 1.5m SFLINT3A.106
+,CHNZ1 ! Neutral drag coef. for heat/moisture @ level1 SFLINT3A.107
&,CDTEMP1 ! Workspace in calc of interpolation coeffs. ARN1F400.102
&,CDTEMP2 ! Workspace in calc of interpolation coeffs. ARN1F400.103
&,CDTEMP3 ! Workspace in calc of interpolation coeffs. ARN1F400.104
SFLINT3A.109
INTEGER SFLINT3A.110
+ I ! Loop counter (horizontal field index). SFLINT3A.111
C* SFLINT3A.112
IF (LTIMER) THEN SFLINT3A.113
CALL TIMER
('SFL_INT ',3) SFLINT3A.114
ENDIF SFLINT3A.115
C SFLINT3A.116
C----------------------------------------------------------------------- SFLINT3A.117
CL 1. This routine uses a generalised formulation of Geleyn (1988) SFLINT3A.118
CL to interpolate to screen and 10m height using surface and bottom SFLINT3A.119
CL level values as well as roughness lengths for heat and momentum. SFLINT3A.120
CL Start of main loop. Set up variables needed later (eg height of SFLINT3A.121
CL envelope orography). SFLINT3A.122
C----------------------------------------------------------------------- SFLINT3A.123
SFLINT3A.124
DO I=1,P_POINTS SFLINT3A.125
SQRTCD(I) = SQRT(CD(I)) SFLINT3A.126
Z1E(I) = Z1(I) + Z0M_EFF(I) SFLINT3A.127
ENDDO SFLINT3A.128
SFLINT3A.129
C----------------------------------------------------------------------- SFLINT3A.130
CL 2. If selected calculate interpolation coefficient for 10m winds. SFLINT3A.131
C----------------------------------------------------------------------- SFLINT3A.132
C SFLINT3A.133
IF(SU10.OR.SV10) THEN SFLINT3A.134
DO I=1,P_POINTS SFLINT3A.135
Z10ME = Z10M + Z0M(I) SFLINT3A.136
CDNZ = LOG( Z10ME / Z0M(I) ) SFLINT3A.137
CDNZ1 = LOG( Z1E(I) / Z0M_EFF(I) ) ARN1F400.105
SQRTCD_K = SQRTCD(I) / VKMAN SFLINT3A.139
Z_OVER_Z1 = Z10M / Z1(I) ARN1F400.106
SFLINT3A.141
IF (RIB(I).GE.0.0) THEN SFLINT3A.142
C ARN1F400.107
C Stable case ARN1F400.108
C ARN1F400.109
CDR10M(I) = Z_OVER_Z1 + SQRTCD_K * SFLINT3A.143
+ (CDNZ - Z_OVER_Z1*CDNZ1) SFLINT3A.144
ELSE SFLINT3A.145
C ARN1F400.110
C Unstable Case ARN1F400.111
C ARN1F400.112
CDTEMP1 = EXP( 1.0 / SQRTCD_K ) ARN1F400.113
CDTEMP2 = Z1E(I) / Z0M_EFF(I) ARN1F400.114
CDTEMP3 = LOG( ( ( Z1E(I) - Z0M(I) ) * CDTEMP1 - ARN1F400.115
& ( Z0M_EFF(I) - Z0M(I) ) * CDTEMP2 ) / ARN1F400.116
& ( ( Z1E(I) - Z10ME ) * CDTEMP1 - ARN1F400.117
& ( Z0M_EFF(I) - Z10ME ) * CDTEMP2 ) ) ARN1F400.118
CDR10M(I) = SQRTCD_K * ( CDNZ + CDTEMP3 ) ARN1F400.119
ENDIF SFLINT3A.150
CDR10M(I) = CDR10M(I) * WIND_PROFILE_FACTOR(I) ARN1F400.120
ENDDO SFLINT3A.151
ENDIF SFLINT3A.152
C SFLINT3A.153
C----------------------------------------------------------------------- SFLINT3A.154
CL 3. If selected calculate interpolation coefficient for 1.5m screen SFLINT3A.155
CL temperature and specific humidity. SFLINT3A.156
C----------------------------------------------------------------------- SFLINT3A.157
C SFLINT3A.158
IF (ST1P5.OR.SQ1P5) THEN SFLINT3A.159
DO I=1,P_POINTS SFLINT3A.160
C SFLINT3A.161
C variables to be used later SFLINT3A.162
Z1S = Z1E(I) - Z0H(I) SFLINT3A.163
Z1P5ME = Z1P5M + Z0H(I) SFLINT3A.164
CHNZ = LOG( Z1P5ME / Z0H(I) ) SFLINT3A.165
CHNZ1 = LOG( Z1E(I) / Z0H(I) ) SFLINT3A.166
SQRTCD(I) = SQRT(CD_STD(I)) ARN1F400.121
SQRTCD_K =0.0 SFLINT3A.167
IF (SQRTCD(I) .GT. 0.0) SQRTCD_K = CH(I) / (SQRTCD(I) * VKMAN) SFLINT3A.168
Z_OVER_Z1 = Z1P5M / Z1S SFLINT3A.169
SFLINT3A.170
C SFLINT3A.171
C Stable case SFLINT3A.172
C ARN1F400.122
IF (RIB(I).GE.0.0) THEN SFLINT3A.173
CHR1P5M(I) = Z_OVER_Z1 + SQRTCD_K * SFLINT3A.174
+ (CHNZ - Z_OVER_Z1 * CHNZ1) SFLINT3A.175
ELSE SFLINT3A.176
C SFLINT3A.177
C Unstable Case SFLINT3A.178
C ARN1F400.123
CDTEMP1 = EXP(1.0 / SQRTCD_K) ARN1F400.124
CDTEMP2 = ( Z1P5M * Z1E(I) + ARN1F400.125
& Z0H(I) * ( Z1S - Z1P5M ) * CDTEMP1 ) / ARN1F400.126
& ( Z1S * Z1P5ME ) ARN1F400.127
CHR1P5M(I) = 1.0 - SQRTCD_K * LOG ( CDTEMP2 ) ARN1F400.128
ENDIF SFLINT3A.183
SFLINT3A.184
C SFLINT3A.185
CER1P5M(I) = ( CHR1P5M(I) - 1.0 ) * RESFT(I) ! P243.123 SFLINT3A.186
ENDDO SFLINT3A.187
ENDIF SFLINT3A.188
C SFLINT3A.189
IF (LTIMER) THEN SFLINT3A.190
CALL TIMER
('SFL_INT ',4) SFLINT3A.191
ENDIF SFLINT3A.192
RETURN SFLINT3A.193
END SFLINT3A.194
*ENDIF SFLINT3A.195