*IF DEF,A03_3A ASJ4F401.5
C ******************************COPYRIGHT****************************** GTS2F400.4429
C (c) CROWN COPYRIGHT 1995, METEOROLOGICAL OFFICE, All Rights Reserved. GTS2F400.4430
C GTS2F400.4431
C Use, duplication or disclosure of this code is subject to the GTS2F400.4432
C restrictions as set forth in the contract. GTS2F400.4433
C GTS2F400.4434
C Meteorological Office GTS2F400.4435
C London Road GTS2F400.4436
C BRACKNELL GTS2F400.4437
C Berkshire UK GTS2F400.4438
C RG12 2SZ GTS2F400.4439
C GTS2F400.4440
C If no contract has been raised with this copy of the code, the use, GTS2F400.4441
C duplication or disclosure of it is strictly prohibited. Permission GTS2F400.4442
C to do so must first be obtained in writing from the Head of Numerical GTS2F400.4443
C Modelling at the above address. GTS2F400.4444
C ******************************COPYRIGHT****************************** GTS2F400.4445
C GTS2F400.4446
C*LL SUBROUTINE IMPL_CAL ---------------------------------------------- IMPLCA2C.3
CLL IMPLCA2C.4
CLL Purpose: Calculate increments for surface temperature and IMPLCA2C.5
CLL atmospheric variables in the boundary layer, using an IMPLCA2C.6
CLL implicit numerical scheme. The tridiagonal matrices are IMPLCA2C.7
CLL inverted using simple Gaussian elimination. See written IMPLCA2C.8
CLL documentation for a very detailed explanation. IMPLCA2C.9
CLL IMPLCA2C.10
CLL Suitable for single column use; activate *IF definition IBM. IMPLCA2C.11
CLL IMPLCA2C.12
CLL Model Modification history from model version 3.0 IMPLCA2C.13
CLL version Date IMPLCA2C.14
CLL 3.1 12/01/93 Alternative, more complete implicit numerical IMPLCA2C.15
CLL scheme made available under *IF DEF,A03_2C. IMPLCA2C.16
CLL 3.3 07/07/93 Only versions 1B & 2C carried forward to UM3.3 RB070793.2
CLL 3.4 06/06/94 DEF TIMER replaced by LOGICAL LTIMER ASJ1F304.360
CLL Argument LTIMER added ASJ1F304.361
CLL S.J.Swarbrick ASJ1F304.362
! 3.5 9/5/95 MPP code: Change updateable area P.Burton APB1F305.329
CLL 4.1 08/05/96 decks A03_2C and A03_3B removed ASJ4F401.6
CLL S D Jackson ASJ4F401.7
CLL 4.2 Oct. 96 T3E migration - *DEF CRAY removed GSS1F402.62
CLL S J Swarbrick GSS1F402.63
!LL 4.3 14/01/97 MPP code : Corrected setting of polar rows GPB1F403.89
!LL P.Burton GPB1F403.90
CLL 4.5 Jul. 98 Kill the IBM specific lines (JCThil) AJC1F405.163
CLL IMPLCA2C.17
CLL ASJ1F304.363
CLL FH, RNBS <- Programmers of some or all of previous code or changes IMPLCA2C.18
CLL IMPLCA2C.19
CLL Programming standard: UM Documentation Paper No 4, Version 2, IMPLCA2C.20
CLL dated 18/1/90 IMPLCA2C.21
CLL IMPLCA2C.22
CLL System component covered: P244 IMPLCA2C.23
CLL IMPLCA2C.24
CLL Project task: P24 IMPLCA2C.25
CLL IMPLCA2C.26
CLL Documentation: UM Documentation Paper No 24. IMPLCA2C.27
CLL IMPLCA2C.28
C*---------------------------------------------------------------------- IMPLCA2C.29
C*L Arguments :- IMPLCA2C.30
SUBROUTINE IMPL_CAL ( 2,8IMPLCA2C.31
+ P_FIELD,U_FIELD,P1,U1, IMPLCA2C.32
+ P_POINTS,U_POINTS,BL_LEVELS,ROW_LENGTH,P_ROWS,U_ROWS IMPLCA2C.33
+,CDR10M,U0,V0,SU10,SV10 IMPLCA2C.34
+,ASOIL_1,DELTA_AK,DELTA_BK,FQW_LEAD,FTL_LEAD IMPLCA2C.35
+,GAMMA_RHOKE,GAMMA_RHOKEA,GAMMA_RHOKES,GAMMA_RHOKESL IMPLCA2C.36
+,GAMMA_RHOKH_RDZ,GAMMA_RHOKM_RDZUV IMPLCA2C.37
+,ICE_FRACT,LAND_MASK,LYING_SNOW,PSTAR,RADNET IMPLCA2C.38
+,SEA_ICE_HTF,SOIL_HT_FLUX,TIMESTEP,TSTAR_NL IMPLCA2C.39
+,EA,ES,ESL,FQW,FTL,QW,GAMMA_RHOKH_1,GAMMA_RHOKM_1,TSTAR,TL,U,V IMPLCA2C.40
+,DQW_1,DQW_RML,DTRDZ,DTRDZ_RML,DTSTAR,E_SEA,H_SEA,TAUX,TAUY IMPLCA2C.41
+,U10M,V10M IMPLCA2C.42
+,NRML IMPLCA2C.43
+,ERROR,LTIMER ASJ1F304.364
+) IMPLCA2C.45
IMPLICIT NONE IMPLCA2C.46
LOGICAL LTIMER ASJ1F304.365
INTEGER IMPLCA2C.47
+ P_FIELD ! IN No. of points in P-grid. IMPLCA2C.48
+,U_FIELD ! IN No. of points in UV-grid. IMPLCA2C.49
+,P1 ! IN First point to be processed in IMPLCA2C.50
C ! P-grid. IMPLCA2C.51
+,U1 ! IN First point to be processed in IMPLCA2C.52
C ! UV-grid. IMPLCA2C.53
+,P_POINTS ! IN Number of P-grid points to be IMPLCA2C.54
C ! processed. IMPLCA2C.55
+,U_POINTS ! IN Number of UV-grid points. IMPLCA2C.59
+,BL_LEVELS ! IN No. of atmospheric levels for IMPLCA2C.63
C ! which boundary layer fluxes are IMPLCA2C.64
C ! calculated. IMPLCA2C.65
+,ROW_LENGTH ! IN No. of points in latitude row. IMPLCA2C.66
+,P_ROWS ! IN No. of P-rows of data to be IMPLCA2C.70
C ! processed. IMPLCA2C.71
+,U_ROWS ! IN No. of UV-rows of data to be IMPLCA2C.75
C ! processed. IMPLCA2C.76
REAL IMPLCA2C.80
+ CDR10M(U_FIELD) ! IN Used in calc of 10m wind - from IMPLCA2C.81
C ! P243 (routine SF_EXCH). First IMPLCA2C.82
C ! and last rows are "missing data" IMPLCA2C.83
C ! and not used. UV-grid. IMPLCA2C.84
+,U0(U_FIELD) ! IN Westerly component of surface IMPLCA2C.85
C ! current (m/s; 0 over land). UVG. IMPLCA2C.86
+,V0(U_FIELD) ! IN Southerly component of surface IMPLCA2C.87
C ! current (m/s; 0 over land). UVG. IMPLCA2C.88
REAL IMPLCA2C.89
+ ASOIL_1(P_FIELD) ! IN Soil thermodynamic coefficient IMPLCA2C.90
C ! from P242 (SOIL_HTF). Sq m K IMPLCA2C.91
C ! per Joule * timestep. IMPLCA2C.92
+,DELTA_AK(BL_LEVELS) ! IN Difference of hybrid 'A' across IMPLCA2C.93
C ! boundary layers (K-1/2 to K+1/2) IMPLCA2C.94
C ! (upper minus lower). IMPLCA2C.95
+,DELTA_BK(BL_LEVELS) ! IN Difference of hybrid 'B' across IMPLCA2C.96
C ! boundary layers (K-1/2 to K+1/2) IMPLCA2C.97
C ! (upper minus lower). IMPLCA2C.98
+,FQW_LEAD(P_FIELD) ! IN "Explicit" surface flux of IMPLCA2C.99
C ! QW (i.e. evaporation), on IMPLCA2C.100
C ! P-grid for leads fraction IMPLCA2C.101
C ! of gridsquare (kg/m2/s). IMPLCA2C.102
C ! Missing data at non-sea-ice pts IMPLCA2C.103
+,FTL_LEAD(P_FIELD) ! IN "Explicit surface flux of TL IMPLCA2C.104
C ! = H/CP (sensible heat/Cp) for IMPLCA2C.105
C ! leads fraction of gridsquare. IMPLCA2C.106
C ! Missing data at non-sea-ice pts IMPLCA2C.107
+,GAMMA_RHOKE(P_FIELD) ! IN Surface exchange coeff. for FQW, IMPLCA2C.108
C ! =GAMMA(1)*RHOKE from SF_EXCH. IMPLCA2C.109
+,GAMMA_RHOKEA(P_FIELD) ! IN Surface exchange coeff. for EA, IMPLCA2C.110
C ! =GAMMA(1)*RHOKEA from SF_EXCH. IMPLCA2C.111
+,GAMMA_RHOKES(P_FIELD) ! IN Surface exchange coeff. for ES, IMPLCA2C.112
C ! =GAMMA(1)*RHOKES from SF_EXCH. IMPLCA2C.113
+,GAMMA_RHOKESL(P_FIELD) ! IN Surface exchange coeff. for ESL, IMPLCA2C.114
C ! =GAMMA(1)*RHOKESL from SF_EXCH. IMPLCA2C.115
+,GAMMA_RHOKH_RDZ(P_FIELD,2:BL_LEVELS) IMPLCA2C.116
C ! IN Exchange coeff for FTL above IMPLCA2C.117
C ! surface, =GAMMA(K)*RHOKH(,K) IMPLCA2C.118
C ! *RDZ(K) for K>=2 (from KMKH). IMPLCA2C.119
+,GAMMA_RHOKM_RDZUV(U_FIELD,2:BL_LEVELS) IMPLCA2C.120
C ! IN Exchange coefficients for IMPLCA2C.121
C ! momentum, on UV-grid with IMPLCA2C.122
C ! first and last rows ignored. IMPLCA2C.123
C ! =GAMMA(K)*RHOKM(,K)*RDZUV(,K) IMPLCA2C.124
C ! for K>=2 (from KMKH). IMPLCA2C.125
REAL ! Split to avoid > 19 continuations. IMPLCA2C.126
+ ICE_FRACT(P_FIELD) ! IN Fraction of grid-box which is IMPLCA2C.127
C ! sea-ice (decimal fraction). IMPLCA2C.128
+,LYING_SNOW(P_FIELD) ! IN Lying snow (kg per sq m ie "mm") IMPLCA2C.129
+,PSTAR(P_FIELD) ! IN Surface pressure (Pa). IMPLCA2C.130
+,RADNET(P_FIELD) ! IN Area weighted ice component IMPLCA2C.131
C ! of surface net radiation. IMPLCA2C.132
C ! (+ve downwards, W per sq m) IMPLCA2C.133
+,SEA_ICE_HTF(P_FIELD) ! IN Heat flux through sea-ice (W per IMPLCA2C.134
C ! sq m, +ve downwards). From P241 IMPLCA2C.135
C ! (routine SICE_HTF). IMPLCA2C.136
+,SOIL_HT_FLUX(P_FIELD) ! IN Heat flux from soil layer 1 to IMPLCA2C.137
C ! soil layer 2 (deep soil layer 1) IMPLCA2C.138
C ! i.e. +ve is downwards (W per sq IMPLCA2C.139
C ! m). From P242 (SOIL_HTF). IMPLCA2C.140
+,TIMESTEP ! IN Timestep in seconds. IMPLCA2C.141
+,TSTAR_NL(P_FIELD) ! IN TSTAR No Leads: surface temp- IMPLCA2C.142
C ! erature of the ice at sea-ice IMPLCA2C.143
C ! points and =TSTAR elsewhere.(K) IMPLCA2C.144
LOGICAL IMPLCA2C.145
+ LAND_MASK(P_FIELD) ! IN T for land, F elsewhere. IMPLCA2C.146
+,SU10 ! IN STASH flag for 10-metre W wind. IMPLCA2C.147
+,SV10 ! IN STASH flag for 10-metre S wind. IMPLCA2C.148
C IMPLCA2C.149
C Next 5 arrays are all IN as "explicit" fluxes from P243 (SF_EXCH and IMPLCA2C.150
C possibly KMKH), and OUT as "implicit" fluxes. IMPLCA2C.151
C IMPLCA2C.152
REAL IMPLCA2C.153
+ EA(P_FIELD) ! INOUT Surface evaporation with only IMPLCA2C.154
C ! aerodynamic resistance (+ve), IMPLCA2C.155
C ! or condensation (-ve), IMPLCA2C.156
C ! averaged over gridbox. IMPLCA2C.157
C ! Kg per sq m per sec. IMPLCA2C.158
+,ES(P_FIELD) ! INOUT Surface evapotranspiration IMPLCA2C.159
C ! (through a not-entirely IMPLCA2C.160
C ! aerodynamic resistance). IMPLCA2C.161
C ! Always non-negative. Kg per IMPLCA2C.162
C ! sq m per sec. IMPLCA2C.163
+,ESL(P_FIELD) ! INOUT ES without fractional IMPLCA2C.164
C ! weighting factor FRACS. "L" IMPLCA2C.165
C ! is for "local". Kg/sq m/sec. IMPLCA2C.166
+,FQW(P_FIELD,BL_LEVELS) ! INOUT Flux of QW (ie., for surface, IMPLCA2C.167
C ! total evaporation). Kg/sq m/s IMPLCA2C.168
+,FTL(P_FIELD,BL_LEVELS) ! INOUT Flux of TL (ie., for surface, IMPLCA2C.169
C ! H/Cp where H is sensible heat IMPLCA2C.170
C ! in W per sq m). IMPLCA2C.171
+,QW(P_FIELD,BL_LEVELS) ! INOUT Total water content (kg per IMPLCA2C.172
C ! kg air). From P243. IMPLCA2C.173
+,GAMMA_RHOKH_1(P_FIELD) ! IN Surface exchange coeffs for IMPLCA2C.174
C ! FTL, =GAMMA(1)*RHOKH(,1) IMPLCA2C.175
C ! from P243 (SF_EXCH). IMPLCA2C.176
C ! OUT =RHOKH_1 to satisfy STASH. IMPLCA2C.177
+,GAMMA_RHOKM_1(U_FIELD) ! IN Surface exchange coeffs for IMPLCA2C.178
C ! momentum, on UV-grid IMPLCA2C.179
C ! with first and last rows IMPLCA2C.180
C ! ignored. =GAMMA(1)*RHOKM(,1) IMPLCA2C.181
C ! from P243 (SF_EXCH). IMPLCA2C.182
C ! OUT =RHOKM_1 to satisfy STASH. IMPLCA2C.183
+,TSTAR(P_FIELD) ! INOUT Gridbox mean surface IMPLCA2C.184
C ! temperature. IMPLCA2C.185
+,TL(P_FIELD,BL_LEVELS) ! INOUT Liquid/frozen water IMPLCA2C.186
C ! temperature (K). From P243. IMPLCA2C.187
+,U(U_FIELD,BL_LEVELS) ! INOUT Westerly wind (m/s). On UV- IMPLCA2C.188
C ! grid, 1st & last rows unused. IMPLCA2C.189
+,V(U_FIELD,BL_LEVELS) ! INOUT Southerly wind (m/s). On UV- IMPLCA2C.190
C ! grid, 1st & last rows unused. IMPLCA2C.191
REAL IMPLCA2C.192
+ DQW_1(P_FIELD) ! OUT Increment for QW(,1) (needed IMPLCA2C.193
C ! again in P245). IMPLCA2C.194
+,DQW_RML(P_FIELD) ! OUT Rapidly mixing layer increment IMPLCA2C.195
C ! to QW (needed again in P245). IMPLCA2C.196
+,DTRDZ(P_FIELD,BL_LEVELS) ! OUT -g.dt/dp for bottom BL_LEVELS IMPLCA2C.197
C ! model layers (needed in P245). IMPLCA2C.198
+,DTRDZ_RML(P_FIELD) ! OUT -g.dt/dp for the rapidly IMPLCA2C.199
C ! mixing layer (needed in P245). IMPLCA2C.200
+,DTSTAR(P_FIELD) ! OUT Surface temperature increment IMPLCA2C.201
C ! (needed again in P245). IMPLCA2C.202
+,E_SEA(P_FIELD) ! OUT Evaporation from sea times IMPLCA2C.203
C ! leads fraction. Zero over land IMPLCA2C.204
C ! (kg per square metre per sec). IMPLCA2C.205
+,H_SEA(P_FIELD) ! OUT Surface sensible heat flux IMPLCA2C.206
C ! over sea times leads fraction. IMPLCA2C.207
C ! Zero over land. ( W/m2) IMPLCA2C.208
+,TAUX(U_FIELD,BL_LEVELS) ! INOUT x-component of turbulent IMPLCA2C.209
C ! stress at levels k-1/2; IMPLCA2C.210
C ! eg. TAUX(,1) is surface stress. IMPLCA2C.211
C ! UV-grid, 1st and last rows set IMPLCA2C.212
C ! to "missing data". (N/sq m) IMPLCA2C.213
+,TAUY(U_FIELD,BL_LEVELS) ! INOUT y-component of turbulent IMPLCA2C.214
C ! stress at levels k-1/2; IMPLCA2C.215
C ! eg. TAUY(,1) is surface stress. IMPLCA2C.216
C ! UV-grid, 1st and last rows set IMPLCA2C.217
C ! to "missing data". (N/sq m) IMPLCA2C.218
+,U10M(U_FIELD) ! OUT Westerly wind at 10m (m/s). IMPLCA2C.219
C ! 1st & last rows "missing data". IMPLCA2C.220
+,V10M(U_FIELD) ! OUT Southerly wind at 10m (m/s). IMPLCA2C.221
C ! 1st & last rows "missing data". IMPLCA2C.222
INTEGER IMPLCA2C.223
+ NRML(P_FIELD) ! IN The number of model layers IMPLCA2C.224
C ! in the unstable rapidly mixing IMPLCA2C.225
C ! layer. Zero if surface layer IMPLCA2C.226
C ! is stable. IMPLCA2C.227
+,ERROR ! OUT 1 if bad arguments, else 0. IMPLCA2C.228
C* IMPLCA2C.229
C*L External references :- IMPLCA2C.230
EXTERNAL QSAT IMPLCA2C.231
*IF -DEF,SCMA AJC1F405.164
EXTERNAL P_TO_UV IMPLCA2C.233
*ENDIF IMPLCA2C.234
EXTERNAL TIMER IMPLCA2C.236
C* IMPLCA2C.238
C*L Local and other symbolic constants :- IMPLCA2C.239
*CALL C_EPSLON
IMPLCA2C.240
*CALL C_G
IMPLCA2C.241
*CALL C_LHEAT
IMPLCA2C.242
*CALL C_R_CP
IMPLCA2C.243
*CALL C_GAMMA
IMPLCA2C.244
*CALL C_SICEHC
! Holds integer AI IMPLCA2C.245
REAL LS IMPLCA2C.246
PARAMETER ( IMPLCA2C.247
+ LS=LC+LF ! Latent heat of sublimation (J per kg). IMPLCA2C.248
+) IMPLCA2C.249
C* IMPLCA2C.250
C*L Workspace :- IMPLCA2C.251
C 6*BL_LEVELS + 2 blocks of real workspace are required. IMPLCA2C.252
REAL IMPLCA2C.254
+ ALPHAS(P_FIELD) ! Partial dQsat/dT* from Clausius- IMPLCA2C.255
C ! Clapeyron relation. IMPLCA2C.256
+,AQ_AM(U_FIELD,BL_LEVELS) ! As AQ: "Q" elements of matrix in IMPLCA2C.257
C ! eqn P244.79 (modified during IMPLCA2C.258
C ! Gaussian elimination process). IMPLCA2C.259
C ! As AM: elements of matrix in eqn IMPLCA2C.260
C ! P244.80 (also get modified). IMPLCA2C.261
+,AQ_RML(U_FIELD) ! Matrix element for humidity in IMPLCA2C.262
C ! rapidly mixing layer. Then briefly IMPLCA2C.263
C ! used for DELTAP on the UV grid. IMPLCA2C.264
+,AT_ATQ(P_FIELD,BL_LEVELS) ! Elements in atmospheric T rows of IMPLCA2C.265
C ! matrix in eqn P244.79 (modified IMPLCA2C.266
C ! during Gaussian elimination). IMPLCA2C.267
+,AT_RML(P_FIELD) ! Matrix element for temperature in IMPLCA2C.268
C ! rapidly mixing layer. IMPLCA2C.269
+,ATSTAR(P_FIELD) ! Element in surface temperature row IMPLCA2C.270
C ! of matrix in eqn P244.79 (modified IMPLCA2C.271
C ! during Gaussian elimination). IMPLCA2C.272
C ! Also used for saturated sp hum, and IMPLCA2C.273
C ! in interpolations to UV-grid. IMPLCA2C.274
+,DELTAP(U_FIELD,BL_LEVELS) ! Vertical pressure difference across IMPLCA2C.275
C ! hybrid layers (upper minus lower) IMPLCA2C.276
C ! (Pa). IMPLCA2C.277
+,DELTAP_RML(U_FIELD) ! Vertical pressure difference across IMPLCA2C.278
C ! the rapidly mixing layer (Pa). IMPLCA2C.279
+,DQW_DU(U_FIELD,BL_LEVELS) ! As DQW: delta QW elements of vector IMPLCA2C.280
C ! on RHS, then LHS, of eqn P244.79. IMPLCA2C.281
C ! As DU: delta U elements of vector IMPLCA2C.282
C ! on RHS, then LHS, of eqn P244.80. IMPLCA2C.283
+,DTL_DV(U_FIELD,BL_LEVELS) ! As DTL: delta TL (for atmosphere) IMPLCA2C.284
C ! elements of vector on RHS, then IMPLCA2C.285
C ! LHS, of eqn P244.79. IMPLCA2C.286
C ! As DV: delta V elements of vector IMPLCA2C.287
C ! on RHS, then LHS, of eqn P244.80. IMPLCA2C.288
+,DTL_RML(U_FIELD) ! Delta TL for rapidly mixing layer. IMPLCA2C.289
+,DTRDZ_UV(U_FIELD,BL_LEVELS) ! -g.dt/dp for model layers IMPLCA2C.290
C ! interpolated to the UV-grid. IMPLCA2C.291
+,FQW_ICE(P_FIELD) ! "Explicit" surface flux of QW for IMPLCA2C.292
C ! sea-ice fraction of gridsquare. IMPLCA2C.293
+,FTL_ICE(P_FIELD) ! "Explicit" surface flux of TL for IMPLCA2C.294
C ! sea-ice fraction of gridsquare. IMPLCA2C.295
C* IMPLCA2C.340
C Local scalars :- IMPLCA2C.341
REAL IMPLCA2C.342
+ CTQ ! Matrix element in P244.??, for local increments to rml IMPLCA2C.343
+,CM ! Matrix element in eqn P244.80. IMPLCA2C.344
+,CQ ! Matrix element in "Q" row in eqn P244.79. IMPLCA2C.345
+,CQ_RML ! As above but for rapidly mixing layer increment. IMPLCA2C.346
+,CT ! Matrix element in "T" row in eqn P244.79. IMPLCA2C.347
+,CT_RML ! As above but for rapidly mixing layer increment. IMPLCA2C.348
+,RBTQ ! Reciprocal of B P244.??, for local increments to rml IMPLCA2C.349
+,RBM ! Reciprocal of BM(') (eqns P244.81, 85, 89). IMPLCA2C.350
+,RBQ ! Reciprocal of BQ(') (eqns P244.98, 101, 104). IMPLCA2C.351
+,RBQ_RML ! As above but for the rapidly mixing layer increment. IMPLCA2C.352
+,RBT ! Reciprocal of BT' (eqns P244.107, 110, 113). IMPLCA2C.353
+,RBT_RML ! As above but for the rapidly mixing layer increment. IMPLCA2C.354
+,DELTDQ ! P244.126 times GAMMA(1) (factor required in eqns IMPLCA2C.355
C ! P244.122-4). IMPLCA2C.356
+,DTRICE ! = DTSTAR/ICE_FRACT at sea-ice points. IMPLCA2C.357
+,DTIMEG ! TIMESTEP * G (used in several places). IMPLCA2C.358
+,LAT_HEAT ! Latent heat (either LC or LS) in section 2.2. IMPLCA2C.359
INTEGER IMPLCA2C.360
+ BLM1 ! BL_LEVELS minus 1. IMPLCA2C.361
+,NRMLP1 ! NRML plus 1. IMPLCA2C.362
+,I ! Loop counter (horizontal field index). IMPLCA2C.363
+,J ! Offset version of I. IMPLCA2C.364
+,K ! Loop counter (vertical index). IMPLCA2C.365
+,KM1 ! K minus 1. IMPLCA2C.366
+,KP1 ! K plus 1. IMPLCA2C.367
*IF DEF,MPP GPB1F403.91
! MPP Common block GPB1F403.92
*CALL PARVARS
GPB1F403.93
*ENDIF GPB1F403.94
C IMPLCA2C.368
C----------------------------------------------------------------------- IMPLCA2C.369
CL 0. Check that the scalars input to define the grid are consistent. IMPLCA2C.370
C See comments to routine SF_EXCH for details. IMPLCA2C.371
C----------------------------------------------------------------------- IMPLCA2C.372
C IMPLCA2C.373
IF (LTIMER) THEN ASJ1F304.366
CALL TIMER
('IMPLCAL ',3) IMPLCA2C.375
ENDIF ASJ1F304.367
ERROR=0 IMPLCA2C.377
*IF DEF,MPP AJC1F405.165
IF ( AJC1F405.166
*ELSEIF DEF,SCMA AJC1F405.167
IF ( U_ROWS .NE. (P_ROWS) .OR. AJC1F405.168
*ELSE AJC1F405.169
IF ( U_ROWS .NE. (P_ROWS+1) .OR. AJC1F405.170
*ENDIF AJC1F405.171
+ U_POINTS .NE. (U_ROWS*ROW_LENGTH) .OR. IMPLCA2C.386
+ P_POINTS .NE. (P_ROWS*ROW_LENGTH) ) THEN IMPLCA2C.387
ERROR=1 IMPLCA2C.388
GOTO999 IMPLCA2C.389
ENDIF IMPLCA2C.390
DTIMEG = TIMESTEP * G IMPLCA2C.392
BLM1 = BL_LEVELS-1 IMPLCA2C.393
C IMPLCA2C.394
CL---------------------------------------------------------------------- IMPLCA2C.395
CL (A) Calculations on P-grid. IMPLCA2C.396
CL---------------------------------------------------------------------- IMPLCA2C.397
C IMPLCA2C.398
C----------------------------------------------------------------------- IMPLCA2C.399
CL 1. Calculate pressure across layer (in hybrid coordinates), DELTAP, IMPLCA2C.400
CL and then -gdt/dP = dt/rho*dz for use throughout section (A) IMPLCA2C.401
C----------------------------------------------------------------------- IMPLCA2C.402
C IMPLCA2C.403
DO 1 K=1,BL_LEVELS IMPLCA2C.404
DO 11 I=P1,P1+P_POINTS-1 IMPLCA2C.405
DELTAP(I,K) = DELTA_AK(K) + PSTAR(I)*DELTA_BK(K) IMPLCA2C.406
DTRDZ(I,K) = -DTIMEG / DELTAP(I,K) IMPLCA2C.407
11 CONTINUE IMPLCA2C.408
1 CONTINUE IMPLCA2C.409
C IMPLCA2C.410
C----------------------------------------------------------------------- IMPLCA2C.411
CL 2. Calculate implicit T and Q increments due to local mixing within IMPLCA2C.412
CL the rapidly mixing layer (where it exists). IMPLCA2C.413
CL The surface fluxes FTL(I,1), FQW(I,1) are used for calculating IMPLCA2C.414
CL the rapidly mixing layer (rml) increments but not here. IMPLCA2C.415
CL Therefore the matrix equation we must solve to find the implicit IMPLCA2C.416
CL T and Q increments due to local mixing within the rml does not IMPLCA2C.417
CL have a "surface" row and we can solve for the T and Q increments IMPLCA2C.418
CL for K = 1 to NRML simultaneously. IMPLCA2C.419
C----------------------------------------------------------------------- IMPLCA2C.420
C IMPLCA2C.421
C----------------------------------------------------------------------- IMPLCA2C.422
CL 2.1 Start 'upward sweep' with lowest model layer, which will be the IMPLCA2C.423
CL bottom of the rapidly mixing layer (rml) if it exists. IMPLCA2C.424
C----------------------------------------------------------------------- IMPLCA2C.425
C IMPLCA2C.426
DO 21 I=P1,P1+P_POINTS-1 IMPLCA2C.427
IF (NRML(I) .GE. 2) THEN IMPLCA2C.428
C IMPLCA2C.429
C "Explicit" increments due to local mixing within the rml. IMPLCA2C.430
C P244.49/31 but surface flux used in rml increment calculations. IMPLCA2C.431
C IMPLCA2C.432
DQW_DU(I,1) = -DTRDZ(I,1) * FQW(I,2) IMPLCA2C.433
DTL_DV(I,1) = -DTRDZ(I,1) * FTL(I,2) IMPLCA2C.434
C IMPLCA2C.435
C Define matrix elements (CTQ always zero for this case). IMPLCA2C.436
C IMPLCA2C.437
AT_ATQ(I,1) = -DTRDZ(I,1) * GAMMA_RHOKH_RDZ(I,2) ! P244.28 IMPLCA2C.438
RBTQ = 1.0 / ( 1.0 - AT_ATQ(I,1) ) ! Reciprocal of P244.110 IMPLCA2C.439
C IMPLCA2C.440
C Now start Gaussian elimination IMPLCA2C.441
C IMPLCA2C.442
DQW_DU(I,1) = RBTQ * DQW_DU(I,1) ! P244.102 IMPLCA2C.443
DTL_DV(I,1) = RBTQ * DTL_DV(I,1) ! P244.111 IMPLCA2C.444
AT_ATQ(I,1) = RBTQ * AT_ATQ(I,1) ! P244.112 IMPLCA2C.445
C IMPLCA2C.446
C Start calculating DELTAP_RML. Mid-level depths added in 2.2 below. IMPLCA2C.447
C IMPLCA2C.448
DELTAP_RML(I) = DELTAP(I,1) IMPLCA2C.449
ELSE ! No rapidly mixing layer calculations IMPLCA2C.450
DTRDZ_RML(I) = 1.E30 IMPLCA2C.451
DQW_RML(I) = 1.E30 IMPLCA2C.452
DTL_RML(I) = 1.E30 IMPLCA2C.453
AQ_RML(I) = 1.E30 IMPLCA2C.454
AT_RML(I) = 1.E30 IMPLCA2C.455
DELTAP_RML(I) = 1.E30 IMPLCA2C.456
ENDIF IMPLCA2C.457
21 CONTINUE IMPLCA2C.458
C IMPLCA2C.459
C----------------------------------------------------------------------- IMPLCA2C.460
CL 2.2 Continue upward sweep through middle of the rapidly mixing layer IMPLCA2C.461
CL (if it exists) and to its top. NB NRML is always <= BLM1. IMPLCA2C.462
C----------------------------------------------------------------------- IMPLCA2C.463
C IMPLCA2C.464
DO 22 K=2,BLM1 IMPLCA2C.465
KP1 = K+1 IMPLCA2C.466
KM1 = K-1 IMPLCA2C.467
DO 221 I=P1,P1+P_POINTS-1 IMPLCA2C.468
C IMPLCA2C.469
C If in the top rapidly mixing layer then do not include flux at its IMPLCA2C.470
C top in the calculation ie FQW(I,NRML+1) and FTL(I,NRML+1) are not IMPLCA2C.471
C included here; they are taken account of in the non-local mixing IMPLCA2C.472
C through the "rapidly mixing layer". IMPLCA2C.473
C IMPLCA2C.474
IF ( K .EQ. NRML(I) ) THEN IMPLCA2C.475
C IMPLCA2C.476
C Add final DELTAP contribution to DELTAP_RML and then calculate IMPLCA2C.477
C DTRDZ_RML. Lower level contributions added in 2.1 and below. IMPLCA2C.478
C IMPLCA2C.479
DELTAP_RML(I) = DELTAP_RML(I) + DELTAP(I,K) IMPLCA2C.480
DTRDZ_RML(I) = -DTIMEG / DELTAP_RML(I) IMPLCA2C.481
C IMPLCA2C.482
C "Explicit" flux divergence across layer giving explicit IMPLCA2C.483
C increment due to the local mixing at the top of rml. IMPLCA2C.484
C IMPLCA2C.485
DQW_DU(I,K) = DTRDZ(I,K) * FQW(I,K) IMPLCA2C.486
DTL_DV(I,K) = DTRDZ(I,K) * FTL(I,K) IMPLCA2C.487
C IMPLCA2C.488
C Define matrix elements (A always zero for this case). IMPLCA2C.489
C IMPLCA2C.490
CTQ = -DTRDZ(I,K) * GAMMA_RHOKH_RDZ(I,K) ! P244.36 IMPLCA2C.491
RBTQ = 1.0 / ( 1.0 - CTQ*( 1.0 + AT_ATQ(I,KM1) ) ) IMPLCA2C.492
C ... Reciprocal of P244.113 IMPLCA2C.493
C Now start Gaussian elimination IMPLCA2C.494
C IMPLCA2C.495
DQW_DU(I,K) = RBTQ * ( DQW_DU(I,K) - CTQ*DQW_DU(I,KM1) ) IMPLCA2C.496
DTL_DV(I,K) = RBTQ * ( DTL_DV(I,K) - CTQ*DTL_DV(I,KM1) ) IMPLCA2C.497
ELSEIF (K .LT. NRML(I)) THEN IMPLCA2C.498
C IMPLCA2C.499
C Add layer depths to form total rml depth. IMPLCA2C.500
C IMPLCA2C.501
DELTAP_RML(I) = DELTAP_RML(I) + DELTAP(I,K) IMPLCA2C.502
C IMPLCA2C.503
C "Explicit" flux divergence across layer giving explicit IMPLCA2C.504
C increment due to the local mixing.P244.54/38 IMPLCA2C.505
C IMPLCA2C.506
DQW_DU(I,K) = -DTRDZ(I,K) * ( FQW(I,KP1) - FQW(I,K) ) IMPLCA2C.507
DTL_DV(I,K) = -DTRDZ(I,K) * ( FTL(I,KP1) - FTL(I,K) ) IMPLCA2C.508
C IMPLCA2C.509
C Define matrix elements. IMPLCA2C.510
C IMPLCA2C.511
AT_ATQ(I,K) = -DTRDZ(I,K) * GAMMA_RHOKH_RDZ(I,KP1) ! P244.35 IMPLCA2C.512
CTQ = -DTRDZ(I,K) * GAMMA_RHOKH_RDZ(I,K) ! P244.36 IMPLCA2C.513
RBTQ = 1.0 / ( 1.0 - AT_ATQ(I,K) IMPLCA2C.514
& - CTQ * ( 1.0 + AT_ATQ(I,KM1) ) ) IMPLCA2C.515
C ... Reciprocal of P244.113 IMPLCA2C.516
C Now start Gaussian elimination IMPLCA2C.517
C IMPLCA2C.518
DQW_DU(I,K) = RBTQ * ( DQW_DU(I,K) - CTQ*DQW_DU(I,KM1) ) IMPLCA2C.519
DTL_DV(I,K) = RBTQ * ( DTL_DV(I,K) - CTQ*DTL_DV(I,KM1) ) IMPLCA2C.520
AT_ATQ(I,K) = RBTQ * AT_ATQ(I,K) ! P244.115 IMPLCA2C.521
ENDIF IMPLCA2C.522
221 CONTINUE IMPLCA2C.523
22 CONTINUE IMPLCA2C.524
C IMPLCA2C.525
C----------------------------------------------------------------------- IMPLCA2C.526
CL 2.3 Downward sweep through matrix. Add implicit increments due to IMPLCA2C.527
CL local mixing within the rapidly mixing layer. Update fluxes of IMPLCA2C.528
CL heat and moisture at the surface and the top-of-rml using IMPLCA2C.529
CL local mixing increments for layers 1 and NRML respectively. IMPLCA2C.530
C----------------------------------------------------------------------- IMPLCA2C.531
C IMPLCA2C.532
DO 23 K=BLM1,1,-1 IMPLCA2C.533
KP1 = K + 1 IMPLCA2C.534
DO 231 I=P1,P1+P_POINTS-1 IMPLCA2C.535
IF ((NRML(I) .GE. 2) .AND. (K .EQ. NRML(I))) THEN IMPLCA2C.536
QW(I,K) = QW(I,K) + DQW_DU(I,K) ! P244.128 IMPLCA2C.537
TL(I,K) = TL(I,K) + DTL_DV(I,K) ! P244.127 IMPLCA2C.538
FQW(I,KP1) = FQW(I,KP1) IMPLCA2C.539
& + GAMMA_RHOKH_RDZ(I,KP1)*DQW_DU(I,K) IMPLCA2C.540
FTL(I,KP1) = FTL(I,KP1) IMPLCA2C.541
& + GAMMA_RHOKH_RDZ(I,KP1)*DTL_DV(I,K) IMPLCA2C.542
ELSEIF ((NRML(I) .GE. 2) .AND. (K .LT. NRML(I))) THEN IMPLCA2C.543
DQW_DU(I,K) = DQW_DU(I,K) IMPLCA2C.544
& - AT_ATQ(I,K)*DQW_DU(I,KP1) ! P244.??? IMPLCA2C.545
DTL_DV(I,K) = DTL_DV(I,K) IMPLCA2C.546
& - AT_ATQ(I,K)*DTL_DV(I,KP1) ! P244.??? IMPLCA2C.547
QW(I,K) = QW(I,K) + DQW_DU(I,K) ! P244.128 IMPLCA2C.548
TL(I,K) = TL(I,K) + DTL_DV(I,K) ! P244.127 IMPLCA2C.549
ENDIF IMPLCA2C.550
IF ((NRML(I) .GE. 2) .AND. (K .EQ. 1)) THEN IMPLCA2C.551
IF (LAND_MASK(I)) THEN IMPLCA2C.552
FTL(I,1) = FTL(I,1) - GAMMA_RHOKH_1(I) * DTL_DV(I,1) IMPLCA2C.553
EA(I) = EA(I) - GAMMA_RHOKEA(I)*DQW_DU(I,1) ! P244.122 IMPLCA2C.554
ESL(I) = ESL(I) - GAMMA_RHOKESL(I)*DQW_DU(I,1) ! P244.124 IMPLCA2C.555
ES(I) = ES(I) - GAMMA_RHOKES(I)*DQW_DU(I,1) ! P244.123 IMPLCA2C.556
FQW(I,1) = EA(I) + ES(I) ! P244.125 IMPLCA2C.557
ELSEIF (ICE_FRACT(I) .GT. 0.0) THEN IMPLCA2C.558
FQW_ICE(I) = FQW(I,1) - FQW_LEAD(I) IMPLCA2C.559
& - GAMMA_RHOKE(I) * ICE_FRACT(I)*DQW_DU(I,1) IMPLCA2C.560
FTL_ICE(I) = FTL(I,1) - FTL_LEAD(I) IMPLCA2C.561
& - GAMMA_RHOKH_1(I)*ICE_FRACT(I)*DTL_DV(I,1) IMPLCA2C.562
FQW_LEAD(I) = FQW_LEAD(I) - GAMMA_RHOKE(I) IMPLCA2C.563
& * ( 1.0 - ICE_FRACT(I) ) * DQW_DU(I,1) IMPLCA2C.564
FTL_LEAD(I) = FTL_LEAD(I) - GAMMA_RHOKH_1(I) IMPLCA2C.565
& * ( 1.0 - ICE_FRACT(I) ) * DTL_DV(I,1) IMPLCA2C.566
FTL(I,1) = FTL_ICE(I) + FTL_LEAD(I) IMPLCA2C.567
FQW(I,1) = FQW_ICE(I) + FQW_LEAD(I) IMPLCA2C.568
ELSE ! ordinary sea point IMPLCA2C.569
FTL(I,1) = FTL(I,1) - GAMMA_RHOKH_1(I)*DTL_DV(I,1) IMPLCA2C.570
FQW(I,1) = FQW(I,1) - GAMMA_RHOKE(I)*DQW_DU(I,1) IMPLCA2C.571
ENDIF IMPLCA2C.572
ENDIF IMPLCA2C.573
231 CONTINUE IMPLCA2C.574
23 CONTINUE IMPLCA2C.575
C IMPLCA2C.576
C----------------------------------------------------------------------- IMPLCA2C.577
CL 3. Calculate those matrix and vector elements on the LHS of eqn IMPLCA2C.578
CL P244.79 which are to do with implicit solution of the moisture IMPLCA2C.579
CL transport problem at the surface, above the rml (if it exists) IMPLCA2C.580
CL and between all levels if it does not. IMPLCA2C.581
CL Begin with "upward sweep" through lower half of matrix). IMPLCA2C.582
C----------------------------------------------------------------------- IMPLCA2C.583
CL 3.1 Row of matrix applying to QW transport into top "boundary" IMPLCA2C.584
CL layer of model atmosphere. IMPLCA2C.585
C----------------------------------------------------------------------- IMPLCA2C.586
C IMPLCA2C.587
DO 31 I=P1,P1+P_POINTS-1 IMPLCA2C.588
DQW_DU(I,BL_LEVELS) = DTRDZ(I,BL_LEVELS) * FQW(I,BL_LEVELS) IMPLCA2C.589
C ... P244.58 IMPLCA2C.590
C IMPLCA2C.591
AQ_AM(I,BL_LEVELS) = -DTRDZ(I,BL_LEVELS) ! P244.56 IMPLCA2C.592
+ * GAMMA_RHOKH_RDZ(I,BL_LEVELS) IMPLCA2C.593
C IMPLCA2C.594
RBQ = 1.0 / ( 1.0 - AQ_AM(I,BL_LEVELS) ) ! Reciprocal P244.98 IMPLCA2C.595
DQW_DU(I,BL_LEVELS) = RBQ * DQW_DU(I,BL_LEVELS) ! P244.99 IMPLCA2C.596
AQ_AM(I,BL_LEVELS) = RBQ * AQ_AM(I,BL_LEVELS) ! P244.100 IMPLCA2C.597
31 CONTINUE IMPLCA2C.598
C IMPLCA2C.599
C----------------------------------------------------------------------- IMPLCA2C.600
CL 3.2 Rows of matrix applying to "middle of boundary layer" model IMPLCA2C.601
CL layers, i.e. all but the topmost and bottom layers. IMPLCA2C.602
C----------------------------------------------------------------------- IMPLCA2C.603
C IMPLCA2C.604
DO 32 K=BLM1,2,-1 IMPLCA2C.605
KP1=K+1 IMPLCA2C.606
DO 321 I=P1,P1+P_POINTS-1 IMPLCA2C.607
C IMPLCA2C.608
C "Explicit" flux divergence across layer giving explicit QW increment. IMPLCA2C.609
C IMPLCA2C.610
IF ( K .GT. NRML(I) ) THEN IMPLCA2C.611
DQW_DU(I,K) = -DTRDZ(I,K) * ( FQW(I,KP1) - FQW(I,K) ) IMPLCA2C.612
C ... P244.54 IMPLCA2C.613
C IMPLCA2C.614
CQ = -DTRDZ(I,K) * GAMMA_RHOKH_RDZ(I,KP1) ! P244.52 IMPLCA2C.615
AQ_AM(I,K) = -DTRDZ(I,K) * GAMMA_RHOKH_RDZ(I,K) ! P244.51 IMPLCA2C.616
RBQ = 1.0 / ( 1.0 - AQ_AM(I,K) - CQ*( 1.0 + AQ_AM(I,KP1) ) ) IMPLCA2C.617
C 1 2 2 1 IMPLCA2C.618
C ... reciprocal of P244.101 IMPLCA2C.619
C IMPLCA2C.620
DQW_DU(I,K) = RBQ * ( DQW_DU(I,K) - CQ*DQW_DU(I,KP1) ) IMPLCA2C.621
C ... P244.102 IMPLCA2C.622
C IMPLCA2C.623
AQ_AM(I,K) = RBQ * AQ_AM(I,K) ! P244.103 IMPLCA2C.624
ENDIF IMPLCA2C.625
321 CONTINUE IMPLCA2C.626
32 CONTINUE IMPLCA2C.627
C IMPLCA2C.628
C----------------------------------------------------------------------- IMPLCA2C.629
CL 3.3 Bottom model layer QW row of matrix equation. IMPLCA2C.630
C----------------------------------------------------------------------- IMPLCA2C.631
C IMPLCA2C.632
DO 33 I=P1,P1+P_POINTS-1 IMPLCA2C.633
IF ( NRML(I) .GE. 2 ) THEN IMPLCA2C.634
C IMPLCA2C.635
C----------------------------------------------------------------------- IMPLCA2C.636
CL 3.3.1 Start calculating rapidly mixing layer increments. IMPLCA2C.637
C----------------------------------------------------------------------- IMPLCA2C.638
C IMPLCA2C.639
NRMLP1 = NRML(I) + 1 IMPLCA2C.640
C IMPLCA2C.641
C "Explicit" QW increment for the rapidly mixing layer. IMPLCA2C.642
C IMPLCA2C.643
DQW_RML(I) = -DTRDZ_RML(I) * ( FQW(I,NRMLP1) - FQW(I,1) ) IMPLCA2C.644
C IMPLCA2C.645
C Define coefficients A,B,C, for implicit calculations. IMPLCA2C.646
C IMPLCA2C.647
AQ_RML(I) = -DTRDZ_RML(I) * GAMMA_RHOKE(I) IMPLCA2C.648
CQ_RML = -DTRDZ_RML(I) * GAMMA_RHOKH_RDZ(I,NRMLP1) IMPLCA2C.649
RBQ_RML = 1.0 / ( 1.0 - AQ_RML(I) IMPLCA2C.650
& - CQ_RML * ( 1.0 + AQ_AM(I,NRMLP1) ) ) IMPLCA2C.651
DQW_RML(I) = RBQ_RML * ( DQW_RML(I) IMPLCA2C.652
& - CQ_RML * DQW_DU(I,NRMLP1) ) IMPLCA2C.653
AQ_RML(I) = RBQ_RML * AQ_RML(I) IMPLCA2C.654
ELSE IMPLCA2C.655
C IMPLCA2C.656
C "Explicit" increment for QW(1) when there is no rapidly mixing IMPLCA2C.657
C layer or when it is only one model layer in depth. IMPLCA2C.658
C IMPLCA2C.659
DQW_DU(I,1) = -DTRDZ(I,1) * ( FQW(I,2) - FQW(I,1) ) ! P244.49 IMPLCA2C.660
AQ_AM(I,1) = -DTRDZ(I,1) * GAMMA_RHOKE(I) ! P244.46 IMPLCA2C.661
CQ = -DTRDZ(I,1) * GAMMA_RHOKH_RDZ(I,2) ! P244.47 IMPLCA2C.662
RBQ = 1.0 / ( 1.0 - AQ_AM(I,1) - CQ*( 1.0 + AQ_AM(I,2) ) ) IMPLCA2C.663
C 1 2 2 1 IMPLCA2C.664
C ... reciprocal of P244.104 IMPLCA2C.665
C IMPLCA2C.666
DQW_DU(I,1) = RBQ * ( DQW_DU(I,1) - CQ*DQW_DU(I,2) )! P244.105 IMPLCA2C.667
AQ_AM(I,1) = RBQ * AQ_AM(I,1) ! P244.106 IMPLCA2C.668
ENDIF IMPLCA2C.669
C IMPLCA2C.670
C----------------------------------------------------------------------- IMPLCA2C.671
CL 4. Calculate those matrix and vector elements on the LHS of eqn IMPLCA2C.672
CL P244.79 which are to do with implicit solution of the heat IMPLCA2C.673
CL transport problem (i.e. for surface temperature and ice/liquid IMPLCA2C.674
CL water temperatures in atmospheric boundary layer), and begin the IMPLCA2C.675
CL solution algorithm (perform "upward sweep" through upper half of IMPLCA2C.676
CL matrix). IMPLCA2C.677
C----------------------------------------------------------------------- IMPLCA2C.678
CL 4.1 Surface temperature row of matrix (different treatments for IMPLCA2C.679
CL different surface types). IMPLCA2C.680
C ALPHAS used as a temporary store for estimate of TSTAR. IMPLCA2C.681
C----------------------------------------------------------------------- IMPLCA2C.682
C IMPLCA2C.683
IF (LAND_MASK(I)) THEN IMPLCA2C.684
IF (LYING_SNOW(I).GT.0.0) THEN IMPLCA2C.685
LAT_HEAT = LS IMPLCA2C.686
ELSE IMPLCA2C.687
LAT_HEAT = LC IMPLCA2C.688
ENDIF IMPLCA2C.689
DTSTAR(I) = ASOIL_1(I) * ( RADNET(I) IMPLCA2C.690
& - CP*FTL(I,1) - LAT_HEAT*FQW(I,1) IMPLCA2C.691
& - SOIL_HT_FLUX(I) ) IMPLCA2C.692
C ... P244.15 (Explicit surface temperature increment; IMPLCA2C.693
C note that ASOIL_1 contains the TIMESTEP factor) IMPLCA2C.694
C IMPLCA2C.695
ALPHAS(I) = TSTAR_NL(I) + DTSTAR(I) IMPLCA2C.696
C ...(Estimate of TSTAR at timelevel n+1 based on the IMPLCA2C.697
C explicit surface temperature increment) IMPLCA2C.698
C IMPLCA2C.699
ELSEIF (ICE_FRACT(I) .GT. 0.0) THEN IMPLCA2C.700
FTL_ICE(I) = FTL(I,1) - FTL_LEAD(I) IMPLCA2C.701
FQW_ICE(I) = FQW(I,1) - FQW_LEAD(I) IMPLCA2C.702
DTSTAR(I) = TIMESTEP * AI * ( RADNET(I) IMPLCA2C.703
& - CP*FTL_ICE(I) - LS*FQW_ICE(I) IMPLCA2C.704
& - SEA_ICE_HTF(I) ) IMPLCA2C.705
C ... P244.19 (Gridbox mean explicit surface temp. increment) IMPLCA2C.706
C IMPLCA2C.707
ALPHAS(I) = TSTAR_NL(I) + DTSTAR(I)/ICE_FRACT(I) IMPLCA2C.708
C ...(Estimate of T*(I), surface temperature of sea-ice, at IMPLCA2C.709
C timelevel n+1 based on the explicit surface IMPLCA2C.710
C temperature increment.) IMPLCA2C.711
C IMPLCA2C.712
ELSE ! Sea without sea-ice. IMPLCA2C.713
DTSTAR(I) = 0.0 ! Surface temperature not updated. IMPLCA2C.714
ALPHAS(I) = TSTAR_NL(I) + DTSTAR(I) IMPLCA2C.715
C ...(Estimate of TSTAR at timelevel n+1 based on the IMPLCA2C.716
C explicit surface temperature increment) IMPLCA2C.717
C IMPLCA2C.718
ENDIF IMPLCA2C.719
33 CONTINUE IMPLCA2C.720
CALL QSAT
(ATSTAR(P1),TSTAR_NL(P1),PSTAR(P1),P_POINTS) IMPLCA2C.721
C ...Qsat for T*n (T*n(I) at sea-ice points) temporarily IMPLCA2C.722
C stored in ATSTAR. IMPLCA2C.723
C IMPLCA2C.724
CALL QSAT
(AT_ATQ(P1,1),ALPHAS(P1),PSTAR(P1),P_POINTS) IMPLCA2C.725
C ...Qsat for (T*n + explicit deltaT*), or (T*n(I) + IMPLCA2C.726
C explicit deltaT*(I) at sea-ice points), temporarily IMPLCA2C.727
C storedin AT_ATQ(*,1). IMPLCA2C.728
C IMPLCA2C.729
DO 41 I = P1,P1+P_POINTS-1 IMPLCA2C.730
IF (LAND_MASK(I)) THEN IMPLCA2C.731
IF (LYING_SNOW(I).GT.0.0) THEN IMPLCA2C.732
LAT_HEAT = LS IMPLCA2C.733
ELSE IMPLCA2C.734
LAT_HEAT = LC IMPLCA2C.735
ENDIF IMPLCA2C.736
IF (DTSTAR(I) .GE. -2.0 .AND. DTSTAR(I) .LE. 2.0) THEN IMPLCA2C.737
ALPHAS(I)=(EPSILON*LAT_HEAT*ATSTAR(I)) IMPLCA2C.738
+ / (R*TSTAR_NL(I)*TSTAR_NL(I)) IMPLCA2C.739
C ...P244.13a(Clausius-Clapeyron formula for partial dQsat/dT* IMPLCA2C.740
C used for small explicit T* increments.) IMPLCA2C.741
C IMPLCA2C.742
ELSE IMPLCA2C.743
ALPHAS(I) = ( AT_ATQ(I,1) - ATSTAR(I) ) / DTSTAR(I) IMPLCA2C.744
C ...P244.13b(partial dQsat/dT* calculated as a finite IMPLCA2C.745
C difference for large explicit T* increments.) IMPLCA2C.746
C IMPLCA2C.747
ENDIF IMPLCA2C.748
ATSTAR(I) = -ASOIL_1(I) * CP * GAMMA_RHOKH_1(I) ! P244.16 IMPLCA2C.749
CT = -GAMMA_RHOKE(I) * ASOIL_1(I) * LAT_HEAT IMPLCA2C.750
C ... P244.17 (Note: ASOIL_1 includes factor TIMESTEP) IMPLCA2C.751
C IMPLCA2C.752
IF ( NRML(I) .GE. 2 ) THEN IMPLCA2C.753
RBT = 1.0 / ( 1.0 - ATSTAR(I) - CT * ALPHAS(I) IMPLCA2C.754
& * ( 1.0 + AQ_RML(I) ) ) IMPLCA2C.755
DTSTAR(I) = RBT * ( DTSTAR(I) - CT * DQW_RML(I) ) IMPLCA2C.756
ELSE IMPLCA2C.757
RBT = 1.0 / ( 1.0 - ATSTAR(I) - CT * ALPHAS(I) IMPLCA2C.758
& * ( 1.0 + AQ_AM(I,1) ) ) IMPLCA2C.759
C 1 2 2 1 IMPLCA2C.760
C ... Reciprocal of P244.107 IMPLCA2C.761
C IMPLCA2C.762
DTSTAR(I) = RBT * ( DTSTAR(I) - CT * DQW_DU(I,1) ) IMPLCA2C.763
C ... P244.108 (giving DTSTAR') IMPLCA2C.764
C IMPLCA2C.765
ENDIF IMPLCA2C.766
ATSTAR(I) = RBT * ATSTAR(I) ! P244.109 IMPLCA2C.767
ELSEIF (ICE_FRACT(I).GT.0.0) THEN IMPLCA2C.768
DTRICE = DTSTAR(I)/ICE_FRACT(I) IMPLCA2C.769
IF (DTRICE .GE. -2.0 .AND. DTRICE .LE. 2.0) THEN IMPLCA2C.770
ALPHAS(I) = (EPSILON*LS*ATSTAR(I)) IMPLCA2C.771
+ / (R*TSTAR_NL(I)*TSTAR_NL(I)) IMPLCA2C.772
C ...P244.13a (Clausius-Clapeyron formula for partial IMPLCA2C.773
C dQsat/dT*(I) used for small explicit T*(I) increments.) IMPLCA2C.774
C IMPLCA2C.775
ELSE IMPLCA2C.776
ALPHAS(I) = ( AT_ATQ(I,1) - ATSTAR(I) ) / DTRICE IMPLCA2C.777
C ...P244.13b (partial dQsat/dT*(I) calculated as a finite IMPLCA2C.778
C difference used for large explicit T*(I) increments.) IMPLCA2C.779
C IMPLCA2C.780
ENDIF IMPLCA2C.781
ATSTAR(I) = -GAMMA_RHOKH_1(I) * TIMESTEP * AI * CP ! P244.20 IMPLCA2C.782
CT = -GAMMA_RHOKE(I) * TIMESTEP * AI * LS ! P244.21 IMPLCA2C.783
IF ( NRML(I) .GE. 2 ) THEN IMPLCA2C.784
RBT = 1.0 / ( 1.0 -ATSTAR(I) - CT * ALPHAS(I) * IMPLCA2C.785
& ( 1.0 + ICE_FRACT(I) * AQ_RML(I) ) ) IMPLCA2C.786
DTSTAR(I) = RBT * ( DTSTAR(I) - ICE_FRACT(I) IMPLCA2C.787
& * CT * DQW_RML(I) ) IMPLCA2C.788
ELSE IMPLCA2C.789
RBT = 1.0 / ( 1.0 - ATSTAR(I) - CT * ALPHAS(I) IMPLCA2C.790
& * ( 1.0 + ICE_FRACT(I)*AQ_AM(I,1) ) ) IMPLCA2C.791
C 1 2 2 1 IMPLCA2C.792
C ... Reciprocal of P2440.14 IMPLCA2C.793
C IMPLCA2C.794
DTSTAR(I) = RBT * ( DTSTAR(I) - ICE_FRACT(I)*CT*DQW_DU(I,1)) IMPLCA2C.795
C ... P2440.15, giving DTSTAR' IMPLCA2C.796
C IMPLCA2C.797
ENDIF IMPLCA2C.798
ATSTAR(I) = RBT * ATSTAR(I) * ICE_FRACT(I) ! P2440.16 IMPLCA2C.799
ELSE ! Sea with no sea-ice. IMPLCA2C.800
ALPHAS(I) = (EPSILON*LC*ATSTAR(I)) IMPLCA2C.801
+ / (R*TSTAR_NL(I)*TSTAR_NL(I)) IMPLCA2C.802
C ...P244.13a(Clausius-Clapeyron formula for partial dQsat/dT*) IMPLCA2C.803
ATSTAR(I) = 0.0 IMPLCA2C.804
ENDIF IMPLCA2C.805
C----------------------------------------------------------------------- IMPLCA2C.806
CL 4.2 Lowest atmospheric layer TL row of matrix. IMPLCA2C.807
C----------------------------------------------------------------------- IMPLCA2C.808
IF (NRML(I) .GE. 2) THEN IMPLCA2C.809
NRMLP1 = NRML(I) + 1 IMPLCA2C.810
C IMPLCA2C.811
C "Explicit" rapidly mixing layer increment for TL. IMPLCA2C.812
C IMPLCA2C.813
DTL_RML(I) = -DTRDZ_RML(I) * ( FTL(I,NRMLP1) - FTL(I,1) ) IMPLCA2C.814
AT_RML(I) = -DTRDZ_RML(I) * GAMMA_RHOKH_RDZ(I,NRMLP1) IMPLCA2C.815
CT_RML = -DTRDZ_RML(I) * GAMMA_RHOKH_1(I) IMPLCA2C.816
RBT_RML = 1.0 / ( 1.0 - AT_RML(I) - CT_RML IMPLCA2C.817
& * ( 1.0 + ATSTAR(I) ) ) IMPLCA2C.818
DTL_RML(I) = RBT_RML * ( DTL_RML(I) IMPLCA2C.819
& - CT_RML * DTSTAR(I) ) IMPLCA2C.820
AT_RML(I) = RBT_RML * AT_RML(I) IMPLCA2C.821
ELSE IMPLCA2C.822
C IMPLCA2C.823
C "Explicit" increment to TL(1) when there is no rapidly mixing layer IMPLCA2C.824
C or it does not extend beyond the bottom model layer. IMPLCA2C.825
C IMPLCA2C.826
DTL_DV(I,1) = -DTRDZ(I,1) * ( FTL(I,2) - FTL(I,1) ) ! P244.31 IMPLCA2C.827
CT = -DTRDZ(I,1) * GAMMA_RHOKH_1(I) ! P244.29 IMPLCA2C.828
AT_ATQ(I,1) = -DTRDZ(I,1) * GAMMA_RHOKH_RDZ(I,2) ! P244.28 IMPLCA2C.829
RBT = 1.0 / ( 1.0 - AT_ATQ(I,1) - CT*( 1.0 + ATSTAR(I) ) ) IMPLCA2C.830
C 1 2 2 1 IMPLCA2C.831
C ... Reciprocal of P244.110 IMPLCA2C.832
C IMPLCA2C.833
DTL_DV(I,1) = RBT * ( DTL_DV(I,1) - CT*DTSTAR(I) ) ! P244.111 IMPLCA2C.834
AT_ATQ(I,1) = RBT * AT_ATQ(I,1) ! P244.112 IMPLCA2C.835
ENDIF IMPLCA2C.836
41 CONTINUE IMPLCA2C.837
C----------------------------------------------------------------------- IMPLCA2C.838
CL 4.3 Rows of matrix applying to TL transport into model layers in the IMPLCA2C.839
CL "middle" of the "boundary" layer, i.e. all but the bottom layer IMPLCA2C.840
CL and the top "boundary" layer. IMPLCA2C.841
C----------------------------------------------------------------------- IMPLCA2C.842
DO 43 K=2,BLM1 IMPLCA2C.843
KP1 = K+1 IMPLCA2C.844
KM1 = K-1 IMPLCA2C.845
DO 431 I=P1,P1+P_POINTS-1 IMPLCA2C.846
C IMPLCA2C.847
C "Explicit" flux divergence across layer giving explicit TL increment IMPLCA2C.848
C due to mixing above rml if it exists or for all levels if it does IMPLCA2C.849
C not. IMPLCA2C.850
C IMPLCA2C.851
NRMLP1 = NRML(I) + 1 IMPLCA2C.852
IF (K .GT. NRML(I)) THEN IMPLCA2C.853
DTL_DV(I,K) = -DTRDZ(I,K) * ( FTL(I,KP1) - FTL(I,K) ) IMPLCA2C.854
C ... P244.38 IMPLCA2C.855
C IMPLCA2C.856
AT_ATQ(I,K) = -DTRDZ(I,K) * GAMMA_RHOKH_RDZ(I,KP1) ! P244.35 IMPLCA2C.857
CT = -DTRDZ(I,K) * GAMMA_RHOKH_RDZ(I,K) ! P244.36 IMPLCA2C.858
IF ((NRML(I) .GE. 2) .AND. (K .EQ. NRMLP1)) THEN IMPLCA2C.859
RBT = 1.0 / ( 1.0 - AT_ATQ(I,K) - CT*( 1.0 + AT_RML(I) ) ) IMPLCA2C.860
DTL_DV(I,K) = RBT * ( DTL_DV(I,K) - CT*DTL_RML(I) ) IMPLCA2C.861
ELSE IMPLCA2C.862
RBT = 1.0 / ( 1.0 - AT_ATQ(I,K) IMPLCA2C.863
& - CT*( 1.0 + AT_ATQ(I,KM1) ) ) IMPLCA2C.864
C 1 2 2 1 IMPLCA2C.865
C ... Reciprocal of P244.113 IMPLCA2C.866
C IMPLCA2C.867
DTL_DV(I,K) = RBT * ( DTL_DV(I,K) - CT*DTL_DV(I,KM1) ) IMPLCA2C.868
C ... P244.114 IMPLCA2C.869
C IMPLCA2C.870
ENDIF IMPLCA2C.871
AT_ATQ(I,K) = RBT * AT_ATQ(I,K) ! P244.115 IMPLCA2C.872
ENDIF IMPLCA2C.873
431 CONTINUE IMPLCA2C.874
43 CONTINUE IMPLCA2C.875
C----------------------------------------------------------------------- IMPLCA2C.876
CL 4.4 Top "boundary" layer TL row of matrix. TL for this layer can IMPLCA2C.877
CL then be, and is, updated. IMPLCA2C.878
C----------------------------------------------------------------------- IMPLCA2C.879
DO 44 I=P1,P1+P_POINTS-1 IMPLCA2C.880
DTL_DV(I,BL_LEVELS) = DTRDZ(I,BL_LEVELS) * FTL(I,BL_LEVELS) IMPLCA2C.881
C ... P244.42 IMPLCA2C.882
C IMPLCA2C.883
CT = -DTRDZ(I,BL_LEVELS) * GAMMA_RHOKH_RDZ(I,BL_LEVELS)! P244.40 IMPLCA2C.884
IF (NRML(I) .EQ. BLM1) THEN IMPLCA2C.885
RBT = 1.0 / ( 1.0 - CT*( 1.0 + AT_RML(I) ) ) IMPLCA2C.886
DTL_DV(I,BL_LEVELS) = RBT * ( DTL_DV(I,BL_LEVELS) IMPLCA2C.887
+ - CT*DTL_RML(I) ) IMPLCA2C.888
ELSE IMPLCA2C.889
RBT = 1.0 / ( 1.0 - CT*( 1.0 + AT_ATQ(I,BLM1) ) ) IMPLCA2C.890
C 1 2 2 1 IMPLCA2C.891
C ... Reciprocal of P244.116 IMPLCA2C.892
C IMPLCA2C.893
DTL_DV(I,BL_LEVELS) = RBT * ( DTL_DV(I,BL_LEVELS) ! P244.117 IMPLCA2C.894
+ - CT*DTL_DV(I,BLM1) ) IMPLCA2C.895
ENDIF IMPLCA2C.896
TL(I,BL_LEVELS) = TL(I,BL_LEVELS) + DTL_DV(I,BL_LEVELS) IMPLCA2C.897
C ... P244.127 IMPLCA2C.898
C IMPLCA2C.899
44 CONTINUE IMPLCA2C.900
C IMPLCA2C.901
C----------------------------------------------------------------------- IMPLCA2C.902
CL 5. "Downward sweep" through whole matrix. TL, QW and TSTAR are IMPLCA2C.903
CL updated when the final implicit increments have been calculated. IMPLCA2C.904
C----------------------------------------------------------------------- IMPLCA2C.905
CL 5.1 Remaining TL rows of matrix and add implicit increments above IMPLCA2C.906
CL the rml or at all levels if it is less than two layers thick. IMPLCA2C.907
C----------------------------------------------------------------------- IMPLCA2C.908
C IMPLCA2C.909
DO 51 K=BLM1,1,-1 IMPLCA2C.910
DO 511 I=P1,P1+P_POINTS-1 IMPLCA2C.911
IF ( (K .GT. NRML(I)) .OR. (NRML(I) .LT. 2) ) THEN IMPLCA2C.912
DTL_DV(I,K) = DTL_DV(I,K) - AT_ATQ(I,K)*DTL_DV(I,K+1) IMPLCA2C.913
! P244.118 IMPLCA2C.914
TL(I,K) = TL(I,K) + DTL_DV(I,K) ! P244.127 IMPLCA2C.915
ENDIF IMPLCA2C.916
511 CONTINUE IMPLCA2C.917
51 CONTINUE IMPLCA2C.918
DO 52 I=P1,P1+P_POINTS-1 IMPLCA2C.919
C IMPLCA2C.920
C----------------------------------------------------------------------- IMPLCA2C.921
CL 5.2 Surface temperature (TSTAR) row of matrix and rapidly mixing IMPLCA2C.922
CL layer increments. IMPLCA2C.923
C----------------------------------------------------------------------- IMPLCA2C.924
IF ( NRML(I) .GE. 2 ) THEN IMPLCA2C.925
NRMLP1 = NRML(I) + 1 IMPLCA2C.926
DTL_RML(I) = DTL_RML(I) - AT_RML(I) * DTL_DV(I,NRMLP1) IMPLCA2C.927
DTSTAR(I) = DTSTAR(I) - ATSTAR(I) * DTL_RML(I) IMPLCA2C.928
DQW_RML(I) = DQW_RML(I) IMPLCA2C.929
& - AQ_RML(I) * ALPHAS(I) * DTSTAR(I) IMPLCA2C.930
TL(I,1) = TL(I,1) + DTL_RML(I) IMPLCA2C.931
QW(I,1) = QW(I,1) + DQW_RML(I) IMPLCA2C.932
ELSE IMPLCA2C.933
DTSTAR(I) = DTSTAR(I) - ATSTAR(I)*DTL_DV(I,1) ! P244.119 IMPLCA2C.934
C IMPLCA2C.935
C----------------------------------------------------------------------- IMPLCA2C.936
CL 5.3 Lowest-level QW row of matrix; local mixing where there is no rml IMPLCA2C.937
C----------------------------------------------------------------------- IMPLCA2C.938
IMPLCA2C.939
DQW_DU(I,1) = DQW_DU(I,1) - AQ_AM(I,1)*ALPHAS(I)*DTSTAR(I) IMPLCA2C.940
C ... P244.120 IMPLCA2C.941
QW(I,1) = QW(I,1) + DQW_DU(I,1) ! P244.128 IMPLCA2C.942
DQW_1(I) = DQW_DU(I,1) IMPLCA2C.943
ENDIF IMPLCA2C.944
TSTAR(I) = TSTAR(I) + DTSTAR(I) IMPLCA2C.945
52 CONTINUE IMPLCA2C.946
C IMPLCA2C.947
C----------------------------------------------------------------------- IMPLCA2C.948
CL 5.4 Remaining QW rows of matrix + updating of QW's. IMPLCA2C.949
CL Add implicit increments due to mixing above rml or at all levels IMPLCA2C.950
CL if there it does not exist. IMPLCA2C.951
C----------------------------------------------------------------------- IMPLCA2C.952
DO 54 K=2,BL_LEVELS IMPLCA2C.953
DO 541 I=P1,P1+P_POINTS-1 IMPLCA2C.954
IF (K .GT. NRML(I)) THEN IMPLCA2C.955
IF ((NRML(I) .GE. 2) .AND. (K-1 .EQ. NRML(I))) THEN IMPLCA2C.956
DQW_DU(I,K) = DQW_DU(I,K) - AQ_AM(I,K)*DQW_RML(I) IMPLCA2C.957
ELSE IMPLCA2C.958
DQW_DU(I,K) = DQW_DU(I,K) - AQ_AM(I,K)*DQW_DU(I,K-1) IMPLCA2C.959
C ... P244.121 IMPLCA2C.960
ENDIF IMPLCA2C.961
QW(I,K) = QW(I,K) + DQW_DU(I,K) ! P244.128 IMPLCA2C.962
ELSE IMPLCA2C.963
C IMPLCA2C.964
C Add the increments due to rapid mixing if in the rapidly mixing layer IMPLCA2C.965
C IMPLCA2C.966
TL(I,K) = TL(I,K) + DTL_RML(I) IMPLCA2C.967
QW(I,K) = QW(I,K) + DQW_RML(I) IMPLCA2C.968
ENDIF IMPLCA2C.969
541 CONTINUE IMPLCA2C.970
54 CONTINUE IMPLCA2C.971
C IMPLCA2C.972
C----------------------------------------------------------------------- IMPLCA2C.973
CL 6. Calculate final implicit fluxes of heat and moisture. IMPLCA2C.974
C----------------------------------------------------------------------- IMPLCA2C.975
CL 6.1 Surface fluxes for the 3 surface types: land, sea-ice, ordinary IMPLCA2C.976
CL sea. Pass out the value of RHOKH(,1) in GAMMA_RHOKH_1 to satisfy IMPLCA2C.977
CL STASH GAMMA_RHOKH_RDZ will contain precisely that on output. IMPLCA2C.978
C----------------------------------------------------------------------- IMPLCA2C.979
C IMPLCA2C.980
DO 61 I=P1,P1+P_POINTS-1 IMPLCA2C.981
IF (LAND_MASK(I)) THEN IMPLCA2C.982
FTL_ICE(I) = 0.0 IMPLCA2C.983
H_SEA(I) = 0.0 IMPLCA2C.984
FQW_ICE(I) = 0.0 IMPLCA2C.985
E_SEA(I) = 0.0 IMPLCA2C.986
IF ( NRML(I) .GE. 2 ) THEN IMPLCA2C.987
FTL(I,1) = FTL(I,1) - GAMMA_RHOKH_1(I) * IMPLCA2C.988
& ( DTL_RML(I) - DTSTAR(I) ) IMPLCA2C.989
DELTDQ = DQW_RML(I) - ALPHAS(I) * DTSTAR(I) IMPLCA2C.990
ELSE IMPLCA2C.991
FTL(I,1) = FTL(I,1) - GAMMA_RHOKH_1(I) * IMPLCA2C.992
& ( DTL_DV(I,1) - DTSTAR(I) ) IMPLCA2C.993
C ... P244.11 (for FTL = H/Cp, rather than for H itself) IMPLCA2C.994
C IMPLCA2C.995
DELTDQ = DQW_DU(I,1) - ALPHAS(I)*DTSTAR(I) ! P244.126 IMPLCA2C.996
ENDIF IMPLCA2C.997
EA(I) = EA(I) - GAMMA_RHOKEA(I)*DELTDQ ! P244.122 IMPLCA2C.998
ES(I) = ES(I) - GAMMA_RHOKES(I)*DELTDQ ! P244.123 IMPLCA2C.999
ESL(I) = ESL(I) - GAMMA_RHOKESL(I)*DELTDQ ! P244.124 IMPLCA2C.1000
FQW(I,1) = EA(I) + ES(I) ! P244.125 IMPLCA2C.1001
ELSEIF (ICE_FRACT(I).GT.0.0) THEN IMPLCA2C.1002
IF ( NRML(I) .GE. 2 ) THEN IMPLCA2C.1003
H_SEA(I) = CP * ( FTL_LEAD(I) - GAMMA_RHOKH_1(I) IMPLCA2C.1004
& * ( 1.0 - ICE_FRACT(I) ) * DTL_RML(I) ) IMPLCA2C.1005
FTL_ICE(I) = FTL_ICE(I) - GAMMA_RHOKH_1(I) IMPLCA2C.1006
& * ( ICE_FRACT(I)*DTL_RML(I) - DTSTAR(I) ) IMPLCA2C.1007
E_SEA(I) = FQW_LEAD(I) - GAMMA_RHOKE(I) IMPLCA2C.1008
& * (1.0 - ICE_FRACT(I)) * DQW_RML(I) IMPLCA2C.1009
DELTDQ = ICE_FRACT(I)*DQW_RML(I) - ALPHAS(I)*DTSTAR(I) IMPLCA2C.1010
ELSE IMPLCA2C.1011
H_SEA(I) = CP * ( FTL_LEAD(I) - GAMMA_RHOKH_1(I) IMPLCA2C.1012
& * ( 1.0 - ICE_FRACT(I) ) * DTL_DV(I,1) ) IMPLCA2C.1013
FTL_ICE(I) = FTL_ICE(I) - GAMMA_RHOKH_1(I) IMPLCA2C.1014
& * ( ICE_FRACT(I)*DTL_DV(I,1) - DTSTAR(I) ) ! P2440.1 IMPLCA2C.1015
E_SEA(I) = FQW_LEAD(I) - GAMMA_RHOKE(I) IMPLCA2C.1016
& * (1.0 - ICE_FRACT(I)) * DQW_DU(I,1) ! P2440.2 IMPLCA2C.1017
DELTDQ = ICE_FRACT(I)*DQW_DU(I,1) - ALPHAS(I)*DTSTAR(I) IMPLCA2C.1018
ENDIF IMPLCA2C.1019
FQW_ICE(I) = FQW_ICE(I) - GAMMA_RHOKE(I) * DELTDQ IMPLCA2C.1020
FTL(I,1) = FTL_ICE(I) + H_SEA(I)/CP IMPLCA2C.1021
FQW(I,1) = FQW_ICE(I) + E_SEA(I) IMPLCA2C.1022
EA(I) = FQW(I,1) IMPLCA2C.1023
ES(I) = 0.0 IMPLCA2C.1024
ESL(I) = 0.0 IMPLCA2C.1025
ELSE ! ordinary sea point IMPLCA2C.1026
IF ( NRML(I) .GE. 2) THEN IMPLCA2C.1027
H_SEA(I) = CP*( FTL(I,1) - GAMMA_RHOKH_1(I)*DTL_RML(I) ) IMPLCA2C.1028
E_SEA(I) = FQW(I,1) - GAMMA_RHOKE(I)*DQW_RML(I) IMPLCA2C.1029
ELSE IMPLCA2C.1030
H_SEA(I) = CP * ( FTL(I,1) - GAMMA_RHOKH_1(I)*DTL_DV(I,1) ) IMPLCA2C.1031
E_SEA(I) = FQW(I,1) - GAMMA_RHOKE(I)*DQW_DU(I,1) IMPLCA2C.1032
ENDIF IMPLCA2C.1033
FTL_ICE(I) = 0.0 IMPLCA2C.1034
FQW_ICE(I) = 0.0 IMPLCA2C.1035
FTL(I,1) = H_SEA(I) / CP IMPLCA2C.1036
FQW(I,1) = E_SEA(I) IMPLCA2C.1037
EA(I) = FQW(I,1) IMPLCA2C.1038
ES(I) = 0.0 IMPLCA2C.1039
ESL(I) = 0.0 IMPLCA2C.1040
ENDIF IMPLCA2C.1041
GAMMA_RHOKH_1(I) = GAMMA_RHOKH_1(I) / GAMMA(1) IMPLCA2C.1042
61 CONTINUE IMPLCA2C.1043
C IMPLCA2C.1044
C----------------------------------------------------------------------- IMPLCA2C.1045
CL 6.2 Fluxes at layer interfaces above the surface. IMPLCA2C.1046
C----------------------------------------------------------------------- IMPLCA2C.1047
DO 62 K=2,BL_LEVELS IMPLCA2C.1048
KM1 = K-1 IMPLCA2C.1049
DO 621 I=P1,P1+P_POINTS-1 IMPLCA2C.1050
C IMPLCA2C.1051
C Calculate and store fluxes due to local mixing. IMPLCA2C.1052
C FTL(local mixing) stored in array AT, IMPLCA2C.1053
C FQW(local mixing) stored in array AQ_AM. IMPLCA2C.1054
C IMPLCA2C.1055
NRMLP1 = NRML(I) + 1 IMPLCA2C.1056
IF ((NRML(I) .GE. 2) .AND. (K .EQ. NRMLP1)) THEN IMPLCA2C.1057
AT_ATQ(I,K) = FTL(I,K) - GAMMA_RHOKH_RDZ(I,K) IMPLCA2C.1058
& * ( DTL_DV(I,K) - DTL_RML(I) ) IMPLCA2C.1059
AQ_AM(I,K) = FQW(I,K) - GAMMA_RHOKH_RDZ(I,K) IMPLCA2C.1060
& * ( DQW_DU(I,K) - DQW_RML(I) ) IMPLCA2C.1061
ELSE IMPLCA2C.1062
AT_ATQ(I,K) = FTL(I,K) - GAMMA_RHOKH_RDZ(I,K) IMPLCA2C.1063
& * ( DTL_DV(I,K) - DTL_DV(I,KM1) ) IMPLCA2C.1064
AQ_AM(I,K) = FQW(I,K) - GAMMA_RHOKH_RDZ(I,K) IMPLCA2C.1065
& * ( DQW_DU(I,K) - DQW_DU(I,KM1) ) IMPLCA2C.1066
ENDIF IMPLCA2C.1067
C IMPLCA2C.1068
C Now calculate the implicit fluxes including both local mixing and IMPLCA2C.1069
C if appropriate also the fluxes due to rapid mixing through layers. IMPLCA2C.1070
C IMPLCA2C.1071
IF ( K .EQ. 2 ) THEN IMPLCA2C.1072
IF ( NRML(I) .GE. 2 ) THEN IMPLCA2C.1073
FTL(I,K) = AT_ATQ(I,K) IMPLCA2C.1074
& + FTL(I,KM1) - DTL_RML(I) / DTRDZ(I,KM1) IMPLCA2C.1075
FQW(I,K) = AQ_AM(I,K) IMPLCA2C.1076
& + FQW(I,KM1) - DQW_RML(I) / DTRDZ(I,KM1) IMPLCA2C.1077
ELSE IMPLCA2C.1078
FTL(I,K) = AT_ATQ(I,K) IMPLCA2C.1079
FQW(I,K) = AQ_AM(I,K) IMPLCA2C.1080
ENDIF IMPLCA2C.1081
ELSEIF ( K .LE. NRML(I) ) THEN IMPLCA2C.1082
FTL(I,K) = AT_ATQ(I,K) - AT_ATQ(I,KM1) IMPLCA2C.1083
& + FTL(I,KM1) - DTL_RML(I) / DTRDZ(I,KM1) IMPLCA2C.1084
FQW(I,K) = AQ_AM(I,K) - AQ_AM(I,KM1) IMPLCA2C.1085
& + FQW(I,KM1) - DQW_RML(I) / DTRDZ(I,KM1) IMPLCA2C.1086
ELSE IMPLCA2C.1087
FTL(I,K) = AT_ATQ(I,K) IMPLCA2C.1088
FQW(I,K) = AQ_AM(I,K) IMPLCA2C.1089
ENDIF IMPLCA2C.1090
621 CONTINUE IMPLCA2C.1091
62 CONTINUE IMPLCA2C.1092
C IMPLCA2C.1093
CL---------------------------------------------------------------------- IMPLCA2C.1094
CL (B) Calculations on UV-grid. IMPLCA2C.1095
CL---------------------------------------------------------------------- IMPLCA2C.1096
C IMPLCA2C.1097
C----------------------------------------------------------------------- IMPLCA2C.1098
CL 7. Solve matrix equation P244.80 for implicit increments to U and V. IMPLCA2C.1099
C----------------------------------------------------------------------- IMPLCA2C.1100
C IMPLCA2C.1101
*IF -DEF,SCMA AJC1F405.172
C----------------------------------------------------------------------- IMPLCA2C.1103
C 7.0 Interpolate to UV-grid the pressure changes across the model IMPLCA2C.1104
C layers IMPLCA2C.1105
C----------------------------------------------------------------------- IMPLCA2C.1106
DO 70 K=1,BL_LEVELS IMPLCA2C.1107
CALL P_TO_UV
(DELTAP(P1,K),AQ_RML(U1+ROW_LENGTH), IMPLCA2C.1108
+ P_POINTS,P_POINTS-ROW_LENGTH,ROW_LENGTH,P_ROWS) IMPLCA2C.1109
DO 701 I=U1+ROW_LENGTH,U1+U_POINTS-ROW_LENGTH-1 IMPLCA2C.1110
DELTAP(I,K) = AQ_RML(I) IMPLCA2C.1111
701 CONTINUE IMPLCA2C.1112
70 CONTINUE IMPLCA2C.1113
C IMPLCA2C.1114
C----------------------------------------------------------------------- IMPLCA2C.1115
CL 7.1 Initial calculations and "upward sweep". IMPLCA2C.1116
CL (a) "Surface" fluxes. IMPLCA2C.1117
C----------------------------------------------------------------------- IMPLCA2C.1118
C IMPLCA2C.1119
DO 71 I=U1+ROW_LENGTH,U1+U_POINTS-ROW_LENGTH-1 IMPLCA2C.1120
*ELSE IMPLCA2C.1121
C IMPLCA2C.1122
C----------------------------------------------------------------------- IMPLCA2C.1123
CL 7.1 Initial calculations and "upward sweep". IMPLCA2C.1124
C----------------------------------------------------------------------- IMPLCA2C.1125
CL (a) "Surface" fluxes. IMPLCA2C.1126
C----------------------------------------------------------------------- IMPLCA2C.1127
DO 71 I=1,U_POINTS IMPLCA2C.1128
*ENDIF IMPLCA2C.1129
DTRDZ_UV(I,1) = -DTIMEG / DELTAP(I,1) IMPLCA2C.1130
DELTAP_RML(I) = 0.0 IMPLCA2C.1131
C IMPLCA2C.1132
C "Explicit" increments to U(1) and V(1) when there is no rapidly IMPLCA2C.1133
C mixing layer or it does not extend beyond the bottom model layer. IMPLCA2C.1134
C IMPLCA2C.1135
DQW_DU(I,1) = DTRDZ_UV(I,1) * ( TAUX(I,2) - TAUX(I,1) ) IMPLCA2C.1136
C ... P244.67 IMPLCA2C.1137
C IMPLCA2C.1138
DTL_DV(I,1) = DTRDZ_UV(I,1) * ( TAUY(I,2) - TAUY(I,1) ) IMPLCA2C.1139
C ... P244.67 IMPLCA2C.1140
C IMPLCA2C.1141
CM = -DTRDZ_UV(I,1) * GAMMA_RHOKM_1(I) ! P244.66 IMPLCA2C.1142
AQ_AM(I,1) = -DTRDZ_UV(I,1) * GAMMA_RHOKM_RDZUV(I,2) ! P244.64 IMPLCA2C.1143
RBM = 1.0 / ( 1.0 - AQ_AM(I,1) - CM ) ! Reciprocal of P244.81 IMPLCA2C.1144
DQW_DU(I,1) = RBM * DQW_DU(I,1) ! P244.82 IMPLCA2C.1145
DTL_DV(I,1) = RBM * DTL_DV(I,1) ! P244.83 IMPLCA2C.1146
AQ_AM(I,1) = RBM * AQ_AM(I,1) ! P244.84 IMPLCA2C.1147
71 CONTINUE IMPLCA2C.1148
C IMPLCA2C.1149
C IMPLCA2C.1150
C----------------------------------------------------------------------- IMPLCA2C.1151
CL (b) Fluxes at (or rows representing) layer interfaces above the IMPLCA2C.1152
CL surface but below the top of the boundary layer. IMPLCA2C.1153
C----------------------------------------------------------------------- IMPLCA2C.1154
DO 72 K=2,BLM1 IMPLCA2C.1155
KP1 = K+1 IMPLCA2C.1156
KM1 = K-1 IMPLCA2C.1157
*IF -DEF,SCMA AJC1F405.173
DO 721 I=U1+ROW_LENGTH,U1+U_POINTS-ROW_LENGTH-1 IMPLCA2C.1159
*ELSE IMPLCA2C.1160
DO 721 I=1,U_POINTS IMPLCA2C.1161
*ENDIF IMPLCA2C.1162
DTRDZ_UV(I,K) = -DTIMEG / DELTAP(I,K) IMPLCA2C.1163
DQW_DU(I,K) = DTRDZ_UV(I,K) * ( TAUX(I,KP1) - TAUX(I,K) ) IMPLCA2C.1164
C ... P244.74 IMPLCA2C.1165
DTL_DV(I,K) = DTRDZ_UV(I,K) * ( TAUY(I,KP1) - TAUY(I,K) ) IMPLCA2C.1166
C ... P244.74 IMPLCA2C.1167
AQ_AM(I,K) = -DTRDZ_UV(I,K) * GAMMA_RHOKM_RDZUV(I,KP1) IMPLCA2C.1168
C ... P244.71 IMPLCA2C.1169
CM = -DTRDZ_UV(I,K) * GAMMA_RHOKM_RDZUV(I,K) ! P244.72 IMPLCA2C.1170
RBM = 1.0 / ( 1.0 - AQ_AM(I,K) -CM*( 1.0 + AQ_AM(I,KM1) ) ) IMPLCA2C.1171
C 1 2 2 1 IMPLCA2C.1172
C ... Reciprocal of P244.85 IMPLCA2C.1173
C IMPLCA2C.1174
DQW_DU(I,K) = RBM * ( DQW_DU(I,K) - CM*DQW_DU(I,KM1) ) IMPLCA2C.1175
C ... P244.86 IMPLCA2C.1176
C IMPLCA2C.1177
DTL_DV(I,K) = RBM * ( DTL_DV(I,K) - CM*DTL_DV(I,KM1) ) IMPLCA2C.1178
C ... P244.87 IMPLCA2C.1179
C IMPLCA2C.1180
AQ_AM(I,K) = RBM * AQ_AM(I,K) ! P244.88 IMPLCA2C.1181
721 CONTINUE IMPLCA2C.1182
72 CONTINUE IMPLCA2C.1183
C----------------------------------------------------------------------- IMPLCA2C.1184
CL (c) Top "boundary" layer; also increment U and V here, as implicit IMPLCA2C.1185
CL flux for this layer is got from "upward sweep" alone. IMPLCA2C.1186
C----------------------------------------------------------------------- IMPLCA2C.1187
*IF -DEF,SCMA AJC1F405.174
DO 73 I=U1+ROW_LENGTH,U1+U_POINTS-ROW_LENGTH-1 IMPLCA2C.1189
*ELSE IMPLCA2C.1190
DO 73 I=1,U_POINTS IMPLCA2C.1191
*ENDIF IMPLCA2C.1192
DTRDZ_UV(I,BL_LEVELS) = -DTIMEG / DELTAP(I,BL_LEVELS) IMPLCA2C.1193
DQW_DU(I,BL_LEVELS) = -DTRDZ_UV(I,BL_LEVELS) * TAUX(I,BL_LEVELS) IMPLCA2C.1194
C ... P244.78 IMPLCA2C.1195
C IMPLCA2C.1196
DTL_DV(I,BL_LEVELS) = -DTRDZ_UV(I,BL_LEVELS) * TAUY(I,BL_LEVELS) IMPLCA2C.1197
C ... P244.78 IMPLCA2C.1198
C IMPLCA2C.1199
CM = -DTRDZ_UV(I,BL_LEVELS) * GAMMA_RHOKM_RDZUV(I,BL_LEVELS) IMPLCA2C.1200
C ... P244.76 IMPLCA2C.1201
RBM = 1.0 / ( 1.0 - CM*( 1.0 + AQ_AM(I,BLM1) ) ) IMPLCA2C.1202
C 1 2 2 1 IMPLCA2C.1203
C ... Reciprocal of P244.89 IMPLCA2C.1204
C IMPLCA2C.1205
DQW_DU(I,BL_LEVELS) = RBM * ( DQW_DU(I,BL_LEVELS) ! P244.90 IMPLCA2C.1206
+ - CM*DQW_DU(I,BLM1) ) IMPLCA2C.1207
DTL_DV(I,BL_LEVELS) = RBM * ( DTL_DV(I,BL_LEVELS) ! P244.91 IMPLCA2C.1208
+ - CM*DTL_DV(I,BLM1) ) IMPLCA2C.1209
U(I,BL_LEVELS) = U(I,BL_LEVELS) + DQW_DU(I,BL_LEVELS) ! P244.94 IMPLCA2C.1210
V(I,BL_LEVELS) = V(I,BL_LEVELS) + DTL_DV(I,BL_LEVELS) ! P244.95 IMPLCA2C.1211
73 CONTINUE IMPLCA2C.1212
IMPLCA2C.1213
C----------------------------------------------------------------------- IMPLCA2C.1214
CL 7.4 Complete solution of matrix equation by performing "downward IMPLCA2C.1215
CL sweep", then update U and V. IMPLCA2C.1216
C----------------------------------------------------------------------- IMPLCA2C.1217
DO 74 K=BLM1,1,-1 IMPLCA2C.1218
KP1 = K+1 IMPLCA2C.1219
*IF -DEF,SCMA AJC1F405.175
DO 741 I=U1+ROW_LENGTH,U1+U_POINTS-ROW_LENGTH-1 IMPLCA2C.1221
*ELSE IMPLCA2C.1222
DO 741 I=1,U_POINTS IMPLCA2C.1223
*ENDIF IMPLCA2C.1224
DQW_DU(I,K) = DQW_DU(I,K) - AQ_AM(I,K)*DQW_DU(I,KP1) ! P244.92 IMPLCA2C.1225
DTL_DV(I,K) = DTL_DV(I,K) - AQ_AM(I,K)*DTL_DV(I,KP1) ! P244.93 IMPLCA2C.1226
U(I,K) = U(I,K) + DQW_DU(I,K) ! P244.94 IMPLCA2C.1227
V(I,K) = V(I,K) + DTL_DV(I,K) ! P244.95 IMPLCA2C.1228
741 CONTINUE IMPLCA2C.1229
74 CONTINUE IMPLCA2C.1230
C IMPLCA2C.1231
C----------------------------------------------------------------------- IMPLCA2C.1232
CL 8. Essentially diagnostic calculations, though some values (i.e. the IMPLCA2C.1233
CL surface wind stresses) are required by the coupled version of the IMPLCA2C.1234
CL model. IMPLCA2C.1235
C----------------------------------------------------------------------- IMPLCA2C.1236
CL 8.1 Surface wind stress components. IMPLCA2C.1237
CL Pass out the value of RHOKM(,1) in GAMMA_RHOKM_1 to satisfy STASH IMPLCA2C.1238
CL GAMMA_RHOKM_RDZ will contain precisely that on output. IMPLCA2C.1239
C----------------------------------------------------------------------- IMPLCA2C.1240
C IMPLCA2C.1241
*IF -DEF,SCMA AJC1F405.176
DO 81 I=U1+ROW_LENGTH,U1+U_POINTS-ROW_LENGTH-1 IMPLCA2C.1243
*ELSE IMPLCA2C.1244
DO 81 I=1,U_POINTS IMPLCA2C.1245
*ENDIF IMPLCA2C.1246
TAUX(I,1) = TAUX(I,1) + GAMMA_RHOKM_1(I)*DQW_DU(I,1) IMPLCA2C.1247
C ... P244.61 IMPLCA2C.1248
C IMPLCA2C.1249
TAUY(I,1) = TAUY(I,1) + GAMMA_RHOKM_1(I)*DTL_DV(I,1) IMPLCA2C.1250
C ... P244.61 IMPLCA2C.1251
C IMPLCA2C.1252
GAMMA_RHOKM_1(I) = GAMMA_RHOKM_1(I) / GAMMA(1) IMPLCA2C.1253
81 CONTINUE IMPLCA2C.1254
*IF -DEF,SCMA AJC1F405.177
C IMPLCA2C.1256
C Set extreme rows to "missing data indicator". IMPLCA2C.1257
C IMPLCA2C.1258
*IF DEF,MPP GPB1F403.95
IF (attop) THEN GPB1F403.96
*ENDIF GPB1F403.97
DO I=U1,U1+ROW_LENGTH-1 GPB1F403.98
TAUX(I,1)=1.E30 GPB1F403.99
TAUY(I,1)=1.E30 GPB1F403.100
ENDDO GPB1F403.101
*IF DEF,MPP GPB1F403.102
ENDIF GPB1F403.103
GPB1F403.104
IF (atbase) THEN GPB1F403.105
*ENDIF GPB1F403.106
DO I= U1 + (U_ROWS-1)*ROW_LENGTH , U1 + U_ROWS*ROW_LENGTH - 1 GPB1F403.107
TAUX(I,1)=1.E30 GPB1F403.108
TAUY(I,1)=1.E30 GPB1F403.109
ENDDO GPB1F403.110
*IF DEF,MPP GPB1F403.111
ENDIF GPB1F403.112
*ENDIF GPB1F403.113
*ENDIF IMPLCA2C.1283
C IMPLCA2C.1284
C----------------------------------------------------------------------- IMPLCA2C.1285
CL 8.2 Wind stress components at layer interfaces above the surface. IMPLCA2C.1286
C----------------------------------------------------------------------- IMPLCA2C.1287
DO 82 K=2,BL_LEVELS IMPLCA2C.1288
KM1 = K-1 IMPLCA2C.1289
*IF -DEF,SCMA AJC1F405.178
DO 821 I=U1+ROW_LENGTH,U1+U_POINTS-ROW_LENGTH-1 IMPLCA2C.1291
*ELSE IMPLCA2C.1292
DO 821 I=1,U_POINTS IMPLCA2C.1293
*ENDIF IMPLCA2C.1294
AQ_AM(I,K) = TAUX(I,K) + GAMMA_RHOKM_RDZUV(I,K) ! P244.61 IMPLCA2C.1295
+ *( DQW_DU(I,K) - DQW_DU(I,KM1) ) IMPLCA2C.1296
AT_ATQ(I,K) = TAUY(I,K) + GAMMA_RHOKM_RDZUV(I,K) ! P244.61 IMPLCA2C.1297
+ *( DTL_DV(I,K) - DTL_DV(I,KM1) ) IMPLCA2C.1298
TAUX(I,K) = AQ_AM(I,K) IMPLCA2C.1299
TAUY(I,K) = AT_ATQ(I,K) IMPLCA2C.1300
821 CONTINUE IMPLCA2C.1301
*IF -DEF,SCMA AJC1F405.179
C IMPLCA2C.1303
C Set extreme rows to "missing data indicator". IMPLCA2C.1304
C IMPLCA2C.1305
*IF DEF,MPP GPB1F403.114
IF (attop) THEN GPB1F403.115
*ENDIF GPB1F403.116
DO I=U1,U1+ROW_LENGTH-1 GPB1F403.117
TAUX(I,K)=1.E30 GPB1F403.118
TAUY(I,K)=1.E30 GPB1F403.119
ENDDO GPB1F403.120
*IF DEF,MPP GPB1F403.121
ENDIF GPB1F403.122
GPB1F403.123
IF (atbase) THEN GPB1F403.124
*ENDIF GPB1F403.125
DO I= U1 + (U_ROWS-1)*ROW_LENGTH , U1 + U_ROWS*ROW_LENGTH - 1 GPB1F403.126
TAUX(I,K)=1.E30 GPB1F403.127
TAUY(I,K)=1.E30 GPB1F403.128
ENDDO GPB1F403.129
*IF DEF,MPP GPB1F403.130
ENDIF GPB1F403.131
*ENDIF GPB1F403.132
*ENDIF IMPLCA2C.1313
82 CONTINUE IMPLCA2C.1314
C IMPLCA2C.1315
C----------------------------------------------------------------------- IMPLCA2C.1316
CL 8.3 Wind components at 10 metres above the surface. IMPLCA2C.1317
C----------------------------------------------------------------------- IMPLCA2C.1318
IF (SU10 .OR. SV10) THEN IMPLCA2C.1319
*IF -DEF,SCMA AJC1F405.180
CMIC$ DO ALL VECTOR SHARED(U_POINTS, ROW_LENGTH, U1, SU10, SV10, U_FIELD IMPLCA2C.1321
CMIC$1 , U0, CDR10M, BL_LEVELS, U, U10M, V0, V, V10M) PRIVATE(I) IMPLCA2C.1322
CDIR$ IVDEP IMPLCA2C.1323
! Fujitsu vectorization directive GRB0F405.361
!OCL NOVREC GRB0F405.362
DO 83 I=U1+ROW_LENGTH,U1+U_POINTS-ROW_LENGTH-1 IMPLCA2C.1324
*ELSE IMPLCA2C.1325
DO 83 I=1,U_POINTS IMPLCA2C.1326
*ENDIF IMPLCA2C.1327
IF (SU10) IMPLCA2C.1328
+ U10M(I) = U0(I) + CDR10M(I)*( U(I,1) - U0(I) ) ! P244.96 IMPLCA2C.1329
IF (SV10) IMPLCA2C.1330
+ V10M(I) = V0(I) + CDR10M(I)*( V(I,1) - V0(I) ) ! P244.97 IMPLCA2C.1331
83 CONTINUE IMPLCA2C.1332
*IF -DEF,SCMA AJC1F405.181
C IMPLCA2C.1334
C Set extreme rows to "missing data indicator". IMPLCA2C.1335
C IMPLCA2C.1336
*IF DEF,MPP GPB1F403.133
IF (attop) THEN GPB1F403.134
*ENDIF GPB1F403.135
DO I=U1,U1+ROW_LENGTH-1 GPB1F403.136
IF (SU10) THEN GPB1F403.137
U10M(I)=1.E30 GPB1F403.138
ENDIF GPB1F403.139
IF (SV10) THEN GPB1F403.140
V10M(I)=1.E30 GPB1F403.141
ENDIF GPB1F403.142
ENDDO GPB1F403.143
*IF DEF,MPP GPB1F403.144
ENDIF GPB1F403.145
GPB1F403.146
IF (atbase) THEN GPB1F403.147
*ENDIF GPB1F403.148
DO I= U1 + (U_ROWS-1)*ROW_LENGTH , U1 + U_ROWS*ROW_LENGTH - 1 GPB1F403.149
IF (SU10) THEN GPB1F403.150
U10M(I)=1.E30 GPB1F403.151
ENDIF GPB1F403.152
IF (SV10) THEN GPB1F403.153
V10M(I)=1.E30 GPB1F403.154
ENDIF GPB1F403.155
ENDDO GPB1F403.156
*IF DEF,MPP GPB1F403.157
ENDIF GPB1F403.158
*ENDIF GPB1F403.159
*ENDIF IMPLCA2C.1369
ENDIF IMPLCA2C.1370
999 CONTINUE ! Branch for error exit. IMPLCA2C.1371
IF (LTIMER) THEN ASJ1F304.368
CALL TIMER
('IMPLCAL ',4) IMPLCA2C.1373
ENDIF ASJ1F304.369
RETURN IMPLCA2C.1375
END IMPLCA2C.1376
*ENDIF IMPLCA2C.1377