New URL for NEMO forge!   http://forge.nemo-ocean.eu

Since March 2022 along with NEMO 4.2 release, the code development moved to a self-hosted GitLab.
This present forge is now archived and remained online for history.
Changeset 3446 for branches/2012/dev_r3438_LOCEAN15_PISLOB/NEMOGCM/NEMO/TOP_SRC/PISCES/P4Z/p4zfechem.F90 – NEMO

Ignore:
Timestamp:
2012-08-10T13:13:55+02:00 (12 years ago)
Author:
cetlod
Message:

branch:2012/dev_r3438_LOCEAN15_PISLOB : 2nd step new PISCES updates from Olivier, see ticket #972

File:
1 edited

Legend:

Unmodified
Added
Removed
  • branches/2012/dev_r3438_LOCEAN15_PISLOB/NEMOGCM/NEMO/TOP_SRC/PISCES/P4Z/p4zfechem.F90

    r3443 r3446  
    3737   REAL(wp), PUBLIC ::  ligand     = 0.6E-9_wp  !: ligand concentration in the ocean  
    3838 
    39    REAL(wp) :: kl1, kl2, kb1, kb2, ks, kpr, spd, con 
     39   REAL(wp) :: kl1, kl2, kb1, kb2, ks, kpr, spd, con, kth 
    4040 
    4141   !!* Substitution 
     
    6565      ! 
    6666      INTEGER  ::   ji, jj, jk, jic 
    67       REAL(wp) ::   zdep, zlam1b, zlamfac 
    68       REAL(wp) ::   zkeq, zfeequi, zfesatur 
    69       REAL(wp) ::   zdenom1, zscave, zaggdfe, zcoag 
     67      REAL(wp) ::   zdep, zlam1a, zlam1b, zlamfac 
     68      REAL(wp) ::   zkeq, zfeequi, zfesatur, zfecoll 
     69      REAL(wp) ::   zdenom1, zscave, zaggdfea, zaggdfeb, zcoag 
    7070      REAL(wp) ::   ztrc, zdust 
    7171#if ! defined key_kriest 
     
    7474      REAL(wp), POINTER, DIMENSION(:,:,:) :: zTL1, zFe3, ztotlig 
    7575      REAL(wp), POINTER, DIMENSION(:,:,:) :: zFeL1, zFeL2, zTL2, zFe2, zFeP 
    76       REAL(wp) :: zkox, zkph1, zkph2, zph, zionic 
     76      REAL(wp) :: zkox, zkph1, zkph2, zph, zionic, ztligand 
    7777      REAL(wp) :: za, zb, zc, zkappa1, zkappa2, za0, za1, za2 
    7878      REAL(wp) :: zxs, zfunc, zp, zq, zd, zr, zphi, zfff, zp3, zq2 
     
    8585      ! 
    8686      ! Allocate temporary workspace 
    87       CALL wrk_alloc( jpi, jpj, jpk, zFe3, zTL1, ztotlig ) 
    88       zFe3   (:,:,:) = 0._wp 
    89       zTL1   (:,:,:) = 0._wp 
    90       ztotlig(:,:,:) = 0._wp 
     87                                          CALL wrk_alloc( jpi, jpj, jpk, zFe3, zTL1, ztotlig ) 
     88      IF( ln_diatrc .AND. lk_iomput )     CALL wrk_alloc( jpi, jpj, jpk, zFeL1 ) 
     89      ! 
    9190      IF( ln_fechem ) THEN 
    92          CALL wrk_alloc( jpi, jpj, jpk, zTL2, zFeP ) 
    93          zTL2(:,:,:) = 0._wp 
    94          zFeP(:,:,:) = 0._wp 
    95          IF( ln_diatrc .AND. lk_iomput ) THEN 
    96             CALL wrk_alloc( jpi, jpj, jpk, zFe2, zFeL1, zFeL2 ) 
    97             zFe2 (:,:,:) = 0._wp 
    98             zFeL1(:,:,:) = 0._wp 
    99             zFeL2(:,:,:) = 0._wp 
    100          ENDIF 
     91                                          CALL wrk_alloc( jpi, jpj, jpk, zTL2, zFeP ) 
     92         IF( ln_diatrc .AND. lk_iomput )  CALL wrk_alloc( jpi, jpj, jpk, zFe2, zFeL2 ) 
    10193      ENDIF 
    10294      ! Total ligand concentration : Ligands can be chosen to be constant or variable 
     
    10496      ! ------------------------------------------------- 
    10597      IF( ln_ligvar ) THEN 
    106          ztotlig(:,:,:) =  0.09 * ( trn(:,:,:,jpdoc) * 1E6 + 40. ) - 3.2 
    107          ztotlig(:,:,:) = MAX( MIN( ztotlig(:,:,:), 10. ), 0.4 ) 
     98         ztotlig(:,:,:) =  0.09 * trn(:,:,:,jpdoc) * 1E6 + ligand * 1E9 
     99         ztotlig(:,:,:) =  MIN( ztotlig(:,:,:), 10. ) 
    108100      ELSE 
    109101         ztotlig(:,:,:) = ligand * 1E9 
     
    120112               DO ji = 1, jpi 
    121113                  ! Calculate ligand concentrations : assume 2/3rd of excess goes to L2 and 1/3rd to L1 
    122                   zTL2(ji,jj,jk) = 0.000001       + 0.67 * ( ztotlig(ji,jj,jk) - 0.4 ) 
    123                   zTL1(ji,jj,jk) = 0.4 - 0.000001 + 0.33 * ( ztotlig(ji,jj,jk) - 0.4 ) 
     114                  ztligand       = ztotlig(ji,jj,jk) - ligand * 1E9 
     115                  zTL2(ji,jj,jk) =                0.000001 + 0.67 * ztligand 
     116                  zTL1(ji,jj,jk) = ligand * 1E9 - 0.000001 + 0.33 * ztligand 
    124117                  ! ionic strength from Millero et al. 1987 
    125118                  zionic = 19.9201 * tsn(ji,jj,jk,jp_sal) / ( 1000. - 1.00488 * tsn(ji,jj,jk,jp_sal) + rtrn ) 
     
    150143                  ! calculate some parameters 
    151144                  za = 1 + ks / kpr 
    152                   zb = 1 + zkph1 / ( zkox + rtrn ) 
     145                  zb = 1 + ( zkph1 + kth ) / ( zkox + rtrn ) 
    153146                  zc = 1 + zkph2 / ( zkox + rtrn ) 
    154                   zkappa1 = ( kb1 + zkph1) / kl1 
    155                   zkappa2 = ( kb2 + zkph2) / kl2 
     147                  zkappa1 = ( kb1 + zkph1 + kth ) / kl1 
     148                  zkappa2 = ( kb2 + zkph2 ) / kl2 
    156149                  za2 = zTL1(ji,jj,jk) * zb / za + zTL2(ji,jj,jk) * zc / za + zkappa1 + zkappa2 - ztfe / za 
    157150                  za1 = zkappa2 * zTL1(ji,jj,jk) * zb / za + zkappa1 * zTL2(ji,jj,jk) * zc / za & 
     
    216209            END DO 
    217210         END DO 
     211         IF( ln_diatrc .AND. lk_iomput )   & 
     212            &  zFeL1(:,:,1:jpkm1) = MAX( 0., trn(:,:,1:jpkm1,jpfer) * 1E9 - zFe3(:,:,1:jpkm1) ) 
    218213         ! 
    219214      ENDIF 
     
    234229               ! Scavenging onto dust is also included as evidenced from the DUNE experiments. 
    235230               ! -------------------------------------------------------------------------------------- 
    236                zfeequi = zFe3(ji,jj,jk) * 1E-9  
    237                IF( ln_fechem ) zfeequi = zfeequi + zFep(ji,jj,jk) * 1E-9 
     231               IF( ln_fechem ) THEN 
     232                  zfeequi = ( zFe3(ji,jj,jk) + zFe2(ji,jj,jk) ) * 1E-9 
     233                  zfecoll = ( 0.3 * zFeL1(ji,jj,jk) + 0.5 * zFeL2(ji,jj,jk) ) * 1E-9 
     234               ELSE 
     235                  zfeequi = zFe3(ji,jj,jk) * 1E-9  
     236                  zfecoll = 0.5 * zFeL1(ji,jj,jk) * 1E-9 
     237               ENDIF 
    238238#if defined key_kriest 
    239239               ztrc   = ( trn(ji,jj,jk,jppoc) + trn(ji,jj,jk,jpcal) + trn(ji,jj,jk,jpgsi) ) * 1.e6  
     
    258258               !  due to increased lithogenic particles and let say it is unknown processes (precipitation, ...) 
    259259               !  ----------------------------------------------------------- 
    260                ztfe = trn(ji,jj,jk,jpfer) * 1.e9 - 1.  
    261                IF( ln_fechem )  ztfe = ztfe +  1. - ztotlig(ji,jj,jk)  
    262                zlam1b  = xlam1 * MAX( 0.e0, ztfe ) 
     260               zlam1b  = xlam1 * MAX( 0.e0, ( trn(ji,jj,jk,jpfer) * 1.e9 - ztotlig(ji,jj,jk) ) ) 
    263261               zcoag   = zfeequi * zlam1b * zstep 
    264262 
     
    269267               zlamfac = MAX( 0.e0, ( gphit(ji,jj) + 55.) / 30. ) 
    270268               zlamfac = MIN( 1.  , zlamfac ) 
    271                zdep    =  MIN(1., 1000. / fsdept(ji,jj,jk) ) 
    272 #if ! defined key_kriest 
    273                zlam1b = (  80.* ( trn(ji,jj,jk,jpdoc) + 35.e-6 )                           & 
    274                   &     + 698.*   trn(ji,jj,jk,jppoc) + 1.05e4 * trn(ji,jj,jk,jpgoc)  )    & 
    275                   &     * xdiss(ji,jj,jk) + 1E-4 * ( 1. - zlamfac ) * zdep 
     269               zdep    = MIN( 1., 1000. / fsdept(ji,jj,jk) ) 
     270               zlam1a  = ( 0.369  * 0.3 * trn(ji,jj,jk,jpdoc) + 102.4  * trn(ji,jj,jk,jppoc) ) * xdiss(ji,jj,jk)    & 
     271                   &   + ( 114.   * 0.3 * trn(ji,jj,jk,jpdoc) + 5.09E3 * trn(ji,jj,jk,jppoc) ) 
     272#if defined key_kriest 
     273               zlam1a   = zlam1a + 1E-4 * ( 1. - zlamfac ) * zdep 
     274               zaggdfea = zlam1a * zstep * zfecoll 
     275               zaggdfeb = 0. 
     276               ! 
     277               tra(ji,jj,jk,jpfer) = tra(ji,jj,jk,jpfer) - zscave - zaggdfea - zaggdfeb - zcoag 
     278               tra(ji,jj,jk,jpsfe) = tra(ji,jj,jk,jpsfe) + zscave * zdenom1 + zaggdfea + zaggdfeb 
    276279#else 
    277                zlam1b = (  80.* (trn(ji,jj,jk,jpdoc) + 35E-6)              & 
    278                   &     + 698.*  trn(ji,jj,jk,jppoc)  )                    & 
    279                   &     * xdiss(ji,jj,jk) + 1E-4 * ( 1. - zlamfac ) * zdep 
    280 #endif 
    281                zaggdfe = zlam1b * zstep * 0.5 * ( trn(ji,jj,jk,jpfer) - zfeequi ) 
    282                tra(ji,jj,jk,jpfer) = tra(ji,jj,jk,jpfer) - zscave - zaggdfe - zcoag 
    283 #if defined key_kriest 
    284                tra(ji,jj,jk,jpsfe) = tra(ji,jj,jk,jpsfe) + zscave * zdenom1 
    285 #else 
    286                tra(ji,jj,jk,jpsfe) = tra(ji,jj,jk,jpsfe) + zscave * zdenom1 
    287                tra(ji,jj,jk,jpbfe) = tra(ji,jj,jk,jpbfe) + zscave * zdenom2 
     280               zlam1b = 3.53E3 *   trn(ji,jj,jk,jpgoc) * xdiss(ji,jj,jk) + 1E-4 * ( 1. - zlamfac ) * zdep  
     281               zaggdfea = zlam1a * zstep * zfecoll 
     282               zaggdfeb = zlam1b * zstep * zfecoll 
     283               ! 
     284               tra(ji,jj,jk,jpfer) = tra(ji,jj,jk,jpfer) - zscave - zaggdfea - zaggdfeb - zcoag 
     285               tra(ji,jj,jk,jpsfe) = tra(ji,jj,jk,jpsfe) + zscave * zdenom1 + zaggdfea 
     286               tra(ji,jj,jk,jpbfe) = tra(ji,jj,jk,jpbfe) + zscave * zdenom2 + zaggdfeb 
    288287#endif 
    289288            END DO 
     
    291290      END DO 
    292291      ! 
     292      !  Define the bioavailable fraction of iron 
     293      !  ---------------------------------------- 
     294      IF( ln_fechem ) THEN 
     295          biron(:,:,:) = MAX( 0., trn(:,:,:,jpfer) - zFeP(:,:,:) * 1E-9 ) 
     296      ELSE 
     297          biron(:,:,:) = trn(:,:,:,jpfer)  
     298      ENDIF 
    293299 
    294300      IF(ln_ctl)   THEN  ! print mean trends (used for debugging) 
     
    302308      IF( ln_diatrc .AND. lk_iomput ) THEN 
    303309         IF( jnt == nrdttrc ) THEN 
    304             CALL iom_put("Fe3"    , zFe3   (:,:,:) * tmask(:,:,:) )   ! Fe3+ 
    305             CALL iom_put("TL1"    , zTL1   (:,:,:) * tmask(:,:,:) )   ! TL1 
    306             CALL iom_put("Totlig" , ztotlig(:,:,:) * tmask(:,:,:) )   ! TL 
     310            CALL iom_put("Fe3"    , zFe3   (:,:,:)       * tmask(:,:,:) )   ! Fe3+ 
     311            CALL iom_put("FeL1"   , zFeL1  (:,:,:)       * tmask(:,:,:) )   ! FeL1 
     312            CALL iom_put("TL1"    , zTL1   (:,:,:)       * tmask(:,:,:) )   ! TL1 
     313            CALL iom_put("Totlig" , ztotlig(:,:,:)       * tmask(:,:,:) )   ! TL 
     314            CALL iom_put("Biron"  , biron  (:,:,:) * 1e9 * tmask(:,:,:) )   ! TL 
    307315            IF( ln_fechem ) THEN 
    308                CALL iom_put("Fe2" , zFe2   (:,:,:) * tmask(:,:,:) )   ! Fe2+ 
    309                CALL iom_put("FeL1", zFeL1  (:,:,:) * tmask(:,:,:) )   ! FeL1 
    310                CALL iom_put("FeL2", zFeL2  (:,:,:) * tmask(:,:,:) )   ! FeL2 
    311                CALL iom_put("FeP" , zFeP   (:,:,:) * tmask(:,:,:) )   ! FeP 
    312                CALL iom_put("TL2" , zTL2   (:,:,:) * tmask(:,:,:) )   ! TL2 
     316               CALL iom_put("Fe2" , zFe2   (:,:,:)       * tmask(:,:,:) )   ! Fe2+ 
     317               CALL iom_put("FeL2", zFeL2  (:,:,:)       * tmask(:,:,:) )   ! FeL2 
     318               CALL iom_put("FeP" , zFeP   (:,:,:)       * tmask(:,:,:) )   ! FeP 
     319               CALL iom_put("TL2" , zTL2   (:,:,:)       * tmask(:,:,:) )   ! TL2 
    313320            ENDIF 
    314321         ENDIF 
     
    321328      ENDIF 
    322329      ! 
    323       CALL wrk_dealloc( jpi, jpj, jpk, zFe3, zTL1, ztotlig ) 
     330                                         CALL wrk_dealloc( jpi, jpj, jpk, zFe3, zTL1, ztotlig ) 
     331      IF( ln_diatrc .AND. lk_iomput )    CALL wrk_dealloc( jpi, jpj, jpk, zFeL1 ) 
    324332      IF( ln_fechem ) THEN 
    325          CALL wrk_dealloc( jpi, jpj, jpk, zTL2, zFeP ) 
    326          IF( ln_diatrc .AND. lk_iomput ) THEN 
    327             CALL wrk_dealloc( jpi, jpj, jpk, zFe2, zFeL1, zFeL2 ) 
    328          ENDIF 
     333                                         CALL wrk_dealloc( jpi, jpj, jpk, zTL2, zFeP ) 
     334         IF( ln_diatrc .AND. lk_iomput ) CALL wrk_dealloc( jpi, jpj, jpk, zFe2, zFeL2 ) 
    329335      ENDIF 
    330336      ! 
     
    365371         ! initialization of some constants used by the complexe chemistry scheme 
    366372         ! ---------------------------------------------------------------------- 
    367          spd = 3600.*24. 
     373         spd = 3600. * 24. 
    368374         con = 1.E9 
    369375         ! LIGAND KINETICS (values from Witter et al. 2000) 
     
    377383         ks  = 0.075 
    378384         kpr = 0.05 
     385         ! thermal reduction of Fe3 
     386         kth = 0.0048 * 24. 
    379387      ENDIF 
    380388      ! 
Note: See TracChangeset for help on using the changeset viewer.