** Warning **
Issuing rollback() due to DESTROY without explicit disconnect() of DBD::mysql::db handle dbname=MITgcm at /usr/local/share/lxr/lib/LXR/Common.pm line 1224.
Last-Modified: Thu, 15 May 2024 05:11:27 GMT
Content-Type: text/html; charset=utf-8
MITgcm/MITgcm/pkg/layers/layers_fluxcalc.F
File indexing completed on 2018-03-09 06:09:17 UTC
view on github raw file Latest commit 84f265d4 on 2018-03-08 18:49:29 UTC
aa668e01f1 Gael* 0001 #include "LAYERS_OPTIONS.h "
0002 #ifdef ALLOW_GMREDI
0003 #include "GMREDI_OPTIONS.h "
0004 #endif
0005
28ff9bcfe6 Jean* 0006
9479dbe556 Mart* 0007
0008
cf336ab6c5 Ryan* 0009
9479dbe556 Mart* 0010
0011
aa668e01f1 Gael* 0012
9479dbe556 Mart* 0013
0014
aa668e01f1 Gael* 0015 SUBROUTINE LAYERS_FLUXCALC (
406891c1a3 Gael* 0016 I uVel ,vVel ,tracer ,iLa ,
50d8304171 Ryan* 0017 O UH ,VH ,Hw ,Hs ,PIw ,PIs ,Uw ,Vs ,
aa668e01f1 Gael* 0018 I myThid )
0019
9479dbe556 Mart* 0020
0021
0022
0023
0024
0025
0026
aa668e01f1 Gael* 0027
0028
0029 IMPLICIT NONE
9479dbe556 Mart* 0030
aa668e01f1 Gael* 0031 #include "SIZE.h "
0032 #include "EEPARAMS.h "
0033 #include "PARAMS.h "
0034 #include "GRID.h "
0035 #include "LAYERS_SIZE.h "
0036 #include "LAYERS.h "
0037 #ifdef ALLOW_GMREDI
0038 # include "GMREDI.h "
0039 #endif
0040
0041
0042
0043
0044
d746ff23b0 Jean* 0045
aa668e01f1 Gael* 0046
0047
0048
0049
14bb2ae926 Ryan* 0050
0051
50d8304171 Ryan* 0052
0053
406891c1a3 Gael* 0054
aa668e01f1 Gael* 0055 INTEGER myThid
28ff9bcfe6 Jean* 0056 _RL uVel (1-OLx :sNx +OLx ,1-OLy :sNy +OLy ,Nr , nSx ,nSy )
0057 _RL vVel (1-OLx :sNx +OLx ,1-OLy :sNy +OLy ,Nr , nSx ,nSy )
0058 _RL tracer (1-OLx :sNx +OLx ,1-OLy :sNy +OLy ,Nr , nSx ,nSy )
6f9adcbc06 Jean* 0059 #ifdef LAYERS_UFLUX
28ff9bcfe6 Jean* 0060 _RL UH (1-OLx :sNx +OLx ,1-OLy :sNy +OLy ,Nlayers ,nSx ,nSy )
6f9adcbc06 Jean* 0061 # ifdef LAYERS_THICKNESS
28ff9bcfe6 Jean* 0062 _RL Hw (1-OLx :sNx +OLx ,1-OLy :sNy +OLy ,Nlayers ,nSx ,nSy )
0063 _RL PIw (1-OLx :sNx +OLx ,1-OLy :sNy +OLy ,Nlayers ,nSx ,nSy )
50d8304171 Ryan* 0064 _RL Uw (1-OLx :sNx +OLx ,1-OLy :sNy +OLy ,Nlayers ,nSx ,nSy )
6f9adcbc06 Jean* 0065 # else
50d8304171 Ryan* 0066 _RL Hw (1), PIw (1), Uw (1)
6f9adcbc06 Jean* 0067 # endif
0068 #else
06e3ddabf4 Jean* 0069 _RL UH (1), Hw (1), PIw (1), Uw (1)
6f9adcbc06 Jean* 0070 #endif
0071 #ifdef LAYERS_VFLUX
0072 _RL VH (1-OLx :sNx +OLx ,1-OLy :sNy +OLy ,Nlayers ,nSx ,nSy )
0073 # ifdef LAYERS_THICKNESS
0074 _RL Hs (1-OLx :sNx +OLx ,1-OLy :sNy +OLy ,Nlayers ,nSx ,nSy )
0075 _RL PIs (1-OLx :sNx +OLx ,1-OLy :sNy +OLy ,Nlayers ,nSx ,nSy )
50d8304171 Ryan* 0076 _RL Vs (1-OLx :sNx +OLx ,1-OLy :sNy +OLy ,Nlayers ,nSx ,nSy )
6f9adcbc06 Jean* 0077 # else
50d8304171 Ryan* 0078 _RL Hs (1), PIs (1), Vs (1)
6f9adcbc06 Jean* 0079 # endif
0080 #else
50d8304171 Ryan* 0081 _RL VH (1), Hs (1), PIs (1), Vs (1)
6f9adcbc06 Jean* 0082 #endif
406891c1a3 Gael* 0083 INTEGER iLa
aa668e01f1 Gael* 0084
0085
0086 #ifdef ALLOW_LAYERS
0087
0088
0089
0090
0091
0092
406891c1a3 Gael* 0093
aa668e01f1 Gael* 0094
0095
9479dbe556 Mart* 0096
0097
6f9adcbc06 Jean* 0098
5b31d86392 Mart* 0099
aa668e01f1 Gael* 0100
06e3ddabf4 Jean* 0101
50d8304171 Ryan* 0102
aa668e01f1 Gael* 0103
0104 INTEGER bi , bj
9479dbe556 Mart* 0105 INTEGER i ,j ,k ,kk ,kg ,kci ,kp1 ,kloc
0106 INTEGER mSteps
aa668e01f1 Gael* 0107 INTEGER kgu (sNx +1,sNy +1), kgv (sNx +1,sNy +1)
5b31d86392 Mart* 0108 _RL TatU (sNx +1,sNy +1), TatV (sNx +1,sNy +1)
f6012f0463 Ryan* 0109 _RL dzfac
d746ff23b0 Jean* 0110 #ifdef ALLOW_GMREDI
aa668e01f1 Gael* 0111 INTEGER kcip1
0112 _RL delPsi , maskp1
0113 #endif
9479dbe556 Mart* 0114 LOGICAL errorFlag
0115 CHARACTER *(MAX_LEN_MBUF ) msgBuf
aa668e01f1 Gael* 0116
0117
0118
28ff9bcfe6 Jean* 0119
9479dbe556 Mart* 0120
0121 mSteps = int(log10(dble(Nlayers ))/log10(2. _d 0))+1
0122
aa668e01f1 Gael* 0123
0124 DO bj =myByLo (myThid ),myByHi (myThid )
0125 DO bi =myBxLo (myThid ),myBxHi (myThid )
0126
0127
0128 DO j = 1,sNy +1
0129 DO i = 1,sNx +1
0130
0131
0132
0133 kgu (i ,j ) = Nlayers
0134 kgv (i ,j ) = Nlayers
0135 ENDDO
0136 ENDDO
0137
0138
0139 DO kg =1,Nlayers
28ff9bcfe6 Jean* 0140 DO j = 1-OLy ,sNy +OLy
0141 DO i = 1-OLx ,sNx +OLx
aa668e01f1 Gael* 0142 #ifdef LAYERS_UFLUX
1377839540 Jean* 0143 UH (i ,j ,kg ,bi ,bj ) = 0. _d 0
aa668e01f1 Gael* 0144 #ifdef LAYERS_THICKNESS
50d8304171 Ryan* 0145 Hw (i ,j ,kg ,bi ,bj ) = 0. _d 0
1377839540 Jean* 0146 PIw (i ,j ,kg ,bi ,bj ) = 0. _d 0
50d8304171 Ryan* 0147 Uw (i ,j ,kg ,bi ,bj ) = 0. _d 0
aa668e01f1 Gael* 0148 #endif /* LAYERS_THICKNESS */
0149 #endif /* UH */
0150 #ifdef LAYERS_VFLUX
1377839540 Jean* 0151 VH (i ,j ,kg ,bi ,bj ) = 0. _d 0
aa668e01f1 Gael* 0152 #ifdef LAYERS_THICKNESS
50d8304171 Ryan* 0153 Hs (i ,j ,kg ,bi ,bj ) = 0. _d 0
1377839540 Jean* 0154 PIs (i ,j ,kg ,bi ,bj ) = 0. _d 0
50d8304171 Ryan* 0155 Vs (i ,j ,kg ,bi ,bj ) = 0. _d 0
aa668e01f1 Gael* 0156 #endif /* LAYERS_THICKNESS */
0157 #endif /* VH */
0158 ENDDO
0159 ENDDO
0160 ENDDO
0161
0162 DO kk =1,NZZ
0163 k = MapIndex (kk )
0164 kci = CellIndex (kk )
d746ff23b0 Jean* 0165 #ifdef ALLOW_GMREDI
0166 kcip1 = MIN(kci +1,Nr )
0167 maskp1 = 1.
0168 IF (kci .GE. Nr ) maskp1 = 0.
0169 #endif /* ALLOW_GMREDI */
5b31d86392 Mart* 0170 #ifdef LAYERS_UFLUX
aa668e01f1 Gael* 0171 DO j = 1,sNy +1
0172 DO i = 1,sNx +1
0173
0174
0175 kp1 =k +1
f6012f0463 Ryan* 0176 IF (maskW (i ,j ,kp1 ,bi ,bj ).EQ. zeroRS ) kp1 =k
0177 TatU (i ,j ) = MapFact (kk ) *
0178 & 0.5 _d 0 * (tracer (i -1,j ,k ,bi ,bj )+tracer (i ,j ,k ,bi ,bj )) +
0179 & (1. _d 0 -MapFact (kk )) *
0180 & 0.5 _d 0 * (tracer (i -1,j ,kp1 ,bi ,bj )+tracer (i ,j ,kp1 ,bi ,bj ))
aa668e01f1 Gael* 0181
5b31d86392 Mart* 0182 ENDDO
0183 ENDDO
aa668e01f1 Gael* 0184
9479dbe556 Mart* 0185
0186 CALL LAYERS_LOCATE (
28ff9bcfe6 Jean* 0187 I layers_bounds (1,iLa ),Nlayers ,mSteps ,sNx ,sNy ,TatU ,
9479dbe556 Mart* 0188 O kgu ,
0189 I myThid )
28ff9bcfe6 Jean* 0190 #ifndef TARGET_NEC_SX
9479dbe556 Mart* 0191
0192 IF ( debugLevel .GE. debLevC ) THEN
0193 errorFlag = .FALSE.
0194 DO j = 1,sNy +1
0195 DO i = 1,sNx +1
0196 IF ( kgu (i ,j ) .LE. 0 ) THEN
0197 WRITE (msgBuf ,'(2A,I3,A,I3,A,1E14.6)' )
0198 & 'S/R LAYERS_LOCATE: Could not find a bin in ' ,
0199 & 'layers_bounds for TatU(' ,i ,',' ,j ,',)=' ,TatU (i ,j )
0200 CALL PRINT_ERROR ( msgBuf , myThid )
0201 errorFlag = .TRUE.
0202 ENDIF
0203 ENDDO
5b31d86392 Mart* 0204 ENDDO
9479dbe556 Mart* 0205 IF ( errorFlag ) STOP 'ABNORMAL END: S/R LAYERS_FLUXCALC'
0206 ENDIF
0207 #endif /* ndef TARGET_NEC_SX */
28ff9bcfe6 Jean* 0208
5b31d86392 Mart* 0209 DO j = 1,sNy +1
0210 DO i = 1,sNx +1
0211
9479dbe556 Mart* 0212 kloc = kgu (i ,j )
f6012f0463 Ryan* 0213 dzfac = dZZf (kk ) * hFacW (i ,j ,kci ,bi ,bj )
0214
aa668e01f1 Gael* 0215
9479dbe556 Mart* 0216 UH (i ,j ,kloc ,bi ,bj ) =
0217 & UH (i ,j ,kloc ,bi ,bj ) +
f6012f0463 Ryan* 0218 & dzfac * uVel (i ,j ,kci ,bi ,bj )
aa668e01f1 Gael* 0219
0220 #ifdef ALLOW_GMREDI
d746ff23b0 Jean* 0221 IF ( layers_bolus (iLa ) ) THEN
0222 IF ( .NOT. GM_AdvForm ) THEN
9ecea62500 Jean* 0223 delPsi = 0.25 _d 0 *(
0224 & ( rA (i -1,j ,bi ,bj )*Kwx (i -1,j ,kcip1 ,bi ,bj )
0225 & +rA ( i ,j ,bi ,bj )*Kwx ( i ,j ,kcip1 ,bi ,bj )
0226 & ) * maskW (i ,j ,kcip1 ,bi ,bj ) * maskp1
0227 & - ( rA (i -1,j ,bi ,bj )*Kwx (i -1,j , kci ,bi ,bj )
0228 & +rA ( i ,j ,bi ,bj )*Kwx ( i ,j , kci ,bi ,bj )
0229 & ) * maskW (i ,j , kci ,bi ,bj )
0230 & ) * recip_rAw (i ,j ,bi ,bj )
d746ff23b0 Jean* 0231 #ifdef GM_BOLUS_ADVEC
0232 ELSE
0233 delPsi = GM_PsiX (i ,j ,kcip1 ,bi ,bj )*maskp1
0234 & - GM_PsiX (i ,j , kci , bi ,bj )
0235 #endif
0236 ENDIF
9479dbe556 Mart* 0237 UH (i ,j ,kloc ,bi ,bj ) = UH (i ,j ,kloc ,bi ,bj )
aa668e01f1 Gael* 0238 & + delPsi *recip_drF (kci )*_recip_hFacW (i ,j ,kci ,bi ,bj )
f6012f0463 Ryan* 0239 & * dzfac
aa668e01f1 Gael* 0240 ENDIF
d746ff23b0 Jean* 0241 #endif /* ALLOW_GMREDI */
aa668e01f1 Gael* 0242
0243 #ifdef LAYERS_THICKNESS
f6012f0463 Ryan* 0244 Hw (i ,j ,kloc ,bi ,bj ) = Hw (i ,j ,kloc ,bi ,bj ) + dzfac
aa668e01f1 Gael* 0245 #endif /* LAYERS_THICKNESS */
0246
5b31d86392 Mart* 0247 ENDDO
0248 ENDDO
aa668e01f1 Gael* 0249 #endif /* LAYERS_UFLUX */
0250
0251 #ifdef LAYERS_VFLUX
5b31d86392 Mart* 0252 DO j = 1,sNy +1
0253 DO i = 1,sNx +1
aa668e01f1 Gael* 0254
0255 kp1 =k +1
f6012f0463 Ryan* 0256 IF (maskS (i ,j ,kp1 ,bi ,bj ).EQ. zeroRS ) kp1 =k
0257 TatV (i ,j ) = MapFact (kk ) *
0258 & 0.5 _d 0 * (tracer (i ,j -1,k ,bi ,bj )+tracer (i ,j ,k ,bi ,bj )) +
0259 & (1. _d 0 -MapFact (kk )) *
0260 & 0.5 _d 0 * (tracer (i ,j -1,kp1 ,bi ,bj )+tracer (i ,j ,kp1 ,bi ,bj ))
aa668e01f1 Gael* 0261
5b31d86392 Mart* 0262 ENDDO
0263 ENDDO
9479dbe556 Mart* 0264
0265
0266 CALL LAYERS_LOCATE (
28ff9bcfe6 Jean* 0267 I layers_bounds (1,iLa ),Nlayers ,mSteps ,sNx ,sNy ,TatV ,
9479dbe556 Mart* 0268 O kgv ,
0269 I myThid )
0270 #ifndef TARGET_NEC_SX
0271 IF ( debugLevel .GE. debLevC ) THEN
0272
0273 errorFlag = .FALSE.
0274 DO j = 1,sNy +1
0275 DO i = 1,sNx +1
0276 IF ( kgv (i ,j ) .LE. 0 ) THEN
0277 WRITE (msgBuf ,'(2A,I3,A,I3,A,1E14.6)' )
0278 & 'S/R LAYERS_LOCATE: Could not find a bin in ' ,
0279 & 'layers_bounds for TatV(' ,i ,',' ,j ,',)=' ,TatV (i ,j )
0280 CALL PRINT_ERROR ( msgBuf , myThid )
0281 errorFlag = .TRUE.
0282 ENDIF
0283 ENDDO
5b31d86392 Mart* 0284 ENDDO
9479dbe556 Mart* 0285 IF ( errorFlag ) STOP 'ABNORMAL END: S/R LAYERS_FLUXCALC'
0286 ENDIF
0287 #endif /* ndef TARGET_NEC_SX */
0288
5b31d86392 Mart* 0289 DO j = 1,sNy +1
0290 DO i = 1,sNx +1
0291
9479dbe556 Mart* 0292 kloc = kgv (i ,j )
50d8304171 Ryan* 0293 dzfac = dZZf (kk ) * hFacS (i ,j ,kci ,bi ,bj )
0294
0295
0296
0297
0298
0299
0300
0301
0302
0303
0304
aa668e01f1 Gael* 0305
9479dbe556 Mart* 0306 VH (i ,j ,kloc ,bi ,bj ) =
50d8304171 Ryan* 0307 & VH (i ,j ,kloc ,bi ,bj ) + dzfac * vVel (i ,j ,kci ,bi ,bj )
aa668e01f1 Gael* 0308
0309 #ifdef ALLOW_GMREDI
d746ff23b0 Jean* 0310 IF ( layers_bolus (iLa ) ) THEN
0311 IF ( .NOT. GM_AdvForm ) THEN
9ecea62500 Jean* 0312 delPsi = 0.25 _d 0 *(
0313 & ( rA (i ,j -1,bi ,bj )*Kwy (i ,j -1,kcip1 ,bi ,bj )
0314 & +rA (i , j ,bi ,bj )*Kwy (i , j ,kcip1 ,bi ,bj )
0315 & ) * maskS (i ,j ,kcip1 ,bi ,bj ) * maskp1
0316 & - ( rA (i ,j -1,bi ,bj )*Kwy (i ,j -1, kci ,bi ,bj )
0317 & +rA (i , j ,bi ,bj )*Kwy (i , j , kci ,bi ,bj )
0318 & ) * maskS (i ,j , kci ,bi ,bj )
0319 & ) * recip_rAs (i ,j ,bi ,bj )
d746ff23b0 Jean* 0320 #ifdef GM_BOLUS_ADVEC
0321 ELSE
0322 delPsi = GM_PsiY (i ,j ,kcip1 ,bi ,bj )*maskp1
0323 & - GM_PsiY (i ,j , kci , bi ,bj )
0324 #endif
0325 ENDIF
9479dbe556 Mart* 0326 VH (i ,j ,kloc ,bi ,bj ) = VH (i ,j ,kloc ,bi ,bj )
aa668e01f1 Gael* 0327 & + delPsi *recip_drF (kci )*_recip_hFacS (i ,j ,kci ,bi ,bj )
50d8304171 Ryan* 0328 & * dzfac
aa668e01f1 Gael* 0329 ENDIF
d746ff23b0 Jean* 0330 #endif /* ALLOW_GMREDI */
aa668e01f1 Gael* 0331
0332 #ifdef LAYERS_THICKNESS
50d8304171 Ryan* 0333 Hs (i ,j ,kloc ,bi ,bj ) = Hs (i ,j ,kloc ,bi ,bj ) + dzfac
aa668e01f1 Gael* 0334 #endif /* LAYERS_THICKNESS */
0335
0336 ENDDO
0337 ENDDO
5b31d86392 Mart* 0338 #endif /* LAYERS_VFLUX */
aa668e01f1 Gael* 0339 ENDDO
d746ff23b0 Jean* 0340
14bb2ae926 Ryan* 0341
0342
0343 #ifdef LAYERS_THICKNESS
0344 DO kg =1,Nlayers
0345 DO j = 1,sNy +1
0346 DO i = 1,sNx +1
0347 #ifdef LAYERS_UFLUX
0348 IF (Hw (i ,j ,kg ,bi ,bj ) .GT. 0.) THEN
0349 PIw (i ,j ,kg ,bi ,bj ) = 1. _d 0
50d8304171 Ryan* 0350 Uw (i ,j ,kg ,bi ,bj ) =
14bb2ae926 Ryan* 0351 & UH (i ,j ,kg ,bi ,bj ) / Hw (i ,j ,kg ,bi ,bj )
0352 ENDIF
0353 #endif /* LAYERS_UFLUX */
0354 #ifdef LAYERS_VFLUX
0355 IF (Hs (i ,j ,kg ,bi ,bj ) .GT. 0.) THEN
0356 PIs (i ,j ,kg ,bi ,bj ) = 1. _d 0
50d8304171 Ryan* 0357 Vs (i ,j ,kg ,bi ,bj ) =
14bb2ae926 Ryan* 0358 & VH (i ,j ,kg ,bi ,bj ) / Hs (i ,j ,kg ,bi ,bj )
0359 ENDIF
0360 #endif /* LAYERS_VFLUX */
0361 ENDDO
0362 ENDDO
0363 ENDDO
0364 #endif /* LAYERS_THICKNESS */
aa668e01f1 Gael* 0365
0366
0367 ENDDO
0368 ENDDO
0369
9479dbe556 Mart* 0370 RETURN
0371 END
28ff9bcfe6 Jean* 0372
0373
cf336ab6c5 Ryan* 0374
0375
0376
0377 SUBROUTINE LAYERS_DIAPYCNAL (
0378 I tracer ,iLa ,
0379 O TtendSurf , TtendDiffh , TtendDiffr ,
50d8304171 Ryan* 0380 O TtendAdvh , TtendAdvr , Ttendtot ,
cf336ab6c5 Ryan* 0381 O StendSurf , StendDiffh , StendDiffr ,
50d8304171 Ryan* 0382 O StendAdvh , StendAdvr , Stendtot ,
cf336ab6c5 Ryan* 0383 O Hc , PIc ,
0384 I myThid )
0385
0386
0387
0388
0389
0390
0391
0392
0393 IMPLICIT NONE
0394 #include "SIZE.h "
0395 #include "EEPARAMS.h "
0396 #include "PARAMS.h "
0397 #include "GRID.h "
0398 #include "LAYERS_SIZE.h "
0399 #include "LAYERS.h "
0400
0401
0402
0403
0404
1377839540 Jean* 0405
cf336ab6c5 Ryan* 0406
0407
50d8304171 Ryan* 0408
0409
06e3ddabf4 Jean* 0410
50d8304171 Ryan* 0411
0412
0413
0414
cf336ab6c5 Ryan* 0415
0416
0417 INTEGER iLa , myThid
0418 _RL tracer (1-OLx :sNx +OLx ,1-OLy :sNy +OLy ,Nr , nSx ,nSy )
0419 _RL Hc (1-OLx :sNx +OLx ,1-OLy :sNy +OLy ,Nlayers ,nSx ,nSy )
0420 _RL PIc (1-OLx :sNx +OLx ,1-OLy :sNy +OLy ,Nlayers ,nSx ,nSy )
50d8304171 Ryan* 0421 _RL TtendSurf (1-OLx :sNx +OLx ,1-OLy :sNy +OLy ,Nlayers -1,nSx ,nSy )
0422 _RL TtendDiffh (1-OLx :sNx +OLx ,1-OLy :sNy +OLy ,Nlayers -1,nSx ,nSy )
06e3ddabf4 Jean* 0423 _RL TtendDiffr (1-OLx :sNx +OLx ,1-OLy :sNy +OLy ,Nlayers -1,nSx ,nSy )
50d8304171 Ryan* 0424 _RL TtendAdvh (1-OLx :sNx +OLx ,1-OLy :sNy +OLy ,Nlayers -1,nSx ,nSy )
06e3ddabf4 Jean* 0425 _RL TtendAdvr (1-OLx :sNx +OLx ,1-OLy :sNy +OLy ,Nlayers -1,nSx ,nSy )
0426 _RL Ttendtot (1-OLx :sNx +OLx ,1-OLy :sNy +OLy ,Nlayers -1,nSx ,nSy )
50d8304171 Ryan* 0427 _RL StendSurf (1-OLx :sNx +OLx ,1-OLy :sNy +OLy ,Nlayers -1,nSx ,nSy )
0428 _RL StendDiffh (1-OLx :sNx +OLx ,1-OLy :sNy +OLy ,Nlayers -1,nSx ,nSy )
06e3ddabf4 Jean* 0429 _RL StendDiffr (1-OLx :sNx +OLx ,1-OLy :sNy +OLy ,Nlayers -1,nSx ,nSy )
50d8304171 Ryan* 0430 _RL StendAdvh (1-OLx :sNx +OLx ,1-OLy :sNy +OLy ,Nlayers -1,nSx ,nSy )
06e3ddabf4 Jean* 0431 _RL StendAdvr (1-OLx :sNx +OLx ,1-OLy :sNy +OLy ,Nlayers -1,nSx ,nSy )
0432 _RL Stendtot (1-OLx :sNx +OLx ,1-OLy :sNy +OLy ,Nlayers -1,nSx ,nSy )
cf336ab6c5 Ryan* 0433
0434
0435 #ifdef LAYERS_THERMODYNAMICS
0436
0437
0438
0439
0440
0441
f6012f0463 Ryan* 0442
cf336ab6c5 Ryan* 0443
0444
0445
0446
0447
0448
da9f56e003 Jean* 0449 _RL Hcw (1-OLx :sNx +OLx ,1-OLy :sNy +OLy ,Nlayers -1,nSx ,nSy )
cf336ab6c5 Ryan* 0450 INTEGER bi , bj
da9f56e003 Jean* 0451 INTEGER i ,j ,k ,kk ,kg ,kci ,kloc
cf336ab6c5 Ryan* 0452 INTEGER mSteps
0453 INTEGER kgc (sNx +1,sNy +1)
06e3ddabf4 Jean* 0454 INTEGER kgcw (sNx +1,sNy +1)
f6012f0463 Ryan* 0455 _RL TatC (sNx +1,sNy +1), dzfac , Tfac , Sfac
cf336ab6c5 Ryan* 0456 LOGICAL errorFlag
0457 CHARACTER *(MAX_LEN_MBUF ) msgBuf
da9f56e003 Jean* 0458 #ifdef LAYERS_FINEGRID_DIAPYCNAL
0459 INTEGER kp1
0460 #endif
cf336ab6c5 Ryan* 0461
50d8304171 Ryan* 0462
0463 Tfac = 1. _d 0
0464 Sfac = 1. _d 0
06e3ddabf4 Jean* 0465
cf336ab6c5 Ryan* 0466
0467
0468 mSteps = int(log10(dble(Nlayers ))/log10(2. _d 0))+1
0469
0470
0471
0472
0473 DO bj =myByLo (myThid ),myByHi (myThid )
0474 DO bi =myBxLo (myThid ),myBxHi (myThid )
0475
0476
0477 DO j = 1,sNy +1
0478 DO i = 1,sNx +1
0479
0480
0481
0482 kgc (i ,j ) = Nlayers
50d8304171 Ryan* 0483 kgcw (i ,j ) = Nlayers -1
cf336ab6c5 Ryan* 0484 ENDDO
0485 ENDDO
0486
0487
50d8304171 Ryan* 0488
0489 DO kg =1,Nlayers -1
cf336ab6c5 Ryan* 0490 DO j = 1-OLy ,sNy +OLy
0491 DO i = 1-OLx ,sNx +OLx
0492 TtendSurf (i ,j ,kg ,bi ,bj ) = 0. _d 0
0493 TtendDiffh (i ,j ,kg ,bi ,bj ) = 0. _d 0
0494 TtendDiffr (i ,j ,kg ,bi ,bj ) = 0. _d 0
50d8304171 Ryan* 0495 TtendAdvh (i ,j ,kg ,bi ,bj ) = 0. _d 0
0496 TtendAdvr (i ,j ,kg ,bi ,bj ) = 0. _d 0
0497 Ttendtot (i ,j ,kg ,bi ,bj ) = 0. _d 0
cf336ab6c5 Ryan* 0498 StendSurf (i ,j ,kg ,bi ,bj ) = 0. _d 0
0499 StendDiffh (i ,j ,kg ,bi ,bj ) = 0. _d 0
06e3ddabf4 Jean* 0500 StendDiffr (i ,j ,kg ,bi ,bj ) = 0. _d 0
50d8304171 Ryan* 0501 StendAdvh (i ,j ,kg ,bi ,bj ) = 0. _d 0
06e3ddabf4 Jean* 0502 StendAdvr (i ,j ,kg ,bi ,bj ) = 0. _d 0
50d8304171 Ryan* 0503 Stendtot (i ,j ,kg ,bi ,bj ) = 0. _d 0
06e3ddabf4 Jean* 0504 Hcw (i ,j ,kg ,bi ,bj ) = 0. _d 0
50d8304171 Ryan* 0505 ENDDO
0506 ENDDO
0507 ENDDO
0508
0509 DO kg =1,Nlayers
0510 DO j = 1-OLy ,sNy +OLy
06e3ddabf4 Jean* 0511 DO i = 1-OLx ,sNx +OLx
cf336ab6c5 Ryan* 0512 Hc (i ,j ,kg ,bi ,bj ) = 0. _d 0
1377839540 Jean* 0513 PIc (i ,j ,kg ,bi ,bj ) = 0. _d 0
cf336ab6c5 Ryan* 0514 ENDDO
0515 ENDDO
0516 ENDDO
0517
cdcbace143 Ryan* 0518 #ifdef LAYERS_FINEGRID_DIAPYCNAL
cf336ab6c5 Ryan* 0519 DO kk =1,NZZ
0520 k = MapIndex (kk )
0521 kci = CellIndex (kk )
0522 DO j = 1,sNy +1
0523 DO i = 1,sNx +1
f6012f0463 Ryan* 0524
0525 kp1 =k +1
0526 IF (maskC (i ,j ,kp1 ,bi ,bj ).EQ. zeroRS ) kp1 =k
0527 TatC (i ,j ) = MapFact (kk ) * tracer (i ,j ,k ,bi ,bj ) +
0528 & (1. _d 0 -MapFact (kk )) * tracer (i ,j ,kp1 ,bi ,bj )
cdcbace143 Ryan* 0529 ENDDO
0530 ENDDO
0531 #else
0532 DO kk =1,Nr
0533 k = kk
0534 kci = kk
0535 DO j = 1,sNy +1
0536 DO i = 1,sNx +1
0537 TatC (i ,j ) = tracer (i ,j ,k ,bi ,bj )
0538 ENDDO
0539 ENDDO
0540 #endif /* LAYERS_FINEGRID_DIAPYCNAL */
da9f56e003 Jean* 0541
f6012f0463 Ryan* 0542
0543
cdcbace143 Ryan* 0544
0545
f6012f0463 Ryan* 0546
0547
0548
0549
0550
0551
0552
0553
0554
0555
0556
cf336ab6c5 Ryan* 0557
50d8304171 Ryan* 0558
cf336ab6c5 Ryan* 0559 CALL LAYERS_LOCATE (
0560 I layers_bounds (1,iLa ),Nlayers ,mSteps ,sNx ,sNy ,TatC ,
0561 O kgc ,
0562 I myThid )
0563 #ifndef TARGET_NEC_SX
0564
0565 IF ( debugLevel .GE. debLevC ) THEN
0566 errorFlag = .FALSE.
0567 DO j = 1,sNy +1
0568 DO i = 1,sNx +1
0569 IF ( kgc (i ,j ) .LE. 0 ) THEN
0570 WRITE (msgBuf ,'(2A,I3,A,I3,A,1E14.6)' )
0571 & 'S/R LAYERS_LOCATE: Could not find a bin in ' ,
0572 & 'layers_bounds for TatC(' ,i ,',' ,j ,',)=' ,TatC (i ,j )
0573 CALL PRINT_ERROR ( msgBuf , myThid )
0574 errorFlag = .TRUE.
0575 ENDIF
0576 ENDDO
0577 ENDDO
0578 IF ( errorFlag ) STOP 'ABNORMAL END: S/R LAYERS_DIAPYCNAL'
0579 ENDIF
0580 #endif /* ndef TARGET_NEC_SX */
50d8304171 Ryan* 0581
0582
0583 CALL LAYERS_LOCATE (
0584 I layers_bounds_w (1,iLa ),Nlayers -1,mSteps ,sNx ,sNy ,TatC ,
0585 O kgcw ,
0586 I myThid )
0587 #ifndef TARGET_NEC_SX
0588
0589 IF ( debugLevel .GE. debLevC ) THEN
0590 errorFlag = .FALSE.
0591 DO j = 1,sNy +1
0592 DO i = 1,sNx +1
0593 IF ( kgcw (i ,j ) .LE. 0 ) THEN
0594 WRITE (msgBuf ,'(2A,I3,A,I3,A,1E14.6)' )
0595 & 'S/R LAYERS_LOCATE: Could not find a bin in ' ,
0596 & 'layers_bounds for TatC(' ,i ,',' ,j ,',)=' ,TatC (i ,j )
0597 CALL PRINT_ERROR ( msgBuf , myThid )
0598 errorFlag = .TRUE.
0599 ENDIF
0600 ENDDO
0601 ENDDO
0602 IF ( errorFlag ) STOP 'ABNORMAL END: S/R LAYERS_DIAPYCNAL'
0603 ENDIF
0604 #endif /* ndef TARGET_NEC_SX */
cf336ab6c5 Ryan* 0605
1377839540 Jean* 0606
cf336ab6c5 Ryan* 0607 DO j = 1,sNy +1
0608 DO i = 1,sNx +1
cdcbace143 Ryan* 0609 #ifdef LAYERS_FINEGRID_DIAPYCNAL
06e3ddabf4 Jean* 0610 dzfac = dZZf (kk ) * hFacC (i ,j ,kci ,bi ,bj )
cdcbace143 Ryan* 0611 #else
0612 dzfac = dRf (kci ) * hFacC (i ,j ,kci ,bi ,bj )
0613 #endif /* LAYERS_FINEGRID_DIAPYCNAL */
50d8304171 Ryan* 0614 kloc = kgcw (i ,j )
0615
0616
0617 Hcw (i ,j ,kloc ,bi ,bj ) = Hcw (i ,j ,kloc ,bi ,bj )
0618 & + dzfac
0619
0620 Hc (i ,j ,kgc (i ,j ),bi ,bj ) = Hc (i ,j ,kgc (i ,j ),bi ,bj )
0621 & + dzfac
06e3ddabf4 Jean* 0622
50d8304171 Ryan* 0623
0624 dzfac = dzfac * layers_recip_delta (kloc ,iLa )
0625
84f265d4b3 dfer* 0626 #ifdef LAYERS_PRHO_REF
50d8304171 Ryan* 0627 IF ( layers_num (iLa ) .EQ. 3 ) THEN
0628 Tfac = layers_alpha (i ,j ,kci ,bi ,bj )
0629 Sfac = layers_beta (i ,j ,kci ,bi ,bj )
0630 ENDIF
84f265d4b3 dfer* 0631 #endif
cf336ab6c5 Ryan* 0632 IF (kci .EQ. 1) THEN
0633
0634 TtendSurf (i ,j ,kloc ,bi ,bj ) =
06e3ddabf4 Jean* 0635 & TtendSurf (i ,j ,kloc ,bi ,bj ) +
50d8304171 Ryan* 0636 & Tfac * dzfac * layers_surfflux (i ,j ,1,1,bi ,bj )
cf336ab6c5 Ryan* 0637 StendSurf (i ,j ,kloc ,bi ,bj ) =
06e3ddabf4 Jean* 0638 & StendSurf (i ,j ,kloc ,bi ,bj ) +
50d8304171 Ryan* 0639 & Sfac * dzfac * layers_surfflux (i ,j ,1,2,bi ,bj )
cf336ab6c5 Ryan* 0640 ENDIF
50d8304171 Ryan* 0641
06e3ddabf4 Jean* 0642 #ifdef SHORTWAVE_HEATING
50d8304171 Ryan* 0643 TtendSurf (i ,j ,kloc ,bi ,bj ) =
0644 & TtendSurf (i ,j ,kloc ,bi ,bj ) +
0645 & Tfac * dzfac * layers_sw (i ,j ,kci ,1,bi ,bj )
0646 #endif /* SHORTWAVE_HEATING */
0647
0648
cf336ab6c5 Ryan* 0649 TtendDiffh (i ,j ,kloc ,bi ,bj ) =
50d8304171 Ryan* 0650 & TtendDiffh (i ,j ,kloc ,bi ,bj ) + dzfac * Tfac *
0651 & (layers_dfx (i ,j ,kci ,1,bi ,bj )+
0652 & layers_dfy (i ,j ,kci ,1,bi ,bj ))
cf336ab6c5 Ryan* 0653 TtendDiffr (i ,j ,kloc ,bi ,bj ) =
06e3ddabf4 Jean* 0654 & TtendDiffr (i ,j ,kloc ,bi ,bj ) +
50d8304171 Ryan* 0655 & dzfac * Tfac * layers_dfr (i ,j ,kci ,1,bi ,bj )
cf336ab6c5 Ryan* 0656 StendDiffh (i ,j ,kloc ,bi ,bj ) =
50d8304171 Ryan* 0657 & StendDiffh (i ,j ,kloc ,bi ,bj ) + dzfac * Sfac *
0658 & (layers_dfx (i ,j ,kci ,2,bi ,bj )+
0659 & layers_dfy (i ,j ,kci ,2,bi ,bj ))
cf336ab6c5 Ryan* 0660 StendDiffr (i ,j ,kloc ,bi ,bj ) =
06e3ddabf4 Jean* 0661 & StendDiffr (i ,j ,kloc ,bi ,bj ) +
50d8304171 Ryan* 0662 & dzfac * Sfac * layers_dfr (i ,j ,kci ,2,bi ,bj )
0663
0664 TtendAdvh (i ,j ,kloc ,bi ,bj ) =
0665 & TtendAdvh (i ,j ,kloc ,bi ,bj ) + dzfac * Tfac *
0666 & (layers_afx (i ,j ,kci ,1,bi ,bj )+
0667 & layers_afy (i ,j ,kci ,1,bi ,bj ))
0668 TtendAdvr (i ,j ,kloc ,bi ,bj ) =
06e3ddabf4 Jean* 0669 & TtendAdvr (i ,j ,kloc ,bi ,bj ) +
50d8304171 Ryan* 0670 & dzfac * Tfac * layers_afr (i ,j ,kci ,1,bi ,bj )
0671 StendAdvh (i ,j ,kloc ,bi ,bj ) =
0672 & StendAdvh (i ,j ,kloc ,bi ,bj ) + dzfac * Sfac *
0673 & (layers_afx (i ,j ,kci ,2,bi ,bj )+
0674 & layers_afy (i ,j ,kci ,2,bi ,bj ))
0675 StendAdvr (i ,j ,kloc ,bi ,bj ) =
06e3ddabf4 Jean* 0676 & StendAdvr (i ,j ,kloc ,bi ,bj ) +
50d8304171 Ryan* 0677 & dzfac * Sfac * layers_afr (i ,j ,kci ,2,bi ,bj )
0678
06e3ddabf4 Jean* 0679 Ttendtot (i ,j ,kloc ,bi ,bj ) =
0680 & Ttendtot (i ,j ,kloc ,bi ,bj ) +
50d8304171 Ryan* 0681 & dzfac * Tfac * layers_tottend (i ,j ,kci ,1,bi ,bj )
06e3ddabf4 Jean* 0682 Stendtot (i ,j ,kloc ,bi ,bj ) =
0683 & Stendtot (i ,j ,kloc ,bi ,bj ) +
50d8304171 Ryan* 0684 & dzfac * Sfac * layers_tottend (i ,j ,kci ,2,bi ,bj )
cf336ab6c5 Ryan* 0685 ENDDO
0686 ENDDO
0687 ENDDO
1377839540 Jean* 0688
cf336ab6c5 Ryan* 0689
0690
0691 DO kg =1,Nlayers
0692 DO j = 1,sNy +1
0693 DO i = 1,sNx +1
0694 IF (Hc (i ,j ,kg ,bi ,bj ) .GT. 0.) THEN
0695 PIc (i ,j ,kg ,bi ,bj ) = 1. _d 0
0696 ENDIF
0697 ENDDO
0698 ENDDO
0699 ENDDO
0700
0701
0702 ENDDO
0703 ENDDO
0704
0705 #endif /* LAYERS_THERMODYNAMICS */
0706
0707 RETURN
0708 END
0709
0710
9479dbe556 Mart* 0711
0712
28ff9bcfe6 Jean* 0713 SUBROUTINE LAYERS_LOCATE (
0714 I xx ,n ,m ,sNx ,sNy ,x ,
0715 O k ,
0716 I myThid )
9479dbe556 Mart* 0717
0718
0719
28ff9bcfe6 Jean* 0720
0721
9479dbe556 Mart* 0722
0723
0724
0725
0726
0727 IMPLICIT NONE
0728
28ff9bcfe6 Jean* 0729
9479dbe556 Mart* 0730
0731
28ff9bcfe6 Jean* 0732
9479dbe556 Mart* 0733
0734
28ff9bcfe6 Jean* 0735
0736 INTEGER n ,m ,sNx ,sNy
9479dbe556 Mart* 0737 _RL xx (1:n +1)
28ff9bcfe6 Jean* 0738 _RL x (snx +1,sny +1)
0739 INTEGER k (snx +1,sny +1)
0740 INTEGER myThid
9479dbe556 Mart* 0741
0742
0743
0744
28ff9bcfe6 Jean* 0745
0746 INTEGER i ,j
9479dbe556 Mart* 0747
0748 #ifdef TARGET_NEC_SX
28ff9bcfe6 Jean* 0749 INTEGER l , km
0750 INTEGER kl (sNx +1,sNy +1), ku (sNx +1,sNy +1)
0751
9479dbe556 Mart* 0752
0753
0754 DO j = 1,sNy +1
0755 DO i = 1,sNx +1
0756 kl (i ,j )=1
0757 ku (i ,j )=n +1
0758 END DO
0759 END DO
0760 DO l = 1,m
0761 DO j = 1,sNy +1
0762 DO i = 1,sNx +1
0763 IF (ku (i ,j )-kl (i ,j ).GT. 1) THEN
0764 km =(ku (i ,j )+kl (i ,j ))/2
0765
0766 IF ( ((xx (n ).GE. xx (1)).AND. (x (i ,j ).GE. xx (km ))).OR.
0767 & ((xx (n ).GE. xx (1)).AND. (x (i ,j ).GE. xx (km ))) ) THEN
0768 kl (i ,j )=km
0769 ELSE
0770 ku (i ,j )=km
0771 END IF
0772 END IF
0773 END DO
0774 END DO
0775 END DO
0776 DO j = 1,sNy +1
0777 DO i = 1,sNx +1
0778 IF ( x (i ,j ).LT. xx (2) ) THEN
0779 k (i ,j )=1
0780 ELSE IF ( x (i ,j ).GE. xx (n ) ) THEN
0781 k (i ,j )=n
0782 ELSE
0783 k (i ,j )=kl (i ,j )
0784 END IF
0785 END DO
0786 END DO
0787 #else
0788
0789 DO j = 1,sNy +1
0790 DO i = 1,sNx +1
0791 IF (x (i ,j ) .GE. xx (n )) THEN
0792
0793 k (i ,j ) = n
0794 ELSE IF (x (i ,j ) .LT. xx (2)) THEN
0795
0796 k (i ,j ) = 1
0797 ELSE IF ( (x (i ,j ) .GE. xx (k (i ,j )))
0798 & .AND. (x (i ,j ) .LT. xx (k (i ,j )+1)) ) THEN
0799
0800 ELSE IF (x (i ,j ) .GE. xx (k (i ,j ))) THEN
0801
0802 DO WHILE (x (i ,j ) .GE. xx (k (i ,j )+1))
0803 k (i ,j ) = k (i ,j ) + 1
0804 ENDDO
0805
0806 ELSE IF (x (i ,j ) .LT. xx (k (i ,j )+1)) THEN
0807
0808 DO WHILE (x (i ,j ) .LT. xx (k (i ,j )))
0809 k (i ,j ) = k (i ,j ) - 1
0810 ENDDO
0811
0812 ELSE
0813
0814 k (i ,j ) = -1
0815 ENDIF
28ff9bcfe6 Jean* 0816
9479dbe556 Mart* 0817 ENDDO
0818 ENDDO
0819 #endif /* TARGET_NEC_SX */
0820
aa668e01f1 Gael* 0821 #endif /* ALLOW_LAYERS */
0822
0823 RETURN
0824 END
cf336ab6c5 Ryan* 0825