*IF DEF,A03_3A ASJ4F401.8
C ******************************COPYRIGHT****************************** GTS2F400.8803
C (c) CROWN COPYRIGHT 1995, METEOROLOGICAL OFFICE, All Rights Reserved. GTS2F400.8804
C GTS2F400.8805
C Use, duplication or disclosure of this code is subject to the GTS2F400.8806
C restrictions as set forth in the contract. GTS2F400.8807
C GTS2F400.8808
C Meteorological Office GTS2F400.8809
C London Road GTS2F400.8810
C BRACKNELL GTS2F400.8811
C Berkshire UK GTS2F400.8812
C RG12 2SZ GTS2F400.8813
C GTS2F400.8814
C If no contract has been raised with this copy of the code, the use, GTS2F400.8815
C duplication or disclosure of it is strictly prohibited. Permission GTS2F400.8816
C to do so must first be obtained in writing from the Head of Numerical GTS2F400.8817
C Modelling at the above address. GTS2F400.8818
C ******************************COPYRIGHT****************************** GTS2F400.8819
C GTS2F400.8820
C*LL SUBROUTINE SF_EVAP------------------------------------------------ SFEVAP2B.3
CLL SFEVAP2B.4
CLL Purpose: Calculate surface evaporation and sublimation amounts SFEVAP2B.5
CLL (without applying them to the surface stores). SFEVAP2B.6
CLL Also calculate heat flux due to sea-ice melting. SFEVAP2B.7
CLL Also calculate 1.5 metre T and Q. SFEVAP2B.8
CLL SFEVAP2B.9
CLL Change to the calculation of soil evaporation when SFEVAP2B.10
CLL lying snow disappears during a timestep. Also two new SFEVAP2B.11
CLL common decks C_HT_M and C_SICEHC SFEVAP2B.12
CLL SFEVAP2B.13
CLL Suitable for single column usage. SFEVAP2B.14
CLL SFEVAP2B.15
CLL Model Modification history: SFEVAP2B.16
CLL version Date SFEVAP2B.17
CLL 3.1 11/1/93 Co-ordinate transformation so that screen SFEVAP2B.18
CLL temperature is calculated at 1.5m above SFEVAP2B.19
CLL surface + 1.5m S.Jackson SFEVAP2B.20
CLL 3.4 06/06/94 DEF TIMER replaced by LOGICAL LTIMER ASJ1F304.382
CLL Argument LTIMER added ASJ1F304.383
CLL S.J.Swarbrick ASJ1F304.384
CLL 4.1 08/05/96 decks A03_2C and A03_3B removed ASJ4F401.9
CLL S D Jackson ASJ4F401.10
CLL 4.2 Oct. 96 T3E migration - *DEF CRAY removed GSS1F402.68
CLL S J Swarbrick GSS1F402.69
CLL 4.5 Jul. 98 Kill the IBM specific lines (JCThil) AJC1F405.127
CLL SFEVAP2B.21
CLL Programming standard: Unified Model Documentation Paper No 4, SFEVAP2B.22
CLL version 2, dated 18/1/90. SFEVAP2B.23
CLL SFEVAP2B.24
CLL Logical component covered: P245. SFEVAP2B.25
CLL SFEVAP2B.26
CLL System task: SFEVAP2B.27
CLL SFEVAP2B.28
CLL Documentation: UMDP 24 SFEVAP2B.29
CLL SFEVAP2B.30
CLL--------------------------------------------------------------------- SFEVAP2B.31
C* SFEVAP2B.32
C*L Arguments :--------------------------------------------------------- SFEVAP2B.33
SUBROUTINE SF_EVAP ( 4,14SFEVAP2B.34
+ P_FIELD,P1,LAND_FIELD,LAND1, SFEVAP2B.35
+ POINTS,BL_LEVELS,LAND_MASK,LAND_PTS,LAND_INDEX SFEVAP2B.39
+,ASOIL_1,CANOPY,CATCH,DTRDZ,DTRDZ_RML,EA,ESL,SOIL_HT_FLUX SFEVAP2B.41
+,E_SEA,ICE_FRACT,LYING_SNOW,RADNET,SMC,TIMESTEP SFEVAP2B.42
+,CER1P5M,CHR1P5M,PSTAR,Z1,Z0M,Z0H,SQ1P5,ST1P5,SMLT SFEVAP2B.43
+,NRML,RHOKH_1,DQW_1,DQW_RML,DTSTAR,FTL,FQW,ES,QW,TL,TSTAR SFEVAP2B.44
+,ECAN,EI SFEVAP2B.45
+,SICE_MLT_HTF,Q1P5M,T1P5M,LTIMER ASJ1F304.385
+) SFEVAP2B.47
IMPLICIT NONE SFEVAP2B.48
LOGICAL LTIMER ASJ1F304.386
INTEGER SFEVAP2B.49
+ P_FIELD ! IN No. of gridpoints in the whole grid. SFEVAP2B.50
+,P1 ! IN 1st P-pt in full field to be processed. SFEVAP2B.51
+,LAND_FIELD ! IN No. of LANDpoints in the whole grid. SFEVAP2B.52
+,LAND1 ! IN 1st L-pt in full field to be processed. SFEVAP2B.53
+,POINTS ! IN No. of gridpoints to be processed. SFEVAP2B.54
+,BL_LEVELS ! IN No. of levels treated by b.l. scheme. SFEVAP2B.55
+,LAND_PTS ! IN No. of land points to be processed. SFEVAP2B.56
LOGICAL SFEVAP2B.57
+ LAND_MASK(P_FIELD) ! IN T for land points, F otherwise. SFEVAP2B.58
INTEGER SFEVAP2B.60
+ LAND_INDEX(P_FIELD) ! IN Index of land points on the P-grid. SFEVAP2B.61
C ! The ith element contains the position SFEVAP2B.62
C ! in whole grid of the ith land point. SFEVAP2B.63
REAL SFEVAP2B.65
+ ASOIL_1(P_FIELD) ! IN Soil coefficient from P242 (sq m K per SFEVAP2B.66
C ! per Joule * timestep). SFEVAP2B.67
+,CANOPY(LAND_FIELD) ! IN Gridbox mean canopy / surface water SFEVAP2B.68
C ! store (kg per sq m). SFEVAP2B.69
+,CATCH(LAND_FIELD) ! IN Canopy / surface water store capacity SFEVAP2B.70
C ! (kg per sq m). SFEVAP2B.71
+,DTRDZ(P_FIELD, ! IN -g.dt/dp for each model layer on p-grid SFEVAP2B.72
+ BL_LEVELS) ! From P244 ((kg/sq m/s)**-1). SFEVAP2B.73
+,DTRDZ_RML(P_FIELD) ! IN -g.dt/dp for the rapidly mixing layer SFEVAP2B.74
C ! (if it exists) on the p-grid from P244. SFEVAP2B.75
+,EA(P_FIELD) ! IN Surface evaporation with only aero- SFEVAP2B.76
C ! dynamic resistance (+ve), or condens- SFEVAP2B.77
C ! ation (-ve), averaged over gridbox. SFEVAP2B.78
C ! From P243 and P244. Kg per sq m per s. SFEVAP2B.79
+,ESL(P_FIELD) ! IN ES (q.v.) without fractional weighting SFEVAP2B.80
C ! factor FRACS ('L' is for 'local'). SFEVAP2B.81
C ! From P243 and P244. Kg per sq m per s. SFEVAP2B.82
+,SOIL_HT_FLUX(P_FIELD) ! IN Heat flux from surface to deep soil SFEVAP2B.83
C ! layer 1 (Watts per sq m). From P242. SFEVAP2B.84
+,E_SEA(P_FIELD) ! IN Evaporation from sea (weighted with SFEVAP2B.85
C ! leads fraction at sea-ice points). SFEVAP2B.86
C ! Diagnostics defined on land and sea. SFEVAP2B.88
+,ICE_FRACT(P_FIELD) ! IN Fraction of gridbox which is covered by SFEVAP2B.90
C ! sea-ice (decimal fraction, but most of SFEVAP2B.91
C ! this sub-component assumes it to be SFEVAP2B.92
C ! either 1.0 or 0.0 precisely). SFEVAP2B.93
+,LYING_SNOW(P_FIELD) ! IN Lying snow (kg per sq m). SFEVAP2B.94
C ! NB Dimension is PFIELD not LAND_FIELD f SFEVAP2B.96
C ! snow on sea-ice in coupled model runs. SFEVAP2B.97
+,RADNET(P_FIELD) ! IN Surface net radiation (Watts per sq m). SFEVAP2B.99
+,SMC(LAND_FIELD) ! IN Soil moisture content (kg per sq m). SFEVAP2B.100
+,TIMESTEP ! IN Timestep (sec). SFEVAP2B.101
LOGICAL SFEVAP2B.102
+ SQ1P5 ! IN STASH flag for 1.5-metre sp humidity. SFEVAP2B.103
+,ST1P5 ! IN STASH flag for 1.5-metre temperature. SFEVAP2B.104
+,SMLT ! IN STASH flag for sea-ice melting ht flux. SFEVAP2B.105
REAL SFEVAP2B.106
+ CER1P5M(P_FIELD) ! IN Transfer coefficient ratio, from P243. SFEVAP2B.107
+,CHR1P5M(P_FIELD) ! IN Transfer coefficient ratio, from P243. SFEVAP2B.108
+,PSTAR(P_FIELD) ! IN Surface pressure (Pa). SFEVAP2B.109
+,Z1(P_FIELD) ! IN Height of lowest atmospheric level SFEVAP2B.110
C ! (i.e. middle of lowest layer). Metres. SFEVAP2B.111
+,Z0M(P_FIELD) ! IN Roughness length for momentum (m) SFEVAP2B.112
+,Z0H(P_FIELD) ! IN Roughness length for heat and moisture SFEVAP2B.113
INTEGER SFEVAP2B.114
& NRML(P_FIELD) ! IN The Number of model layers in the SFEVAP2B.115
C ! Rapidly Mixing Layer. SFEVAP2B.116
REAL SFEVAP2B.117
+ RHOKH_1(P_FIELD) ! IN Turbulent surface exchange SFEVAP2B.118
C ! coefficient for sensible heat. SFEVAP2B.119
+,DQW_1(P_FIELD) ! IN Increment for lowest-level total SFEVAP2B.120
C ! water content. From P244. SFEVAP2B.121
+,DQW_RML(P_FIELD) ! IN Increment to QW due to rapid mixing. SFEVAP2B.122
+,DTSTAR(P_FIELD) ! IN Increment for surface temperature. SFEVAP2B.123
C ! From P244. SFEVAP2B.124
+,FTL(P_FIELD, ! INOUT Sensible heat flux from layer k-1 to SFEVAP2B.125
+ BL_LEVELS) ! layer k (W/sq m). From P243 and SFEVAP2B.126
C ! P244, units changed in P24 top level. SFEVAP2B.127
+,FQW(P_FIELD, ! INOUT Turbulent moisture flux from level SFEVAP2B.128
+ BL_LEVELS) ! k-1 to k (kg/sq m/s). From P243/4. SFEVAP2B.129
+,ES(P_FIELD) ! INOUT Surface evapotranspiration (through SFEVAP2B.130
C ! a resistance which is not entirely SFEVAP2B.131
C ! aerodynamic). Always non-negative. SFEVAP2B.132
C ! Kg per sq m per sec. From P243/4. SFEVAP2B.133
C ! Diagnostics defined on land and sea. SFEVAP2B.135
+,QW(P_FIELD,BL_LEVELS)! INOUT Total water content (kg(water)/ SFEVAP2B.137
C ! kg(air)). From P243/4. SFEVAP2B.138
+,TL(P_FIELD,BL_LEVELS)! INOUT Liquid/frozen water temperature (K). SFEVAP2B.139
+,TSTAR(P_FIELD) ! INOUT Surface temperature (K). SFEVAP2B.140
REAL SFEVAP2B.141
+ ECAN(P_FIELD) ! OUT Gridbox mean evaporation from canopy/ SFEVAP2B.142
C ! surface store (kg per sq m per s). SFEVAP2B.143
C ! Zero over sea and sea-ice. SFEVAP2B.144
C ! Diagnostics defined on land and sea. SFEVAP2B.146
+,EI(P_FIELD) ! OUT Sublimation from lying snow or sea- SFEVAP2B.148
C ! ice (kg per sq m per s). SFEVAP2B.149
REAL SFEVAP2B.150
+ SICE_MLT_HTF(P_FIELD)! OUT Heat flux due to melting of sea-ice SFEVAP2B.151
C ! (Watts per square metre). SFEVAP2B.152
+,Q1P5M(P_FIELD) ! OUT Specific humidity at screen height of SFEVAP2B.153
C ! 1.5 metres (kg water per kg air). SFEVAP2B.154
+,T1P5M(P_FIELD) ! OUT Temperature at 1.5 metres above the SFEVAP2B.155
C ! surface (K). SFEVAP2B.156
C* SFEVAP2B.157
C*L External subprogram(s) required :- SFEVAP2B.158
EXTERNAL QSAT SFEVAP2B.159
EXTERNAL TIMER SFEVAP2B.161
C* SFEVAP2B.163
C*L Local and other symbolic constants used :- SFEVAP2B.164
*CALL C_0_DG_C
SFEVAP2B.165
*CALL C_LHEAT
SFEVAP2B.166
*CALL C_G
SFEVAP2B.167
*CALL C_HT_M
! Contains Z1P5M SFEVAP2B.168
*CALL C_R_CP
SFEVAP2B.169
*CALL C_SICEHC
! Contains AI SFEVAP2B.170
REAL GRCP SFEVAP2B.171
PARAMETER ( SFEVAP2B.172
+ GRCP=G/CP ! Accn due to gravity / standard heat capacity of SFEVAP2B.173
C ! air at const pressure. Used in diagnosis of 1.5 SFEVAP2B.174
C ! metre temperature. SFEVAP2B.175
+) SFEVAP2B.176
C* SFEVAP2B.177
REAL SFEVAP2B.179
+ QS(P_FIELD) ! Used for saturated specific humidity SFEVAP2B.180
C ! at surface, in Q1P5M calculation. SFEVAP2B.181
REAL SFEVAP2B.182
+ DTL_RML(P_FIELD) ! Adjustment increment to TL in the SFEVAP2B.183
C ! rapidly mixing layer. SFEVAP2B.184
+,DFTL(P_FIELD, ! Adjustment increment to the sensible SFEVAP2B.185
+ BL_LEVELS) ! heat flux. SFEVAP2B.186
+,DFQW(P_FIELD, ! Adjustment increment to the flux of SFEVAP2B.187
+ BL_LEVELS) ! total water. SFEVAP2B.188
C Local scalars SFEVAP2B.202
REAL SFEVAP2B.203
+ EADT ! EA (q.v.) integrated over timestep. SFEVAP2B.204
+,ECANDT ! ECAN (q.v.) integrated over timestep. SFEVAP2B.205
+,EDT ! E=FQW(,1) (q.v.) integrated over timestep. SFEVAP2B.206
+,EIDT ! EI (q.v.) integrated over timestep. SFEVAP2B.207
+,ESDT ! ES (q.v.) integrated over timestep. SFEVAP2B.208
+,ESLDT ! ESL (q.v.) integrated over timestep. SFEVAP2B.209
+,EW ! Total surface flux of water, excluding SFEVAP2B.210
C ! sublimation/frost deposition, over land. SFEVAP2B.211
+,FRACA ! Fraction of gridbox at which moisture flux SFEVAP2B.212
C ! is going at the potential rate, i.e. with SFEVAP2B.213
C ! only aerodynamic resistance. SFEVAP2B.214
+,FRACS ! Fraction of gridbox at which moisture flux SFEVAP2B.215
C ! is additionally impeded by a surface and/or SFEVAP2B.216
C ! stomatal resistance. SFEVAP2B.217
+,CT_1 ! Parameter in the implicit adjustment of SFEVAP2B.218
C ! TSTAR and TL(,1). SFEVAP2B.219
+,CT_RML ! Parameter in the implicit adjustment of SFEVAP2B.220
C ! TSTAR and TL in the rapidly mixing layer. SFEVAP2B.221
+,TSTARMAX ! Maximum gridbox mean surface temperature in SFEVAP2B.222
C ! sea-ice melting calculation. See comments SFEVAP2B.223
C ! to section 1.6.3(b) or 2.5.3(b) below. SFEVAP2B.224
INTEGER SFEVAP2B.225
+ I ! Loop counter - full horizontal field index. SFEVAP2B.226
+,L ! Loop counter - land field index. SFEVAP2B.227
+,K ! Loop counter in the vertical. SFEVAP2B.228
+,KM1 ! K - 1 SFEVAP2B.229
+,NRMLP1 ! NRML + 1 SFEVAP2B.230
C SFEVAP2B.231
C----------------------------------------------------------------------- SFEVAP2B.232
CL 1. Initialise some output variables to zero. SFEVAP2B.233
C----------------------------------------------------------------------- SFEVAP2B.234
C SFEVAP2B.235
IF (LTIMER) THEN ASJ1F304.387
CALL TIMER
('SFEVAP ',3) SFEVAP2B.237
ENDIF ASJ1F304.388
DO 1 I=P1,P1+POINTS-1 SFEVAP2B.239
ECAN(I) = 0.0 SFEVAP2B.240
EI(I) = 0.0 SFEVAP2B.241
IF (SMLT) SFEVAP2B.242
+ SICE_MLT_HTF(I) = 0.0 SFEVAP2B.243
1 CONTINUE SFEVAP2B.244
C SFEVAP2B.245
C--------------------------------------------------------------------- SFEVAP2B.246
CL 2. Do calculations for land points. SFEVAP2B.247
C--------------------------------------------------------------------- SFEVAP2B.248
C SFEVAP2B.249
CMIC$ DO ALL VECTOR SHARED(P_FIELD, LAND_FIELD, BL_LEVELS, LAND1, SFEVAP2B.250
CMIC$1 LAND_PTS, LAND_INDEX, ESL, TIMESTEP, ES, LYING_SNOW, ECAN, SFEVAP2B.251
CMIC$2 EA, CATCH, CANOPY, SMC, EI, TSTAR, FQW) PRIVATE(I, L, ESLDT, SFEVAP2B.252
CMIC$3 ESDT, EADT, EDT, ECANDT, FRACA, FRACS, EW, EIDT) SFEVAP2B.253
CDIR$ IVDEP SFEVAP2B.254
! Fujitsu vectorization directive GRB0F405.443
!OCL NOVREC GRB0F405.444
DO 2 L=LAND1,LAND1+LAND_PTS-1 SFEVAP2B.260
I = LAND_INDEX(L) SFEVAP2B.261
C SFEVAP2B.263
C----------------------------------------------------------------------- SFEVAP2B.264
CL 2.1 Calculate fluxes integrated over timestep. SFEVAP2B.265
C----------------------------------------------------------------------- SFEVAP2B.266
C SFEVAP2B.267
ESLDT = ESL(I) * TIMESTEP SFEVAP2B.268
EADT = EA(I) * TIMESTEP SFEVAP2B.269
ESDT = ES(I) * TIMESTEP SFEVAP2B.270
EDT = EADT + ESDT SFEVAP2B.271
C SFEVAP2B.272
C----------------------------------------------------------------------- SFEVAP2B.273
CL 2.2 Do calculations for snow-free land. Canopy processes operate. SFEVAP2B.274
CL LYING_SNOW is defined on sea and land points for snow on sea-ice SFEVAP2B.276
CL in coupled model runs. SFEVAP2B.277
C----------------------------------------------------------------------- SFEVAP2B.279
C SFEVAP2B.280
IF (LYING_SNOW(I).LE.0.0) THEN SFEVAP2B.281
IF (EDT.GE.0.0) THEN SFEVAP2B.282
C SFEVAP2B.283
C----------------------------------------------------------------------- SFEVAP2B.284
CL 2.2.1 Non-negative moisture flux over snow-free land. SFEVAP2B.285
C----------------------------------------------------------------------- SFEVAP2B.286
C SFEVAP2B.287
C- - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - SFEVAP2B.288
CL (a) Water in canopy and soil is assumed to be liquid, so all SFEVAP2B.289
CL positive moisture flux over snow-free land is evaporation SFEVAP2B.290
CL rather than sublimation, even if TSTAR is less than or equal SFEVAP2B.291
CL to TM. SFEVAP2B.292
C- - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - SFEVAP2B.293
C SFEVAP2B.294
ECAN(I) = EA(I) SFEVAP2B.295
ECANDT = EADT SFEVAP2B.296
C SFEVAP2B.297
C If EDT is non-negative, then ECANDT must be non-negative. SFEVAP2B.298
C SFEVAP2B.299
FRACA = 0.0 SFEVAP2B.300
IF (CATCH(L).GT.0.0) SFEVAP2B.301
+ FRACA = CANOPY(L) / CATCH(L) SFEVAP2B.302
IF (CANOPY(L).LT.ECANDT) THEN SFEVAP2B.303
C SFEVAP2B.304
C- - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - SFEVAP2B.305
CL (b) It is assumed that any 'canopy' moisture flux in excess of the SFEVAP2B.306
CL current canopy water amount is in fact soil evaporation. SFEVAP2B.307
C- - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - SFEVAP2B.308
C SFEVAP2B.309
C This situation is highly improbable - it will occur at, at SFEVAP2B.310
C most, a few gridpoints in any given timestep. SFEVAP2B.311
C SFEVAP2B.312
FRACS = 1.0 - FRACA*( CANOPY(L) / ECANDT ) SFEVAP2B.313
ESDT = ESLDT * FRACS SFEVAP2B.314
ECANDT = CANOPY(L) SFEVAP2B.315
ECAN(I) = ECANDT / TIMESTEP SFEVAP2B.316
ES(I) = ESDT / TIMESTEP SFEVAP2B.317
ENDIF SFEVAP2B.318
C SFEVAP2B.319
C (The canopy store is depleted by evaporation in P252, and not here, SFEVAP2B.320
C according to the formula: CANOPY=CANOPY-ECANDT) SFEVAP2B.321
C SFEVAP2B.322
C- - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - SFEVAP2B.323
CL (c) Adjustments to evaporation from soil as calculated so far :- SFEVAP2B.324
C- - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - SFEVAP2B.325
C SFEVAP2B.326
IF (SMC(L).LE.0.0) THEN SFEVAP2B.327
C SFEVAP2B.328
CL (i) If there is currently no soil moisture, there must be no SFEVAP2B.329
CL evaporation of soil moisture, so this flux is set to zero. SFEVAP2B.330
C SFEVAP2B.331
ESDT = 0.0 SFEVAP2B.332
ES(I) = 0.0 SFEVAP2B.333
ELSEIF (SMC(L).LT.ESDT) THEN SFEVAP2B.334
C SFEVAP2B.335
CL (ii) Ensure that the soil evaporation is not greater than the SFEVAP2B.336
CL current soil moisture store. SFEVAP2B.337
C This situation is extremely unlikely at any given gridpoint SFEVAP2B.338
C at any given timestep. SFEVAP2B.339
C SFEVAP2B.340
ESDT = SMC(L) SFEVAP2B.341
ES(I) = ESDT / TIMESTEP SFEVAP2B.342
ENDIF SFEVAP2B.343
C SFEVAP2B.344
C (The soil moisture store is depleted by evaporation in P253, and not SFEVAP2B.345
C here, using the formula: SMC=SMC-ESDT) SFEVAP2B.346
C SFEVAP2B.347
EW = ECAN(I) + ES(I) SFEVAP2B.348
EI(I) = 0.0 SFEVAP2B.349
C SFEVAP2B.350
C----------------------------------------------------------------------- SFEVAP2B.351
CL 2.2.2 Negative moisture flux onto snow-free land above freezing SFEVAP2B.352
C----------------------------------------------------------------------- SFEVAP2B.353
CL (i.e. condensation onto snow-free land). The whole flux is SFEVAP2B.354
CL into the surface/canopy store. SFEVAP2B.355
C SFEVAP2B.356
ELSEIF (TSTAR(I).GT.TM) THEN ! ELSE of evaporation / SFEVAP2B.357
C ! condensation block. SFEVAP2B.358
C SFEVAP2B.359
C Condensation implies ES=0, so ECAN=EA=EW=E (=FQW(,1)) SFEVAP2B.360
C SFEVAP2B.361
ECAN(I) = FQW(I,1) SFEVAP2B.362
ES(I) = 0.0 SFEVAP2B.363
EW = ECAN(I) SFEVAP2B.364
EI(I) = 0.0 SFEVAP2B.365
C SFEVAP2B.366
C (The canopy store is augmented by interception of condensation at SFEVAP2B.367
C P252, and not here.) SFEVAP2B.368
C SFEVAP2B.369
C----------------------------------------------------------------------- SFEVAP2B.370
CL 2.2.3 Negative moisture flux onto snow-free land below freezing SFEVAP2B.371
CL (i.e. deposition of frost). SFEVAP2B.372
C----------------------------------------------------------------------- SFEVAP2B.373
C SFEVAP2B.374
ELSE ! ELSE of condensation / frost deposition block. SFEVAP2B.375
EI(I) = FQW(I,1) SFEVAP2B.376
ES(I) = 0.0 SFEVAP2B.377
EW = 0.0 SFEVAP2B.378
C SFEVAP2B.379
C (Negative EI is used to increment the snowdepth store - there is SFEVAP2B.380
C no separate "frost" store. This incrementing is done in P251, SFEVAP2B.381
C according to: LYING_SNOW = LYING_SNOW - EI*TIMESTEP) SFEVAP2B.382
C SFEVAP2B.383
ENDIF ! End of evaporation/condensation/deposition block. SFEVAP2B.384
C SFEVAP2B.385
C----------------------------------------------------------------------- SFEVAP2B.386
CL 2.3 Do calculations for snow-covered land. SFEVAP2B.387
C----------------------------------------------------------------------- SFEVAP2B.388
C SFEVAP2B.389
ELSEIF (LYING_SNOW(I).LE.EDT) THEN ! ELSEIF of no-snow. SFEVAP2B.390
C SFEVAP2B.391
C----------------------------------------------------------------------- SFEVAP2B.392
CL 2.3.1 Shallow snow (lying snow or frost which is being exhausted SFEVAP2B.393
CL by evaporation). All the snow is sublimated, the remaining SFEVAP2B.394
CL moisture flux being taken from the canopy and soil, with all SFEVAP2B.395
CL the palaver of section 1.2.1 above. SFEVAP2B.396
C----------------------------------------------------------------------- SFEVAP2B.397
C SFEVAP2B.398
C This is extremely unlikely at more than one or two gridpoints SFEVAP2B.399
C at any given timestep, yet the complicated logic probably SFEVAP2B.400
C slows down the routine considerably - this section is a SFEVAP2B.401
C suitable candidate for further consideration as regards SFEVAP2B.402
C making the model optimally efficient. SFEVAP2B.403
C SFEVAP2B.404
EI(I) = LYING_SNOW(I) / TIMESTEP SFEVAP2B.405
EIDT = LYING_SNOW(I) SFEVAP2B.406
C SFEVAP2B.407
C Set EDT = ( E - SNOSUB ) * TIMESTEP. This is the moisture in kg per SFEVAP2B.408
C square metre left over to be evaporated from the canopy and soil. SFEVAP2B.409
C N.B. E=FQW(,1) SFEVAP2B.410
C SFEVAP2B.411
EDT = EDT - EIDT SFEVAP2B.412
C SFEVAP2B.413
C (Snowdepth is decreased using EI at P251, and not here. The formula SFEVAP2B.414
C used is simply: LYING_SNOW = LYING_SNOW - EI*TIMESTEP.) SFEVAP2B.415
C SFEVAP2B.416
C Now that all the snow has sublimed, canopy processes come into SFEVAP2B.417
C operation (FRACA no longer necessarily equal to 1). SFEVAP2B.418
C SFEVAP2B.419
FRACA = 0.0 SFEVAP2B.420
IF (CATCH(L).GT.0.0) SFEVAP2B.421
+ FRACA = CANOPY(L) / CATCH(L) SFEVAP2B.422
ECANDT = EDT * FRACA SFEVAP2B.423
IF (CANOPY(L).LT.ECANDT) THEN SFEVAP2B.424
C SFEVAP2B.425
C Dry out the canopy completely and assume the remaining moisture flux SFEVAP2B.426
C is soil evaporation. SFEVAP2B.427
C SFEVAP2B.428
FRACS = 1.0 - FRACA*( CANOPY(L) / ECANDT ) SFEVAP2B.429
ESDT = EDT * FRACS SFEVAP2B.430
ECANDT = CANOPY(L) SFEVAP2B.431
ELSE SFEVAP2B.432
C SFEVAP2B.433
C Calculate soil evaporation. SFEVAP2B.434
C SFEVAP2B.435
FRACS = 1.0 - FRACA SFEVAP2B.436
ESDT = EDT * FRACS SFEVAP2B.437
ENDIF SFEVAP2B.438
ECAN(I) = ECANDT / TIMESTEP SFEVAP2B.439
ES(I) = ESDT / TIMESTEP SFEVAP2B.440
C SFEVAP2B.441
C (ECAN is used to deplete the canopy store at P252, and not here. The SFEVAP2B.442
C formula used is simply: CANOPY = CANOPY - ECAN*TIMESTEP.) SFEVAP2B.443
C SFEVAP2B.444
C Evaporation from soil. SFEVAP2B.445
C SFEVAP2B.446
IF (SMC(L).LE.0.0) THEN SFEVAP2B.447
C SFEVAP2B.448
C No evaporation from soil possible when there is no soil moisture. SFEVAP2B.449
C SFEVAP2B.450
ESDT = 0.0 SFEVAP2B.451
ES(I) = 0.0 SFEVAP2B.452
ELSEIF (SMC(L).LT.ESDT) THEN SFEVAP2B.453
C SFEVAP2B.454
C Limit evaporation of soil moisture in the extremely unlikely event SFEVAP2B.455
C that soil moisture is exhausted by the evaporation left over from SFEVAP2B.456
C sublimation which exhausted the snow store. SFEVAP2B.457
C SFEVAP2B.458
ESDT = SMC(L) SFEVAP2B.459
ES(I) = ESDT / TIMESTEP SFEVAP2B.460
ENDIF SFEVAP2B.461
C SFEVAP2B.462
C (ES is used to deplete the soil moisture store at P253, and not here, SFEVAP2B.463
C according to the formula: SMC = SMC - ES*TIMESTEP.) SFEVAP2B.464
C SFEVAP2B.465
EW = ECAN(I) + ES(I) SFEVAP2B.466
C SFEVAP2B.467
C----------------------------------------------------------------------- SFEVAP2B.468
CL 2.3.2 Deep snow (i.e. not being exhausted by evaporation). This SFEVAP2B.469
CL covers two cases: (a) sublimation from deep snow (if total SFEVAP2B.470
CL moisture flux over the timestep is non-negative but less than SFEVAP2B.471
CL the lying snow amount), and (b) deposition onto an already SFEVAP2B.472
CL snowy surface (if the total moisture flux is negative and SFEVAP2B.473
CL the lying snow amount is positive). SFEVAP2B.474
C----------------------------------------------------------------------- SFEVAP2B.475
C SFEVAP2B.476
ELSE ! ELSE of shallow snow / deep snow block. SFEVAP2B.477
EI(I) = FQW(I,1) SFEVAP2B.478
EW = 0.0 SFEVAP2B.479
C SFEVAP2B.480
C (EI is used to increase or decrease the snowdepth at P251, and not SFEVAP2B.481
C here, according to the formula: SFEVAP2B.482
C LYING_SNOW = LYING_SNOW - EI*TIMESTEP . ) SFEVAP2B.483
C SFEVAP2B.484
ENDIF ! End of no snow/shallow snow/deep snow block. SFEVAP2B.485
FQW(I,1) = EW + EI(I) SFEVAP2B.486
2 CONTINUE SFEVAP2B.488
C SFEVAP2B.489
C Split loop 2 here so that it will vectorise. SFEVAP2B.490
C SFEVAP2B.491
CDIR$ IVDEP SFEVAP2B.492
! Fujitsu vectorization directive GRB0F405.445
!OCL NOVREC GRB0F405.446
DO 24 L=LAND1,LAND1+LAND_PTS-1 SFEVAP2B.493
I = LAND_INDEX(L) SFEVAP2B.494
C SFEVAP2B.496
C----------------------------------------------------------------------- SFEVAP2B.497
CL 2.4 Do calculations which are the same at all land points. Calculate SFEVAP2B.498
CL and apply increments to (a) total atmospheric water content QW, SFEVAP2B.499
CL (b) surface temperature TSTAR, (c) TL SFEVAP2B.500
C----------------------------------------------------------------------- SFEVAP2B.501
C SFEVAP2B.502
DTSTAR(I) = ASOIL_1(I) * ( RADNET(I) - FTL(I,1) SFEVAP2B.503
& - LC*FQW(I,1) - LF*EI(I) - SOIL_HT_FLUX(I) ) SFEVAP2B.504
& - DTSTAR(I) SFEVAP2B.505
IF ( NRML(I) .GE. 2 ) THEN SFEVAP2B.506
NRMLP1 = NRML(I) + 1 SFEVAP2B.507
DQW_RML(I) = -DTRDZ_RML(I) * ( FQW(I,NRMLP1) - FQW(I,1) ) SFEVAP2B.508
& -DQW_RML(I) SFEVAP2B.509
DFQW(I,1) = DQW_RML(I) / DTRDZ_RML(I) SFEVAP2B.510
CT_RML = -DTRDZ_RML(I) * RHOKH_1(I) SFEVAP2B.511
DTSTAR(I) = DTSTAR(I) / SFEVAP2B.512
& ( 1.0 - CT_RML + ASOIL_1(I)*CP*RHOKH_1(I) ) SFEVAP2B.513
DFTL(I,1) = CP * RHOKH_1(I) * DTSTAR(I) SFEVAP2B.514
FTL(I,1) = FTL(I,1) + DFTL(I,1) SFEVAP2B.515
DTL_RML(I) = -CT_RML * DTSTAR(I) SFEVAP2B.516
DTSTAR(I) = DTSTAR(I) + DTL_RML(I) SFEVAP2B.517
ELSE SFEVAP2B.518
DQW_1(I) = -DTRDZ(I,1) * ( FQW(I,2) - FQW(I,1) ) SFEVAP2B.519
& - DQW_1(I) SFEVAP2B.520
QW(I,1) = QW(I,1) + DQW_1(I) SFEVAP2B.521
CT_1 = -DTRDZ(I,1) * RHOKH_1(I) SFEVAP2B.522
DTSTAR(I) = DTSTAR(I) / SFEVAP2B.523
& ( 1.0 - CT_1 + ASOIL_1(I)*CP*RHOKH_1(I) ) SFEVAP2B.524
FTL(I,1) = FTL(I,1) + CP * RHOKH_1(I) * DTSTAR(I) SFEVAP2B.525
TL(I,1) = TL(I,1) - CT_1 * DTSTAR(I) SFEVAP2B.526
DTSTAR(I) = (1.0 - CT_1) * DTSTAR(I) SFEVAP2B.527
ENDIF SFEVAP2B.528
TSTAR(I) = TSTAR(I) + DTSTAR(I) SFEVAP2B.529
24 CONTINUE SFEVAP2B.533
C SFEVAP2B.535
C----------------------------------------------------------------------- SFEVAP2B.536
CL 2.5 Do calculations for sea points. SFEVAP2B.537
C----------------------------------------------------------------------- SFEVAP2B.538
C SFEVAP2B.539
CMIC$ DO ALL VECTOR SHARED(P_FIELD, BL_LEVELS, P1, POINTS, SFEVAP2B.541
CMIC$1 LAND_MASK, ES, ECAN, EI, ICE_FRACT, FQW, E_SEA, SFEVAP2B.542
CMIC$2 TSTAR, SMLT, SICE_MLT_HTF, TIMESTEP) PRIVATE(I, TSTARMAX) SFEVAP2B.543
CDIR$ IVDEP SFEVAP2B.544
! Fujitsu vectorization directive GRB0F405.447
!OCL NOVREC GRB0F405.448
DO 25 I=P1,P1+POINTS-1 SFEVAP2B.545
IF (.NOT.LAND_MASK(I)) THEN SFEVAP2B.546
C SFEVAP2B.548
C----------------------------------------------------------------------- SFEVAP2B.549
CL 2.5.1 Set soil and canopy evaporation amounts to zero, and set SFEVAP2B.550
CL sublimation to zero for liquid sea points. SFEVAP2B.551
C----------------------------------------------------------------------- SFEVAP2B.552
C SFEVAP2B.553
ES(I) = 0.0 SFEVAP2B.554
ECAN(I) = 0.0 SFEVAP2B.555
EI(I) = 0.0 SFEVAP2B.556
C----------------------------------------------------------------------- SFEVAP2B.557
CL 2.5.3 For sea-ice points :- SFEVAP2B.558
C----------------------------------------------------------------------- SFEVAP2B.559
C SFEVAP2B.560
IF (ICE_FRACT(I).GT.0.0) THEN SFEVAP2B.561
C SFEVAP2B.562
CL (a) Set sublimation to total moisture flux less evaporation SFEVAP2B.563
CL from sea. SFEVAP2B.564
C SFEVAP2B.565
EI(I) = FQW(I,1) - E_SEA(I) SFEVAP2B.566
C SFEVAP2B.567
CL (b) If the surface temperature is above the melting point of SFEVAP2B.568
CL sea-ice, it is reset to be equal to the melting point. SFEVAP2B.569
C (Strictly speaking, the gridbox mean surface temperature is SFEVAP2B.570
C adjusted so that the temperature of the icy fraction is not SFEVAP2B.571
C above the melting point; unlike the rest of this brick, this SFEVAP2B.572
C sub-section is valid for any ice fraction between 1 and 0 SFEVAP2B.573
C inclusive, and not just for fractions of precisely 1 or 0.) SFEVAP2B.574
CL This corresponds to a melting of some of the sea-ice. If SFEVAP2B.575
CL requested, the heat flux associated with this melting is SFEVAP2B.576
CL calculated. This code is indifferent to whether the melting SFEVAP2B.577
CL reduces the thickness or areal extent of the sea-ice - neither SFEVAP2B.578
CL is altered here. SFEVAP2B.579
C (This too works correctly for a general ice fraction, by SFEVAP2B.580
C weighting the heat flux with the ice fraction, to give the SFEVAP2B.581
C gridbox mean heat flux.) SFEVAP2B.582
C SFEVAP2B.583
TSTARMAX = ICE_FRACT(I)*TM + (1.0-ICE_FRACT(I))*TFS SFEVAP2B.584
IF (TSTAR(I).GT.TSTARMAX) THEN SFEVAP2B.585
IF (SMLT) SFEVAP2B.586
+ SICE_MLT_HTF(I) = (TSTAR(I)-TSTARMAX) / (AI*TIMESTEP) SFEVAP2B.587
TSTAR(I) = TSTARMAX SFEVAP2B.588
ENDIF ! End of TSTAR adjustment block. SFEVAP2B.589
ENDIF ! End of liquid sea/sea-ice block. SFEVAP2B.590
ENDIF ! End of sea point calculations. SFEVAP2B.595
25 CONTINUE SFEVAP2B.596
C SFEVAP2B.598
C----------------------------------------------------------------------- SFEVAP2B.599
CL 3. Update heat and moisture fluxes over land due to a rapidly SFEVAP2B.600
CL layer. SFEVAP2B.601
C----------------------------------------------------------------------- SFEVAP2B.602
C SFEVAP2B.603
DO 3 K = 1,BL_LEVELS-1 SFEVAP2B.604
KM1 = K - 1 SFEVAP2B.605
CDIR$ IVDEP SFEVAP2B.611
! Fujitsu vectorization directive GRB0F405.449
!OCL NOVREC GRB0F405.450
DO L = LAND1,LAND1+LAND_PTS-1 SFEVAP2B.612
I = LAND_INDEX(L) SFEVAP2B.613
IF ( (K .LE. NRML(I)) .AND. (NRML(I) .GE. 2) ) THEN SFEVAP2B.615
TL(I,K) = TL(I,K) + DTL_RML(I) SFEVAP2B.616
QW(I,K) = QW(I,K) + DQW_RML(I) SFEVAP2B.617
IF ( K .GE. 2 ) THEN SFEVAP2B.618
DFTL(I,K) = DFTL(I,KM1) - CP * DTL_RML(I) / DTRDZ(I,KM1) SFEVAP2B.619
DFQW(I,K) = DFQW(I,KM1) - DQW_RML(I) / DTRDZ(I,KM1) SFEVAP2B.620
FTL(I,K) = FTL(I,K) + DFTL(I,K) SFEVAP2B.621
FQW(I,K) = FQW(I,K) + DFQW(I,K) SFEVAP2B.622
ENDIF ! K .GE. 2 SFEVAP2B.623
ENDIF ! Rapidly mixing layer exists and is more than one SFEVAP2B.624
C ! model layer deep and current layer is in it. SFEVAP2B.625
ENDDO ! Loop over points SFEVAP2B.629
3 CONTINUE! Loop over levels SFEVAP2B.630
C SFEVAP2B.631
C----------------------------------------------------------------------- SFEVAP2B.632
CL 4. Diagnose temperature and/or specific humidity at screen height SFEVAP2B.633
CL (1.5 metres), as requested via the STASH flags. SFEVAP2B.634
C----------------------------------------------------------------------- SFEVAP2B.635
C SFEVAP2B.636
IF (SQ1P5 .OR. ST1P5) THEN SFEVAP2B.637
IF (SQ1P5) CALL QSAT
(QS(P1),TSTAR(P1),PSTAR(P1),POINTS) SFEVAP2B.638
DO 4 I=P1,P1+POINTS-1 SFEVAP2B.639
IF (ST1P5) T1P5M(I) = TSTAR(I) - GRCP*Z1P5M + CHR1P5M(I) * SFEVAP2B.640
+ ( TL(I,1) - TSTAR(I) + GRCP*(Z1(I)+Z0M(I)-Z0H(I)) ) SFEVAP2B.641
IF (SQ1P5) Q1P5M(I) = QW(I,1) + CER1P5M(I)*( QW(I,1) - QS(I) ) SFEVAP2B.642
4 CONTINUE SFEVAP2B.643
ENDIF SFEVAP2B.644
IF (LTIMER) THEN ASJ1F304.389
CALL TIMER
('SFEVAP ',4) SFEVAP2B.646
ENDIF ASJ1F304.390
RETURN SFEVAP2B.648
END SFEVAP2B.649
*ENDIF SFEVAP2B.650