Changeset 216 for trunk/NEMO/OPA_SRC/TRD
- Timestamp:
- 2005-03-17T15:02:38+01:00 (19 years ago)
- Location:
- trunk/NEMO/OPA_SRC/TRD
- Files:
-
- 2 edited
-
trdmld.F90 (modified) (20 diffs)
-
trdvor.F90 (modified) (18 diffs)
Legend:
- Unmodified
- Added
- Removed
-
trunk/NEMO/OPA_SRC/TRD/trdmld.F90
r101 r216 8 8 !! 'key_trdmld' mixed layer trend diagnostics 9 9 !!---------------------------------------------------------------------- 10 !! trd_mld : T and S trends averaged over the mixed layer 10 !! trd_mld : T and S cumulated trends averaged over the mixed layer 11 !! trd_mld_zint : T and S trends vertical integration 12 !! trd_mld_init : initialization step 11 13 !!---------------------------------------------------------------------- 12 14 !! * Modules used 13 15 USE oce ! ocean dynamics and tracers variables 14 16 USE dom_oce ! ocean space and time domain variables 15 USE ldftra_oce ! ocean active tracers: lateral physics16 USE trdtra_oce ! ocean active tracer trend variables17 USE trdmod_oce ! ocean variables trends 18 USE ldftra_oce ! ocean active tracers lateral physics 17 19 USE zdf_oce ! ocean vertical physics 18 20 USE in_out_manager ! I/O manager 19 20 21 USE phycst ! Define parameters for the routines 21 22 USE daymod ! calendar 22 23 USE dianam ! build the name of file (routine) 23 USE ldfslp ! iso-neutral slopes 24 USE zdfmxl 25 USE zdfddm 26 USE ioipsl 27 USE lbclnk 28 #if defined key_dimgout 29 USE diawri, ONLY : dia_wri_dimg 30 #endif 24 USE ldfslp ! iso-neutral slopes 25 USE zdfmxl ! mixed layer depth 26 USE zdfddm ! ocean vertical physics: double diffusion 27 USE ioipsl ! NetCDF library 28 USE lbclnk ! ocean lateral boundary conditions (or mpp link) 29 USE diadimg ! dimg direct access file format output 31 30 32 31 IMPLICIT NONE … … 34 33 35 34 !! * Accessibility 36 PUBLIC trd_mld ! routine called by step.F90 35 PUBLIC trd_mld ! routine called by step.F90 36 PUBLIC trd_mld_init ! routine called by opa.F90 37 PUBLIC trd_mld_zint ! routine called by tracers routines 37 38 38 39 !! * Shared module variables … … 40 41 41 42 !! * Module variables 43 INTEGER :: & 44 nh_t, nmoymltrd, & ! ??? 45 nidtrd, & 46 ndextrd1(jpi*jpj), & 47 ndimtrd1 48 INTEGER, SAVE :: & 49 ionce, icount, & 50 idebug ! (0/1) set it to 1 in case of problem to have more print 51 42 52 INTEGER, DIMENSION(jpi,jpj) :: & 43 nmld, & ! mixed layer depth 44 nbol 45 46 INTEGER :: & 47 nh_t, nmoymltrd, & ! ??? 48 nidtrd,nhoridtrd, & 49 ndextrd1(jpi*jpj),ndimtrd1 53 nmld, & ! mixed layer depth 54 nbol 50 55 51 56 REAL(wp), DIMENSION(jpi,jpj) :: & … … 59 64 tmltrdm, smltrdm ! 60 65 61 ! REAL(wp), DIMENSION(jpi,jpj,jpltrd) :: & ! Must be jpk for mpp lbc_lnk62 ! TO BE FIXED ???63 66 REAL(wp), DIMENSION(jpi,jpj,jpk) :: & 64 67 tmltrd , & ! total cumulative trends of temperature and 65 smltrd 66 67 INTEGER :: iyear,imon,iday 68 CHARACTER(LEN=80) :: clname , cltext, clmode68 smltrd , & ! salinity over nwrite-1 time steps 69 wkx 70 71 CHARACTER(LEN=80) :: clname 69 72 70 73 !! * Substitutions … … 78 81 CONTAINS 79 82 80 SUBROUTINE trd_mld( kt)83 SUBROUTINE trd_mld_zint( pttrdmld, pstrdmld, ktrd, ctype ) 81 84 !!---------------------------------------------------------------------- 82 !! *** ROUTINE trd_mld ***85 !! *** ROUTINE trd_mld_zint *** 83 86 !! 84 87 !! ** Purpose : computation of vertically integrated T and S budgets 85 !! from ocean surface down to control surface (NetCDF output)88 !! from ocean surface down to control surface 86 89 !! 87 90 !! ** Method/usage : … … 134 137 !! ! 99-09 (E. Guilyardi) Re-writing + netCDF output 135 138 !! 8.5 ! 02-06 (G. Madec) F90: Free form and module 139 !! 9.0 ! 04-08 (C. Talandier) New trends organization 136 140 !!---------------------------------------------------------------------- 137 141 !! * Arguments 138 INTEGER, INTENT( in ) :: kt ! ocean time-step index 142 INTEGER, INTENT( in ) :: ktrd ! ocean trend index 143 144 CHARACTER(len=2), INTENT( in ) :: & 145 ctype ! surface/bottom (2D arrays) or 146 ! interior (3D arrays) physics 147 148 REAL(wp), DIMENSION(jpi,jpj,jpk), INTENT( in ) :: & 149 pttrdmld, & ! Temperature trend 150 pstrdmld ! Salinity trend 139 151 140 152 !! * Local declarations 141 INTEGER :: ilseq 142 INTEGER :: ji, jj, jk, jl, ik, ikb, idebug, isum, it 143 INTEGER, DIMENSION(jpi,jpj) :: zvlmsk 144 145 REAL(wp) :: zmean, zavt, zjulian, zsto, zout 146 REAL(wp) ,DIMENSION(jpi,jpj,jpktrd) :: zwkx 147 REAL(wp) ,DIMENSION(jpi,jpj) :: & 148 & ztmltot, ztmlres, z2d, zsmltot, zsmlres 149 150 CHARACTER (len=21) :: & 151 clold ='OLD' , & ! open specifier (direct access files) 152 clunf ='UNFORMATTED', & ! open specifier (direct access files) 153 clseq ='SEQUENTIAL' ! open specifier (direct access files) 154 CHARACTER (len=80) :: clname 155 CHARACTER (len=40) :: clhstnam 156 CHARACTER (len=40) :: clop 157 CHARACTER (len=12) :: clmxl 158 159 NAMELIST/namtrd/ ntrd, nctls 153 INTEGER :: ji, jj, jk, isum 154 # if defined key_trabbl_dif 155 INTEGER :: ikb 156 # endif 157 158 REAL(wp), DIMENSION(jpi,jpj) :: & 159 zvlmsk 160 160 !!---------------------------------------------------------------------- 161 161 162 ! =================== 163 ! 0. initialization 164 ! =================== 165 166 ! Open specifier 167 ilseq = 1 168 idebug = 0 ! set it to 1 in case of problem to have more print 169 170 IF( kt == nit000 ) THEN 171 ! namelist namtrd : trend diagnostic 172 REWIND( numnam ) 173 READ ( numnam, namtrd ) 174 175 IF(lwp) THEN 176 WRITE(numout,*) 'namtrd' 177 WRITE(numout,*) ' ' 178 WRITE(numout,*) ' time step frequency trend ntrd = ',ntrd 179 WRITE(numout,*) ' control surface for trends nctls = ',nctls 180 WRITE(numout,*) ' ' 181 ENDIF 182 183 ! cumulated trends array init 184 nmoymltrd = 0 185 tmltrdm(:,:) = 0. 186 smltrdm(:,:) = 0. 187 ENDIF 188 189 ! set before values of vertically average T and S 190 191 IF( kt > nit000 ) THEN 192 tmlb(:,:) = tml(:,:) 193 smlb(:,:) = sml(:,:) 194 ENDIF 195 196 ! read control surface from file ctlsurf_idx 197 198 IF( kt == nit000 .and. nctls == - 1 ) THEN 199 clname ='ctlsurf_idx' 200 CALL ctlopn(numbol,clname,clold,clunf,clseq, & 201 ilseq,numout,lwp,1) 202 REWIND (numbol) 203 READ(numbol) nbol 204 ENDIF 205 206 IF( idebug /= 0 ) THEN 207 WRITE(numout,*) ' debuging trd_mld: 0. done ' 208 CALL FLUSH(numout) 209 ENDIF 210 211 212 ! ======================================================== 213 ! I. definition of control surface and associated fields 214 ! ======================================================== 215 216 ! I.1 set nmld(ji,jj) = index of first T point below control surface 217 ! ------------------- or outside mixed-layer 218 219 ! clmxl = legend root for netCDF output 220 221 IF( nctls == 0 ) THEN 222 ! control surface = mixed-layer with density criterion 223 ! (array nmln computed in zdfmxl.F90) 224 nmld(:,:) = nmln(:,:) 225 clmxl = 'Mixed Layer ' 226 ELSE IF( nctls == 1 ) THEN 227 ! control surface = read index from file 228 nmld(:,:) = nbol(:,:) 229 clmxl = ' Bowl ' 230 ELSE IF( nctls >= 2 ) THEN 231 ! control surface = model level 232 nctls = MIN( nctls, jpktrd - 1 ) 233 nmld(:,:) = nctls + 1 234 WRITE(clmxl,'(A9,I2,1X)') 'Levels 1-', nctls 235 ENDIF 236 237 ! Check of validity : nmld(ji,jj) =< jpktrd 238 isum = 0 239 240 IF( jpktrd < jpk ) THEN 241 DO jj = 1, jpj 242 DO ji = 1, jpi 243 IF( nmld(ji,jj) <= jpktrd ) THEN 244 zvlmsk(ji,jj) = tmask(ji,jj,1) 245 ELSE 246 isum = isum + 1 247 zvlmsk(ji,jj) = 0. 248 ENDIF 162 IF( icount == 1 ) THEN 163 164 zvlmsk(:,:) = 0.e0 165 tmltrd(:,:,:) = 0.e0 166 smltrd(:,:,:) = 0.e0 167 168 ! This computation should be done only once per time step 169 170 ! ======================================================== 171 ! I. definition of control surface and associated fields 172 ! ======================================================== 173 174 ! I.1 set nmld(ji,jj) = index of first T point below control surface 175 ! ------------------- or outside mixed-layer 176 177 IF( nctls == 0 ) THEN 178 ! control surface = mixed-layer with density criterion 179 ! (array nmln computed in zdfmxl.F90) 180 nmld(:,:) = nmln(:,:) 181 ELSE IF( nctls == 1 ) THEN 182 ! control surface = read index from file 183 nmld(:,:) = nbol(:,:) 184 ELSE IF( nctls >= 2 ) THEN 185 ! control surface = model level 186 nctls = MIN( nctls, jpktrd - 1 ) 187 nmld(:,:) = nctls + 1 188 ENDIF 189 190 IF( ionce == 1 ) THEN ! compute ndextrd1 and ndimtrd1 only once 191 ! Check of validity : nmld(ji,jj) =< jpktrd 192 isum = 0 193 194 IF( jpktrd < jpk ) THEN 195 DO jj = 1, jpj 196 DO ji = 1, jpi 197 IF( nmld(ji,jj) <= jpktrd ) THEN 198 zvlmsk(ji,jj) = tmask(ji,jj,1) 199 ELSE 200 isum = isum + 1 201 zvlmsk(ji,jj) = 0. 202 ENDIF 203 END DO 204 END DO 205 ENDIF 206 207 ! Index of ocean points (2D only) 208 IF( isum > 0 ) THEN 209 WRITE(numout,*)' Number of invalid points nmld > jpktrd', isum 210 CALL wheneq( jpi*jpj, zvlmsk(:,:) , 1, 1., ndextrd1, ndimtrd1 ) ! surface 211 ELSE 212 CALL wheneq( jpi*jpj, tmask(:,:,1), 1, 1., ndextrd1, ndimtrd1 ) ! surface 213 ENDIF 214 215 ! no more pass here 216 ionce = 0 217 218 ENDIF 219 220 IF( idebug /= 0 ) THEN 221 ! CALL prihre (zvlmsk,jpi,jpj,1,jpi,2,1,jpj,2,3,numout) 222 WRITE(numout,*) ' debuging trd_mld_zint: I.1 done ' 223 CALL FLUSH(numout) 224 ENDIF 225 226 227 ! I.2 probability density function of presence in mixed-layer 228 ! -------------------------------- 229 ! (i.e. weight of each grid point in vertical integration : wkx(ji,jj,jk) 230 231 232 ! initialize wkx with vertical scale factor in mixed-layer 233 234 wkx(:,:,:) = 0.e0 235 DO jk = 1, jpktrd 236 DO jj = 1,jpj 237 DO ji = 1,jpi 238 IF( jk - nmld(ji,jj) < 0. ) wkx(ji,jj,jk) = fse3t(ji,jj,jk) * tmask(ji,jj,jk) 239 END DO 249 240 END DO 250 241 END DO 251 ENDIF 252 242 243 ! compute mixed-layer depth : rmld 244 245 rmld(:,:) = 0. 246 DO jk = 1, jpktrd 247 rmld(:,:) = rmld(:,:) + wkx(:,:,jk) 248 END DO 249 250 ! compute PDF 251 252 DO jk = 1, jpktrd 253 wkx(:,:,jk) = wkx(:,:,jk) / MAX( 1., rmld(:,:) ) 254 END DO 255 256 IF( idebug /= 0 ) THEN 257 WRITE(numout,*) ' debuging trd_mld_zint: I.2 done ' 258 CALL FLUSH(numout) 259 ENDIF 260 261 ! Set counter icount to 0 to avoid this part at each time step 262 icount = 0 263 264 ENDIF 265 266 267 ! ==================================================== 268 ! II. vertical integration of trends in mixed-layer 269 ! ==================================================== 270 271 ! II.1 vertical integration of 3D and 2D trends 272 ! --------------------------------------------- 273 274 SELECT CASE (ctype) 275 276 CASE ('3D') ! 3D treatment 277 278 ! trends terms in the mixed-layer 279 DO jk = 1, jpktrd 280 ! Temperature 281 tmltrd(:,:,ktrd) = tmltrd(:,:,ktrd) + pttrdmld(:,:,jk) * wkx(:,:,jk) 282 283 ! Salinity 284 smltrd(:,:,ktrd) = smltrd(:,:,ktrd) + pstrdmld(:,:,jk) * wkx(:,:,jk) 285 ENDDO 286 287 CASE ('2D') ! 2D treatment 288 289 SELECT CASE (ktrd) 290 291 CASE (jpmldldf) 292 293 # if defined key_trabbl_dif 294 ! trends terms from Beckman over-flow parameterization 295 DO jj = 1,jpj 296 DO ji = 1,jpi 297 ikb = MAX( mbathy(ji,jj)-1, 1 ) 298 ! beckmann component -> horiz. part of lateral diffusion 299 tmltrd(ji,jj,ktrd) = tmltrd(ji,jj,ktrd) + pttrdmld(ji,jj,1) * wkx(ji,jj,ikb) 300 smltrd(ji,jj,ktrd) = smltrd(ji,jj,ktrd) + pstrdmld(ji,jj,1) * wkx(ji,jj,ikb) 301 END DO 302 END DO 303 # endif 304 305 CASE DEFAULT 306 307 ! trends terms at upper boundary of mixed-layer 308 309 ! forcing term (non penetrative) 310 ! Temperature 311 tmltrd(:,:,ktrd) = tmltrd(:,:,ktrd) + pttrdmld(:,:,1) * wkx(:,:,1) 312 313 ! forcing term 314 ! Salinity 315 smltrd(:,:,ktrd) = smltrd(:,:,ktrd) + pstrdmld(:,:,1) * wkx(:,:,1) 316 317 END SELECT 318 319 END SELECT 320 253 321 IF( idebug /= 0 ) THEN 254 ! CALL prihre (zvlmsk,jpi,jpj,1,jpi,2,1,jpj,2,3,numout) 255 WRITE(numout,*) ' debuging trd_mld: I.1 done ' 322 IF(lwp) WRITE(numout,*) ' debuging trd_mld_zint: II.1 done' 256 323 CALL FLUSH(numout) 257 324 ENDIF 258 325 259 260 ! I.2 probability density function of presence in mixed-layer 261 ! -------------------------------- 262 ! (i.e. weight of each grid point in vertical integration : zwkx(ji,jj,jk) 263 264 265 ! initialize zwkx with vertical scale factor in mixed-layer 266 267 zwkx(:,:,:) = 0.e0 268 DO jk = 1, jpktrd 269 DO jj = 1,jpj 270 DO ji = 1,jpi 271 IF( jk - nmld(ji,jj) < 0. ) zwkx(ji,jj,jk) = fse3t(ji,jj,jk) * tmask(ji,jj,jk) 272 END DO 273 END DO 274 END DO 275 276 ! compute mixed-layer depth : rmld 277 278 rmld(:,:) = 0. 279 DO jk = 1, jpktrd 280 rmld(:,:) = rmld(:,:) + zwkx(:,:,jk) 281 END DO 282 283 ! compute PDF 284 285 DO jk = 1, jpktrd 286 zwkx(:,:,jk) = zwkx(:,:,jk) / MAX( 1., rmld(:,:) ) 287 END DO 288 289 IF( idebug /= 0 ) THEN 290 WRITE(numout,*) ' debuging trd_mld: I.2 done ' 291 CALL FLUSH(numout) 292 ENDIF 293 294 295 ! I.3 vertically integrated T and S 296 ! --------------------------------- 297 298 tml(:,:) = 0. 299 sml(:,:) = 0. 300 301 DO jk = 1, jpktrd - 1 302 tml(:,:) = tml(:,:) + zwkx(:,:,jk) * tn(:,:,jk) 303 sml(:,:) = sml(:,:) + zwkx(:,:,jk) * sn(:,:,jk) 304 END DO 305 306 IF(idebug /= 0) THEN 307 WRITE(numout,*) ' debuging trd_mld: I.3 done' 308 CALL FLUSH(numout) 309 ENDIF 310 311 312 ! =================================== 313 ! II. netCDF output initialization 314 ! =================================== 315 316 #if defined key_dimgout 317 318 #else 319 #include "trdmld_ncinit.h90" 326 END SUBROUTINE trd_mld_zint 327 328 329 330 SUBROUTINE trd_mld( kt ) 331 !!---------------------------------------------------------------------- 332 !! *** ROUTINE trd_mld *** 333 !! 334 !! ** Purpose : computation of cumulated trends over analysis period 335 !! and make outputs (NetCDF or DIMG format) 336 !! 337 !! ** Method/usage : 338 !! 339 !! History : 340 !! 9.0 ! 04-08 (C. Talandier) New trends organization 341 !!---------------------------------------------------------------------- 342 !! * Arguments 343 INTEGER, INTENT( in ) :: kt ! ocean time-step index 344 345 !! * Local declarations 346 INTEGER :: ji, jj, jk, jl, ik, it 347 348 REAL(wp) :: zmean, zavt 349 350 REAL(wp) ,DIMENSION(jpi,jpj) :: & 351 ztmltot, ztmlres, & 352 zsmltot, zsmlres, & 353 z2d 354 355 #if defined key_dimgout 356 INTEGER :: iyear,imon,iday 357 CHARACTER(LEN=80) :: cltext, clmode 320 358 #endif 321 322 IF( idebug /= 0 ) THEN 323 WRITE(numout,*) ' debuging trd_mld: II. done' 324 CALL FLUSH(numout) 325 ENDIF 326 327 ! ==================================================== 328 ! III. vertical integration of trends in mixed-layer 329 ! ==================================================== 330 331 332 ! III.0 initializations 333 ! --------------------- 334 335 tmltrd(:,:,:) = 0.e0 336 smltrd(:,:,:) = 0.e0 337 338 IF( idebug /= 0 ) THEN 339 WRITE(numout,*) ' debuging trd_mld: III.0 done' 340 CALL FLUSH(numout) 341 ENDIF 342 343 344 ! III.1 vertical integration of 3D trends 345 ! --------------------------------------- 346 347 DO jk = 1,jpktrd 348 349 ! Temperature 350 tmltrd(:,:,1) = tmltrd(:,:,1) + ttrdh(:,:,jk,1) * zwkx(:,:,jk) ! zonal advection 351 tmltrd(:,:,2) = tmltrd(:,:,2) + ttrdh(:,:,jk,2) * zwkx(:,:,jk) ! meridional advection 352 tmltrd(:,:,3) = tmltrd(:,:,3) + ttrd (:,:,jk,2) * zwkx(:,:,jk) ! vertical advection 353 tmltrd(:,:,4) = tmltrd(:,:,4) + ttrd (:,:,jk,3) * zwkx(:,:,jk) ! lateral diffusion (hor. part) 354 tmltrd(:,:,5) = tmltrd(:,:,5) + ttrd (:,:,jk,7) * zwkx(:,:,jk) ! forcing (penetrative) 355 IF( l_traldf_iso ) THEN 356 tmltrd(:,:,7) = tmltrd(:,:,7) + ttrd (:,:,jk,4) * zwkx(:,:,jk) ! lateral diffusion (explicit 357 ! ! vert. part (isopycnal diff.) 358 ENDIF 359 !#if defined key_traldf_eiv 360 ! tmltrd(:,:,8 ) = tmltrdg(:,:,8) + ttrdh(:,:,jk,3) * zwkx(:,:,jk) ! eddy induced zonal advection 361 ! tmltrd(:,:,9 ) = tmltrdg(:,:,9) + ttrdh(:,:,jk,4) * zwkx(:,:,jk) ! eddy induced merid. advection 362 ! tmltrd(:,:,10) = tmltrdg(:,:,10) + ttrd(:,:,jk,6) * zwkx(:,:,jk) ! eddy induced vert. advection 363 !#endif 364 365 ! Salinity 366 smltrd(:,:,1) = smltrd(:,:,1) + strdh(:,:,jk,1) * zwkx(:,:,jk) ! zonal advection 367 smltrd(:,:,2) = smltrd(:,:,2) + strdh(:,:,jk,2) * zwkx(:,:,jk) ! meridional advection 368 smltrd(:,:,3) = smltrd(:,:,3) + strd (:,:,jk,2) * zwkx(:,:,jk) ! vertical advection 369 smltrd(:,:,4) = smltrd(:,:,4) + strd (:,:,jk,3) * zwkx(:,:,jk) ! lateral diffusion (hor. part) 370 IF( l_traldf_iso ) THEN 371 smltrd(:,:,7) = smltrd(:,:,7) + strd (:,:,jk,4) * zwkx(:,:,jk) ! lateral diffusion (explicit 372 ! ! vert. part (isopycnal diff.) 373 ENDIF 374 !#if defined key_traldf_eiv 375 ! smltrd(:,:,8) = smltrdg(:,:,8) + strdh(:,:,jk,3) * zwkx(:,:,jk) ! eddy induced zonal advection 376 ! smltrd(:,:,9) = smltrdg(:,:,9) + strdh(:,:,jk,4) * zwkx(:,:,jk) ! eddy induced merid. advection 377 ! smltrd(:,:,10) = smltrdg(:,:,10) + strd(:,:,jk,6) * zwkx(:,:,jk) ! eddy induced vert. advection 378 !#endif 379 380 END DO 381 382 383 IF( idebug /= 0 ) THEN 384 IF(lwp) WRITE(numout,*) ' debuging trd_mld: III.1 done' 385 CALL FLUSH(numout) 386 ENDIF 387 388 389 ! III.2 trends terms at upper and lower boundaries of mixed-layer 390 ! --------------------------------------------------------------- 359 !!---------------------------------------------------------------------- 360 361 ! I. trends terms at lower boundary of mixed-layer 362 ! ------------------------------------------------ 391 363 392 364 DO jj = 1,jpj … … 396 368 397 369 ! Temperature 398 399 ! forcing (non penetrative)400 tmltrd(ji,jj,5) = tmltrd(ji,jj,5) + flxtrd(ji,jj,1) * zwkx(ji,jj,1)401 370 ! entrainment due to vertical diffusion 402 371 ! - due to vertical mixing scheme (TKE) 403 372 zavt = avt(ji,jj,ik) 404 tmltrd(ji,jj,6) = - 1. * zavt / fse3w(ji,jj,ik) * tmask(ji,jj,ik) & 405 * ( tn(ji,jj,ik-1) - tn(ji,jj,ik) ) & 406 / MAX( 1., rmld(ji,jj) ) * tmask(ji,jj,1) 407 373 tmltrd(ji,jj,jpmldevd) = - 1. * zavt / fse3w(ji,jj,ik) * tmask(ji,jj,ik) & 374 & * ( tn(ji,jj,ik-1) - tn(ji,jj,ik) ) & 375 & / MAX( 1., rmld(ji,jj) ) * tmask(ji,jj,1) 408 376 ! Salinity 409 410 ! forcing411 smltrd(ji,jj,5) = flxtrd(ji,jj,2) * zwkx(ji,jj,1)412 377 ! entrainment due to vertical diffusion 413 378 ! - due to vertical mixing scheme (TKE) 414 379 zavt = fsavs(ji,jj,ik) 415 smltrd(ji,jj, 6) = -1. * zavt / fse3w(ji,jj,ik) * tmask(ji,jj,ik) &380 smltrd(ji,jj,jpmldevd) = -1. * zavt / fse3w(ji,jj,ik) * tmask(ji,jj,ik) & 416 381 & * ( sn(ji,jj,ik-1) - sn(ji,jj,ik) ) & 417 382 & / MAX( 1., rmld(ji,jj) ) * tmask(ji,jj,1) … … 420 385 421 386 IF( l_traldf_iso ) THEN 422 !!Clem On retire de la diffusion verticale TOTALE calculee par tmltrd(:,:,7) 423 !!Clem ds trazdf.isopycnal et implicit, la partie verticale due au Kz afin de ne garder 424 !!Clem effectivement que la diffusion verticale isopycnale (ie composante de la 425 !!Clem diff isopycnale sur la verticale):426 tmltrd(:,:,7) = tmltrd(:,:,7) - tmltrd(:,:,6) ! - due to isopycnal mixing scheme (implicit part)427 smltrd(:,:,7) = smltrd(:,:,7) - smltrd(:,:,6) ! - due to isopycnal mixing scheme (implicit part)387 ! We substract to the TOTAL vertical diffusion tmltrd(:,:,jpmldzdf) 388 ! computed in subroutines trazdf_iso.F90 or trazdf_imp.F90 389 ! the vertical part du to the Kz in order to keep only the vertical 390 ! isopycnal diffusion (i.e the isopycnal diffusion componant on the vertical): 391 tmltrd(:,:,jpmldzdf) = tmltrd(:,:,jpmldzdf) - tmltrd(:,:,jpmldevd) ! - due to isopycnal mixing scheme (implicit part) 392 smltrd(:,:,jpmldzdf) = smltrd(:,:,jpmldzdf) - smltrd(:,:,jpmldevd) ! - due to isopycnal mixing scheme (implicit part) 428 393 ENDIF 429 394 … … 435 400 436 401 IF( idebug /= 0 ) THEN 437 WRITE(numout,*) ' debuging trd_mld: I II.2done'402 WRITE(numout,*) ' debuging trd_mld: I. done' 438 403 CALL FLUSH(numout) 439 404 ENDIF 440 405 441 #if defined key_trabbl_dif 442 ! III.3 trends terms from beckman over-flow parameterization 443 ! ---------------------------------------------------------- 444 445 DO jj = 1,jpj 446 DO ji = 1,jpi 447 ikb = MAX( mbathy(ji,jj)-1, 1 ) 448 ! beckmann component -> horiz. part of lateral diffusion 449 tmltrd(ji,jj,4) = tmltrd(ji,jj,4) + bbltrd(ji,jj,1) * zwkx(ji,jj,ikb) 450 smltrd(ji,jj,4) = smltrd(ji,jj,4) + bbltrd(ji,jj,2) * zwkx(ji,jj,ikb) 451 END DO 406 ! ================================= 407 ! II. Cumulated trends 408 ! ================================= 409 410 ! II.1 set before values of vertically average T and S 411 ! --------------------------------------------------- 412 413 IF( kt > nit000 ) THEN 414 tmlb(:,:) = tml(:,:) 415 smlb(:,:) = sml(:,:) 416 ENDIF 417 418 ! II.2 vertically integrated T and S 419 ! --------------------------------- 420 421 tml(:,:) = 0. 422 sml(:,:) = 0. 423 424 DO jk = 1, jpktrd - 1 425 tml(:,:) = tml(:,:) + wkx(:,:,jk) * tn(:,:,jk) 426 sml(:,:) = sml(:,:) + wkx(:,:,jk) * sn(:,:,jk) 452 427 END DO 453 454 #endif 455 456 IF( idebug /= 0 ) THEN 457 WRITE(numout,*) ' debuging trd_mld: III.3 done' 428 429 IF(idebug /= 0) THEN 430 WRITE(numout,*) ' debuging trd_mld: II.2 done' 458 431 CALL FLUSH(numout) 459 432 ENDIF 460 433 461 462 ! ================================= 463 ! IV. Cumulated trends 464 ! ================================= 465 466 467 468 ! IV.1 set `before' mixed layer values for kt = nit000+1 434 ! II.3 set `before' mixed layer values for kt = nit000+1 469 435 ! -------------------------------------------------------- 470 436 … … 477 443 478 444 IF( idebug /= 0 ) THEN 479 WRITE(numout,*) ' debuging trd_mld: I V.1done'445 WRITE(numout,*) ' debuging trd_mld: II.3 done' 480 446 CALL FLUSH(numout) 481 447 ENDIF 482 448 483 484 ! IV.2 cumulated trends over analysis period (kt=2 to nwrite) 485 ! ---------------------- 449 ! II.4 cumulated trends over analysis period (kt=2 to nwrite) 450 ! ----------------------------------------------------------- 486 451 487 452 ! trends cumulated over nwrite-2 time steps … … 496 461 497 462 IF( idebug /= 0 ) THEN 498 WRITE(numout,*) ' debuging trd_mld: I V.2done'463 WRITE(numout,*) ' debuging trd_mld: II.4 done' 499 464 CALL FLUSH(numout) 500 465 ENDIF 501 466 502 503 467 ! ============================================= 504 ! V. Output in netCDF + residual computation468 ! III. Output in netCDF + residual computation 505 469 ! ============================================= 506 470 … … 512 476 IF( MOD( kt - nit000+1, nwrite ) == 0 ) THEN 513 477 514 ! V.1 compute total trend478 ! III.1 compute total trend 515 479 ! ------------------------ 516 480 … … 522 486 IF(idebug /= 0) THEN 523 487 WRITE(numout,*) ' zmean = ',zmean 524 WRITE(numout,*) ' debuging trd_mld: V.1 done'488 WRITE(numout,*) ' debuging trd_mld: III.1 done' 525 489 CALL FLUSH(numout) 526 490 ENDIF 527 491 528 492 529 ! V.2 compute residual493 ! III.2 compute residual 530 494 ! --------------------- 531 495 … … 542 506 543 507 IF( idebug /= 0 ) THEN 544 WRITE(numout,*) ' debuging trd_mld: V.2 done'508 WRITE(numout,*) ' debuging trd_mld: III.2 done' 545 509 CALL FLUSH(numout) 546 510 ENDIF 547 511 548 512 549 ! V.3 time evolution array swap513 ! III.3 time evolution array swap 550 514 ! ------------------------------ 551 515 … … 556 520 557 521 IF( idebug /= 0 ) THEN 558 WRITE(numout,*) ' debuging trd_mld: V.3 done'522 WRITE(numout,*) ' debuging trd_mld: III.3 done' 559 523 CALL FLUSH(numout) 560 524 ENDIF 561 525 562 526 563 ! V.4 zero cumulative array527 ! III.4 zero cumulative array 564 528 ! --------------------------- 565 529 … … 570 534 571 535 IF(idebug /= 0) THEN 572 WRITE(numout,*) ' debuging trd_mld: I V.4 done'536 WRITE(numout,*) ' debuging trd_mld: III.4 done' 573 537 CALL FLUSH(numout) 574 538 ENDIF … … 576 540 ENDIF 577 541 578 ! I V.5 write trends to output542 ! III.5 write trends to output 579 543 ! --------------------------- 580 544 … … 595 559 IF( kt >= nit000+1 ) THEN 596 560 597 #include "trdmld_ncwrite.h90" 561 ! define time axis 562 it= kt-nit000+1 563 IF( lwp .AND. MOD( kt, nwrite ) == 0 ) THEN 564 WRITE(numout,*) ' trd_mld : write NetCDF fields' 565 ENDIF 566 567 CALL histwrite( nidtrd,"somlttml",it,rmld ,ndimtrd1,ndextrd1) ! Mixed-layer depth 568 569 ! Temperature trends 570 ! ------------------ 571 CALL histwrite( nidtrd,"somltemp",it,tml ,ndimtrd1,ndextrd1) ! Mixed-layer temperature 572 CALL histwrite( nidtrd,"somlttto",it,ztmltot ,ndimtrd1,ndextrd1) ! total 573 CALL histwrite( nidtrd,"somlttax",it,tmltrd(:,:, 1),ndimtrd1,ndextrd1) ! i- adv. 574 CALL histwrite( nidtrd,"somlttay",it,tmltrd(:,:, 2),ndimtrd1,ndextrd1) ! j- adv. 575 CALL histwrite( nidtrd,"somlttaz",it,tmltrd(:,:, 3),ndimtrd1,ndextrd1) ! vertical adv. 576 CALL histwrite( nidtrd,"somlttdh",it,tmltrd(:,:, 4),ndimtrd1,ndextrd1) ! hor. lateral diff. 577 CALL histwrite( nidtrd,"somlttfo",it,tmltrd(:,:, 5),ndimtrd1,ndextrd1) ! forcing 578 579 CALL histwrite( nidtrd,"somlbtdz",it,tmltrd(:,:, 6),ndimtrd1,ndextrd1) ! vert. diffusion 580 CALL histwrite( nidtrd,"somlbtdt",it,ztmlres ,ndimtrd1,ndextrd1) ! dh/dt entrainment (residual) 581 IF( l_traldf_iso ) THEN 582 CALL histwrite( nidtrd,"somlbtdv",it,tmltrd(:,:, 7),ndimtrd1,ndextrd1) ! vert. lateral diff. 583 ENDIF 584 #if defined key_traldf_eiv 585 CALL histwrite( nidtrd,"somlgtax",it,tmltrd(:,:, 8),ndimtrd1,ndextrd1) ! i- adv. (eiv) 586 CALL histwrite( nidtrd,"somlgtay",it,tmltrd(:,:, 9),ndimtrd1,ndextrd1) ! j- adv. (eiv) 587 CALL histwrite( nidtrd,"somlgtaz",it,tmltrd(:,:,10),ndimtrd1,ndextrd1) ! vert. adv. (eiv) 588 z2d(:,:) = tmltrd(:,:,8) + tmltrd(:,:,9) + tmltrd(:,:,10) 589 CALL histwrite( nidtrd,"somlgtat",it,z2d ,ndimtrd1,ndextrd1) ! total adv. (eiv) 590 #endif 591 592 ! Salinity trends 593 ! --------------- 594 CALL histwrite( nidtrd,"somlsalt",it,sml ,ndimtrd1,ndextrd1) ! Mixed-layer salinity 595 CALL histwrite( nidtrd,"somltsto",it,zsmltot ,ndimtrd1,ndextrd1) ! total 596 CALL histwrite( nidtrd,"somltsax",it,smltrd(:,:, 1),ndimtrd1,ndextrd1) ! i- adv. 597 CALL histwrite( nidtrd,"somltsay",it,smltrd(:,:, 2),ndimtrd1,ndextrd1) ! j- adv. 598 CALL histwrite( nidtrd,"somltsaz",it,smltrd(:,:, 3),ndimtrd1,ndextrd1) ! vert. adv. 599 CALL histwrite( nidtrd,"somltsdh",it,smltrd(:,:, 4),ndimtrd1,ndextrd1) ! hor. lateral diff. 600 CALL histwrite( nidtrd,"somltsfo",it,smltrd(:,:, 5),ndimtrd1,ndextrd1) ! forcing 601 CALL histwrite( nidtrd,"somlbsdz",it,smltrd(:,:, 6),ndimtrd1,ndextrd1) ! vert. diff. 602 CALL histwrite( nidtrd,"somlbsdt",it,zsmlres ,ndimtrd1,ndextrd1) ! dh/dt entrainment (residual) 603 IF( l_traldf_iso ) THEN 604 CALL histwrite( nidtrd,"somlbsdv",it,smltrd(:,:, 7),ndimtrd1,ndextrd1) ! vert. lateral diff; 605 ENDIF 606 #if defined key_traldf_eiv 607 CALL histwrite( nidtrd,"somlgsax",it,smltrd(:,:, 8),ndimtrd1,ndextrd1) ! i-adv. (eiv) 608 CALL histwrite( nidtrd,"somlgsay",it,smltrd(:,:, 9),ndimtrd1,ndextrd1) ! j-adv. (eiv) 609 CALL histwrite( nidtrd,"somlgsaz",it,smltrd(:,:,10),ndimtrd1,ndextrd1) ! vert. adv. (eiv) 610 z2d(:,:) = smltrd(:,:,8) + smltrd(:,:,9) + smltrd(:,:,10) 611 CALL histwrite( nidtrd,"somlgsat",it,z2d ,ndimtrd1,ndextrd1) ! total adv. (eiv) 612 #endif 598 613 599 614 IF( idebug /= 0 ) THEN 600 WRITE(numout,*) ' debuging trd_mld: I V.5 done'615 WRITE(numout,*) ' debuging trd_mld: III.5 done' 601 616 CALL FLUSH(numout) 602 617 ENDIF 603 618 604 ENDIF 619 ! set counter icount to one to allow the calculation 620 ! of the surface control in the next time step in the trd_mld_zint subroutine 621 icount = 1 622 623 ENDIF 624 625 ! At the end of the 1st time step, set icount to 1 to be 626 ! able to compute the surface control at the beginning of 627 ! the second time step 628 IF( kt == nit000 ) icount = 1 605 629 606 630 IF( kt == nitend ) CALL histclo( nidtrd ) … … 608 632 609 633 END SUBROUTINE trd_mld 634 635 636 637 SUBROUTINE trd_mld_init 638 !!---------------------------------------------------------------------- 639 !! *** ROUTINE trd_mld_init *** 640 !! 641 !! ** Purpose : computation of vertically integrated T and S budgets 642 !! from ocean surface down to control surface (NetCDF output) 643 !! 644 !! ** Method/usage : 645 !! 646 !! History : 647 !! ! 95-04 (J. Vialard) Original code 648 !! ! 97-02 (E. Guilyardi) Adaptation global + base cmo 649 !! ! 99-09 (E. Guilyardi) Re-writing + netCDF output 650 !! 8.5 ! 02-06 (G. Madec) F90: Free form and module 651 !! 9.0 ! 04-08 (C. Talandier) New trends organization 652 !!---------------------------------------------------------------------- 653 !! * Local declarations 654 INTEGER :: ilseq 655 656 REAL(wp) :: zjulian, zsto, zout 657 658 CHARACTER (LEN=21) :: & 659 clold ='OLD' , & ! open specifier (direct access files) 660 clunf ='UNFORMATTED', & ! open specifier (direct access files) 661 clseq ='SEQUENTIAL' ! open specifier (direct access files) 662 CHARACTER (LEN=40) :: clhstnam 663 CHARACTER (LEN=40) :: clop 664 CHARACTER (LEN=12) :: clmxl 665 666 NAMELIST/namtrd/ ntrd, nctls 667 !!---------------------------------------------------------------------- 668 669 ! =================== 670 ! I. initialization 671 ! =================== 672 673 ! Open specifier 674 ilseq = 1 675 idebug = 0 ! set it to 1 in case of problem to have more print 676 icount = 1 677 ionce = 1 678 679 ! namelist namtrd : trend diagnostic 680 REWIND( numnam ) 681 READ ( numnam, namtrd ) 682 683 IF(lwp) THEN 684 WRITE(numout,*) ' ' 685 WRITE(numout,*) 'trd_mld_init: mixed layer heat & freshwater budget trends' 686 WRITE(numout,*) '~~~~~~~~~~~~~' 687 WRITE(numout,*) ' ' 688 WRITE(numout,*) ' Namelist namtrd : ' 689 WRITE(numout,*) ' control surface for trends nctls = ',nctls 690 WRITE(numout,*) ' ' 691 ENDIF 692 693 ! cumulated trends array init 694 nmoymltrd = 0 695 tmltrdm(:,:) = 0.e0 696 smltrdm(:,:) = 0.e0 697 698 ! read control surface from file ctlsurf_idx 699 700 IF( nctls == 1 ) THEN 701 clname ='ctlsurf_idx' 702 CALL ctlopn(numbol,clname,clold,clunf,clseq, & 703 ilseq,numout,lwp,1) 704 REWIND (numbol) 705 READ(numbol) nbol 706 ENDIF 707 708 709 IF( idebug /= 0 ) THEN 710 WRITE(numout,*) ' debuging trd_mld_init: 0. done ' 711 CALL FLUSH(numout) 712 ENDIF 713 714 ! =================================== 715 ! II. netCDF output initialization 716 ! =================================== 717 718 #if defined key_dimgout 719 720 #else 721 ! clmxl = legend root for netCDF output 722 IF( nctls == 0 ) THEN 723 ! control surface = mixed-layer with density criterion 724 ! (array nmln computed in zdfmxl.F90) 725 clmxl = 'Mixed Layer ' 726 ELSE IF( nctls == 1 ) THEN 727 ! control surface = read index from file 728 clmxl = ' Bowl ' 729 ELSE IF( nctls >= 2 ) THEN 730 ! control surface = model level 731 WRITE(clmxl,'(A9,I2,1X)') 'Levels 1-', nctls 732 ENDIF 733 734 !----------------------------------------- 735 ! II.1 Define frequency of output and means 736 ! ----------------------------------------- 737 738 #if defined key_diainstant 739 zsto = nwrite*rdt 740 clop ="inst(x)" 741 #else 742 zsto = rdt 743 clop ="ave(x)" 744 #endif 745 zout = nwrite*rdt 746 747 IF(lwp) WRITE (numout,*) ' trdmld_ncinit: netCDF initialization' 748 749 ! II.2 Compute julian date from starting date of the run 750 ! ------------------------ 751 752 CALL ymds2ju( nyear, nmonth, nday, 0.e0, zjulian ) 753 IF (lwp) WRITE(numout,*)' ' 754 IF (lwp) WRITE(numout,*)' Date 0 used :',nit000 & 755 ,' YEAR ', nyear,' MONTH ', nmonth,' DAY ', nday & 756 ,'Julian day : ', zjulian 757 758 759 ! II.3 Define the T grid trend file (nidtrd) 760 ! --------------------------------- 761 762 CALL dia_nam( clhstnam, nwrite, 'trends' ) ! filename 763 IF(lwp) WRITE(numout,*) ' Name of NETCDF file ', clhstnam 764 CALL histbeg( clhstnam, jpi, glamt, jpj, gphit,1, jpi, & ! Horizontal grid : glamt and gphit 765 1, jpj, 0, zjulian, rdt, nh_t, nidtrd) 766 767 ! Declare output fields as netCDF variables 768 769 ! Mixed layer Depth 770 CALL histdef( nidtrd, "somlttml", clmxl//"Depth" , "m" , & ! hmlp 771 & jpi, jpj, nh_t, 1 , 1, 1 , -99 , 32, clop, zsto, zout ) 772 773 ! Temperature 774 CALL histdef( nidtrd, "somltemp", clmxl//"Temperature" , "C" , & ! ??? 775 & jpi, jpj, nh_t, 1 , 1, 1 , -99 , 32, clop, zsto, zout ) 776 ! Temperature trends 777 CALL histdef( nidtrd, "somlttto", clmxl//"T Total" , "C/s", & ! total 778 & jpi, jpj, nh_t, 1 , 1, 1 , -99 , 32, clop, zout, zout ) 779 CALL histdef( nidtrd, "somlttax", clmxl//"T Zonal Advection", "C/s", & ! i-adv. 780 & jpi, jpj, nh_t, 1 , 1, 1 , -99 , 32, clop, zsto, zout ) 781 CALL histdef( nidtrd, "somlttay", clmxl//"T Meridional Advection", "C/s", & ! j-adv. 782 & jpi, jpj, nh_t, 1 , 1, 1 , -99 , 32, clop, zsto, zout ) 783 CALL histdef( nidtrd, "somlttaz", clmxl//"T Vertical Advection", "C/s", & ! vert. adv. 784 & jpi, jpj, nh_t, 1 , 1, 1 , -99 , 32, clop, zsto, zout ) 785 CALL histdef( nidtrd, "somlttdh", clmxl//"T Horizontal Diffusion ", "C/s", & ! hor. lateral diff. 786 & jpi, jpj, nh_t, 1 , 1, 1 , -99 , 32, clop, zsto, zout ) 787 CALL histdef( nidtrd, "somlttfo", clmxl//"T Forcing", "C/s", & ! forcing 788 & jpi, jpj, nh_t, 1 , 1, 1 , -99 , 32, clop, zsto, zout ) 789 CALL histdef( nidtrd, "somlbtdz", clmxl//"T Vertical Diffusion", "C/s", & ! vert. diff. 790 & jpi, jpj, nh_t, 1 , 1, 1 , -99 , 32, clop, zsto, zout ) 791 CALL histdef( nidtrd, "somlbtdt", clmxl//"T dh/dt Entrainment (Residual)", "C/s", & ! T * dh/dt 792 & jpi, jpj, nh_t, 1 , 1, 1 , -99 , 32, clop, zout, zout ) 793 IF( l_traldf_iso ) THEN 794 CALL histdef( nidtrd, "somlbtdv", clmxl//"T Vert. lateral Diffusion","C/s", & ! vertical diffusion entrainment (ISO) 795 & jpi, jpj, nh_t, 1 , 1, 1 , -99 , 32, clop, zsto, zout ) 796 ENDIF 797 #if defined key_traldf_eiv 798 CALL histdef( nidtrd, "somlgtax", clmxl//"T Zonal EIV Advection", "C/s", & ! i-adv. (eiv) 799 & jpi, jpj, nh_t, 1 , 1, 1 , -99 , 32, clop, zsto, zout ) 800 CALL histdef( nidtrd, "somlgtay", clmxl//"T Meridional EIV Advection", "C/s", & ! j-adv. (eiv) 801 & jpi, jpj, nh_t, 1 , 1, 1 , -99 , 32, clop, zsto, zout ) 802 CALL histdef( nidtrd, "somlgtaz", clmxl//"T Vertical EIV Advection", "C/s", & ! vert. adv. (eiv) 803 & jpi, jpj, nh_t, 1 , 1, 1 , -99 , 32, clop, zsto, zout ) 804 CALL histdef( nidtrd, "somlgtat", clmxl//"T Total EIV Advection", "C/s", & ! total advection (eiv) 805 & jpi, jpj, nh_t, 1 , 1, 1 , -99 , 32, clop, zsto, zout ) 806 #endif 807 ! Salinity 808 CALL histdef( nidtrd, "somlsalt", clmxl//"Salinity", "PSU", & ! Mixed-layer salinity 809 & jpi, jpj, nh_t, 1 , 1, 1 , -99 , 32, clop, zsto, zout ) 810 ! Salinity trends 811 CALL histdef( nidtrd, "somltsto", clmxl//"S Total", "PSU/s", & ! total 812 & jpi, jpj, nh_t, 1 , 1, 1 , -99 , 32, clop, zsto, zout ) 813 CALL histdef( nidtrd, "somltsax", clmxl//"S Zonal Advection", "PSU/s", & ! i-advection 814 & jpi, jpj, nh_t, 1 , 1, 1 , -99 , 32, clop, zsto, zout ) 815 CALL histdef( nidtrd, "somltsay", clmxl//"S Meridional Advection", "PSU/s", & ! j-advection 816 & jpi, jpj, nh_t, 1 , 1, 1 , -99 , 32, clop, zsto, zout ) 817 CALL histdef( nidtrd, "somltsaz", clmxl//"S Vertical Advection", "PSU/s", & ! vertical advection 818 & jpi, jpj, nh_t, 1 , 1, 1 , -99 , 32, clop, zsto, zout ) 819 CALL histdef( nidtrd, "somltsdh", clmxl//"S Horizontal Diffusion ", "PSU/s", & ! hor. lat. diff. 820 & jpi, jpj, nh_t, 1 , 1, 1 , -99 , 32, clop, zsto, zout ) 821 CALL histdef( nidtrd, "somltsfo", clmxl//"S Forcing", "PSU/s", & ! forcing 822 & jpi, jpj, nh_t, 1 , 1, 1 , -99 , 32, clop, zsto, zout ) 823 824 CALL histdef( nidtrd, "somlbsdz", clmxl//"S Vertical Diffusion", "PSU/s", & ! vert. diff. 825 & jpi, jpj, nh_t, 1 , 1, 1 , -99 , 32, clop, zsto, zout ) 826 CALL histdef( nidtrd, "somlbsdt", clmxl//"S dh/dt Entrainment (Residual)", "PSU/s", & ! S * dh/dt 827 & jpi, jpj, nh_t, 1 , 1, 1 , -99 , 32, clop, zsto, zout ) 828 IF( l_traldf_iso ) THEN 829 ! vertical diffusion entrainment (ISO) 830 CALL histdef( nidtrd, "somlbsdv", clmxl//"S Vertical lateral Diffusion", "PSU/s", & ! vert. lat. diff. 831 & jpi, jpj, nh_t, 1 , 1, 1 , -99 , 32, clop, zsto, zout ) 832 ENDIF 833 #if defined key_traldf_eiv 834 CALL histdef( nidtrd, "somlgsax", clmxl//"S Zonal EIV Advection", "PSU/s", & ! i-advection (eiv) 835 & jpi, jpj, nh_t, 1 , 1, 1 , -99 , 32, clop, zsto, zout ) 836 CALL histdef( nidtrd, "somlgsay", clmxl//"S Meridional EIV Advection", "PSU/s", & ! j-advection (eiv) 837 & jpi, jpj, nh_t, 1 , 1, 1 , -99 , 32, clop, zsto, zout ) 838 CALL histdef( nidtrd, "somlgsaz", clmxl//"S Vertical EIV Advection", "PSU/s", & ! vert. adv. (eiv) 839 & jpi, jpj, nh_t, 1 , 1, 1 , -99 , 32, clop, zsto, zout ) 840 CALL histdef( nidtrd, "somlgsat", clmxl//"S Total EIV Advection", "PSU/s", & ! total adv. (eiv) 841 & jpi, jpj, nh_t, 1 , 1, 1 , -99 , 32, clop, zsto, zout ) 842 #endif 843 CALL histend( nidtrd ) 844 #endif 845 846 IF( idebug /= 0 ) THEN 847 WRITE(numout,*) ' debuging trd_mld_init: II. done' 848 CALL FLUSH(numout) 849 ENDIF 850 851 852 END SUBROUTINE trd_mld_init 610 853 611 854 #else … … 615 858 LOGICAL, PUBLIC, PARAMETER :: lk_trdmld = .FALSE. !: momentum trend flag 616 859 CONTAINS 617 SUBROUTINE trd_mld( kt ) ! Empty routine 860 SUBROUTINE trd_mld( kt ) ! Empty routine 861 INTEGER, INTENT( in) :: kt 618 862 WRITE(*,*) 'trd_mld: You should not have seen this print! error?', kt 619 863 END SUBROUTINE trd_mld 864 SUBROUTINE trd_mld_zint( pttrdmld, pstrdmld, ktrd, ctype ) 865 REAL, DIMENSION(:,:,:), INTENT( in ) :: & 866 pttrdmld, pstrdmld ! Temperature and Salinity trends 867 INTEGER, INTENT( in ) :: ktrd ! ocean trend index 868 CHARACTER(len=2), INTENT( in ) :: & 869 ctype ! surface/bottom (2D arrays) or 870 ! ! interior (3D arrays) physics 871 WRITE(*,*) 'trd_mld_zint: You should not have seen this print! error?', pttrdmld(1,1,1) 872 WRITE(*,*) ' " " : You should not have seen this print! error?', pstrdmld(1,1,1) 873 WRITE(*,*) ' " " : You should not have seen this print! error?', ctype 874 WRITE(*,*) ' " " : You should not have seen this print! error?', ktrd 875 END SUBROUTINE trd_mld_zint 876 SUBROUTINE trd_mld_init ! Empty routine 877 WRITE(*,*) 'trd_mld_init: You should not have seen this print! error?' 878 END SUBROUTINE trd_mld_init 620 879 #endif 621 880 -
trunk/NEMO/OPA_SRC/TRD/trdvor.F90
r129 r216 5 5 !!===================================================================== 6 6 7 #if defined key_trd _vor || defined key_esopa7 #if defined key_trdvor || defined key_esopa 8 8 !!---------------------------------------------------------------------- 9 !! 'key_trd _vor' : momentum trend diagnostics9 !! 'key_trdvor' : momentum trend diagnostics 10 10 !!---------------------------------------------------------------------- 11 11 !! trd_vor : momentum trends averaged over the depth 12 !! trd_vor_zint : vorticity vertical integration 13 !! trd_vor_init : initialization step 12 14 !!---------------------------------------------------------------------- 13 15 !! * Modules used 14 16 USE oce ! ocean dynamics and tracers variables 15 17 USE dom_oce ! ocean space and time domain variables 16 USE trd dyn_oce ! ocean active tracer trend variables18 USE trdmod_oce ! ocean variables trends 17 19 USE zdf_oce ! ocean vertical physics 18 20 USE in_out_manager ! I/O manager 19 20 21 USE phycst ! Define parameters for the routines 21 22 USE ldfdyn_oce ! ocean active tracers: lateral physics … … 23 24 USE dianam ! build the name of file (routine) 24 25 USE ldfslp ! iso-neutral slopes 25 USE zdfmxl 26 USE ioipsl 27 USE lbclnk 26 USE zdfmxl ! mixed layer depth 27 USE ioipsl ! NetCDF library 28 USE lbclnk ! ocean lateral boundary conditions (or mpp link) 29 28 30 29 31 IMPLICIT NONE 30 32 PRIVATE 31 33 34 !! * Interfaces 35 INTERFACE trd_vor_zint 36 MODULE PROCEDURE trd_vor_zint_2d, trd_vor_zint_3d 37 END INTERFACE 32 38 33 39 !! * Accessibility 34 PUBLIC trd_vor ! routine called by step.F90 40 PUBLIC trd_vor ! routine called by step.F90 41 PUBLIC trd_vor_zint ! routine called by dynamics routines 42 PUBLIC trd_vor_init ! routine called by opa.F90 35 43 36 44 !! * Shared module variables … … 38 46 39 47 !! * Module variables 40 INTEGER,PARAMETER :: jplvor=1141 48 INTEGER :: & 42 49 nh_t, nmoydpvor , & 43 50 nidvor, nhoridvor, & 44 51 ndexvor1(jpi*jpj), & 45 ndimvor1 52 ndimvor1, icount, & 53 idebug ! (0/1) set it to 1 in case of problem to have more print 46 54 47 55 REAL(wp), DIMENSION(jpi,jpj) :: & … … 49 57 vor_avrb , & ! before vorticity (kt-1) 50 58 vor_avrbb , & ! vorticity at begining of the nwrite-1 timestep averaging period 51 vor_avrbn, & ! after vorticity at time step after the52 rotot, & ! begining of the NWRITE-1 timesteps53 udpvor , & ! total cumulative trends54 vdpvor ! " " "59 vor_avrbn , & ! after vorticity at time step after the 60 rotot , & ! begining of the NWRITE-1 timesteps 61 vor_avrtot , & 62 vor_avrres 55 63 56 64 REAL(wp), DIMENSION(jpi,jpj,jplvor):: & !: curl of trends 57 65 vortrd 66 67 CHARACTER(len=12) :: cvort 58 68 59 69 !! * Substitutions … … 68 78 CONTAINS 69 79 70 SUBROUTINE trd_vor ( kt)80 SUBROUTINE trd_vor_zint_2d( putrdvor, pvtrdvor, ktrd ) 71 81 !!---------------------------------------------------------------------------- 72 !! *** ROUTINE trd_vor ***82 !! *** ROUTINE trd_vor_zint *** 73 83 !! 74 84 !! ** Purpose : computation of vertically integrated vorticity budgets … … 80 90 !! 81 91 !! ** Action : 82 !! /com mld/ :92 !! /comvor/ : 83 93 !! vor_avr average 84 94 !! vor_avrb vorticity at kt-1 … … 107 117 !! 108 118 !! trends output in netCDF format using ioipsl 119 !! 120 !! History : 121 !! 9.0 ! 04-06 (L. Brunier, A-M. Treguier) Original code 122 !! ! 04-08 (C. Talandier) New trends organization 109 123 !!---------------------------------------------------------------------- 110 124 !! * Arguments 111 INTEGER, INTENT( in ) :: kt ! ocean time-step index 125 INTEGER, INTENT( in ) :: ktrd ! ocean trend index 126 127 REAL(wp), DIMENSION(jpi,jpj), INTENT( inout ) :: & 128 putrdvor, & ! u vorticity trend 129 pvtrdvor ! v vorticity trend 112 130 113 131 !! * Local declarations 114 INTEGER ilseq 115 INTEGER ji, jj, jk, jl, idebug, it 132 INTEGER :: ji, jj 133 INTEGER :: ikbu, ikbum1, ikbv, ikbvm1 134 REAL(wp), DIMENSION(jpi,jpj) :: & 135 zudpvor, & ! total cmulative trends 136 zvdpvor ! " " " 137 !!---------------------------------------------------------------------- 138 139 ! Initialization 140 zudpvor(:,:) = 0.e0 141 zvdpvor(:,:) = 0.e0 142 143 CALL lbc_lnk( putrdvor, 'U' , -1. ) 144 CALL lbc_lnk( pvtrdvor, 'V' , -1. ) 145 146 ! ===================================== 147 ! I vertical integration of 2D trends 148 ! ===================================== 149 150 SELECT CASE (ktrd) 151 152 CASE (jpvorbfr) ! bottom friction 153 154 DO jj = 2, jpjm1 155 DO ji = fs_2, fs_jpim1 156 ikbu = min( mbathy(ji+1,jj), mbathy(ji,jj) ) 157 ikbum1 = max( ikbu-1, 1 ) 158 ikbv = min( mbathy(ji,jj+1), mbathy(ji,jj) ) 159 ikbvm1 = max( ikbv-1, 1 ) 160 161 zudpvor(ji,jj) = putrdvor(ji,jj) * fse3u(ji,jj,ikbum1) * e1u(ji,jj) * umask(ji,jj,ikbum1) 162 zvdpvor(ji,jj) = pvtrdvor(ji,jj) * fse3v(ji,jj,ikbvm1) * e2v(ji,jj) * vmask(ji,jj,ikbvm1) 163 END DO 164 END DO 165 166 CASE (jpvorswf) ! wind stress 167 168 zudpvor(:,:) = putrdvor(:,:) * fse3u(:,:,1) * e1u(:,:) * umask(:,:,1) 169 zvdpvor(:,:) = pvtrdvor(:,:) * fse3v(:,:,1) * e2v(:,:) * vmask(:,:,1) 170 171 END SELECT 172 173 ! Average except for Beta.V 174 zudpvor(:,:) = zudpvor(:,:) * hur(:,:) 175 zvdpvor(:,:) = zvdpvor(:,:) * hvr(:,:) 176 177 ! Curl 178 DO ji=1,jpim1 179 DO jj=1,jpjm1 180 vortrd(ji,jj,ktrd) = ( zvdpvor(ji+1,jj) - zvdpvor(ji,jj) & 181 & - ( zudpvor(ji,jj+1) - zudpvor(ji,jj) ) ) & 182 & / ( e1f(ji,jj) * e2f(ji,jj) ) 183 END DO 184 END DO 185 186 ! Surface mask 187 vortrd(:,:,ktrd) = vortrd(:,:,ktrd) * fmask(:,:,1) 188 189 IF( idebug /= 0 ) THEN 190 IF(lwp) WRITE(numout,*) ' debuging trd_vor_zint: I done' 191 CALL FLUSH(numout) 192 ENDIF 193 194 END SUBROUTINE trd_vor_zint_2d 195 196 197 198 SUBROUTINE trd_vor_zint_3d( putrdvor, pvtrdvor, ktrd ) 199 !!---------------------------------------------------------------------------- 200 !! *** ROUTINE trd_vor_zint *** 201 !! 202 !! ** Purpose : computation of vertically integrated vorticity budgets 203 !! from ocean surface down to control surface (NetCDF output) 204 !! 205 !! ** Method/usage : 206 !! integration done over nwrite-1 time steps 207 !! 208 !! 209 !! ** Action : 210 !! /comvor/ : 211 !! vor_avr average 212 !! vor_avrb vorticity at kt-1 213 !! vor_avrbb vorticity at begining of the NWRITE-1 214 !! time steps averaging period 215 !! vor_avrbn vorticity at time step after the 216 !! begining of the NWRITE-1 time 217 !! steps averaging period 218 !! 219 !! trends : 220 !! 221 !! vortrd (,,1) = Pressure Gradient Trend 222 !! vortrd (,,2) = KE Gradient Trend 223 !! vortrd (,,3) = Relative Vorticity Trend 224 !! vortrd (,,4) = Coriolis Term Trend 225 !! vortrd (,,5) = Horizontal Diffusion Trend 226 !! vortrd (,,6) = Vertical Advection Trend 227 !! vortrd (,,7) = Vertical Diffusion Trend 228 !! vortrd (,,8) = Surface Pressure Grad. Trend 229 !! vortrd (,,9) = Beta V 230 !! vortrd (,,10) = forcing term 231 !! vortrd (,,11) = bottom friction term 232 !! rotot(,) : total cumulative trends over nwrite-1 time steps 233 !! vor_avrtot(,) : first membre of vrticity equation 234 !! vor_avrres(,) : residual = dh/dt entrainment 235 !! 236 !! trends output in netCDF format using ioipsl 237 !! 238 !! History : 239 !! 9.0 ! 04-06 (L. Brunier, A-M. Treguier) Original code 240 !! ! 04-08 (C. Talandier) New trends organization 241 !!---------------------------------------------------------------------- 242 !! * Arguments 243 INTEGER, INTENT( in ) :: ktrd ! ocean trend index 244 245 REAL(wp), DIMENSION(jpi,jpj,jpk), INTENT( inout ) :: & 246 putrdvor, & ! u vorticity trend 247 pvtrdvor ! v vorticity trend 248 249 !! * Local declarations 250 INTEGER :: ji, jj, jk 251 252 REAL(wp), DIMENSION(jpi,jpj) :: & 253 zubet, & ! u Beta.V case 254 zvbet, & ! v Beta.V case 255 zudpvor, & ! total cmulative trends 256 zvdpvor ! " " " 257 !!---------------------------------------------------------------------- 258 259 ! Initialization 260 zubet(:,:) = 0.e0 261 zvbet(:,:) = 0.e0 262 zudpvor(:,:) = 0.e0 263 zvdpvor(:,:) = 0.e0 264 265 ! ===================================== 266 ! I vertical integration of 3D trends 267 ! ===================================== 268 269 CALL lbc_lnk( putrdvor, 'U' , -1. ) 270 CALL lbc_lnk( pvtrdvor, 'V' , -1. ) 271 272 ! putrdvor and pvtrdvor terms 273 DO jk = 1,jpk 274 zudpvor(:,:) = zudpvor(:,:) + putrdvor(:,:,jk) * fse3u(:,:,jk) * e1u(:,:) * umask(:,:,jk) 275 zvdpvor(:,:) = zvdpvor(:,:) + pvtrdvor(:,:,jk) * fse3v(:,:,jk) * e2v(:,:) * vmask(:,:,jk) 276 END DO 277 278 ! Save Beta.V term to avoid average before Curl 279 ! Beta.V : intergration, no average 280 IF( ktrd == jpvorbev ) THEN 281 zubet(:,:) = zudpvor(:,:) 282 zvbet(:,:) = zvdpvor(:,:) 283 ENDIF 284 285 ! Average except for Beta.V 286 zudpvor(:,:) = zudpvor(:,:) * hur(:,:) 287 zvdpvor(:,:) = zvdpvor(:,:) * hvr(:,:) 288 289 ! Curl 290 DO ji=1,jpim1 291 DO jj=1,jpjm1 292 vortrd(ji,jj,ktrd) = ( zvdpvor(ji+1,jj) - zvdpvor(ji,jj) - & 293 & ( zudpvor(ji,jj+1) - zudpvor(ji,jj) ) ) & 294 & / ( e1f(ji,jj) * e2f(ji,jj) ) 295 END DO 296 END DO 297 298 ! Surface mask 299 vortrd(:,:,ktrd) = vortrd(:,:,ktrd) * fmask(:,:,1) 300 301 ! Special treatement for the Beta.V term 302 ! Compute the Curl of the Beta.V term which is not averaged 303 IF( ktrd == jpvorbev ) THEN 304 DO ji=1,jpim1 305 DO jj=1,jpjm1 306 vortrd(ji,jj,jpvorbev) = ( zvbet(ji+1,jj) - zvbet(ji,jj) - & 307 & ( zubet(ji,jj+1) - zubet(ji,jj) ) ) & 308 & / ( e1f(ji,jj) * e2f(ji,jj) ) 309 END DO 310 END DO 311 312 ! Average on the Curl 313 vortrd(:,:,jpvorbev) = vortrd(:,:,jpvorbev) * hur(:,:) 314 315 ! Surface mask 316 vortrd(:,:,jpvorbev) = vortrd(:,:,jpvorbev) * fmask(:,:,1) 317 ENDIF 318 319 IF( idebug /= 0 ) THEN 320 IF(lwp) WRITE(numout,*) ' debuging trd_vor_zint: I done' 321 CALL FLUSH(numout) 322 ENDIF 323 324 END SUBROUTINE trd_vor_zint_3d 325 326 327 328 SUBROUTINE trd_vor( kt ) 329 !!---------------------------------------------------------------------- 330 !! *** ROUTINE trd_vor *** 331 !! 332 !! ** Purpose : computation of cumulated trends over analysis period 333 !! and make outputs (NetCDF or DIMG format) 334 !! 335 !! ** Method/usage : 336 !! 337 !! History : 338 !! 9.0 ! 04-06 (L. Brunier, A-M. Treguier) Original code 339 !! ! 04-08 (C. Talandier) New trends organization 340 !!---------------------------------------------------------------------- 341 !! * Arguments 342 INTEGER, INTENT( in ) :: kt ! ocean time-step index 343 344 !! * Local declarations 345 INTEGER :: ji, jj, jk, jl, it 116 346 117 347 REAL(wp) :: zmean 118 REAL(wp) :: zun(jpi,jpj), zvn(jpi,jpj) 119 REAL(wp) :: zjulian, zsto, zout 120 REAL(wp) :: vor_avrtot(jpi,jpj), vor_avrres(jpi,jpj) 121 INTEGER(wp) :: ikbu,ikbum1,ikbv,ikbvm1 122 CHARACTER (len=12) :: cvort 123 CHARACTER (len=40) :: clhstnam 124 CHARACTER (len=40) :: clop 125 126 NAMELIST/namtrd/ ntrd,nctls 127 !!---------------------------------------------------------------------- 348 349 REAL(wp) ,DIMENSION(jpi,jpj) :: & 350 zun, zvn 351 !!---------------------------------------------------------------------- 352 353 ! ================= 354 ! I. Initialization 355 ! ================= 128 356 129 ! =================== 130 ! 0. initialization 131 ! =================== 132 133 cvort='averaged-vor' 134 135 ! Open specifier 136 ilseq = 1 137 idebug = 0 ! set it to 1 in case of problem to have more Print 138 139 IF( kt == nit000 ) THEN 140 141 ! namelist namtrd : trend diagnostic 142 REWIND( numnam ) 143 READ ( numnam, namtrd ) 144 145 IF(lwp) THEN 146 WRITE(numout,*) 'namtrd' 147 WRITE(numout,*) ' ' 148 WRITE(numout,*) ' time step frequency trend ntrd = ',ntrd 149 WRITE(numout,*) ' ' 150 ENDIF 151 152 ! cumulated trends array init 153 nmoydpvor = 0 154 rotot(:,:)=0 155 vor_avrtot(:,:)=0 156 vor_avrres(:,:)=0 157 ENDIF 158 159 ! set before values of vertically average u and v 357 358 ! I.1 set before values of vertically average u and v 359 ! --------------------------------------------------- 160 360 161 361 IF( kt > nit000 ) THEN … … 164 364 165 365 IF( idebug /= 0 ) THEN 166 WRITE(numout,*) ' debuging trd_vor: 0.done '366 WRITE(numout,*) ' debuging trd_vor: I.1 done ' 167 367 CALL FLUSH(numout) 168 368 ENDIF 169 369 170 ! ================================= 171 ! I. vertically integrated vorticity 172 ! ================================= 370 ! I.2 vertically integrated vorticity 371 ! ---------------------------------- 173 372 174 373 vor_avr(:,:) = 0. … … 178 377 vor_avrres(:,:)=0 179 378 180 ! vertically averaged velocity379 ! Vertically averaged velocity 181 380 DO jk = 1, jpk - 1 182 381 zun(:,:)=zun(:,:) + e1u(:,:)*un(:,:,jk)*fse3u(:,:,jk) 183 382 zvn(:,:)=zvn(:,:) + e2v(:,:)*vn(:,:,jk)*fse3v(:,:,jk) 184 383 END DO 185 186 384 187 385 zun(:,:)=zun(:,:)*hur(:,:) 188 386 zvn(:,:)=zvn(:,:)*hvr(:,:) 189 387 190 ! Curl388 ! Curl 191 389 DO ji=1,jpim1 192 390 DO jj=1,jpjm1 … … 198 396 END DO 199 397 200 201 398 IF(idebug /= 0) THEN 202 WRITE(numout,*) ' debuging trd_vor: I done'399 WRITE(numout,*) ' debuging trd_vor: I.2 done' 203 400 CALL FLUSH(numout) 204 401 ENDIF 205 402 206 403 ! ================================= 207 ! II. netCDF output initialization404 ! II. Cumulated trends 208 405 ! ================================= 209 406 210 # include "trdvor_ncinit.h90" 211 212 IF( idebug /= 0 ) THEN 213 WRITE(numout,*) ' debuging trd_vor: II. done' 214 CALL FLUSH(numout) 215 ENDIF 216 217 ! ===================================== 218 ! III vertical integration of 3D trends 219 ! ===================================== 220 ! Beta.V : intergration, no average 221 utrd(:,:,:,9)=utrd(:,:,:,4) 222 vtrd(:,:,:,9)=vtrd(:,:,:,4) 223 224 DO jl=1,jplvor 225 226 udpvor(:,:)=0 227 vdpvor(:,:)=0 228 229 !bottom friction 230 IF( jl == jplvor ) THEN 231 232 CALL lbc_lnk( tautrd(:,:,3), 'U' , -1. ) 233 CALL lbc_lnk( tautrd(:,:,4), 'V' , -1. ) 234 235 DO jj = 2, jpjm1 236 DO ji = fs_2, fs_jpim1 237 ikbu = min( mbathy(ji+1,jj), mbathy(ji,jj) ) 238 ikbum1 = max( ikbu-1, 1 ) 239 ikbv = min( mbathy(ji,jj+1), mbathy(ji,jj) ) 240 ikbvm1 = max( ikbv-1, 1 ) 241 242 udpvor(ji,jj)=tautrd(ji,jj,3)*fse3u(ji,jj,ikbum1)*e1u(ji,jj)*umask(ji,jj,ikbum1) 243 vdpvor(ji,jj)=tautrd(ji,jj,4)*fse3v(ji,jj,ikbvm1)*e2v(ji,jj)*vmask(ji,jj,ikbvm1) 244 END DO 245 END DO 246 247 !wind stress 248 ELSE IF( jl == (jplvor-1) ) THEN 249 250 CALL lbc_lnk( tautrd(:,:,1), 'U' , -1. ) 251 CALL lbc_lnk( tautrd(:,:,2), 'V' , -1. ) 252 253 udpvor(:,:)=tautrd(:,:,1)*fse3u(:,:,1)*e1u(:,:)*umask(:,:,1) 254 vdpvor(:,:)=tautrd(:,:,2)*fse3v(:,:,1)*e2v(:,:)*vmask(:,:,1) 255 256 ELSE 257 258 CALL lbc_lnk( utrd(:,:,:,jl), 'U' , -1. ) 259 CALL lbc_lnk( vtrd(:,:,:,jl), 'V' , -1. ) 260 261 !utrd and vtrd terms 262 DO jk = 1,jpk 263 udpvor(:,:)=udpvor(:,:)+utrd(:,:,jk,jl)*fse3u(:,:,jk)*e1u(:,:)*umask(:,:,jk) 264 vdpvor(:,:)=vdpvor(:,:)+vtrd(:,:,jk,jl)*fse3v(:,:,jk)*e2v(:,:)*vmask(:,:,jk) 265 END DO 266 267 ENDIF 268 269 !average except for Beta.V 270 IF (jl/=9) THEN 271 udpvor(:,:) = udpvor(:,:) * hur(:,:) 272 vdpvor(:,:) = vdpvor(:,:) * hvr(:,:) 273 ENDIF 274 275 !Curl 276 DO ji=1,jpim1 277 DO jj=1,jpjm1 278 vortrd(ji,jj,jl)=( vdpvor(ji+1,jj)-vdpvor(ji,jj) & 279 - ( udpvor(ji,jj+1)-udpvor(ji,jj) ) ) & 280 / ( e1f(ji,jj) * e2f(ji,jj) ) 281 END DO 282 END DO 283 284 vortrd(:,:,9)=vortrd(:,:,9)*hur(:,:) 285 286 !surface mask 287 DO ji=1,jpi 288 DO jj=1,jpj 289 vortrd(ji,jj,jl)=vortrd(ji,jj,jl)*fmask(ji,jj,1) !surface mask 290 END DO 291 END DO 292 293 END DO 294 295 IF( idebug /= 0 ) THEN 296 IF(lwp) WRITE(numout,*) ' debuging trd_vor: III done' 297 CALL FLUSH(numout) 298 ENDIF 299 300 ! ================================= 301 ! IV. Cumulated trends 302 ! ================================= 303 304 ! IV.1 set `before' mixed layer values for kt = nit000+1 305 ! -------------------------------------------------------- 407 ! II.1 set `before' mixed layer values for kt = nit000+1 408 ! ------------------------------------------------------ 306 409 IF( kt == nit000+1 ) THEN 307 410 vor_avrbb(:,:) = vor_avrb(:,:) … … 310 413 311 414 IF( idebug /= 0 ) THEN 312 WRITE(numout,*) ' debuging trd_vor: I V.1 done'415 WRITE(numout,*) ' debuging trd_vor: I1.1 done' 313 416 CALL FLUSH(numout) 314 417 ENDIF 315 418 316 ! I V.2 cumulated trends over analysis period (kt=2 to nwrite)419 ! II.2 cumulated trends over analysis period (kt=2 to nwrite) 317 420 ! ---------------------- 318 421 ! trends cumulated over nwrite-2 time steps … … 328 431 329 432 IF( idebug /= 0 ) THEN 330 WRITE(numout,*) ' debuging trd_vor: I V.2 done'433 WRITE(numout,*) ' debuging trd_vor: II.2 done' 331 434 CALL FLUSH(numout) 332 435 ENDIF 333 436 334 437 ! ============================================= 335 ! V. Output in netCDF + residual computation438 ! III. Output in netCDF + residual computation 336 439 ! ============================================= 440 337 441 IF( MOD( kt - nit000+1, ntrd ) == 0 ) THEN 338 442 339 ! V.1 compute total trend443 ! III.1 compute total trend 340 444 ! ------------------------ 341 445 zmean = float(nmoydpvor) … … 346 450 IF( idebug /= 0 ) THEN 347 451 WRITE(numout,*) ' zmean = ',zmean 348 WRITE(numout,*) ' debuging trd_vor: V.1 done'452 WRITE(numout,*) ' debuging trd_vor: III.1 done' 349 453 CALL FLUSH(numout) 350 454 ENDIF 351 455 352 ! V.2 compute residual456 ! III.2 compute residual 353 457 ! --------------------- 354 458 vor_avrres(:,:) = vor_avrtot(:,:) - rotot(:,:) / zmean … … 359 463 360 464 IF( idebug /= 0 ) THEN 361 WRITE(numout,*) ' debuging trd_vor: V.2 done'465 WRITE(numout,*) ' debuging trd_vor: III.2 done' 362 466 CALL FLUSH(numout) 363 467 ENDIF 364 468 365 ! V.3 time evolution array swap469 ! III.3 time evolution array swap 366 470 ! ------------------------------ 367 471 vor_avrbb(:,:) = vor_avrb(:,:) … … 369 473 370 474 IF( idebug /= 0 ) THEN 371 WRITE(numout,*) ' debuging trd_vor: V.3 done'475 WRITE(numout,*) ' debuging trd_vor: III.3 done' 372 476 CALL FLUSH(numout) 373 477 ENDIF … … 377 481 ENDIF 378 482 379 ! V.5write trends to output483 ! III.4 write trends to output 380 484 ! --------------------------- 485 381 486 IF( kt >= nit000+1 ) THEN 382 487 383 #include "trdvor_ncwrite.h90" 488 ! define time axis 489 it= kt-nit000+1 490 IF( lwp .AND. MOD( kt, ntrd ) == 0 ) THEN 491 WRITE(numout,*) ' trdvor_ncwrite : write NetCDF fields' 492 ENDIF 493 494 CALL histwrite( nidvor,"sovortPh",it,vortrd(:,:,1),ndimvor1,ndexvor1) ! grad Ph 495 CALL histwrite( nidvor,"sovortEk",it,vortrd(:,:,2),ndimvor1,ndexvor1) ! Energy 496 CALL histwrite( nidvor,"sovozeta",it,vortrd(:,:,3),ndimvor1,ndexvor1) ! rel vorticity 497 CALL histwrite( nidvor,"sovortif",it,vortrd(:,:,4),ndimvor1,ndexvor1) ! coriolis 498 CALL histwrite( nidvor,"sovodifl",it,vortrd(:,:,5),ndimvor1,ndexvor1) ! lat diff 499 CALL histwrite( nidvor,"sovoadvv",it,vortrd(:,:,6),ndimvor1,ndexvor1) ! vert adv 500 CALL histwrite( nidvor,"sovodifv",it,vortrd(:,:,7),ndimvor1,ndexvor1) ! vert diff 501 CALL histwrite( nidvor,"sovortPs",it,vortrd(:,:,8),ndimvor1,ndexvor1) ! grad Ps 502 CALL histwrite( nidvor,"sovortbv",it,vortrd(:,:,9),ndimvor1,ndexvor1) ! beta.V 503 CALL histwrite( nidvor,"sovowind",it,vortrd(:,:,10),ndimvor1,ndexvor1) ! wind stress 504 CALL histwrite( nidvor,"sovobfri",it,vortrd(:,:,11),ndimvor1,ndexvor1) ! bottom friction 505 CALL histwrite( nidvor,"1st_mbre",it,vor_avrtot ,ndimvor1,ndexvor1) ! First membre 506 CALL histwrite( nidvor,"sovorgap",it,vor_avrres ,ndimvor1,ndexvor1) ! gap between 1st and 2 nd mbre 384 507 385 508 IF( idebug /= 0 ) THEN 386 WRITE(numout,*) ' debuging trd_vor: I V.5done'509 WRITE(numout,*) ' debuging trd_vor: III.4 done' 387 510 CALL FLUSH(numout) 388 511 ENDIF … … 390 513 ENDIF 391 514 392 IF( MOD( kt - nit000+1, ntrd ) == 0 ) THEN 393 rotot(:,:)=0 394 ENDIF 395 396 IF( kt == nitend ) THEN 397 CALL histclo( nidvor ) 398 ENDIF 515 IF( MOD( kt - nit000+1, ntrd ) == 0 ) rotot(:,:)=0 516 517 IF( kt == nitend ) CALL histclo( nidvor ) 399 518 400 519 END SUBROUTINE trd_vor 520 521 522 523 SUBROUTINE trd_vor_init 524 !!---------------------------------------------------------------------- 525 !! *** ROUTINE trd_vor_init *** 526 !! 527 !! ** Purpose : computation of vertically integrated T and S budgets 528 !! from ocean surface down to control surface (NetCDF output) 529 !! 530 !! ** Method/usage : 531 !! 532 !! History : 533 !! 9.0 ! 04-06 (L. Brunier, A-M. Treguier) Original code 534 !! ! 04-08 (C. Talandier) New trends organization 535 !!---------------------------------------------------------------------- 536 !! * Local declarations 537 REAL(wp) :: zjulian, zsto, zout 538 539 CHARACTER (len=40) :: clhstnam 540 CHARACTER (len=40) :: clop 541 542 NAMELIST/namtrd/ ntrd,nctls 543 !!---------------------------------------------------------------------- 544 545 ! =================== 546 ! I. initialization 547 ! =================== 548 549 cvort='averaged-vor' 550 551 ! Open specifier 552 idebug = 0 ! set it to 1 in case of problem to have more Print 553 554 ! namelist namtrd : trend diagnostic 555 REWIND( numnam ) 556 READ ( numnam, namtrd ) 557 558 IF(lwp) THEN 559 WRITE(numout,*) ' ' 560 WRITE(numout,*) 'trd_vor_init: vorticity trends' 561 WRITE(numout,*) '~~~~~~~~~~~~~' 562 WRITE(numout,*) ' ' 563 WRITE(numout,*) ' Namelist namtrd : ' 564 WRITE(numout,*) ' time step frequency trend ntrd = ',ntrd 565 WRITE(numout,*) ' ' 566 WRITE(numout,*) '##########################################################################' 567 WRITE(numout,*) ' CAUTION: The interpretation of the vorticity trends is' 568 WRITE(numout,*) ' not obvious, please contact Anne-Marie TREGUIER at: treguier@ifremer.fr ' 569 WRITE(numout,*) '##########################################################################' 570 WRITE(numout,*) ' ' 571 ENDIF 572 573 ! cumulated trends array init 574 nmoydpvor = 0 575 rotot(:,:)=0 576 vor_avrtot(:,:)=0 577 vor_avrres(:,:)=0 578 579 IF( idebug /= 0 ) THEN 580 WRITE(numout,*) ' debuging trd_vor_init: I. done' 581 CALL FLUSH(numout) 582 ENDIF 583 584 ! ================================= 585 ! II. netCDF output initialization 586 ! ================================= 587 588 !----------------------------------------- 589 ! II.1 Define frequency of output and means 590 ! ----------------------------------------- 591 #if defined key_diainstant 592 zsto = nwrite*rdt 593 clop ="inst(x)" 594 #else 595 zsto = rdt 596 clop ="ave(x)" 597 #endif 598 zout = ntrd*rdt 599 600 IF(lwp) WRITE (numout,*) ' trdvor_ncinit: netCDF initialization' 601 602 ! II.2 Compute julian date from starting date of the run 603 ! ------------------------ 604 CALL ymds2ju( nyear, nmonth, nday, 0.e0, zjulian ) 605 IF (lwp) WRITE(numout,*)' ' 606 IF (lwp) WRITE(numout,*)' Date 0 used :',nit000 & 607 ,' YEAR ', nyear,' MONTH ', nmonth,' DAY ', nday & 608 ,'Julian day : ', zjulian 609 610 ! II.3 Define the T grid trend file (nidvor) 611 ! --------------------------------- 612 CALL dia_nam( clhstnam, ntrd, 'vort' ) ! filename 613 IF(lwp) WRITE(numout,*) ' Name of NETCDF file ', clhstnam 614 CALL histbeg( clhstnam, jpi, glamf, jpj, gphif,1, jpi, & ! Horizontal grid : glamt and gphit 615 & 1, jpj, 0, zjulian, rdt, nh_t, nidvor) 616 CALL wheneq( jpi*jpj, fmask, 1, 1., ndexvor1, ndimvor1 ) ! surface 617 618 ! Declare output fields as netCDF variables 619 CALL histdef( nidvor, "sovortPh", cvort//"grad Ph" , "s-2", & ! grad Ph 620 & jpi,jpj,nh_t,1,1,1,-99,32,clop,zsto,zout) 621 CALL histdef( nidvor, "sovortEk", cvort//"Energy", "s-2", & ! Energy 622 & jpi,jpj,nh_t,1,1,1,-99,32,clop,zsto,zout) 623 CALL histdef( nidvor, "sovozeta", cvort//"rel vorticity", "s-2", & ! rel vorticity 624 & jpi,jpj,nh_t,1,1,1,-99,32,clop,zsto,zout) 625 CALL histdef( nidvor, "sovortif", cvort//"coriolis", "s-2", & ! coriolis 626 & jpi,jpj,nh_t,1,1,1,-99,32,clop,zsto,zout) 627 CALL histdef( nidvor, "sovodifl", cvort//"lat diff ", "s-2", & ! lat diff 628 & jpi,jpj,nh_t,1,1,1,-99,32,clop,zsto,zout) 629 CALL histdef( nidvor, "sovoadvv", cvort//"vert adv", "s-2", & ! vert adv 630 & jpi,jpj,nh_t,1,1,1,-99,32,clop,zsto,zout) 631 CALL histdef( nidvor, "sovodifv", cvort//"vert diff" , "s-2", & ! vert diff 632 & jpi,jpj,nh_t,1,1,1,-99,32,clop,zsto,zout) 633 CALL histdef( nidvor, "sovortPs", cvort//"grad Ps", "s-2", & ! grad Ps 634 & jpi,jpj,nh_t,1,1,1,-99,32,clop,zsto,zout) 635 CALL histdef( nidvor, "sovortbv", cvort//"Beta V", "s-2", & ! beta.V 636 & jpi,jpj,nh_t,1,1,1,-99,32,clop,zsto,zout) 637 CALL histdef( nidvor, "sovowind", cvort//"wind stress", "s-2", & ! wind stress 638 & jpi,jpj,nh_t,1,1,1,-99,32,clop,zsto,zout) 639 CALL histdef( nidvor, "sovobfri", cvort//"bottom friction", "s-2", & ! bottom friction 640 & jpi,jpj,nh_t,1,1,1,-99,32,clop,zsto,zout) 641 CALL histdef( nidvor, "1st_mbre", cvort//"1st mbre", "s-2", & ! First membre 642 & jpi,jpj,nh_t,1,1,1,-99,32,clop,zsto,zout) 643 CALL histdef( nidvor, "sovorgap", cvort//"gap", "s-2", & ! gap between 1st and 2 nd mbre 644 & jpi,jpj,nh_t,1,1,1,-99,32,clop,zsto,zout) 645 CALL histend( nidvor ) 646 647 IF( idebug /= 0 ) THEN 648 WRITE(numout,*) ' debuging trd_vor_init: II. done' 649 CALL FLUSH(numout) 650 ENDIF 651 652 END SUBROUTINE trd_vor_init 401 653 402 654 #else … … 405 657 !!---------------------------------------------------------------------- 406 658 LOGICAL, PUBLIC :: lk_trdvor = .FALSE. ! momentum trend flag 659 660 !! * Interfaces 661 INTERFACE trd_vor_zint 662 MODULE PROCEDURE trd_vor_zint_2d, trd_vor_zint_3d 663 END INTERFACE 664 407 665 CONTAINS 408 666 SUBROUTINE trd_vor( kt ) ! Empty routine 409 667 WRITE(*,*) 'trd_vor: You should not have seen this print! error?', kt 410 668 END SUBROUTINE trd_vor 669 SUBROUTINE trd_vor_zint_2d( putrdvor, pvtrdvor, ktrd ) 670 REAL, DIMENSION(:,:), INTENT( inout ) :: & 671 putrdvor, pvtrdvor ! U and V momentum trends 672 INTEGER, INTENT( in ) :: ktrd ! ocean trend index 673 WRITE(*,*) 'trd_vor_zint_2d: You should not have seen this print! error?', putrdvor(1,1) 674 WRITE(*,*) ' " " : You should not have seen this print! error?', pvtrdvor(1,1) 675 WRITE(*,*) ' " " : You should not have seen this print! error?', ktrd 676 END SUBROUTINE trd_vor_zint_2d 677 SUBROUTINE trd_vor_zint_3d( putrdvor, pvtrdvor, ktrd ) 678 REAL, DIMENSION(:,:,:), INTENT( inout ) :: & 679 putrdvor, pvtrdvor ! U and V momentum trends 680 INTEGER, INTENT( in ) :: ktrd ! ocean trend index 681 WRITE(*,*) 'trd_vor_zint_3d: You should not have seen this print! error?', putrdvor(1,1,1) 682 WRITE(*,*) ' " " : You should not have seen this print! error?', pvtrdvor(1,1,1) 683 WRITE(*,*) ' " " : You should not have seen this print! error?', ktrd 684 END SUBROUTINE trd_vor_zint_3d 685 SUBROUTINE trd_vor_init ! Empty routine 686 WRITE(*,*) 'trd_vor_init: You should not have seen this print! error?' 687 END SUBROUTINE trd_vor_init 411 688 #endif 412 413 689 !!====================================================================== 414 690 END MODULE trdvor
Note: See TracChangeset
for help on using the changeset viewer.
