*IF DEF,A03_7A SFFLUX7A.2
C *****************************COPYRIGHT****************************** SFFLUX7A.3
C (c) CROWN COPYRIGHT 1997, METEOROLOGICAL OFFICE, All Rights Reserved. SFFLUX7A.4
C SFFLUX7A.5
C Use, duplication or disclosure of this code is subject to the SFFLUX7A.6
C restrictions as set forth in the contract. SFFLUX7A.7
C SFFLUX7A.8
C Meteorological Office SFFLUX7A.9
C London Road SFFLUX7A.10
C BRACKNELL SFFLUX7A.11
C Berkshire UK SFFLUX7A.12
C RG12 2SZ SFFLUX7A.13
C SFFLUX7A.14
C If no contract has been raised with this copy of the code, the use, SFFLUX7A.15
C duplication or disclosure of it is strictly prohibited. Permission SFFLUX7A.16
C to do so must first be obtained in writing from the Head of Numerical SFFLUX7A.17
C Modelling at the above address. SFFLUX7A.18
C ******************************COPYRIGHT****************************** SFFLUX7A.19
SFFLUX7A.20
!----------------------------------------------------------------------- SFFLUX7A.21
! SFFLUX7A.22
! Subroutines SF_FLUX_LAND and SF_FLUX_SEA to calculate explicit surface SFFLUX7A.23
! fluxes of heat and moisture SFFLUX7A.24
! SFFLUX7A.25
!----------------------------------------------------------------------- SFFLUX7A.26
SFFLUX7A.27
! SUBROUTINE SF_FLUX_LAND------------------------------------------- SFFLUX7A.28
! SFFLUX7A.29
! Calculate explicit surface fluxes of heat and moisture over SFFLUX7A.30
! land tiles SFFLUX7A.31
! SFFLUX7A.32
! ------------------------------------------------------------------ SFFLUX7A.33
SUBROUTINE SF_FLUX_LAND ( 2,2SFFLUX7A.34
& P_FIELD,LAND_FIELD,TILE_PTS,LAND_INDEX,TILE_INDEX, SFFLUX7A.35
& ASHTF,LH,QS1,QSTAR,QW_1,RADNET,RESFT,RHOKH_1,TILE_FRAC, SFFLUX7A.36
& TL_1,TS1,TSTAR,Z0H,Z0M_EFF,Z1_TQ, SFFLUX7A.37
& FQW_1_GB,FTL_1_GB, SFFLUX7A.38
& ALPHA1,FQW_1,FTL_1,RHOKPM,LTIMER SFFLUX7A.39
& ) SFFLUX7A.40
SFFLUX7A.41
IMPLICIT NONE SFFLUX7A.42
SFFLUX7A.43
INTEGER SFFLUX7A.44
& P_FIELD ! IN Total number of P-grid points. SFFLUX7A.45
&,LAND_FIELD ! IN Total number of land points. SFFLUX7A.46
&,TILE_PTS ! IN Number of tile points. SFFLUX7A.47
&,LAND_INDEX(P_FIELD) ! IN Index of land points. SFFLUX7A.48
&,TILE_INDEX(LAND_FIELD)! IN Index of tile points. SFFLUX7A.49
SFFLUX7A.50
LOGICAL SFFLUX7A.51
& LTIMER ! IN Logical for TIMER SFFLUX7A.52
SFFLUX7A.53
REAL SFFLUX7A.54
& ASHTF(P_FIELD) ! IN Coefficient to calculate surface SFFLUX7A.55
! ! heat flux into soil (W/m2/K). SFFLUX7A.56
&,LH ! IN Latent heat (J/K/kg). SFFLUX7A.57
&,QS1(P_FIELD) ! IN Sat. specific humidity qsat(TL_1,PSTAR) SFFLUX7A.58
&,QSTAR(LAND_FIELD) ! IN Surface qsat. SFFLUX7A.59
&,QW_1(P_FIELD) ! IN Total water content of lowest SFFLUX7A.60
! ! atmospheric layer (kg per kg air). SFFLUX7A.61
&,RADNET(P_FIELD) ! IN Net surface radiation (W/m2) positive SFFLUX7A.62
! ! downwards SFFLUX7A.63
&,RESFT(LAND_FIELD) ! IN Total resistance factor. SFFLUX7A.64
&,RHOKH_1(LAND_FIELD) ! IN Surface exchange coefficient. SFFLUX7A.65
&,TILE_FRAC(LAND_FIELD) SFFLUX7A.66
! ! IN Tile fraction. SFFLUX7A.67
&,TL_1(P_FIELD) ! IN Liquid/frozen water temperature for SFFLUX7A.68
! ! lowest atmospheric layer (K). SFFLUX7A.69
&,TS1(LAND_FIELD) ! IN Temperature of surface layer (K). SFFLUX7A.70
&,TSTAR(LAND_FIELD) ! IN Surface temperature (K). SFFLUX7A.71
&,Z0H(LAND_FIELD) ! IN Roughness length for heat and moisture SFFLUX7A.72
&,Z0M_EFF(LAND_FIELD) ! IN Effective roughness length for momentum SFFLUX7A.73
&,Z1_TQ(P_FIELD) ! IN Height of lowest atmospheric level (m). SFFLUX7A.74
SFFLUX7A.75
REAL SFFLUX7A.76
& FQW_1_GB(P_FIELD) ! INOUT GBM surface flux of QW (kg/m2/s). SFFLUX7A.77
&,FTL_1_GB(P_FIELD) ! INOUT GBM surface flux of TL. SFFLUX7A.78
SFFLUX7A.79
REAL SFFLUX7A.80
& ALPHA1(LAND_FIELD) ! OUT Gradient of saturated specific humidity SFFLUX7A.81
! ! with respect to temperature between the SFFLUX7A.82
! ! bottom model layer and the surface. SFFLUX7A.83
&,FQW_1(LAND_FIELD) ! OUT Local surface flux of QW (kg/m2/s). SFFLUX7A.84
&,FTL_1(LAND_FIELD) ! OUT Local surface flux of TL. SFFLUX7A.85
&,RHOKPM(LAND_FIELD) ! OUT Modified surface exchange coefficient. SFFLUX7A.86
SFFLUX7A.87
*CALL C_EPSLON 
SFFLUX7A.88
*CALL C_0_DG_C SFFLUX7A.89
*CALL C_G SFFLUX7A.90
*CALL C_LHEAT SFFLUX7A.91
*CALL C_R_CP SFFLUX7A.92
SFFLUX7A.93
! Derived local parameters SFFLUX7A.94
! Derived local parameters SFFLUX7A.95
REAL GRCP,LS SFFLUX7A.96
PARAMETER ( SFFLUX7A.97
& GRCP=G/CP SFFLUX7A.98
&,LS=LF+LC ! Latent heat of sublimation. SFFLUX7A.99
& ) SFFLUX7A.100
SFFLUX7A.101
! Scalars SFFLUX7A.102
INTEGER SFFLUX7A.103
& I ! Horizontal field index. SFFLUX7A.104
&,J ! Tile field index. SFFLUX7A.105
&,L ! Land point field index. SFFLUX7A.106
SFFLUX7A.107
REAL SFFLUX7A.108
& DQ1 ! (qsat(TL_1,PSTAR)-QW_1) + g/cp*alpha1*Z1 SFFLUX7A.109
&,D_T ! Temporary in calculation of alpha1. SFFLUX7A.110
&,RAD_REDUC ! Radiation term required for surface flux SFFLUX7A.111
SFFLUX7A.112
EXTERNAL TIMER SFFLUX7A.113
SFFLUX7A.114
IF (LTIMER) THEN SFFLUX7A.115
CALL TIMER('SF_FLUX ',3) SFFLUX7A.116
ENDIF SFFLUX7A.117
SFFLUX7A.118
!----------------------------------------------------------------------- SFFLUX7A.119
!! 1 Calculate gradient of saturated specific humidity for use in SFFLUX7A.120
!! calculation of surface fluxes SFFLUX7A.121
!----------------------------------------------------------------------- SFFLUX7A.122
DO J=1,TILE_PTS SFFLUX7A.123
L = TILE_INDEX(J) SFFLUX7A.124
I = LAND_INDEX(L) SFFLUX7A.125
D_T = TSTAR(L) - TL_1(I) SFFLUX7A.126
IF (D_T .GT. 0.05 .OR. D_T .LT. -0.05) THEN SFFLUX7A.127
ALPHA1(L) = (QSTAR(L) - QS1(I)) / D_T SFFLUX7A.128
ELSEIF (TL_1(I) .GT. TM) THEN SFFLUX7A.129
ALPHA1(L) = EPSILON*LC*QS1(I)*(1. + C_VIRTUAL*QS1(I)) SFFLUX7A.130
& / ( R*TL_1(I)*TL_1(I)) SFFLUX7A.131
ELSE SFFLUX7A.132
ALPHA1(L) = EPSILON*LS*QS1(I)*(1. + C_VIRTUAL*QS1(I)) SFFLUX7A.133
& / ( R*TL_1(I)*TL_1(I)) SFFLUX7A.134
ENDIF SFFLUX7A.135
ENDDO SFFLUX7A.136
SFFLUX7A.137
DO J=1,TILE_PTS SFFLUX7A.138
L = TILE_INDEX(J) SFFLUX7A.139
I = LAND_INDEX(L) SFFLUX7A.140
SFFLUX7A.141
RHOKPM(L) = RHOKH_1(L) / ( ASHTF(I) + SFFLUX7A.142
& RHOKH_1(L)*(LH*ALPHA1(L)*RESFT(L) + CP) ) SFFLUX7A.143
RAD_REDUC = RADNET(I) - ASHTF(I) * ( TL_1(I) - TS1(L) SFFLUX7A.144
& + GRCP*(Z1_TQ(I) + Z0M_EFF(L) - Z0H(L)) ) SFFLUX7A.145
DQ1 = QS1(I) - QW_1(I) + SFFLUX7A.146
& GRCP*ALPHA1(L)*(Z1_TQ(I) + Z0M_EFF(L) - Z0H(L)) SFFLUX7A.147
FQW_1(L) = RESFT(L)*RHOKPM(L)*( ALPHA1(L)*RAD_REDUC SFFLUX7A.148
& + (CP*RHOKH_1(L) + ASHTF(I))*DQ1 ) SFFLUX7A.149
FTL_1(L) = RHOKPM(L)*(RAD_REDUC - LH*RESFT(L)*RHOKH_1(L)*DQ1) SFFLUX7A.150
SFFLUX7A.151
FTL_1_GB(I) = FTL_1_GB(I) + TILE_FRAC(L)*FTL_1(L) SFFLUX7A.152
FQW_1_GB(I) = FQW_1_GB(I) + TILE_FRAC(L)*FQW_1(L) SFFLUX7A.153
SFFLUX7A.154
ENDDO SFFLUX7A.155
SFFLUX7A.156
IF (LTIMER) THEN SFFLUX7A.157
CALL TIMER('SF_FLUX ',4) SFFLUX7A.158
ENDIF SFFLUX7A.159
SFFLUX7A.160
RETURN SFFLUX7A.161
END SFFLUX7A.162
SFFLUX7A.163
! SUBROUTINE SF_FLUX_SEA-------------------------------------------- SFFLUX7A.164
! SFFLUX7A.165
! Calculate explicit surface fluxes of heat and moisture over sea SFFLUX7A.166
! and sea-ice SFFLUX7A.167
! SFFLUX7A.168
! ------------------------------------------------------------------ SFFLUX7A.169
SUBROUTINE SF_FLUX_SEA ( 1,2SFFLUX7A.170
& P_POINTS,P_FIELD,P1,NSICE,SICE_INDEX,LAND_MASK, SFFLUX7A.171
& ASHTF,ICE_FRACT,QS1,QSTAR_ICE,QSTAR_SEA,QW_1,RADNET,RHOKH_1,TI, SFFLUX7A.172
& TL_1,TSTAR_ICE,TSTAR_SEA,Z0H_ICE,Z0M_ICE,Z0H_SEA,Z0M_SEA,Z1_TQ, SFFLUX7A.173
& ALPHA1,E_SEA,FQW_ICE,FQW_1,FTL_ICE,FTL_1,H_SEA,RHOKPM,LTIMER SFFLUX7A.174
& ) SFFLUX7A.175
SFFLUX7A.176
IMPLICIT NONE SFFLUX7A.177
SFFLUX7A.178
INTEGER SFFLUX7A.179
& P_POINTS ! IN Number of P-grid points to be processed. SFFLUX7A.180
&,P_FIELD ! IN Total number of P-grid points. SFFLUX7A.181
&,P1 ! IN First P-point to be processed. SFFLUX7A.182
&,NSICE ! IN Number of sea-ice points. SFFLUX7A.183
&,SICE_INDEX(P_FIELD) ! IN Index of sea-ice points SFFLUX7A.184
SFFLUX7A.185
LOGICAL SFFLUX7A.186
& LTIMER ! IN Logical for TIMER SFFLUX7A.187
&,LAND_MASK(P_FIELD) ! IN .TRUE. for land, .FALSE. elsewhere. SFFLUX7A.188
SFFLUX7A.189
REAL SFFLUX7A.190
& ASHTF(P_FIELD) ! IN Coefficient to calculate surface SFFLUX7A.191
! ! heat flux into sea-ice (W/m2/K). SFFLUX7A.192
&,ICE_FRACT(P_FIELD) ! IN Fraction of gridbox which is sea-ice. SFFLUX7A.193
&,QS1(P_FIELD) ! IN Sat. specific humidity qsat(TL_1,PSTAR) SFFLUX7A.194
&,QSTAR_ICE(P_FIELD) ! IN Surface qsat for sea-ice. SFFLUX7A.195
&,QSTAR_SEA(P_FIELD) ! IN Surface qsat for sea or sea-ice leads. SFFLUX7A.196
&,QW_1(P_FIELD) ! IN Total water content of lowest SFFLUX7A.197
! ! atmospheric layer (kg per kg air). SFFLUX7A.198
&,RADNET(P_FIELD) ! IN Net surface radiation (W/m2) positive SFFLUX7A.199
! ! downwards SFFLUX7A.200
&,RHOKH_1(P_FIELD) ! IN Surface exchange coefficient. SFFLUX7A.201
&,TI(P_FIELD) ! IN Temperature of sea-ice surface layer (K) SFFLUX7A.202
&,TL_1(P_FIELD) ! IN Liquid/frozen water temperature for SFFLUX7A.203
! ! lowest atmospheric layer (K). SFFLUX7A.204
&,TSTAR_ICE(P_FIELD) ! IN Sea-ice surface temperature (K). SFFLUX7A.205
&,TSTAR_SEA(P_FIELD) ! IN Sea surface temperature (K). SFFLUX7A.206
&,Z0H_ICE(P_FIELD) ! IN Sea-ice heat and moisture roughness SFFLUX7A.207
! ! length (m). SFFLUX7A.208
&,Z0M_ICE(P_FIELD) ! IN Sea-ice momentum roughness length (m). SFFLUX7A.209
&,Z0H_SEA(P_FIELD) ! IN Sea and lead heat and moisture roughness SFFLUX7A.210
! ! length (m). SFFLUX7A.211
&,Z0M_SEA(P_FIELD) ! IN Sea and lead momentum roughness length. SFFLUX7A.212
&,Z1_TQ(P_FIELD) ! IN Height of lowest atmospheric level (m). SFFLUX7A.213
SFFLUX7A.214
REAL SFFLUX7A.215
& ALPHA1(P_FIELD) ! OUT Gradient of saturated specific humidity SFFLUX7A.216
! ! with respect to temperature between the SFFLUX7A.217
! ! bottom model layer and the surface. SFFLUX7A.218
&,E_SEA(P_FIELD) ! OUT Evaporation from sea times leads SFFLUX7A.219
! ! fraction (kg/m2/s). SFFLUX7A.220
&,FQW_ICE(P_FIELD) ! OUT Surface flux of QW for sea-ice. SFFLUX7A.221
&,FQW_1(P_FIELD) ! OUT GBM surface flux of QW (kg/m2/s). SFFLUX7A.222
&,FTL_ICE(P_FIELD) ! OUT Surface flux of TL for sea-ice. SFFLUX7A.223
&,FTL_1(P_FIELD) ! OUT GBM surface flux of TL. SFFLUX7A.224
&,H_SEA(P_FIELD) ! OUT Surface sensible heat flux over sea SFFLUX7A.225
! ! times leads fraction (W/m2). SFFLUX7A.226
&,RHOKPM(P_FIELD) ! OUT Modified surface exchange coefficient. SFFLUX7A.227
SFFLUX7A.228
*CALL C_EPSLON SFFLUX7A.229
*CALL C_0_DG_C SFFLUX7A.230
*CALL C_G SFFLUX7A.231
*CALL C_LHEAT SFFLUX7A.232
*CALL C_R_CP SFFLUX7A.233
SFFLUX7A.234
! Derived local parameters. SFFLUX7A.235
REAL GRCP,LS SFFLUX7A.236
PARAMETER ( SFFLUX7A.237
& GRCP=G/CP SFFLUX7A.238
&,LS=LF+LC ! Latent heat of sublimation. SFFLUX7A.239
& ) SFFLUX7A.240
SFFLUX7A.241
! Scalars SFFLUX7A.242
INTEGER SFFLUX7A.243
& I ! Horizontal field index. SFFLUX7A.244
&,J ! Sea-ice field index. SFFLUX7A.245
REAL SFFLUX7A.246
& DQ1 ! (qsat(TL_1,PSTAR)-QW_1) + g/cp*alpha1*Z1 SFFLUX7A.247
&,D_T ! Temporary in calculation of alpha1. SFFLUX7A.248
&,RAD_REDUC ! Radiation term required for surface flux SFFLUX7A.249
! ! calcs. SFFLUX7A.250
SFFLUX7A.251
EXTERNAL TIMER SFFLUX7A.252
SFFLUX7A.253
IF (LTIMER) THEN SFFLUX7A.254
CALL TIMER('SF_FLUX ',3) SFFLUX7A.255
ENDIF SFFLUX7A.256
SFFLUX7A.257
DO I=P1,P1+P_POINTS-1 SFFLUX7A.258
ALPHA1(I) = 0. SFFLUX7A.259
E_SEA(I) = 0. SFFLUX7A.260
H_SEA(I) = 0. SFFLUX7A.261
FQW_ICE(I) = 0. SFFLUX7A.262
FTL_ICE(I) = 0. SFFLUX7A.263
RHOKPM(I) = 0. SFFLUX7A.264
ENDDO SFFLUX7A.265
SFFLUX7A.266
!---------------------------------------------------------------------- SFFLUX7A.267
!! 1 Calculate gradient of saturated specific humidity for use in SFFLUX7A.268
!! calculation of surface fluxes - only required for sea-ice points SFFLUX7A.269
!---------------------------------------------------------------------- SFFLUX7A.270
DO J=1,NSICE SFFLUX7A.271
I = SICE_INDEX(J) SFFLUX7A.272
D_T = TSTAR_ICE(I) - TL_1(I) SFFLUX7A.273
IF (D_T .GT. 0.05 .OR. D_T .LT. -0.05) THEN SFFLUX7A.274
ALPHA1(I) = (QSTAR_ICE(I) - QS1(I)) / D_T SFFLUX7A.275
ELSEIF (TL_1(I) .GT. TM) THEN SFFLUX7A.276
ALPHA1(I) = EPSILON*LC*QS1(I)*( 1.0+C_VIRTUAL*QS1(I) ) SFFLUX7A.277
& /(R*TL_1(I)*TL_1(I)) SFFLUX7A.278
ELSE SFFLUX7A.279
ALPHA1(I) = EPSILON*LS*QS1(I)*( 1.0+C_VIRTUAL*QS1(I) ) SFFLUX7A.280
& /(R*TL_1(I)*TL_1(I)) SFFLUX7A.281
ENDIF SFFLUX7A.282
ENDDO SFFLUX7A.283
SFFLUX7A.284
DO I=P1,P1+P_POINTS-1 SFFLUX7A.285
IF ( .NOT. LAND_MASK(I) ) THEN SFFLUX7A.286
SFFLUX7A.287
E_SEA(I) = - (1. - ICE_FRACT(I)) * SFFLUX7A.288
& RHOKH_1(I)*(QW_1(I) - QSTAR_SEA(I)) SFFLUX7A.289
H_SEA(I) = - (1. - ICE_FRACT(I))*CP*RHOKH_1(I) * SFFLUX7A.290
& ( TL_1(I) - TSTAR_SEA(I) SFFLUX7A.291
& + GRCP*(Z1_TQ(I) + Z0M_SEA(I) - Z0H_SEA(I)) ) SFFLUX7A.292
SFFLUX7A.293
IF ( ICE_FRACT(I) .GT. 0. ) THEN SFFLUX7A.294
! Sea-ice SFFLUX7A.295
RHOKPM(I) = RHOKH_1(I) / ( ASHTF(I) + SFFLUX7A.296
& RHOKH_1(I)*(LS*ALPHA1(I) + CP) ) SFFLUX7A.297
RAD_REDUC = RADNET(I) - ICE_FRACT(I) * ASHTF(I) * SFFLUX7A.298
& ( TL_1(I) - TI(I) + SFFLUX7A.299
& GRCP*(Z1_TQ(I) + Z0M_ICE(I) - Z0H_ICE(I)) ) SFFLUX7A.300
DQ1 = QS1(I) - QW_1(I) + SFFLUX7A.301
& GRCP*ALPHA1(I)*(Z1_TQ(I) + Z0M_ICE(I) - Z0H_ICE(I)) SFFLUX7A.302
FQW_ICE(I) = RHOKPM(I) * ( ALPHA1(I)*RAD_REDUC + SFFLUX7A.303
& (CP*RHOKH_1(I) + ASHTF(I))*DQ1*ICE_FRACT(I) ) SFFLUX7A.304
FTL_ICE(I) = RHOKPM(I) * ( RAD_REDUC - SFFLUX7A.305
& ICE_FRACT(I)*LS*RHOKH_1(I)*DQ1 ) SFFLUX7A.306
SFFLUX7A.307
ENDIF SFFLUX7A.308
SFFLUX7A.309
FTL_1(I) = FTL_ICE(I) + H_SEA(I) / CP SFFLUX7A.310
FQW_1(I) = FQW_ICE(I) + E_SEA(I) SFFLUX7A.311
SFFLUX7A.312
ENDIF SFFLUX7A.313
ENDDO SFFLUX7A.314
SFFLUX7A.315
IF (LTIMER) THEN SFFLUX7A.316
CALL TIMER('SF_FLUX ',4) SFFLUX7A.317
ENDIF SFFLUX7A.318
SFFLUX7A.319
RETURN SFFLUX7A.320
END SFFLUX7A.321
SFFLUX7A.322
*ENDIF SFFLUX7A.323