*IF DEF,OCEAN @DYALLOC.4679
C ******************************COPYRIGHT****************************** GTS2F400.11503
C (c) CROWN COPYRIGHT 1995, METEOROLOGICAL OFFICE, All Rights Reserved. GTS2F400.11504
C GTS2F400.11505
C Use, duplication or disclosure of this code is subject to the GTS2F400.11506
C restrictions as set forth in the contract. GTS2F400.11507
C GTS2F400.11508
C Meteorological Office GTS2F400.11509
C London Road GTS2F400.11510
C BRACKNELL GTS2F400.11511
C Berkshire UK GTS2F400.11512
C RG12 2SZ GTS2F400.11513
C GTS2F400.11514
C If no contract has been raised with this copy of the code, the use, GTS2F400.11515
C duplication or disclosure of it is strictly prohibited. Permission GTS2F400.11516
C to do so must first be obtained in writing from the Head of Numerical GTS2F400.11517
C Modelling at the above address. GTS2F400.11518
C ******************************COPYRIGHT****************************** GTS2F400.11519
C GTS2F400.11520
C*LL Subroutine VDIFCALC VDIFCALC.3
CLL VDIFCALC.4
CLL Can run on any compiler accepting lower case variables. VDIFCALC.5
CLL VDIFCALC.6
CLL The code must be precompiled by the UPDOC system. VDIFCALC.7
CLL Option L selects the code for Isopycnal diffusion. VDIFCALC.8
CLL VDIFCALC.9
CLL Author: D.J.Carrington VDIFCALC.10
CLL Date: 12 December 1990 VDIFCALC.11
CLL Reviewer: R.A.Wood VDIFCALC.12
CLL Date: 19 December 1990 VDIFCALC.13
CLL Version 1.0 VDIFCALC.14
CLL VDIFCALC.15
CLL External documentation: VDIFCALC.16
CLL Unified Model Documentation Paper No. 51. VDIFCALC.17
CLL VDIFCALC.18
CLL Naming convention of variables: Cox naming convention is used, VDIFCALC.19
CLL with the addition that lower-case variables are VDIFCALC.20
CLL introduced by the Isopycnal Diffusion Scheme. VDIFCALC.21
CLL VDIFCALC.22
CLL Purpose of Subroutine: VDIFCALC.23
CLL Solves the vertical diffusion equation for momentum. VDIFCALC.24
CLL This separate treatment of the vertical diffusion allows a VDIFCALC.25
CLL greater time-step than would otherwise be possible. VDIFCALC.26
CLL Note that the surface fluxes of momentum (wind stresses) ARE VDIFCALC.27
CLL introduced in this routine. VDIFCALC.28
CLL VDIFCALC.29
CLL List of subroutines required for implementation of Isopycnal VDIFCALC.30
CLL Diffusion Scheme (in order of being called): VDIFCALC.31
CLL VERTCOFC * VDIFCALC.32
CLL VDIFCALC VDIFCALC.33
CLL VERTCOFT * VDIFCALC.34
CLL IPDCOFCL VDIFCALC.35
CLL IPDFLXCL VDIFCALC.36
CLL VDIFCALT * K-theory mixing scheme VDIFCALC.37
CLL VDIFCALC.38
CLLEND------------------------------------------------------------------ VDIFCALC.39
C* VDIFCALC.40
C*L---- Arguments ------------------------------------------------------ VDIFCALC.41
C VDIFCALC.42
SUBROUTINE VDIFCALC 2VDIFCALC.43
& ( J,IMT,IMTM1,KM,KMP1,KMM1,NT, VDIFCALC.44
& UA,UB,VA,VB, VDIFCALC.45
& DZ,DZZ2RQ,DZ2RQ,C2DTUV, VDIFCALC.46
& WSX,WSY,GM,gnu VDIFCALC.47
& ) VDIFCALC.48
C VDIFCALC.49
C VDIFCALC.50
IMPLICIT NONE VDIFCALC.51
C VDIFCALC.52
INTEGER VDIFCALC.53
& I,J,K,M,IMT,IMTM1,KM,KMP1,KMM1,NT VDIFCALC.54
C VDIFCALC.55
REAL VDIFCALC.56
& UA(IMT,KM),VA(IMT,KM),UB(IMT,KM),VB(IMT,KM), VDIFCALC.57
& DZ(KM),DZ2RQ(IMT,KM),DZZ2RQ(IMT,KM),C2DTUV, VDIFCALC.58
& GM(IMT,KM),WSX(IMT),WSY(IMT), VDIFCALC.59
& gnu(IMT,KM) ! IN Vert. diffusivity of momentum at TOP of VDIFCALC.60
C box (calculated by VERTCOFC) (cm2/s) VDIFCALC.61
C VDIFCALC.62
C VDIFCALC.63
C Decalre local variables and arrays VDIFCALC.64
C VDIFCALC.65
REAL VDIFCALC.66
& aa(IMT,KM), ! c(1)A(k) } Eq.3.2 VDIFCALC.67
& cc(IMT,KM), ! c(2)A(k+1) } VDIFCALC.68
& bb(IMT,KM), ! aa+cc+1 VDIFCALC.69
& ee(IMT,KMP1), ! x(k-1) } Eq.1.18 VDIFCALC.70
& ff(IMT,KMP1,NT), ! y(k-1) } VDIFCALC.71
& efdr(IMT), ! 1/(bb-cc*ee) VDIFCALC.72
& tempa(IMT,KMP1), ! workspace VDIFCALC.73
& tempb(IMT,KMP1), ! workspace VDIFCALC.74
& fxa, ! } VDIFCALC.75
& fxb, ! } local csts. VDIFCALC.76
& fxc, ! } VDIFCALC.77
& conv ! SI-cgs conversion factor for windstress VDIFCALC.78
C* VDIFCALC.79
C -------------------------------------------------------------- VDIFCALC.80
CL The arrays are set up here, VDIFCALC.81
CL in preparation for the recurrence method that follows. VDIFCALC.82
C -------------------------------------------------------------- VDIFCALC.83
C VDIFCALC.84
fxa=4. VDIFCALC.85
fxb=1. VDIFCALC.86
fxc=0. VDIFCALC.87
C VDIFCALC.88
CL Certain quantities for the surface & bottom are set to zero. VDIFCALC.89
C VDIFCALC.90
DO 210 I=1,IMT VDIFCALC.91
aa(I, 1)=fxc VDIFCALC.92
cc(I, KM)=fxc VDIFCALC.93
ee(I,KM+1)=fxc VDIFCALC.94
210 CONTINUE VDIFCALC.95
DO 220 M=1,2 VDIFCALC.96
DO 220 I=1,IMT VDIFCALC.97
ff(I,KM+1,M)=fxc VDIFCALC.98
220 CONTINUE VDIFCALC.99
C VDIFCALC.100
DO 230 K=1,KM VDIFCALC.101
DO 230 I=1,IMT VDIFCALC.102
tempa(I,K)=C2DTUV*DZ2RQ(I,K)*fxa*GM(I,K) VDIFCALC.103
tempb(I,K)=gnu(I,K)*DZZ2RQ(I,K) VDIFCALC.104
230 CONTINUE VDIFCALC.105
DO 240 K=2,KM VDIFCALC.106
DO 240 I=1,IMT VDIFCALC.107
aa(I,K)=tempa(I,K)*tempb(I,K) VDIFCALC.108
240 CONTINUE VDIFCALC.109
DO 250 K=1,KMM1 VDIFCALC.110
DO 250 I=1,IMT VDIFCALC.111
cc(I,K)=tempa(I,K)*tempb(I,K+1)*GM(I,K+1) VDIFCALC.112
250 CONTINUE VDIFCALC.113
DO 260 K=1,KM VDIFCALC.114
DO 260 I=1,IMT VDIFCALC.115
bb(I,K)=aa(I,K)+cc(I,K)+fxb VDIFCALC.116
260 CONTINUE VDIFCALC.117
C VDIFCALC.118
C -------------------------------------------------------------- VDIFCALC.119
CL The terms x(k-1) (contained in ee) AND y(k-1) (contained in ff) VDIFCALC.120
CL in Equation (1.18) - and hence in Equation (3.2) - are calculated. VDIFCALC.121
C -------------------------------------------------------------- VDIFCALC.122
C VDIFCALC.123
DO 300 K=KM,1,-1 VDIFCALC.124
DO 270 I=2,IMTM1 VDIFCALC.125
efdr(I)=fxb/(bb(I,K)-cc(I,K)*ee(I,K+1)) VDIFCALC.126
270 CONTINUE VDIFCALC.127
DO 280 I=2,IMTM1 VDIFCALC.128
ee(I,K)=aa(I,K)*efdr(I) VDIFCALC.129
280 CONTINUE VDIFCALC.130
DO 290 I=2,IMTM1 VDIFCALC.131
C VDIFCALC.132
C The velocities at the start of the timestep are computed from VDIFCALC.133
C UB (velocity) and UA (acceleration). VDIFCALC.134
C VDIFCALC.135
ff(I,K,1)=(UB(I,K)+UA(I,K)*C2DTUV+cc(I,K)*ff(I,K+1,1))*efdr(I) VDIFCALC.136
ff(I,K,2)=(VB(I,K)+VA(I,K)*C2DTUV+cc(I,K)*ff(I,K+1,2))*efdr(I) VDIFCALC.137
290 CONTINUE VDIFCALC.138
300 CONTINUE VDIFCALC.139
C VDIFCALC.140
C -------------------------------------------------------------- VDIFCALC.141
CL The surface fluxes of momentum are introduced here, and the new VDIFCALC.142
CL velocities for the top model level are calculated using VDIFCALC.143
CL Equation (3.3). This is done by assuming zero diffusivity VDIFCALC.144
CL at the surface, but adding in the momentum flux explicitly VDIFCALC.145
CL to UB (U(k,n-1) in eqs. (3.2), (3.3)). Note that rho0=1.0 is VDIFCALC.146
CL assumed in this code (OK in cgs units!). VDIFCALC.147
CL VDIFCALC.148
CL Vble. 'conv' converts the wind stress from N/m2 to dyn/cm2. VDIFCALC.149
CL Variable bb is re-used, this time with dimensions of velocity. VDIFCALC.150
C -------------------------------------------------------------- VDIFCALC.151
C VDIFCALC.152
conv=10.0 VDIFCALC.153
DO 360 M=1,2 VDIFCALC.154
IF(M.EQ.1) THEN VDIFCALC.155
DO 310 I=2,IMTM1 VDIFCALC.156
bb(I,1)=(WSX(I)*conv*C2DTUV/DZ(1)+UB(I,1)+UA(I,1)*C2DTUV VDIFCALC.157
* +cc(I,1)*ff(I,2,1))*efdr(I) VDIFCALC.158
310 CONTINUE VDIFCALC.159
ELSE VDIFCALC.160
DO 320 I=2,IMTM1 VDIFCALC.161
bb(I,1)=(WSY(I)*conv*C2DTUV/DZ(1)+VB(I,1)+VA(I,1)*C2DTUV VDIFCALC.162
* +cc(I,1)*ff(I,2,2))*efdr(I) VDIFCALC.163
320 CONTINUE VDIFCALC.164
ENDIF VDIFCALC.165
C VDIFCALC.166
C -------------------------------------------------------------- VDIFCALC.167
CL The new velocities are calculated for the rest of the model VDIFCALC.168
CL levels, using Equation (3.2). VDIFCALC.169
C -------------------------------------------------------------- VDIFCALC.170
C VDIFCALC.171
DO 330 K=2,KM VDIFCALC.172
DO 330 I=2,IMTM1 VDIFCALC.173
bb(I,K)=ee(I,K)*bb(I,K-1)+ff(I,K,M) VDIFCALC.174
330 CONTINUE VDIFCALC.175
C VDIFCALC.176
C -------------------------------------------------------------- VDIFCALC.177
CL The acceleration terms are updated using these velocities. VDIFCALC.178
C -------------------------------------------------------------- VDIFCALC.179
C VDIFCALC.180
IF(M.EQ.1) THEN VDIFCALC.181
DO 340 K=1,KM VDIFCALC.182
DO 340 I=1,IMT VDIFCALC.183
UA(I,K)=(bb(I,K)*GM(I,K)-UB(I,K))/C2DTUV VDIFCALC.184
340 CONTINUE VDIFCALC.185
ELSE VDIFCALC.186
DO 350 K=1,KM VDIFCALC.187
DO 350 I=1,IMT VDIFCALC.188
VA(I,K)=(bb(I,K)*GM(I,K)-VB(I,K))/C2DTUV VDIFCALC.189
350 CONTINUE VDIFCALC.190
ENDIF VDIFCALC.191
360 CONTINUE VDIFCALC.192
C VDIFCALC.193
RETURN VDIFCALC.194
END VDIFCALC.195
*ENDIF @DYALLOC.4680