*IF DEF,OCEAN CALCDIFF.2
C ******************************COPYRIGHT****************************** CALCDIFF.3
C (c) CROWN COPYRIGHT 1998, METEOROLOGICAL OFFICE, All Rights Reserved. CALCDIFF.4
C CALCDIFF.5
C Use, duplication or disclosure of this code is subject to the CALCDIFF.6
C restrictions as set forth in the contract. CALCDIFF.7
C CALCDIFF.8
C Meteorological Office CALCDIFF.9
C London Road CALCDIFF.10
C BRACKNELL CALCDIFF.11
C Berkshire UK CALCDIFF.12
C RG12 2SZ CALCDIFF.13
C CALCDIFF.14
C If no contract has been raised with this copy of the code, the use, CALCDIFF.15
C duplication or disclosure of it is strictly prohibited. Permission CALCDIFF.16
C to do so must first be obtained in writing from the Head of Numerical CALCDIFF.17
C Modelling at the above address. CALCDIFF.18
C ******************************COPYRIGHT****************************** CALCDIFF.19
C CALCDIFF.20
CLL Subroutine CALCDIFF CALCDIFF.21
CLL CALCDIFF.22
CLL Author: M J Roberts CALCDIFF.23
CLL CALCDIFF.24
CLL Description: Mirrors subroutine CALCESAV for the Griffies CALCDIFF.25
CLL isopycnal diffusion scheme - calls isopycnal CALCDIFF.26
CLL diffusion routines to initialise variables for the CALCDIFF.27
CLL first row of a block. CALCDIFF.28
CLL CALCDIFF.29
CLL External documentation: CALCDIFF.30
CLL CALCDIFF.31
CLL Implemented at UM vn 4.5 26 June 1998 CALCDIFF.32
CLL CALCDIFF.33
CLL Modification History : CALCDIFF.34
CLL Version Date Comment & Name CALCDIFF.35
CLL ------- -------- -------------------------------------------- CALCDIFF.36
CLL CALCDIFF.37
CLL Subroutine dependencies. CALCDIFF.38
CLL Called by BLOKINIT CALCDIFF.39
CLL Calls subroutines VERTCOFT, ISOPYC_M, ISOFLUX, STATED, STATEC, CALCDIFF.40
CLL CALC_MLD_LARGEG or CALC_MLD_LARGEB. CALCDIFF.41
CLL CALCDIFF.42
CLLEND------------------------------------------------------------------ CALCDIFF.43
C CALCDIFF.44
SUBROUTINE CALCDIFF( 1,9CALCDIFF.45
*CALL ARGSIZE 
CALCDIFF.46
*CALL COCAWRKA CALCDIFF.47
& ,j,cstr,dyur,dxur,dz2r,dzz,dzz2r,athkdftu,athkdftv,ahi,athkdf, CALCDIFF.48
& dz,dyu,dxu,cs,dxt4r,dyt4r,dxtr,dytr,cst,csjm,dyurjm,j_1, CALCDIFF.49
& KMP,NERGY,FKAPB_SI,CSR,ITT, CALCDIFF.50
& KAPPA_B_SI,J_OFFSET, CALCDIFF.51
& WSXM,WSYM,OCEANHEATFLUX,CARYHEAT,FLXTOICE, CALCDIFF.52
& max_Large_levels,no_layers_in_lev, CALCDIFF.53
& waterflux_ice,L_OBULKRI,L_OWINDMIX,L_OBULKMAXMLD, CALCDIFF.54
& lambda_Large,phit, CALCDIFF.55
& ISX,ISY,WSX_LEADS,WSY_LEADS, CALCDIFF.56
& ISXM,ISYM,WSX_LEADSM,WSY_LEADSM, CALCDIFF.57
& ZDZZ,ZDZ, CALCDIFF.58
& adv_vetiso, CALCDIFF.59
& adv_vbtiso,adv_fbiso) CALCDIFF.60
C CALCDIFF.61
IMPLICIT NONE CALCDIFF.62
C CALCDIFF.63
*CALL OARRYSIZ CALCDIFF.64
*CALL TYPSIZE CALCDIFF.65
*CALL CNTLOCN CALCDIFF.66
*CALL COCTWRKA CALCDIFF.67
*CALL UMSCALAR CALCDIFF.68
C CALCDIFF.69
INTEGER CALCDIFF.70
& I,J,K,M,N,L, CALCDIFF.71
& KMP, CALCDIFF.72
& NERGY,ITT,KS,J_OFFSET,j_1 CALCDIFF.73
C CALCDIFF.74
REAL CALCDIFF.75
& CSR(JMT),CSTR(JMT),DXUR(IMT),DYUR(JMT),DZ2R(km),dzz(km+1), CALCDIFF.76
& dzz2r(km+1), CALCDIFF.77
& csjm,dyurjm, ! IN cs and dyur at j-1 CALCDIFF.78
& athkdftu(imt_vis,jmt_vis), ! IN thickness diffusion coeff CALCDIFF.79
& athkdftv(imt_vis,jmt_vis), ! IN thickness diffusion coeff CALCDIFF.80
c ! where u and v points are on a C-grid CALCDIFF.81
& ahi(km), ! isopycnal diffusion CALCDIFF.82
& athkdf(km), ! thickness diffusion CALCDIFF.83
CALCDIFF.84
& dz(km),dyu(jmt),dxu(imt),cs(jmt),DXT4R(IMT),DYT4R(JMT), CALCDIFF.85
& dxtr(imt),dytr(jmt),cst(jmt),ZDZ(KM),ZDZZ(KM+1), CALCDIFF.86
& gnu(IMT,KM), ! IN coeff. of vertical diffusion at TOP of CALCDIFF.87
C box. Av in Eq.1.10 (cm2/s) CALCDIFF.88
& rhosrn(IMT,KM),RHOSRNA(IMT,KM+1),RHOSRNB(IMT,KM+1), ! CALCDIFF.89
& KAPPA_B_SI(KM),CGFLUX(IMT,KM,NT) CALCDIFF.90
&, mld_Large(imt) ! mixed layer depth on T grid, row j, (cm) CALCDIFF.91
&, waterflux_ice(imt) ! IN Water flux from ice (Kg/m2/s) CALCDIFF.92
! on row j CALCDIFF.93
&, lambda_Large ! IN Value used in calculating minimum mld CALCDIFF.94
&, phit ! IN latitude of row j on tracer grid (radians) CALCDIFF.95
C CALCDIFF.96
INTEGER CALCDIFF.97
& max_Large_levels ! IN Max no of levs for Large scheme CALCDIFF.98
&, no_layers_in_lev ! IN No of layers in a level to calc mld_Large CALCDIFF.99
LOGICAL L_OBULKRI,L_OWINDMIX,L_OBULKMAXMLD CALCDIFF.100
CALCDIFF.101
REAL CALCDIFF.102
& OCEANHEATFLUX(IMT) !HTN:NON-PEN SURF HEATFLUX INTO OCEAN BUDGET CALCDIFF.103
&,CARYHEAT(IMT) !MISCELLANEOUS HEATFLUX FROM ICE CALCDIFF.104
&,FLXTOICE(IMT) !OCEAN TO ICE HEATFLUX CALCDIFF.105
CALCDIFF.106
REAL WSXM(IMT),WSYM(IMT) ! wind stress on row j-1 CALCDIFF.107
&,RiT(imt_gnu,km_gnu-1) ! Richardson number for tracer grid CALCDIFF.108
&,hT(imt_qlarge) ! Maximum depth, quadratic Large scheme CALCDIFF.109
&,L_T(imt) ! Monin Obukhov Length, Large scheme (tracer) CALCDIFF.110
&,Rimldcalc(imt,kmm1) ! Richardson no from mld calculation CALCDIFF.111
&,ISX(IMT_idr) ! IN Stress under sea ice fraction, row j CALCDIFF.112
&,ISY(IMT_idr) ! IN Stress under sea ice fraction, row j CALCDIFF.113
&,WSX_LEADS(IMT_idr) ! IN Stress under leads fraction, row j CALCDIFF.114
&,WSY_LEADS(IMT_idr) ! IN Stress under leads fraction, row j CALCDIFF.115
&,ISXM(IMT_idr) ! IN Stress under sea ice fraction, row j-1 CALCDIFF.116
&,ISYM(IMT_idr) ! IN Stress under sea ice fraction, row j-1 CALCDIFF.117
&,WSX_LEADSM(IMT_idr) ! IN Stress under leads fraction, row j-1 CALCDIFF.118
&,WSY_LEADSM(IMT_idr) ! IN Stress under leads fraction, row j-1 CALCDIFF.119
CALCDIFF.120
REAL FXA CALCDIFF.121
CALCDIFF.122
REAL Ai_ez(imt_iso,km_iso,0:1,0:1),Ai_bx(imt_iso,km_iso,0:1,0:1), CALCDIFF.123
& Ai_by(imt_iso,km_iso,0:1,0:1),Ai_nz(imt_iso,km_iso,0:1,0:1), CALCDIFF.124
& K11(imt_iso,km_iso),K22(imt_iso,km_iso),K33(imt_iso,km_iso) CALCDIFF.125
CALCDIFF.126
C Declare local variables CALCDIFF.127
C CALCDIFF.128
REAL CALCDIFF.129
& worka(IMT,KM), ! workspace CALCDIFF.130
& workb(IMT,KM), ! workspace CALCDIFF.131
& workc(IMT,KM) ! workspace CALCDIFF.132
CALCDIFF.133
REAL CALCDIFF.134
& XSTRESS_ICE(IMT_idr) !IN Total stress under sea ice row j CALCDIFF.135
&,YSTRESS_ICE(IMT_idr) ! (Wind stress at ice free points) CALCDIFF.136
&,XSTRESS_ICEM(IMT_idr) !IN Total stress under sea ice row j-1 CALCDIFF.137
&,YSTRESS_ICEM(IMT_idr) ! (Wind stress at ice free points) CALCDIFF.138
C CALCDIFF.139
C OISOGM variables CALCDIFF.140
REAL adv_vetiso(imt_gmm,km_gmm), CALCDIFF.141
& adv_vbtiso(imt_gmm,0:km_gmm),adv_fbiso(imt_gmm,0:km_gmm) CALCDIFF.142
REAL srho(imt,km) CALCDIFF.143
C CALCDIFF.144
C CALCDIFF.145
C* CALCDIFF.146
C*L---- External subroutines called ------------------------------------ CALCDIFF.147
C VERTCOFT CALCDIFF.148
C ISOPYC_M CALCDIFF.149
C ISOFLUX CALCDIFF.150
C STATED CALCDIFF.151
C ---------------------------------------------------------------- CALCDIFF.152
C First we must calculate gnu CALCDIFF.153
C ---------------------------------------------------------------- CALCDIFF.154
C CALCDIFF.155
IF (L_ORICHARD) THEN CALCDIFF.156
C CALCDIFF.157
IF(L_OSTATEC)THEN CALCDIFF.158
C ---------------------------------------------------------------- CALCDIFF.159
C COMPUTE RHOSRNA/RHOSRNB FOR CALL TO VERTCOFT USING STATEC CALCDIFF.160
C ---------------------------------------------------------------- CALCDIFF.161
CALL STATEC(TB(1,1,1),TB(1,1,2),RHOSRNA,WORKA,WORKB,1, CALCDIFF.162
& IMT,KM,J,JMT) CALCDIFF.163
CALL STATEC(TB(1,1,1),TB(1,1,2),RHOSRNB,WORKA,WORKB,2, CALCDIFF.164
& IMT,KM,J,JMT) CALCDIFF.165
DO I=1,IMT CALCDIFF.166
RHOSRNA(I,KM+1)=RHOSRNA(I,KM) CALCDIFF.167
RHOSRNB(I,KM+1)=RHOSRNB(I,KM) CALCDIFF.168
ENDDO CALCDIFF.169
C CALCDIFF.170
ELSE CALCDIFF.171
C CALCDIFF.172
C ---------------------------------------------------------------- CALCDIFF.173
C Compute rhosrn for call to VERTCOFT USING STATED CALCDIFF.174
C ---------------------------------------------------------------- CALCDIFF.175
CALL STATED(TB(1,1,1),TB(1,1,2),rhosrn,worka,workb,IMT,KM CALCDIFF.176
& ,J,KM,JMT CALCDIFF.177
& ) CALCDIFF.178
C CALCDIFF.179
ENDIF ! L_OSTATEC CALCDIFF.180
C CALCDIFF.181
C ---------------------------------------------------------------- CALCDIFF.182
C Calculate mixed layer depth CALCDIFF.183
C ---------------------------------------------------------------- CALCDIFF.184
IF (L_OFULARGE) THEN CALCDIFF.185
IF (L_OBULKRI) THEN CALCDIFF.186
CALL CALC_MLD_LARGEB(J,KM,IMT,NT CALCDIFF.187
&, RHO_WATER_SI,GRAV_SI,SPECIFIC_HEAT_SI,CRIT_RI_FL CALCDIFF.188
&, ZDZ,DZ,ZDZZ,DZZ,L_OCYCLIC,L_SEAICE CALCDIFF.189
&, UB,VB,UBM,VBM,TB,RHOSRN,HTN CALCDIFF.190
&, PME,WATERFLUX_ICE,SOL,WME,OCEANHEATFLUX,CARYHEAT,FLXTOICE CALCDIFF.191
&, FM,MLD_LARGE,RIMLDCALC CALCDIFF.192
& ) CALCDIFF.193
ELSE CALCDIFF.194
CALL CALC_MLD_LARGEG(J,KM,IMT CALCDIFF.195
&, RHO_WATER_SI,GRAV_SI,CRIT_RI_FL,NO_LAYERS_IN_LEV CALCDIFF.196
&, ZDZ,DZ,DZZ,L_OCYCLIC CALCDIFF.197
&, UB,VB,UBM,VBM CALCDIFF.198
&, FM,RHOSRN,MLD_LARGE,RIMLDCALC CALCDIFF.199
& ) CALCDIFF.200
ENDIF CALCDIFF.201
ENDIF ! L_OFULARGE CALCDIFF.202
C CALCDIFF.203
C ---------------------------------------------------------------- CALCDIFF.204
C Call subroutine to calculate vertical coefficient of diffusion CALCDIFF.205
C ---------------------------------------------------------------- CALCDIFF.206
C CALCDIFF.207
IF (L_SEAICE.AND.L_ICEFREEDR) THEN CALCDIFF.208
do i=1,imt CALCDIFF.209
XSTRESS_ICE(I) = WSX_LEADS(I)+ISX(I) CALCDIFF.210
YSTRESS_ICE(I) = WSY_LEADS(I)+ISY(I) CALCDIFF.211
XSTRESS_ICEM(I) = WSX_LEADSM(I)+ISXM(I) CALCDIFF.212
YSTRESS_ICEM(I) = WSY_LEADSM(I)+ISYM(I) CALCDIFF.213
enddo CALCDIFF.214
C CALCDIFF.215
CALL VERTCOFT CALCDIFF.216
& (J,IMT,KM,KMM1,NT, CALCDIFF.217
& IMT_QLARGE,L_OQLARGE,L_OCYCLIC, CALCDIFF.218
& IMT_FULARGE,L_OFULARGE,L_OUSTARWME,L_OWINDMIX,L_OBULKMAXMLD, CALCDIFF.219
& L_OROTATE,L_OSTATEC,L_SEAICE,L_OPANDP,PHIT,TB, CALCDIFF.220
& UB,VB,UBM,VBM, CALCDIFF.221
& XSTRESS_ICE,YSTRESS_ICE,XSTRESS_ICEM,YSTRESS_ICEM, CALCDIFF.222
& ZDZZ,ZDZ,DZ,DZZ, CALCDIFF.223
& DZZ2RQ,DZ2RQ, CALCDIFF.224
& NERGY,GRAV_SI,FM, CALCDIFF.225
& RHOSRN,RHOSRNA,RHOSRNB, CALCDIFF.226
& HT,RIT,MAX_QLARGE_DEPTH,CRIT_RI, CALCDIFF.227
& MLD_LARGE,MAX_LARGE_DEPTH,MAX_LARGE_LEVELS,RHO_WATER_SI, CALCDIFF.228
& HTN,PME,WATERFLUX_ICE,SOL,WME,L_T,LAMBDA_LARGE, CALCDIFF.229
& SPECIFIC_HEAT_SI, CALCDIFF.230
& GNU(1,1),FNU0_SI,FNUB_SI,STABLM_SI,FKAPB_SI,GNUMINT_SI, CALCDIFF.231
& KAPPA_B_SI,OCEANHEATFLUX,CARYHEAT,FLXTOICE,CGFLUX ) CALCDIFF.232
ELSE CALCDIFF.233
CALL VERTCOFT CALCDIFF.234
& (J,IMT,KM,KMM1,NT, CALCDIFF.235
& IMT_QLARGE,L_OQLARGE,L_OCYCLIC, CALCDIFF.236
& IMT_FULARGE,L_OFULARGE,L_OUSTARWME,L_OWINDMIX,L_OBULKMAXMLD, CALCDIFF.237
& L_OROTATE,L_OSTATEC,L_SEAICE,L_OPANDP,PHIT,TB, CALCDIFF.238
& UB,VB,UBM,VBM, CALCDIFF.239
& WSX, WSY,WSXM,WSYM, CALCDIFF.240
& ZDZZ,ZDZ,DZ,DZZ, CALCDIFF.241
& DZZ2RQ,DZ2RQ, CALCDIFF.242
& NERGY,GRAV_SI,FM, CALCDIFF.243
& RHOSRN,RHOSRNA,RHOSRNB, CALCDIFF.244
& HT,RIT,MAX_QLARGE_DEPTH,CRIT_RI, CALCDIFF.245
& MLD_LARGE,MAX_LARGE_DEPTH,MAX_LARGE_LEVELS,RHO_WATER_SI, CALCDIFF.246
& HTN,PME,WATERFLUX_ICE,SOL,WME,L_T,LAMBDA_LARGE, CALCDIFF.247
& SPECIFIC_HEAT_SI, CALCDIFF.248
& GNU(1,1),FNU0_SI,FNUB_SI,STABLM_SI,FKAPB_SI,GNUMINT_SI, CALCDIFF.249
& KAPPA_B_SI,OCEANHEATFLUX,CARYHEAT,FLXTOICE,CGFLUX ) CALCDIFF.250
ENDIF CALCDIFF.251
C CALCDIFF.252
ELSE CALCDIFF.253
DO K=1,KM CALCDIFF.254
DO I=1,IMT CALCDIFF.255
gnu(I,K)=FKPH CALCDIFF.256
ENDDO ! Over I CALCDIFF.257
ENDDO ! Over K CALCDIFF.258
ENDIF CALCDIFF.259
C CALCDIFF.260
do k=1,km CALCDIFF.261
do i=1,imt CALCDIFF.262
srho(i,k)=0. CALCDIFF.263
enddo CALCDIFF.264
enddo CALCDIFF.265
C ------------------------------------------------------------------ CALCDIFF.266
C Having obtained required values of gnu, CALCDIFF.267
C call subroutine to calculate isopycnal diffusion tensor CALCDIFF.268
C ------------------------------------------------------------------ CALCDIFF.269
C CALCDIFF.270
c CALCDIFF.271
CALL ISOPYC_M( CALCDIFF.272
*CALL ARGSIZE CALCDIFF.273
*CALL COCAWRKA CALCDIFF.274
& ,j,cstr,dyur,dxur,dz2r,dzz,dzz2r,athkdftu,athkdftv,ahi,athkdf, CALCDIFF.275
& slope_max,dz,dyu,dxu,cs,dxt4r,dyt4r,slopec,dslope,dxtr,dytr, CALCDIFF.276
& cst,csr,Ai_ez,Ai_nz,Ai_bx,Ai_by,K11,K22,K33,gnu(1,1), CALCDIFF.277
& adv_vetiso,adv_vbtiso,adv_fbiso,J_OFFSET,athkdf_bi, CALCDIFF.278
& csjm,dyurjm,j_1,srho) CALCDIFF.279
C CALCDIFF.280
C CALCDIFF.281
C -------------------------------------------------------------- CALCDIFF.282
C Calculate diff_fn (equivalent to esav) ready for calculations CALCDIFF.283
C involving the first row in this parallel block. CALCDIFF.284
C (Tracer values are unavoidably populated with spurious figures CALCDIFF.285
C here but since they are kept local to this subroutine, the CALCDIFF.286
C only problem should be one of performance.) CALCDIFF.287
C -------------------------------------------------------------- CALCDIFF.288
! Initialise TA and diff_fn CALCDIFF.289
DO M = 1, NT CALCDIFF.290
DO K = 1, KM CALCDIFF.291
DO I = 1, IMT CALCDIFF.292
TA(I,K,M) = 0.0 CALCDIFF.293
diff_fn(I,K,M,0) = 0.0 CALCDIFF.294
diff_fn(I,K,M,1) = 0.0 CALCDIFF.295
ENDDO ! over I CALCDIFF.296
ENDDO ! over K CALCDIFF.297
ENDDO ! over M CALCDIFF.298
C CALCDIFF.299
c CALCDIFF.300
DO M=1,NT CALCDIFF.301
C CALCDIFF.302
CALL ISOFLUX( CALCDIFF.303
*CALL ARGSIZE CALCDIFF.304
*CALL COCAWRKA CALCDIFF.305
& ,m,j,cstr,dyur,dxur,dz2r,cs,csjm, CALCDIFF.306
& dytr,dxtr,cst,dxt4r,dyt4r,Ai_ez,Ai_nz,Ai_bx,Ai_by,K11,K22,K33, CALCDIFF.307
& itt,adv_vbtiso,adv_fbiso,worka,workb,workc,j_1) CALCDIFF.308
C CALCDIFF.309
C CALCDIFF.310
ENDDO CALCDIFF.311
CALCDIFF.312
RETURN CALCDIFF.313
END CALCDIFF.314
*ENDIF CALCDIFF.315