Changeset 6225 for branches/2014/dev_r4704_NOC5_MPP_BDY_UPDATE/NEMOGCM/NEMO/TOP_SRC/PISCES/P4Z/p4zsink.F90
- Timestamp:
- 2016-01-08T10:35:19+01:00 (8 years ago)
- File:
-
- 1 edited
Legend:
- Unmodified
- Added
- Removed
-
branches/2014/dev_r4704_NOC5_MPP_BDY_UPDATE/NEMOGCM/NEMO/TOP_SRC/PISCES/P4Z/p4zsink.F90
r4624 r6225 41 41 #endif 42 42 43 INTEGER :: ik sed = 1043 INTEGER :: ik100 44 44 45 45 #if defined key_kriest … … 65 65 #endif 66 66 67 !!* Substitution68 # include "top_substitute.h90"69 67 !!---------------------------------------------------------------------- 70 68 !! NEMO/TOP 3.3 , NEMO Consortium (2010) … … 79 77 !!---------------------------------------------------------------------- 80 78 81 SUBROUTINE p4z_sink ( kt, jnt )79 SUBROUTINE p4z_sink ( kt, knt ) 82 80 !!--------------------------------------------------------------------- 83 81 !! *** ROUTINE p4z_sink *** … … 88 86 !! ** Method : - ??? 89 87 !!--------------------------------------------------------------------- 90 INTEGER, INTENT(in) :: kt, jnt88 INTEGER, INTENT(in) :: kt, knt 91 89 INTEGER :: ji, jj, jk, jit 92 90 INTEGER :: iiter1, iiter2 … … 94 92 REAL(wp) :: zagg , zaggfe, zaggdoc, zaggdoc2, zaggdoc3 95 93 REAL(wp) :: zfact, zwsmax, zmax, zstep 96 REAL(wp) :: zrfact297 INTEGER :: ik198 94 CHARACTER (len=25) :: charout 95 REAL(wp), POINTER, DIMENSION(:,:,:) :: zw3d 96 REAL(wp), POINTER, DIMENSION(:,: ) :: zw2d 99 97 !!--------------------------------------------------------------------- 100 98 ! … … 108 106 DO ji = 1,jpi 109 107 zmax = MAX( heup(ji,jj), hmld(ji,jj) ) 110 zfact = MAX( 0., fsdepw(ji,jj,jk+1) - zmax ) / 5000._wp108 zfact = MAX( 0., gdepw_n(ji,jj,jk+1) - zmax ) / 5000._wp 111 109 wsbio4(ji,jj,jk) = wsbio2 + ( 200.- wsbio2 ) * zfact 112 110 END DO … … 137 135 DO ji = 1, jpi 138 136 IF( tmask(ji,jj,jk) == 1) THEN 139 zwsmax = 0.5 * fse3t(ji,jj,jk) / xstep137 zwsmax = 0.5 * e3t_n(ji,jj,jk) / xstep 140 138 iiter1 = MAX( iiter1, INT( wsbio3(ji,jj,jk) / zwsmax ) ) 141 139 iiter2 = MAX( iiter2, INT( wsbio4(ji,jj,jk) / zwsmax ) ) … … 156 154 DO ji = 1, jpi 157 155 IF( tmask(ji,jj,jk) == 1 ) THEN 158 zwsmax = 0.5 * fse3t(ji,jj,jk) / xstep156 zwsmax = 0.5 * e3t_n(ji,jj,jk) / xstep 159 157 wsbio3(ji,jj,jk) = MIN( wsbio3(ji,jj,jk), zwsmax * FLOAT( iiter1 ) ) 160 158 wsbio4(ji,jj,jk) = MIN( wsbio4(ji,jj,jk), zwsmax * FLOAT( iiter2 ) ) … … 199 197 zfact = zstep * xdiss(ji,jj,jk) 200 198 ! Part I : Coagulation dependent on turbulence 201 zagg1 = 25.9 * zfact * tr n(ji,jj,jk,jppoc) * trn(ji,jj,jk,jppoc)202 zagg2 = 4452. * zfact * tr n(ji,jj,jk,jppoc) * trn(ji,jj,jk,jpgoc)199 zagg1 = 25.9 * zfact * trb(ji,jj,jk,jppoc) * trb(ji,jj,jk,jppoc) 200 zagg2 = 4452. * zfact * trb(ji,jj,jk,jppoc) * trb(ji,jj,jk,jpgoc) 203 201 204 202 ! Part II : Differential settling 205 203 206 204 ! Aggregation of small into large particles 207 zagg3 = 47.1 * zstep * tr n(ji,jj,jk,jppoc) * trn(ji,jj,jk,jpgoc)208 zagg4 = 3.3 * zstep * tr n(ji,jj,jk,jppoc) * trn(ji,jj,jk,jppoc)205 zagg3 = 47.1 * zstep * trb(ji,jj,jk,jppoc) * trb(ji,jj,jk,jpgoc) 206 zagg4 = 3.3 * zstep * trb(ji,jj,jk,jppoc) * trb(ji,jj,jk,jppoc) 209 207 210 208 zagg = zagg1 + zagg2 + zagg3 + zagg4 211 zaggfe = zagg * tr n(ji,jj,jk,jpsfe) / ( trn(ji,jj,jk,jppoc) + rtrn )209 zaggfe = zagg * trb(ji,jj,jk,jpsfe) / ( trb(ji,jj,jk,jppoc) + rtrn ) 212 210 213 211 ! Aggregation of DOC to POC : … … 215 213 ! 2nd term is shear aggregation of DOC-POC 216 214 ! 3rd term is differential settling of DOC-POC 217 zaggdoc = ( ( 0.369 * 0.3 * tr n(ji,jj,jk,jpdoc) + 102.4 * trn(ji,jj,jk,jppoc) ) * zfact &218 & + 2.4 * zstep * tr n(ji,jj,jk,jppoc) ) * 0.3 * trn(ji,jj,jk,jpdoc)215 zaggdoc = ( ( 0.369 * 0.3 * trb(ji,jj,jk,jpdoc) + 102.4 * trb(ji,jj,jk,jppoc) ) * zfact & 216 & + 2.4 * zstep * trb(ji,jj,jk,jppoc) ) * 0.3 * trb(ji,jj,jk,jpdoc) 219 217 ! transfer of DOC to GOC : 220 218 ! 1st term is shear aggregation 221 219 ! 2nd term is differential settling 222 zaggdoc2 = ( 3.53E3 * zfact + 0.1 * zstep ) * tr n(ji,jj,jk,jpgoc) * 0.3 * trn(ji,jj,jk,jpdoc)220 zaggdoc2 = ( 3.53E3 * zfact + 0.1 * zstep ) * trb(ji,jj,jk,jpgoc) * 0.3 * trb(ji,jj,jk,jpdoc) 223 221 ! tranfer of DOC to POC due to brownian motion 224 zaggdoc3 = ( 5095. * tr n(ji,jj,jk,jppoc) + 114. * 0.3 * trn(ji,jj,jk,jpdoc) ) *zstep * 0.3 * trn(ji,jj,jk,jpdoc)222 zaggdoc3 = ( 5095. * trb(ji,jj,jk,jppoc) + 114. * 0.3 * trb(ji,jj,jk,jpdoc) ) *zstep * 0.3 * trb(ji,jj,jk,jpdoc) 225 223 226 224 ! Update the trends … … 235 233 END DO 236 234 237 ! Total primary production per year 238 t_oce_co2_exp = t_oce_co2_exp + glob_sum( ( sinking(:,:,iksed+1) + sinking2(:,:,iksed+1) ) * e1e2t(:,:) * tmask(:,:,1) ) 235 236 ! Total carbon export per year 237 IF( iom_use( "tcexp" ) .OR. ( ln_check_mass .AND. kt == nitend .AND. knt == nrdttrc ) ) & 238 & t_oce_co2_exp = glob_sum( ( sinking(:,:,ik100) + sinking2(:,:,ik100) ) * e1e2t(:,:) * tmask(:,:,1) ) 239 239 ! 240 IF( ln_diatrc ) THEN 241 zrfact2 = 1.e3 * rfact2r 242 ik1 = iksed + 1 243 IF( lk_iomput ) THEN 244 IF( jnt == nrdttrc ) THEN 245 CALL iom_put( "EPC100" , ( sinking(:,:,ik1) + sinking2(:,:,ik1) ) * zrfact2 * tmask(:,:,1) ) ! Export of carbon at 100m 246 CALL iom_put( "EPFE100" , ( sinkfer(:,:,ik1) + sinkfer2(:,:,ik1) ) * zrfact2 * tmask(:,:,1) ) ! Export of iron at 100m 247 CALL iom_put( "EPCAL100", sinkcal(:,:,ik1) * zrfact2 * tmask(:,:,1) ) ! Export of calcite at 100m 248 CALL iom_put( "EPSI100" , sinksil(:,:,ik1) * zrfact2 * tmask(:,:,1) ) ! Export of biogenic silica at 100m 249 ENDIF 250 ELSE 251 trc2d(:,:,jp_pcs0_2d + 4) = sinking (:,:,ik1) * zrfact2 * tmask(:,:,1) 252 trc2d(:,:,jp_pcs0_2d + 5) = sinking2(:,:,ik1) * zrfact2 * tmask(:,:,1) 253 trc2d(:,:,jp_pcs0_2d + 6) = sinkfer (:,:,ik1) * zrfact2 * tmask(:,:,1) 254 trc2d(:,:,jp_pcs0_2d + 7) = sinkfer2(:,:,ik1) * zrfact2 * tmask(:,:,1) 255 trc2d(:,:,jp_pcs0_2d + 8) = sinksil (:,:,ik1) * zrfact2 * tmask(:,:,1) 256 trc2d(:,:,jp_pcs0_2d + 9) = sinkcal (:,:,ik1) * zrfact2 * tmask(:,:,1) 240 IF( lk_iomput ) THEN 241 IF( knt == nrdttrc ) THEN 242 CALL wrk_alloc( jpi, jpj, zw2d ) 243 CALL wrk_alloc( jpi, jpj, jpk, zw3d ) 244 zfact = 1.e+3 * rfact2r ! conversion from mol/l/kt to mol/m3/s 245 ! 246 IF( iom_use( "EPC100" ) ) THEN 247 zw2d(:,:) = ( sinking(:,:,ik100) + sinking2(:,:,ik100) ) * zfact * tmask(:,:,1) ! Export of carbon at 100m 248 CALL iom_put( "EPC100" , zw2d ) 249 ENDIF 250 IF( iom_use( "EPFE100" ) ) THEN 251 zw2d(:,:) = ( sinkfer(:,:,ik100) + sinkfer2(:,:,ik100) ) * zfact * tmask(:,:,1) ! Export of iron at 100m 252 CALL iom_put( "EPFE100" , zw2d ) 253 ENDIF 254 IF( iom_use( "EPCAL100" ) ) THEN 255 zw2d(:,:) = sinkcal(:,:,ik100) * zfact * tmask(:,:,1) ! Export of calcite at 100m 256 CALL iom_put( "EPCAL100" , zw2d ) 257 ENDIF 258 IF( iom_use( "EPSI100" ) ) THEN 259 zw2d(:,:) = sinksil(:,:,ik100) * zfact * tmask(:,:,1) ! Export of bigenic silica at 100m 260 CALL iom_put( "EPSI100" , zw2d ) 261 ENDIF 262 IF( iom_use( "EXPC" ) ) THEN 263 zw3d(:,:,:) = ( sinking(:,:,:) + sinking2(:,:,:) ) * zfact * tmask(:,:,:) ! Export of carbon in the water column 264 CALL iom_put( "EXPC" , zw3d ) 265 ENDIF 266 IF( iom_use( "EXPFE" ) ) THEN 267 zw3d(:,:,:) = ( sinkfer(:,:,:) + sinkfer2(:,:,:) ) * zfact * tmask(:,:,:) ! Export of iron 268 CALL iom_put( "EXPFE" , zw3d ) 269 ENDIF 270 IF( iom_use( "EXPCAL" ) ) THEN 271 zw3d(:,:,:) = sinkcal(:,:,:) * zfact * tmask(:,:,:) ! Export of calcite 272 CALL iom_put( "EXPCAL" , zw3d ) 273 ENDIF 274 IF( iom_use( "EXPSI" ) ) THEN 275 zw3d(:,:,:) = sinksil(:,:,:) * zfact * tmask(:,:,:) ! Export of bigenic silica 276 CALL iom_put( "EXPSI" , zw3d ) 277 ENDIF 278 IF( iom_use( "tcexp" ) ) CALL iom_put( "tcexp" , t_oce_co2_exp * zfact ) ! molC/s 279 ! 280 CALL wrk_dealloc( jpi, jpj, zw2d ) 281 CALL wrk_dealloc( jpi, jpj, jpk, zw3d ) 282 ENDIF 283 ELSE 284 IF( ln_diatrc ) THEN 285 zfact = 1.e3 * rfact2r 286 trc2d(:,:,jp_pcs0_2d + 4) = sinking (:,:,ik100) * zfact * tmask(:,:,1) 287 trc2d(:,:,jp_pcs0_2d + 5) = sinking2(:,:,ik100) * zfact * tmask(:,:,1) 288 trc2d(:,:,jp_pcs0_2d + 6) = sinkfer (:,:,ik100) * zfact * tmask(:,:,1) 289 trc2d(:,:,jp_pcs0_2d + 7) = sinkfer2(:,:,ik100) * zfact * tmask(:,:,1) 290 trc2d(:,:,jp_pcs0_2d + 8) = sinksil (:,:,ik100) * zfact * tmask(:,:,1) 291 trc2d(:,:,jp_pcs0_2d + 9) = sinkcal (:,:,ik100) * zfact * tmask(:,:,1) 257 292 ENDIF 258 293 ENDIF … … 272 307 !! *** ROUTINE p4z_sink_init *** 273 308 !!---------------------------------------------------------------------- 274 309 INTEGER :: jk 310 311 ik100 = 10 ! last level where depth less than 100 m 312 DO jk = jpkm1, 1, -1 313 IF( gdept_1d(jk) > 100. ) ik100 = jk - 1 314 END DO 315 IF (lwp) WRITE(numout,*) 316 IF (lwp) WRITE(numout,*) ' Level corresponding to 100m depth ', ik100 + 1 317 IF (lwp) WRITE(numout,*) 318 ! 275 319 t_oce_co2_exp = 0._wp 276 320 ! … … 282 326 !!---------------------------------------------------------------------- 283 327 284 SUBROUTINE p4z_sink ( kt, jnt )328 SUBROUTINE p4z_sink ( kt, knt ) 285 329 !!--------------------------------------------------------------------- 286 330 !! *** ROUTINE p4z_sink *** … … 292 336 !!--------------------------------------------------------------------- 293 337 ! 294 INTEGER, INTENT(in) :: kt, jnt338 INTEGER, INTENT(in) :: kt, knt 295 339 ! 296 340 INTEGER :: ji, jj, jk, jit, niter1, niter2 … … 300 344 REAL(wp) :: zdiv , zdiv1, zdiv2, zdiv3, zdiv4, zdiv5 301 345 REAL(wp) :: zval1, zval2, zval3, zval4 302 REAL(wp) :: z rfact2346 REAL(wp) :: zfact 303 347 INTEGER :: ik1 304 348 CHARACTER (len=25) :: charout 305 349 REAL(wp), POINTER, DIMENSION(:,:,:) :: znum3d 350 REAL(wp), POINTER, DIMENSION(:,:,:) :: zw3d 351 REAL(wp), POINTER, DIMENSION(:,: ) :: zw2d 306 352 !!--------------------------------------------------------------------- 307 353 ! … … 325 371 DO ji = 1, jpi 326 372 IF( tmask(ji,jj,jk) /= 0.e0 ) THEN 327 znum = tr n(ji,jj,jk,jppoc) / ( trn(ji,jj,jk,jpnum) + rtrn ) / xkr_massp373 znum = trb(ji,jj,jk,jppoc) / ( trb(ji,jj,jk,jpnum) + rtrn ) / xkr_massp 328 374 ! -------------- To avoid sinking speed over 50 m/day ------- 329 375 znum = MIN( xnumm(jk), znum ) … … 387 433 IF( tmask(ji,jj,jk) /= 0.e0 ) THEN 388 434 389 znum = tr n(ji,jj,jk,jppoc)/(trn(ji,jj,jk,jpnum)+rtrn) / xkr_massp435 znum = trb(ji,jj,jk,jppoc)/(trb(ji,jj,jk,jpnum)+rtrn) / xkr_massp 390 436 !-------------- To avoid sinking speed over 50 m/day ------- 391 437 znum = min(xnumm(jk),znum) … … 405 451 ! ---------------------------------------------- 406 452 407 zagg1 = 0.163 * tr n(ji,jj,jk,jpnum)**2 &453 zagg1 = 0.163 * trb(ji,jj,jk,jpnum)**2 & 408 454 & * 2.*( (zfm-1.)*(zfm*xkr_mass_max**3-xkr_mass_min**3) & 409 455 & * (zeps-1)/zdiv1 + 3.*(zfm*xkr_mass_max-xkr_mass_min) & 410 456 & * (zfm*xkr_mass_max**2-xkr_mass_min**2) & 411 457 & * (zeps-1.)**2/(zdiv2*zdiv3)) 412 zagg2 = 2*0.163*tr n(ji,jj,jk,jpnum)**2*zfm* &458 zagg2 = 2*0.163*trb(ji,jj,jk,jpnum)**2*zfm* & 413 459 & ((xkr_mass_max**3+3.*(xkr_mass_max**2 & 414 460 & *xkr_mass_min*(zeps-1.)/zdiv2 & … … 418 464 & (zeps-2.)+(zeps-1.)/zdiv3)+(zeps-1.)/zdiv1)) 419 465 420 zagg3 = 0.163*tr n(ji,jj,jk,jpnum)**2*zfm**2*8. * xkr_mass_max**3466 zagg3 = 0.163*trb(ji,jj,jk,jpnum)**2*zfm**2*8. * xkr_mass_max**3 421 467 422 468 ! Aggregation of small into large particles … … 424 470 ! ---------------------------------------------- 425 471 426 zagg4 = 2.*3.141*0.125*tr n(ji,jj,jk,jpnum)**2* &472 zagg4 = 2.*3.141*0.125*trb(ji,jj,jk,jpnum)**2* & 427 473 & xkr_wsbio_min*(zeps-1.)**2 & 428 474 & *(xkr_mass_min**2*((1.-zsm*zfm)/(zdiv3*zdiv4) & … … 431 477 & *xkr_eta)/(zdiv*zdiv3*zdiv5) ) 432 478 433 zagg5 = 2.*3.141*0.125*tr n(ji,jj,jk,jpnum)**2 &479 zagg5 = 2.*3.141*0.125*trb(ji,jj,jk,jpnum)**2 & 434 480 & *(zeps-1.)*zfm*xkr_wsbio_min & 435 481 & *(zsm*(xkr_mass_min**2-zfm*xkr_mass_max**2) & … … 441 487 ! ------------------------------------ 442 488 443 zfract = 2.*3.141*0.125*tr n(ji,jj,jk,jpmes)*12./0.12/0.06**3*trn(ji,jj,jk,jpnum) &489 zfract = 2.*3.141*0.125*trb(ji,jj,jk,jpmes)*12./0.12/0.06**3*trb(ji,jj,jk,jpnum) & 444 490 & * (0.01/xkr_mass_min)**(1.-zeps)*0.1**2 & 445 491 & * 10000.*xstep … … 448 494 ! -------------------------------------- 449 495 450 zaggdoc = 0.83 * tr n(ji,jj,jk,jpdoc) * xstep * xdiss(ji,jj,jk) * trn(ji,jj,jk,jpdoc) &451 & + 0.005 * 231. * tr n(ji,jj,jk,jpdoc) * xstep * trn(ji,jj,jk,jpdoc)452 zaggdoc1 = 271. * tr n(ji,jj,jk,jppoc) * xstep * xdiss(ji,jj,jk) * trn(ji,jj,jk,jpdoc) &453 & + 0.02 * 16706. * tr n(ji,jj,jk,jppoc) * xstep * trn(ji,jj,jk,jpdoc)496 zaggdoc = 0.83 * trb(ji,jj,jk,jpdoc) * xstep * xdiss(ji,jj,jk) * trb(ji,jj,jk,jpdoc) & 497 & + 0.005 * 231. * trb(ji,jj,jk,jpdoc) * xstep * trb(ji,jj,jk,jpdoc) 498 zaggdoc1 = 271. * trb(ji,jj,jk,jppoc) * xstep * xdiss(ji,jj,jk) * trb(ji,jj,jk,jpdoc) & 499 & + 0.02 * 16706. * trb(ji,jj,jk,jppoc) * xstep * trb(ji,jj,jk,jpdoc) 454 500 455 501 # if defined key_degrad … … 466 512 zagg = 0.5 * xkr_stick * ( zaggsh + zaggsi ) 467 513 ! 468 znumdoc = tr n(ji,jj,jk,jpnum) / ( trn(ji,jj,jk,jppoc) + rtrn )514 znumdoc = trb(ji,jj,jk,jpnum) / ( trb(ji,jj,jk,jppoc) + rtrn ) 469 515 tra(ji,jj,jk,jppoc) = tra(ji,jj,jk,jppoc) + zaggdoc + zaggdoc1 470 516 tra(ji,jj,jk,jpnum) = tra(ji,jj,jk,jpnum) + zfract + zaggdoc / xkr_massp - zagg … … 477 523 478 524 ! Total primary production per year 479 t_oce_co2_exp = t_oce_co2_exp + glob_sum( ( sinking(:,:, :) ) * cvol(:,:,:) )525 t_oce_co2_exp = t_oce_co2_exp + glob_sum( ( sinking(:,:,ik100) * e1e2t(:,:) * tmask(:,:,1) ) 480 526 ! 481 IF( ln_diatrc ) THEN 482 ! 483 ik1 = iksed + 1 484 zrfact2 = 1.e3 * rfact2r 485 IF( jnt == nrdttrc ) THEN 486 CALL iom_put( "POCFlx" , sinking (:,:,:) * zrfact2 * tmask(:,:,:) ) ! POC export 487 CALL iom_put( "NumFlx" , sinking2 (:,:,:) * zrfact2 * tmask(:,:,:) ) ! Num export 488 CALL iom_put( "SiFlx" , sinksil (:,:,:) * zrfact2 * tmask(:,:,:) ) ! Silica export 489 CALL iom_put( "CaCO3Flx", sinkcal (:,:,:) * zrfact2 * tmask(:,:,:) ) ! Calcite export 490 CALL iom_put( "xnum" , znum3d (:,:,:) * tmask(:,:,:) ) ! Number of particles in aggregats 491 CALL iom_put( "W1" , wsbio3 (:,:,:) * tmask(:,:,:) ) ! sinking speed of POC 492 CALL iom_put( "W2" , wsbio4 (:,:,:) * tmask(:,:,:) ) ! sinking speed of aggregats 527 IF( lk_iomput ) THEN 528 IF( knt == nrdttrc ) THEN 529 CALL wrk_alloc( jpi, jpj, zw2d ) 530 CALL wrk_alloc( jpi, jpj, jpk, zw3d ) 531 zfact = 1.e+3 * rfact2r ! conversion from mol/l/kt to mol/m3/s 532 ! 533 IF( iom_use( "EPC100" ) ) THEN 534 zw2d(:,:) = sinking(:,:,ik100) * zfact * tmask(:,:,1) ! Export of carbon at 100m 535 CALL iom_put( "EPC100" , zw2d ) 536 ENDIF 537 IF( iom_use( "EPN100" ) ) THEN 538 zw2d(:,:) = sinking2(:,:,ik100) * zfact * tmask(:,:,1) ! Export of number of aggregates ? 539 CALL iom_put( "EPN100" , zw2d ) 540 ENDIF 541 IF( iom_use( "EPCAL100" ) ) THEN 542 zw2d(:,:) = sinkcal(:,:,ik100) * zfact * tmask(:,:,1) ! Export of calcite at 100m 543 CALL iom_put( "EPCAL100" , zw2d ) 544 ENDIF 545 IF( iom_use( "EPSI100" ) ) THEN 546 zw2d(:,:) = sinksil(:,:,ik100) * zfact * tmask(:,:,1) ! Export of bigenic silica at 100m 547 CALL iom_put( "EPSI100" , zw2d ) 548 ENDIF 549 IF( iom_use( "EXPC" ) ) THEN 550 zw3d(:,:,:) = sinking(:,:,:) * zfact * tmask(:,:,:) ! Export of carbon in the water column 551 CALL iom_put( "EXPC" , zw3d ) 552 ENDIF 553 IF( iom_use( "EXPN" ) ) THEN 554 zw3d(:,:,:) = sinking(:,:,:) * zfact * tmask(:,:,:) ! Export of carbon in the water column 555 CALL iom_put( "EXPN" , zw3d ) 556 ENDIF 557 IF( iom_use( "EXPCAL" ) ) THEN 558 zw3d(:,:,:) = sinkcal(:,:,:) * zfact * tmask(:,:,:) ! Export of calcite 559 CALL iom_put( "EXPCAL" , zw3d ) 560 ENDIF 561 IF( iom_use( "EXPSI" ) ) THEN 562 zw3d(:,:,:) = sinksil(:,:,:) * zfact * tmask(:,:,:) ! Export of bigenic silica 563 CALL iom_put( "EXPSI" , zw3d ) 564 ENDIF 565 IF( iom_use( "XNUM" ) ) THEN 566 zw3d(:,:,:) = znum3d(:,:,:) * tmask(:,:,:) ! Number of particles on aggregats 567 CALL iom_put( "XNUM" , zw3d ) 568 ENDIF 569 IF( iom_use( "WSC" ) ) THEN 570 zw3d(:,:,:) = wsbio3(:,:,:) * tmask(:,:,:) ! Sinking speed of carbon particles 571 CALL iom_put( "WSC" , zw3d ) 572 ENDIF 573 IF( iom_use( "WSN" ) ) THEN 574 zw3d(:,:,:) = wsbio4(:,:,:) * tmask(:,:,:) ! Sinking speed of particles number 575 CALL iom_put( "WSN" , zw3d ) 576 ENDIF 577 ! 578 CALL wrk_dealloc( jpi, jpj, zw2d ) 579 CALL wrk_dealloc( jpi, jpj, jpk, zw3d ) 580 ELSE 581 IF( ln_diatrc ) THEN 582 zfact = 1.e3 * rfact2r 583 trc2d(:,: ,jp_pcs0_2d + 4) = sinking (:,:,ik100) * zfact * tmask(:,:,1) 584 trc2d(:,: ,jp_pcs0_2d + 5) = sinking2(:,:,ik100) * zfact * tmask(:,:,1) 585 trc2d(:,: ,jp_pcs0_2d + 6) = sinkfer (:,:,ik100) * zfact * tmask(:,:,1) 586 trc2d(:,: ,jp_pcs0_2d + 7) = sinksil (:,:,ik100) * zfact * tmask(:,:,1) 587 trc2d(:,: ,jp_pcs0_2d + 8) = sinkcal (:,:,ik100) * zfact * tmask(:,:,1) 588 trc3d(:,:,:,jp_pcs0_3d + 11) = sinking (:,:,:) * zfact * tmask(:,:,:) 589 trc3d(:,:,:,jp_pcs0_3d + 12) = sinking2(:,:,:) * zfact * tmask(:,:,:) 590 trc3d(:,:,:,jp_pcs0_3d + 13) = sinksil (:,:,:) * zfact * tmask(:,:,:) 591 trc3d(:,:,:,jp_pcs0_3d + 14) = sinkcal (:,:,:) * zfact * tmask(:,:,:) 592 trc3d(:,:,:,jp_pcs0_3d + 15) = znum3d (:,:,:) * tmask(:,:,:) 593 trc3d(:,:,:,jp_pcs0_3d + 16) = wsbio3 (:,:,:) * tmask(:,:,:) 594 trc3d(:,:,:,jp_pcs0_3d + 17) = wsbio4 (:,:,:) * tmask(:,:,:) 493 595 ENDIF 494 # if ! defined key_iomput495 trc2d(:,: ,jp_pcs0_2d + 4) = sinking (:,:,ik1) * zrfact2 * tmask(:,:,1)496 trc2d(:,: ,jp_pcs0_2d + 5) = sinking2(:,:,ik1) * zrfact2 * tmask(:,:,1)497 trc2d(:,: ,jp_pcs0_2d + 6) = sinkfer (:,:,ik1) * zrfact2 * tmask(:,:,1)498 trc2d(:,: ,jp_pcs0_2d + 7) = sinksil (:,:,ik1) * zrfact2 * tmask(:,:,1)499 trc2d(:,: ,jp_pcs0_2d + 8) = sinkcal (:,:,ik1) * zrfact2 * tmask(:,:,1)500 trc3d(:,:,:,jp_pcs0_3d + 11) = sinking (:,:,:) * zrfact2 * tmask(:,:,:)501 trc3d(:,:,:,jp_pcs0_3d + 12) = sinking2(:,:,:) * zrfact2 * tmask(:,:,:)502 trc3d(:,:,:,jp_pcs0_3d + 13) = sinksil (:,:,:) * zrfact2 * tmask(:,:,:)503 trc3d(:,:,:,jp_pcs0_3d + 14) = sinkcal (:,:,:) * zrfact2 * tmask(:,:,:)504 trc3d(:,:,:,jp_pcs0_3d + 15) = znum3d (:,:,:) * tmask(:,:,:)505 trc3d(:,:,:,jp_pcs0_3d + 16) = wsbio3 (:,:,:) * tmask(:,:,:)506 trc3d(:,:,:,jp_pcs0_3d + 17) = wsbio4 (:,:,:) * tmask(:,:,:)507 # endif508 !509 596 ENDIF 597 510 598 ! 511 599 IF(ln_ctl) THEN ! print mean trends (used for debugging) … … 610 698 zl = zmin 611 699 zr = zmax 612 wmax = 0.5 * fse3t(1,1,jk) * rday * float(niter1max) / rfact2700 wmax = 0.5 * e3t_n(1,1,jk) * rday * float(niter1max) / rfact2 613 701 zdiv = xkr_zeta + xkr_eta - xkr_eta * zl 614 702 znum = zl - 1. … … 663 751 END DO 664 752 ! 753 ik100 = 10 ! last level where depth less than 100 m 754 DO jk = jpkm1, 1, -1 755 IF( gdept_1d(jk) > 100. ) iksed = jk - 1 756 END DO 757 IF (lwp) WRITE(numout,*) 758 IF (lwp) WRITE(numout,*) ' Level corresponding to 100m depth ', ik100 + 1 759 IF (lwp) WRITE(numout,*) 760 ! 665 761 t_oce_co2_exp = 0._wp 666 762 ! … … 702 798 ztraz(:,:,:) = 0.e0 703 799 zakz (:,:,:) = 0.e0 704 ztrb (:,:,:) = tr n(:,:,:,jp_tra)800 ztrb (:,:,:) = trb(:,:,:,jp_tra) 705 801 706 802 DO jk = 1, jpkm1 … … 717 813 ! first guess of the slopes interior values 718 814 DO jk = 2, jpkm1 719 ztraz(:,:,jk) = ( tr n(:,:,jk-1,jp_tra) - trn(:,:,jk,jp_tra) ) * tmask(:,:,jk)815 ztraz(:,:,jk) = ( trb(:,:,jk-1,jp_tra) - trb(:,:,jk,jp_tra) ) * tmask(:,:,jk) 720 816 END DO 721 817 ztraz(:,:,1 ) = 0.0 … … 746 842 DO jj = 1, jpj 747 843 DO ji = 1, jpi 748 zigma = zwsink2(ji,jj,jk+1) * zstep / fse3w(ji,jj,jk+1)844 zigma = zwsink2(ji,jj,jk+1) * zstep / e3w_n(ji,jj,jk+1) 749 845 zew = zwsink2(ji,jj,jk+1) 750 psinkflx(ji,jj,jk+1) = -zew * ( tr n(ji,jj,jk,jp_tra) - 0.5 * ( 1 + zigma ) * zakz(ji,jj,jk) ) * zstep846 psinkflx(ji,jj,jk+1) = -zew * ( trb(ji,jj,jk,jp_tra) - 0.5 * ( 1 + zigma ) * zakz(ji,jj,jk) ) * zstep 751 847 END DO 752 848 END DO … … 760 856 DO jj = 1,jpj 761 857 DO ji = 1, jpi 762 zflx = ( psinkflx(ji,jj,jk) - psinkflx(ji,jj,jk+1) ) / fse3t(ji,jj,jk)763 tr n(ji,jj,jk,jp_tra) = trn(ji,jj,jk,jp_tra) + zflx858 zflx = ( psinkflx(ji,jj,jk) - psinkflx(ji,jj,jk+1) ) / e3t_n(ji,jj,jk) 859 trb(ji,jj,jk,jp_tra) = trb(ji,jj,jk,jp_tra) + zflx 764 860 END DO 765 861 END DO … … 771 867 DO jj = 1,jpj 772 868 DO ji = 1, jpi 773 zflx = ( psinkflx(ji,jj,jk) - psinkflx(ji,jj,jk+1) ) / fse3t(ji,jj,jk)869 zflx = ( psinkflx(ji,jj,jk) - psinkflx(ji,jj,jk+1) ) / e3t_n(ji,jj,jk) 774 870 ztrb(ji,jj,jk) = ztrb(ji,jj,jk) + 2. * zflx 775 871 END DO … … 777 873 END DO 778 874 779 tr n(:,:,:,jp_tra) = ztrb(:,:,:)875 trb(:,:,:,jp_tra) = ztrb(:,:,:) 780 876 psinkflx(:,:,:) = 2. * psinkflx(:,:,:) 781 877 ! … … 815 911 816 912 !!====================================================================== 817 END MODULE 913 END MODULE p4zsink
Note: See TracChangeset
for help on using the changeset viewer.
