Changeset 2528 for trunk/NEMOGCM/NEMO/TOP_SRC/PISCES/p4zsink.F90
- Timestamp:
- 2010-12-27T18:33:53+01:00 (13 years ago)
- File:
-
- 1 edited
Legend:
- Unmodified
- Added
- Removed
-
trunk/NEMOGCM/NEMO/TOP_SRC/PISCES/p4zsink.F90
- Property svn:executable deleted
r1836 r2528 19 19 PRIVATE 20 20 21 PUBLIC p4z_sink ! called in p4zbio.F90 21 PUBLIC p4z_sink ! called in p4zbio.F90 22 PUBLIC p4z_sink_init ! called in trcsms_pisces.F90 22 23 23 24 !! * Shared module variables … … 31 32 sinkcal, sinksil, & !: CaCO3 and BSi sinking fluxes 32 33 sinkfer !: Small BFe sinking flux 33 34 REAL(wp) :: &35 xstep , xstep2 !: Time step duration for biology36 34 37 35 INTEGER :: & … … 71 69 # include "top_substitute.h90" 72 70 !!---------------------------------------------------------------------- 73 !! NEMO/TOP 2.0 , LOCEAN-IPSL (2007)71 !! NEMO/TOP 3.3 , NEMO Consortium (2010) 74 72 !! $Id$ 75 !! Software governed by the CeCILL licence ( modipsl/doc/NEMO_CeCILL.txt)73 !! Software governed by the CeCILL licence (NEMOGCM/NEMO_CeCILL.txt) 76 74 !!---------------------------------------------------------------------- 77 75 … … 97 95 REAL(wp) :: zdiv , zdiv1, zdiv2, zdiv3, zdiv4, zdiv5 98 96 REAL(wp) :: zval1, zval2, zval3, zval4 99 #if defined key_ trc_diaadd97 #if defined key_diatrc 100 98 REAL(wp) :: zrfact2 101 99 INTEGER :: ik1 … … 106 104 !!--------------------------------------------------------------------- 107 105 108 IF( ( kt * jnt ) == nittrc000 ) THEN 109 CALL p4z_sink_init ! Initialization (first time-step only) 110 xstep = rfact2 / rday ! Time step duration for biology 111 xstep2 = rfact2 / 2. 112 ENDIF 113 114 ! Initialisation of variables used to compute Sinking Speed 115 ! --------------------------------------------------------- 106 ! Initialisation of variables used to compute Sinking Speed 107 ! --------------------------------------------------------- 116 108 117 109 znum3d(:,:,:) = 0.e0 … … 120 112 zval3 = 1. + xkr_eta 121 113 122 ! Computation of the vertical sinking speed : Kriest et Evans, 2000123 ! -----------------------------------------------------------------114 ! Computation of the vertical sinking speed : Kriest et Evans, 2000 115 ! ----------------------------------------------------------------- 124 116 125 117 DO jk = 1, jpkm1 … … 128 120 IF( tmask(ji,jj,jk) /= 0.e0 ) THEN 129 121 znum = trn(ji,jj,jk,jppoc) / ( trn(ji,jj,jk,jpnum) + rtrn ) / xkr_massp 130 ! -------------- To avoid sinking speed over 50 m/day -------122 ! -------------- To avoid sinking speed over 50 m/day ------- 131 123 znum = MIN( xnumm(jk), znum ) 132 124 znum = MAX( 1.1 , znum ) 133 125 znum3d(ji,jj,jk) = znum 134 !------------------------------------------------------------126 !------------------------------------------------------------ 135 127 zeps = ( zval1 * znum - 1. )/ ( znum - 1. ) 136 128 zfm = xkr_frac**( 1. - zeps ) … … 150 142 wscal(:,:,:) = MAX( wsbio3(:,:,:), 50. ) 151 143 152 153 ! INITIALIZE TO ZERO ALL THE SINKING ARRAYS 154 ! ----------------------------------------- 144 ! INITIALIZE TO ZERO ALL THE SINKING ARRAYS 145 ! ----------------------------------------- 155 146 156 147 sinking (:,:,:) = 0.e0 … … 160 151 sinksil (:,:,:) = 0.e0 161 152 162 ! Compute the sedimentation term using p4zsink2 for all 163 ! the sinking particles 164 ! ----------------------------------------------------- 153 ! Compute the sedimentation term using p4zsink2 for all the sinking particles 154 ! ----------------------------------------------------- 165 155 166 156 CALL p4z_sink2( wsbio3, sinking , jppoc ) … … 170 160 CALL p4z_sink2( wscal , sinkcal , jpcal ) 171 161 172 ! Exchange between organic matter compartments due to 173 ! coagulation/disaggregation 174 ! --------------------------------------------------- 162 ! Exchange between organic matter compartments due to coagulation/disaggregation 163 ! --------------------------------------------------- 175 164 176 165 zval1 = 1. + xkr_zeta … … 185 174 186 175 znum = trn(ji,jj,jk,jppoc)/(trn(ji,jj,jk,jpnum)+rtrn) / xkr_massp 187 !-------------- To avoid sinking speed over 50 m/day -------176 !-------------- To avoid sinking speed over 50 m/day ------- 188 177 znum = min(xnumm(jk),znum) 189 178 znum = MAX( 1.1,znum) 190 !------------------------------------------------------------179 !------------------------------------------------------------ 191 180 zeps = ( zval1 * znum - 1.) / ( znum - 1.) 192 181 zdiv = MAX( 1.e-4, ABS( zeps - zval3) ) * SIGN( 1., zeps - zval3 ) … … 199 188 zsm = xkr_frac**xkr_eta 200 189 201 ! Part I : Coagulation dependant on turbulence202 ! ----------------------------------------------190 ! Part I : Coagulation dependant on turbulence 191 ! ---------------------------------------------- 203 192 204 193 zagg1 = ( 0.163 * trn(ji,jj,jk,jpnum)**2 & … … 207 196 & * (zfm*xkr_mass_max**2-xkr_mass_min**2) & 208 197 & * (zeps-1.)**2/(zdiv2*zdiv3)) & 209 # if defined key_ off_degrad198 # if defined key_degrad 210 199 & *facvol(ji,jj,jk) & 211 200 # endif … … 219 208 & -zfm*xkr_mass_max**3*(1.+3.*((zeps-1.)/ & 220 209 & (zeps-2.)+(zeps-1.)/zdiv3)+(zeps-1.)/zdiv1)) & 221 # if defined key_ off_degrad210 # if defined key_degrad 222 211 & *facvol(ji,jj,jk) & 223 212 # endif … … 225 214 226 215 zagg3 = ( 0.163*trn(ji,jj,jk,jpnum)**2*zfm**2*8. * xkr_mass_max**3 & 227 # if defined key_ off_degrad216 # if defined key_degrad 228 217 & *facvol(ji,jj,jk) & 229 218 # endif … … 232 221 zaggsh = ( zagg1 + zagg2 + zagg3 ) * rfact2 * xdiss(ji,jj,jk) / 1000. 233 222 234 ! Aggregation of small into large particles235 ! Part II : Differential settling236 ! ----------------------------------------------223 ! Aggregation of small into large particles 224 ! Part II : Differential settling 225 ! ---------------------------------------------- 237 226 238 227 zagg4 = ( 2.*3.141*0.125*trn(ji,jj,jk,jpnum)**2* & … … 242 231 & ((zfm*zfm*xkr_mass_max**2*zsm-xkr_mass_min**2) & 243 232 & *xkr_eta)/(zdiv*zdiv3*zdiv5) ) & 244 # if defined key_ off_degrad233 # if defined key_degrad 245 234 & *facvol(ji,jj,jk) & 246 235 # endif … … 252 241 & /zdiv3-(xkr_mass_min**2-zfm*zsm*xkr_mass_max**2) & 253 242 & /zdiv) & 254 # if defined key_ off_degrad243 # if defined key_degrad 255 244 & *facvol(ji,jj,jk) & 256 245 # endif … … 261 250 zagg = 0.5 * xkr_stick * ( zaggsh + zaggsi ) 262 251 263 ! Aggregation of DOC to small particles264 ! --------------------------------------252 ! Aggregation of DOC to small particles 253 ! -------------------------------------- 265 254 266 255 zaggdoc = ( 0.4 * trn(ji,jj,jk,jpdoc) & 267 256 & + 1018. * trn(ji,jj,jk,jppoc) ) * xstep & 268 # if defined key_ off_degrad257 # if defined key_degrad 269 258 & * facvol(ji,jj,jk) & 270 259 # endif … … 281 270 END DO 282 271 283 #if defined key_ trc_diaadd272 #if defined key_diatrc 284 273 zrfact2 = 1.e3 * rfact2r 285 274 ik1 = iksed + 1 … … 332 321 !! 333 322 !! ** Method : Read the nampiskrs namelist and check the parameters 334 !! called at the first timestep (nittrc000)323 !! called at the first timestep 335 324 !! 336 325 !! ** input : Namelist nampiskrs … … 473 462 REAL(wp) :: zagg1, zagg2, zagg3, zagg4 474 463 REAL(wp) :: zagg , zaggfe, zaggdoc, zaggdoc2 475 REAL(wp) :: zfact, zwsmax 476 #if defined key_ trc_dia3d464 REAL(wp) :: zfact, zwsmax, zstep 465 #if defined key_diatrc 477 466 REAL(wp) :: zrfact2 478 467 INTEGER :: ik1 … … 481 470 !!--------------------------------------------------------------------- 482 471 483 IF( ( kt * jnt ) == nittrc000 ) THEN 484 xstep = rfact2 / rday ! Timestep duration for biology 485 xstep2 = rfact2 / 2. 486 ENDIF 487 488 ! Sinking speeds of detritus is increased with depth as shown 489 ! by data and from the coagulation theory 490 ! ----------------------------------------------------------- 472 ! Sinking speeds of detritus is increased with depth as shown 473 ! by data and from the coagulation theory 474 ! ----------------------------------------------------------- 491 475 DO jk = 1, jpkm1 492 476 DO jj = 1, jpj 493 477 DO ji=1,jpi 494 zfact = MAX( 0., fsdepw(ji,jj,jk+1) -hmld(ji,jj) ) / 4000.478 zfact = MAX( 0., fsdepw(ji,jj,jk+1) - hmld(ji,jj) ) / 4000. 495 479 wsbio4(ji,jj,jk) = wsbio2 + ( 200.- wsbio2 ) * zfact 496 480 END DO … … 498 482 END DO 499 483 500 ! LIMIT THE VALUES OF THE SINKING SPEEDS 501 ! TO AVOID NUMERICAL INSTABILITIES 502 484 ! limit the values of the sinking speeds to avoid numerical instabilities 503 485 wsbio3(:,:,:) = wsbio 504 !505 ! OA Below, this is garbage. the ideal would be to find a time-splitting 506 ! OA algorithm that does not increase the computing cost by too much507 ! OA In ROMS, I have included a time-splitting procedure. But it is508 ! OA too expensive as the loop is computed globally. Thus, a small e3t509 ! OA at one place determines the number of subtimesteps globally510 ! OA AWFULLY EXPENSIVE !! Not able to find a better approach. Damned !!486 ! 487 ! OA Below, this is garbage. the ideal would be to find a time-splitting 488 ! OA algorithm that does not increase the computing cost by too much 489 ! OA In ROMS, I have included a time-splitting procedure. But it is 490 ! OA too expensive as the loop is computed globally. Thus, a small e3t 491 ! OA at one place determines the number of subtimesteps globally 492 ! OA AWFULLY EXPENSIVE !! Not able to find a better approach. Damned !! 511 493 512 494 DO jk = 1,jpkm1 … … 522 504 wscal(:,:,:) = wsbio4(:,:,:) 523 505 524 ! INITIALIZE TO ZERO ALL THE SINKING ARRAYS 525 ! -----------------------------------------506 ! Initializa to zero all the sinking arrays 507 ! ----------------------------------------- 526 508 527 509 sinking (:,:,:) = 0.e0 … … 532 514 sinkfer2(:,:,:) = 0.e0 533 515 534 ! Compute the sedimentation term using p4zsink2 for all 535 ! the sinking particles 536 ! ----------------------------------------------------- 516 ! Compute the sedimentation term using p4zsink2 for all the sinking particles 517 ! ----------------------------------------------------- 537 518 538 519 CALL p4z_sink2( wsbio3, sinking , jppoc ) … … 543 524 CALL p4z_sink2( wscal , sinkcal , jpcal ) 544 525 545 ! Exchange between organic matter compartments due to 546 ! coagulation/disaggregation 547 ! --------------------------------------------------- 526 ! Exchange between organic matter compartments due to coagulation/disaggregation 527 ! --------------------------------------------------- 548 528 549 529 DO jk = 1, jpkm1 550 530 DO jj = 1, jpj 551 531 DO ji = 1, jpi 552 zfact = xstep * xdiss(ji,jj,jk) 532 # if defined key_degrad 533 zstep = xstep * facvol(ji,jj,jk) 534 # else 535 zstep = xstep 536 # endif 537 zfact = zstep * xdiss(ji,jj,jk) 553 538 ! Part I : Coagulation dependent on turbulence 554 # if defined key_off_degrad555 zagg1 = 940.* zfact * trn(ji,jj,jk,jppoc) * trn(ji,jj,jk,jppoc) * facvol(ji,jj,jk)556 zagg2 = 1.054e4 * zfact * trn(ji,jj,jk,jppoc) * trn(ji,jj,jk,jpgoc) * facvol(ji,jj,jk)557 # else558 539 zagg1 = 940.* zfact * trn(ji,jj,jk,jppoc) * trn(ji,jj,jk,jppoc) 559 540 zagg2 = 1.054e4 * zfact * trn(ji,jj,jk,jppoc) * trn(ji,jj,jk,jpgoc) 560 # endif561 541 562 542 ! Part II : Differential settling 563 543 564 544 ! Aggregation of small into large particles 565 # if defined key_off_degrad 566 zagg3 = 0.66 * xstep * trn(ji,jj,jk,jppoc) * trn(ji,jj,jk,jpgoc) * facvol(ji,jj,jk) 567 zagg4 = 0.e0 * xstep * trn(ji,jj,jk,jppoc) * trn(ji,jj,jk,jppoc) * facvol(ji,jj,jk) 568 # else 569 zagg3 = 0.66 * xstep * trn(ji,jj,jk,jppoc) * trn(ji,jj,jk,jpgoc) 570 zagg4 = 0.e0 * xstep * trn(ji,jj,jk,jppoc) * trn(ji,jj,jk,jppoc) 571 # endif 545 zagg3 = 0.66 * zstep * trn(ji,jj,jk,jppoc) * trn(ji,jj,jk,jpgoc) 546 zagg4 = 0.e0 * zstep * trn(ji,jj,jk,jppoc) * trn(ji,jj,jk,jppoc) 547 572 548 zagg = zagg1 + zagg2 + zagg3 + zagg4 573 549 zaggfe = zagg * trn(ji,jj,jk,jpsfe) / ( trn(ji,jj,jk,jppoc) + rtrn ) 574 550 575 551 ! Aggregation of DOC to small particles 576 #if defined key_off_degrad 577 zaggdoc = ( 80.* trn(ji,jj,jk,jpdoc) + 698. * trn(ji,jj,jk,jppoc) ) & 578 & * facvol(ji,jj,jk) * zfact * trn(ji,jj,jk,jpdoc) 579 zaggdoc2 = 1.05e4 * zfact * trn(ji,jj,jk,jpgoc) & 580 & * facvol(ji,jj,jk) * trn(ji,jj,jk,jpdoc) 581 #else 582 zaggdoc = ( 80.* trn(ji,jj,jk,jpdoc) + 698. * trn(ji,jj,jk,jppoc) ) & 583 & * zfact * trn(ji,jj,jk,jpdoc) 552 zaggdoc = ( 80.* trn(ji,jj,jk,jpdoc) + 698. * trn(ji,jj,jk,jppoc) ) * zfact * trn(ji,jj,jk,jpdoc) 584 553 zaggdoc2 = 1.05e4 * zfact * trn(ji,jj,jk,jpgoc) * trn(ji,jj,jk,jpdoc) 585 #endif 554 586 555 ! Update the trends 587 556 tra(ji,jj,jk,jppoc) = tra(ji,jj,jk,jppoc) - zagg + zaggdoc … … 595 564 END DO 596 565 597 #if defined key_ trc_diaadd566 #if defined key_diatrc 598 567 zrfact2 = 1.e3 * rfact2r 599 568 ik1 = iksed + 1 … … 623 592 END SUBROUTINE p4z_sink 624 593 594 SUBROUTINE p4z_sink_init 595 !!---------------------------------------------------------------------- 596 !! *** ROUTINE p4z_sink_init *** 597 !!---------------------------------------------------------------------- 598 END SUBROUTINE p4z_sink_init 599 625 600 #endif 626 601 … … 641 616 !! 642 617 INTEGER :: ji, jj, jk, jn 643 REAL(wp) :: zigma,zew,zign, zflx 618 REAL(wp) :: zigma,zew,zign, zflx, zstep 644 619 REAL(wp), DIMENSION(jpi,jpj,jpk) :: ztraz, zakz 645 620 REAL(wp), DIMENSION(jpi,jpj,jpk) :: zwsink2 … … 647 622 648 623 624 zstep = rfact2 / 2. 625 649 626 ztraz(:,:,:) = 0.e0 650 627 zakz (:,:,:) = 0.e0 651 628 652 629 DO jk = 1, jpkm1 653 # if defined key_ off_degrad630 # if defined key_degrad 654 631 zwsink2(:,:,jk+1) = -pwsink(:,:,jk) / rday * tmask(:,:,jk+1) * facvol(:,:,jk) 655 632 # else … … 693 670 DO jj = 1, jpj 694 671 DO ji = 1, jpi 695 zigma = zwsink2(ji,jj,jk+1) * xstep2/ fse3w(ji,jj,jk+1)672 zigma = zwsink2(ji,jj,jk+1) * zstep / fse3w(ji,jj,jk+1) 696 673 zew = zwsink2(ji,jj,jk+1) 697 psinkflx(ji,jj,jk+1) = -zew * ( trn(ji,jj,jk,jp_tra) - 0.5 * ( 1 + zigma ) * zakz(ji,jj,jk) ) * xstep2674 psinkflx(ji,jj,jk+1) = -zew * ( trn(ji,jj,jk,jp_tra) - 0.5 * ( 1 + zigma ) * zakz(ji,jj,jk) ) * zstep 698 675 END DO 699 676 END DO
Note: See TracChangeset
for help on using the changeset viewer.
