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 10368 for NEMO/branches/2018/dev_r10164_HPC09_ESIWACE_PREP_MERGE/src – NEMO

Ignore:
Timestamp:
2018-12-03T12:45:01+01:00 (5 years ago)
Author:
smasson
Message:

dev_r10164_HPC09_ESIWACE_PREP_MERGE: merge with trunk@10365, see #2133

Location:
NEMO/branches/2018/dev_r10164_HPC09_ESIWACE_PREP_MERGE/src
Files:
43 edited

Legend:

Unmodified
Added
Removed
  • NEMO/branches/2018/dev_r10164_HPC09_ESIWACE_PREP_MERGE/src/OCE/ASM/asminc.F90

    r10170 r10368  
    911911            ! Nudge sea ice depth to bring it up to a required minimum depth 
    912912            WHERE( zseaicendg(:,:) > 0.0_wp .AND. hm_i(:,:) < zhicifmin )  
    913                zhicifinc(:,:) = (zhicifmin - hm_i(:,:)) * zincwgt     
     913               zhicifinc(:,:) = zhicifmin - hm_i(:,:) 
    914914            ELSEWHERE 
    915915               zhicifinc(:,:) = 0.0_wp 
  • NEMO/branches/2018/dev_r10164_HPC09_ESIWACE_PREP_MERGE/src/OCE/C1D/dyndmp.F90

    r10297 r10368  
    182182            DO jj = 2, jpjm1 
    183183               DO ji = fs_2, fs_jpim1   ! vector opt. 
    184                   IF( avt(ji,jj,jk) <= 5.e-4_wp ) THEN 
     184                  IF( avt(ji,jj,jk) <= avt_c ) THEN 
    185185                     zua = resto_uv(ji,jj,jk) * ( zuv_dta(ji,jj,jk,1) - ub(ji,jj,jk) ) 
    186186                     zva = resto_uv(ji,jj,jk) * ( zuv_dta(ji,jj,jk,2) - vb(ji,jj,jk) ) 
  • NEMO/branches/2018/dev_r10164_HPC09_ESIWACE_PREP_MERGE/src/OCE/DYN/dynnxt.F90

    r10170 r10368  
    246246                  DO jj = 1, jpj 
    247247                     DO ji = 1, jpi 
    248                         IF( misfkt(ji,jj) <= jk .AND. jk <= misfkb(ji,jj) ) THEN 
    249                            e3t_b(ji,jj,jk) =   e3t_b(ji,jj,jk) - zcoef * ( fwfisf_b(ji,jj) - fwfisf(ji,jj) ) & 
     248                        IF( misfkt(ji,jj) <=jk .and. jk < misfkb(ji,jj) ) THEN 
     249                           e3t_b(ji,jj,jk) = e3t_b(ji,jj,jk) - zcoef * ( fwfisf_b(ji,jj) - fwfisf(ji,jj) ) & 
    250250                                &          * ( e3t_n(ji,jj,jk) * r1_hisf_tbl(ji,jj) ) * tmask(ji,jj,jk) 
     251                        ELSEIF ( jk==misfkb(ji,jj) ) THEN 
     252                           e3t_b(ji,jj,jk) = e3t_b(ji,jj,jk) - zcoef * ( fwfisf_b(ji,jj) - fwfisf(ji,jj) ) & 
     253                                &          * ( e3t_n(ji,jj,jk) * r1_hisf_tbl(ji,jj) ) * ralpha(ji,jj) * tmask(ji,jj,jk) 
    251254                        ENDIF 
    252255                     END DO 
  • NEMO/branches/2018/dev_r10164_HPC09_ESIWACE_PREP_MERGE/src/OCE/DYN/dynzdf.F90

    r10068 r10368  
    7171      REAL(wp) ::   zzwi, ze3ua, zdt   ! local scalars 
    7272      REAL(wp) ::   zzws, ze3va        !   -      - 
     73      REAL(wp) ::   z1_e3ua, z1_e3va   !   -      - 
     74      REAL(wp) ::   zWu , zWv          !   -      - 
     75      REAL(wp) ::   zWui, zWvi         !   -      - 
     76      REAL(wp) ::   zWus, zWvs         !   -      - 
    7377      REAL(wp), DIMENSION(jpi,jpj,jpk)        ::  zwi, zwd, zws   ! 3D workspace  
    7478      REAL(wp), DIMENSION(:,:,:), ALLOCATABLE ::   ztrdu, ztrdv   !  -      - 
     
    155159      !                    !* Matrix construction 
    156160      zdt = r2dt * 0.5 
    157       SELECT CASE( nldf_dyn ) 
    158       CASE( np_lap_i )           ! rotated lateral mixing: add its vertical mixing (akzu) 
    159          DO jk = 1, jpkm1 
    160             DO jj = 2, jpjm1  
    161                DO ji = fs_2, fs_jpim1   ! vector opt. 
    162                   ze3ua =  ( 1._wp - r_vvl ) * e3u_n(ji,jj,jk) + r_vvl * e3u_a(ji,jj,jk)   ! after scale factor at T-point 
    163                   zzwi = - zdt * ( avm(ji+1,jj,jk  ) + avm(ji,jj,jk  ) + akzu(ji,jj,jk  ) )   & 
    164                      &         / ( ze3ua * e3uw_n(ji,jj,jk  ) ) * wumask(ji,jj,jk  ) 
    165                   zzws = - zdt * ( avm(ji+1,jj,jk+1) + avm(ji,jj,jk+1) + akzu(ji,jj,jk+1) )   & 
    166                      &         / ( ze3ua * e3uw_n(ji,jj,jk+1) ) * wumask(ji,jj,jk+1) 
    167                   zwi(ji,jj,jk) = zzwi 
    168                   zws(ji,jj,jk) = zzws 
    169                   zwd(ji,jj,jk) = 1._wp - zzwi - zzws 
    170                END DO 
    171             END DO 
    172          END DO 
    173       CASE DEFAULT               ! iso-level lateral mixing 
    174          DO jk = 1, jpkm1 
    175             DO jj = 2, jpjm1  
    176                DO ji = fs_2, fs_jpim1   ! vector opt. 
    177                   ze3ua =  ( 1._wp - r_vvl ) * e3u_n(ji,jj,jk) + r_vvl * e3u_a(ji,jj,jk)   ! after scale factor at T-point 
    178                   zzwi = - zdt * ( avm(ji+1,jj,jk  ) + avm(ji,jj,jk  ) ) / ( ze3ua * e3uw_n(ji,jj,jk  ) ) * wumask(ji,jj,jk  ) 
    179                   zzws = - zdt * ( avm(ji+1,jj,jk+1) + avm(ji,jj,jk+1) ) / ( ze3ua * e3uw_n(ji,jj,jk+1) ) * wumask(ji,jj,jk+1) 
    180                   zwi(ji,jj,jk) = zzwi 
    181                   zws(ji,jj,jk) = zzws 
    182                   zwd(ji,jj,jk) = 1._wp - zzwi - zzws 
    183                END DO 
    184             END DO 
    185          END DO 
    186       END SELECT 
    187       ! 
    188       DO jj = 2, jpjm1     !* Surface boundary conditions 
    189          DO ji = fs_2, fs_jpim1   ! vector opt. 
    190             zwi(ji,jj,1) = 0._wp 
    191             zwd(ji,jj,1) = 1._wp - zws(ji,jj,1) 
    192          END DO 
    193       END DO 
     161      IF( ln_zad_Aimp ) THEN   !! 
     162         SELECT CASE( nldf_dyn ) 
     163         CASE( np_lap_i )           ! rotated lateral mixing: add its vertical mixing (akzu) 
     164            DO jk = 1, jpkm1 
     165               DO jj = 2, jpjm1  
     166                  DO ji = fs_2, fs_jpim1   ! vector opt. 
     167                     ze3ua =  ( 1._wp - r_vvl ) * e3u_n(ji,jj,jk) + r_vvl * e3u_a(ji,jj,jk)   ! after scale factor at U-point 
     168                     zzwi = - zdt * ( avm(ji+1,jj,jk  ) + avm(ji,jj,jk  ) + akzu(ji,jj,jk  ) )   & 
     169                        &         / ( ze3ua * e3uw_n(ji,jj,jk  ) ) * wumask(ji,jj,jk  ) 
     170                     zzws = - zdt * ( avm(ji+1,jj,jk+1) + avm(ji,jj,jk+1) + akzu(ji,jj,jk+1) )   & 
     171                        &         / ( ze3ua * e3uw_n(ji,jj,jk+1) ) * wumask(ji,jj,jk+1) 
     172                     zWui = 0.5_wp * ( wi(ji,jj,jk  ) + wi(ji+1,jj,jk  ) ) 
     173                     zWus = 0.5_wp * ( wi(ji,jj,jk+1) + wi(ji+1,jj,jk+1) ) 
     174                     zwi(ji,jj,jk) = zzwi + zdt * MIN( zWui, 0._wp )  
     175                     zws(ji,jj,jk) = zzws - zdt * MAX( zWus, 0._wp ) 
     176                     zwd(ji,jj,jk) = 1._wp - zzwi - zzws + zdt * ( MAX( zWui, 0._wp ) - MIN( zWus, 0._wp ) ) 
     177                  END DO 
     178               END DO 
     179            END DO 
     180         CASE DEFAULT               ! iso-level lateral mixing 
     181            DO jk = 1, jpkm1 
     182               DO jj = 2, jpjm1  
     183                  DO ji = fs_2, fs_jpim1   ! vector opt. 
     184                     ze3ua =  ( 1._wp - r_vvl ) * e3u_n(ji,jj,jk) + r_vvl * e3u_a(ji,jj,jk)   ! after scale factor at U-point 
     185                     zzwi = - zdt * ( avm(ji+1,jj,jk  ) + avm(ji,jj,jk  ) ) / ( ze3ua * e3uw_n(ji,jj,jk  ) ) * wumask(ji,jj,jk  ) 
     186                     zzws = - zdt * ( avm(ji+1,jj,jk+1) + avm(ji,jj,jk+1) ) / ( ze3ua * e3uw_n(ji,jj,jk+1) ) * wumask(ji,jj,jk+1) 
     187                     zWui = 0.5_wp * ( wi(ji,jj,jk  ) + wi(ji+1,jj,jk  ) ) 
     188                     zWus = 0.5_wp * ( wi(ji,jj,jk+1) + wi(ji+1,jj,jk+1) ) 
     189                     zwi(ji,jj,jk) = zzwi + zdt * MIN( zWui, 0._wp ) 
     190                     zws(ji,jj,jk) = zzws - zdt * MAX( zWus, 0._wp ) 
     191                     zwd(ji,jj,jk) = 1._wp - zzwi - zzws + zdt * ( MAX( zWui, 0._wp ) - MIN( zWus, 0._wp ) ) 
     192                  END DO 
     193               END DO 
     194            END DO 
     195         END SELECT 
     196         DO jj = 2, jpjm1     !* Surface boundary conditions 
     197            DO ji = fs_2, fs_jpim1   ! vector opt. 
     198               zwi(ji,jj,1) = 0._wp 
     199               ze3ua =  ( 1._wp - r_vvl ) * e3u_n(ji,jj,1) + r_vvl * e3u_a(ji,jj,1) 
     200               zzws = - zdt * ( avm(ji+1,jj,2) + avm(ji  ,jj,2) ) / ( ze3ua * e3uw_n(ji,jj,2) ) * wumask(ji,jj,2) 
     201               zWus = 0.5_wp * ( wi(ji  ,jj,2) +  wi(ji+1,jj,2) ) 
     202               zws(ji,jj,1 ) = zzws - zdt * MAX( zWus, 0._wp ) 
     203               zwd(ji,jj,1 ) = 1._wp - zzws - zdt * ( MIN( zWus, 0._wp ) ) 
     204            END DO 
     205         END DO 
     206      ELSE 
     207         SELECT CASE( nldf_dyn ) 
     208         CASE( np_lap_i )           ! rotated lateral mixing: add its vertical mixing (akzu) 
     209            DO jk = 1, jpkm1 
     210               DO jj = 2, jpjm1  
     211                  DO ji = fs_2, fs_jpim1   ! vector opt. 
     212                     ze3ua =  ( 1._wp - r_vvl ) * e3u_n(ji,jj,jk) + r_vvl * e3u_a(ji,jj,jk)   ! after scale factor at U-point 
     213                     zzwi = - zdt * ( avm(ji+1,jj,jk  ) + avm(ji,jj,jk  ) + akzu(ji,jj,jk  ) )   & 
     214                        &         / ( ze3ua * e3uw_n(ji,jj,jk  ) ) * wumask(ji,jj,jk  ) 
     215                     zzws = - zdt * ( avm(ji+1,jj,jk+1) + avm(ji,jj,jk+1) + akzu(ji,jj,jk+1) )   & 
     216                        &         / ( ze3ua * e3uw_n(ji,jj,jk+1) ) * wumask(ji,jj,jk+1) 
     217                     zwi(ji,jj,jk) = zzwi 
     218                     zws(ji,jj,jk) = zzws 
     219                     zwd(ji,jj,jk) = 1._wp - zzwi - zzws 
     220                  END DO 
     221               END DO 
     222            END DO 
     223         CASE DEFAULT               ! iso-level lateral mixing 
     224            DO jk = 1, jpkm1 
     225               DO jj = 2, jpjm1  
     226                  DO ji = fs_2, fs_jpim1   ! vector opt. 
     227                     ze3ua =  ( 1._wp - r_vvl ) * e3u_n(ji,jj,jk) + r_vvl * e3u_a(ji,jj,jk)   ! after scale factor at U-point 
     228                     zzwi = - zdt * ( avm(ji+1,jj,jk  ) + avm(ji,jj,jk  ) ) / ( ze3ua * e3uw_n(ji,jj,jk  ) ) * wumask(ji,jj,jk  ) 
     229                     zzws = - zdt * ( avm(ji+1,jj,jk+1) + avm(ji,jj,jk+1) ) / ( ze3ua * e3uw_n(ji,jj,jk+1) ) * wumask(ji,jj,jk+1) 
     230                     zwi(ji,jj,jk) = zzwi 
     231                     zws(ji,jj,jk) = zzws 
     232                     zwd(ji,jj,jk) = 1._wp - zzwi - zzws 
     233                  END DO 
     234               END DO 
     235            END DO 
     236         END SELECT 
     237         DO jj = 2, jpjm1     !* Surface boundary conditions 
     238            DO ji = fs_2, fs_jpim1   ! vector opt. 
     239               zwi(ji,jj,1) = 0._wp 
     240               zwd(ji,jj,1) = 1._wp - zws(ji,jj,1) 
     241            END DO 
     242         END DO 
     243      ENDIF 
     244      ! 
    194245      ! 
    195246      !              !==  Apply semi-implicit bottom friction  ==! 
     
    274325      !                       !* Matrix construction 
    275326      zdt = r2dt * 0.5 
    276       SELECT CASE( nldf_dyn ) 
    277       CASE( np_lap_i )           ! rotated lateral mixing: add its vertical mixing (akzu) 
    278          DO jk = 1, jpkm1 
    279             DO jj = 2, jpjm1    
    280                DO ji = fs_2, fs_jpim1   ! vector opt. 
    281                   ze3va =  ( 1._wp - r_vvl ) * e3v_n(ji,jj,jk) + r_vvl * e3v_a(ji,jj,jk)   ! after scale factor at T-point 
    282                   zzwi = - zdt * ( avm(ji,jj+1,jk  ) + avm(ji,jj,jk  ) + akzv(ji,jj,jk  ) )   & 
    283                      &         / ( ze3va * e3vw_n(ji,jj,jk  ) ) * wvmask(ji,jj,jk  ) 
    284                   zzws = - zdt * ( avm(ji,jj+1,jk+1) + avm(ji,jj,jk+1) + akzv(ji,jj,jk+1) )   & 
    285                      &         / ( ze3va * e3vw_n(ji,jj,jk+1) ) * wvmask(ji,jj,jk+1) 
    286                   zwi(ji,jj,jk) = zzwi * wvmask(ji,jj,jk  ) 
    287                   zws(ji,jj,jk) = zzws * wvmask(ji,jj,jk+1) 
    288                   zwd(ji,jj,jk) = 1._wp - zzwi - zzws 
    289                END DO 
    290             END DO 
    291          END DO 
    292       CASE DEFAULT               ! iso-level lateral mixing 
    293          DO jk = 1, jpkm1 
    294             DO jj = 2, jpjm1    
    295                DO ji = fs_2, fs_jpim1   ! vector opt. 
    296                   ze3va =  ( 1._wp - r_vvl ) * e3v_n(ji,jj,jk) + r_vvl * e3v_a(ji,jj,jk)   ! after scale factor at T-point 
    297                   zzwi = - zdt * ( avm(ji,jj+1,jk  ) + avm(ji,jj,jk  ) ) / ( ze3va * e3vw_n(ji,jj,jk  ) ) * wvmask(ji,jj,jk  ) 
    298                   zzws = - zdt * ( avm(ji,jj+1,jk+1) + avm(ji,jj,jk+1) ) / ( ze3va * e3vw_n(ji,jj,jk+1) ) * wvmask(ji,jj,jk+1) 
    299                   zwi(ji,jj,jk) = zzwi * wvmask(ji,jj,jk  ) 
    300                   zws(ji,jj,jk) = zzws * wvmask(ji,jj,jk+1) 
    301                   zwd(ji,jj,jk) = 1._wp - zzwi - zzws 
    302                END DO 
    303             END DO 
    304          END DO 
    305       END SELECT 
    306       ! 
    307       DO jj = 2, jpjm1        !* Surface boundary conditions 
    308          DO ji = fs_2, fs_jpim1   ! vector opt. 
    309             zwi(ji,jj,1) = 0._wp 
    310             zwd(ji,jj,1) = 1._wp - zws(ji,jj,1) 
    311          END DO 
    312       END DO 
     327      IF( ln_zad_Aimp ) THEN   !! 
     328         SELECT CASE( nldf_dyn ) 
     329         CASE( np_lap_i )           ! rotated lateral mixing: add its vertical mixing (akzv) 
     330            DO jk = 1, jpkm1 
     331               DO jj = 2, jpjm1  
     332                  DO ji = fs_2, fs_jpim1   ! vector opt. 
     333                     ze3va =  ( 1._wp - r_vvl ) * e3v_n(ji,jj,jk) + r_vvl * e3v_a(ji,jj,jk)   ! after scale factor at V-point 
     334                     zzwi = - zdt * ( avm(ji,jj+1,jk  ) + avm(ji,jj,jk  ) + akzv(ji,jj,jk  ) )   & 
     335                        &         / ( ze3va * e3vw_n(ji,jj,jk  ) ) * wvmask(ji,jj,jk  ) 
     336                     zzws = - zdt * ( avm(ji,jj+1,jk+1) + avm(ji,jj,jk+1) + akzv(ji,jj,jk+1) )   & 
     337                        &         / ( ze3va * e3vw_n(ji,jj,jk+1) ) * wvmask(ji,jj,jk+1) 
     338                     zWvi = 0.5_wp * ( wi(ji,jj,jk  ) + wi(ji,jj+1,jk  ) ) * wvmask(ji,jj,jk  ) 
     339                     zWvs = 0.5_wp * ( wi(ji,jj,jk+1) + wi(ji,jj+1,jk+1) ) * wvmask(ji,jj,jk+1) 
     340                     zwi(ji,jj,jk) = zzwi + zdt * MIN( zWvi, 0._wp ) 
     341                     zws(ji,jj,jk) = zzws - zdt * MAX( zWvs, 0._wp ) 
     342                     zwd(ji,jj,jk) = 1._wp - zzwi - zzws - zdt * ( - MAX( zWvi, 0._wp ) + MIN( zWvs, 0._wp ) ) 
     343                  END DO 
     344               END DO 
     345            END DO 
     346         CASE DEFAULT               ! iso-level lateral mixing 
     347            DO jk = 1, jpkm1 
     348               DO jj = 2, jpjm1  
     349                  DO ji = fs_2, fs_jpim1   ! vector opt. 
     350                     ze3va =  ( 1._wp - r_vvl ) * e3v_n(ji,jj,jk) + r_vvl * e3v_a(ji,jj,jk)   ! after scale factor at V-point 
     351                     zzwi = - zdt * ( avm(ji,jj+1,jk  ) + avm(ji,jj,jk  ) ) / ( ze3va * e3vw_n(ji,jj,jk  ) ) * wvmask(ji,jj,jk  ) 
     352                     zzws = - zdt * ( avm(ji,jj+1,jk+1) + avm(ji,jj,jk+1) ) / ( ze3va * e3vw_n(ji,jj,jk+1) ) * wvmask(ji,jj,jk+1) 
     353                     zWvi = 0.5_wp * ( wi(ji,jj,jk  ) + wi(ji,jj+1,jk  ) ) * wvmask(ji,jj,jk  ) 
     354                     zWvs = 0.5_wp * ( wi(ji,jj,jk+1) + wi(ji,jj+1,jk+1) ) * wvmask(ji,jj,jk+1) 
     355                     zwi(ji,jj,jk) = zzwi  + zdt * MIN( zWvi, 0._wp ) 
     356                     zws(ji,jj,jk) = zzws  - zdt * MAX( zWvs, 0._wp ) 
     357                     zwd(ji,jj,jk) = 1._wp - zzwi - zzws - zdt * ( - MAX( zWvi, 0._wp ) + MIN( zWvs, 0._wp ) ) 
     358                  END DO 
     359               END DO 
     360            END DO 
     361         END SELECT 
     362         DO jj = 2, jpjm1     !* Surface boundary conditions 
     363            DO ji = fs_2, fs_jpim1   ! vector opt. 
     364               zwi(ji,jj,1) = 0._wp 
     365               ze3va =  ( 1._wp - r_vvl ) * e3v_n(ji,jj,1) + r_vvl * e3v_a(ji,jj,1) 
     366               zzws = - zdt * ( avm(ji,jj+1,2) + avm(ji,jj,2) ) / ( ze3va * e3vw_n(ji,jj,2) ) * wvmask(ji,jj,2) 
     367               zWvs = 0.5_wp * ( wi(ji,jj  ,2) +  wi(ji,jj+1,2) ) 
     368               zws(ji,jj,1 ) = zzws - zdt * MAX( zWvs, 0._wp ) 
     369               zwd(ji,jj,1 ) = 1._wp - zzws - zdt * ( MIN( zWvs, 0._wp ) ) 
     370            END DO 
     371         END DO 
     372      ELSE 
     373         SELECT CASE( nldf_dyn ) 
     374         CASE( np_lap_i )           ! rotated lateral mixing: add its vertical mixing (akzu) 
     375            DO jk = 1, jpkm1 
     376               DO jj = 2, jpjm1    
     377                  DO ji = fs_2, fs_jpim1   ! vector opt. 
     378                     ze3va =  ( 1._wp - r_vvl ) * e3v_n(ji,jj,jk) + r_vvl * e3v_a(ji,jj,jk)   ! after scale factor at V-point 
     379                     zzwi = - zdt * ( avm(ji,jj+1,jk  ) + avm(ji,jj,jk  ) + akzv(ji,jj,jk  ) )   & 
     380                        &         / ( ze3va * e3vw_n(ji,jj,jk  ) ) * wvmask(ji,jj,jk  ) 
     381                     zzws = - zdt * ( avm(ji,jj+1,jk+1) + avm(ji,jj,jk+1) + akzv(ji,jj,jk+1) )   & 
     382                        &         / ( ze3va * e3vw_n(ji,jj,jk+1) ) * wvmask(ji,jj,jk+1) 
     383                     zwi(ji,jj,jk) = zzwi 
     384                     zws(ji,jj,jk) = zzws 
     385                     zwd(ji,jj,jk) = 1._wp - zzwi - zzws 
     386                  END DO 
     387               END DO 
     388            END DO 
     389         CASE DEFAULT               ! iso-level lateral mixing 
     390            DO jk = 1, jpkm1 
     391               DO jj = 2, jpjm1    
     392                  DO ji = fs_2, fs_jpim1   ! vector opt. 
     393                     ze3va =  ( 1._wp - r_vvl ) * e3v_n(ji,jj,jk) + r_vvl * e3v_a(ji,jj,jk)   ! after scale factor at V-point 
     394                     zzwi = - zdt * ( avm(ji,jj+1,jk  ) + avm(ji,jj,jk  ) ) / ( ze3va * e3vw_n(ji,jj,jk  ) ) * wvmask(ji,jj,jk  ) 
     395                     zzws = - zdt * ( avm(ji,jj+1,jk+1) + avm(ji,jj,jk+1) ) / ( ze3va * e3vw_n(ji,jj,jk+1) ) * wvmask(ji,jj,jk+1) 
     396                     zwi(ji,jj,jk) = zzwi 
     397                     zws(ji,jj,jk) = zzws 
     398                     zwd(ji,jj,jk) = 1._wp - zzwi - zzws 
     399                  END DO 
     400               END DO 
     401            END DO 
     402         END SELECT 
     403         DO jj = 2, jpjm1        !* Surface boundary conditions 
     404            DO ji = fs_2, fs_jpim1   ! vector opt. 
     405               zwi(ji,jj,1) = 0._wp 
     406               zwd(ji,jj,1) = 1._wp - zws(ji,jj,1) 
     407            END DO 
     408         END DO 
     409      ENDIF 
     410      ! 
    313411      !              !==  Apply semi-implicit top/bottom friction  ==! 
    314412      ! 
  • NEMO/branches/2018/dev_r10164_HPC09_ESIWACE_PREP_MERGE/src/OCE/DYN/sshwzv.F90

    r10170 r10368  
    2828#endif 
    2929   ! 
     30   USE iom  
    3031   USE in_out_manager ! I/O manager 
    3132   USE restart        ! only for lrst_oce 
     
    4142   PUBLIC   ssh_nxt    ! called by step.F90 
    4243   PUBLIC   wzv        ! called by step.F90 
     44   PUBLIC   wAimp      ! called by step.F90 
    4345   PUBLIC   ssh_swp    ! called by step.F90 
    4446 
     
    265267   END SUBROUTINE ssh_swp 
    266268 
     269   SUBROUTINE wAimp( kt ) 
     270      !!---------------------------------------------------------------------- 
     271      !!                ***  ROUTINE wAimp  *** 
     272      !!                    
     273      !! ** Purpose :   compute the Courant number and partition vertical velocity 
     274      !!                if a proportion needs to be treated implicitly 
     275      !! 
     276      !! ** Method  : -  
     277      !! 
     278      !! ** action  :   wn      : now vertical velocity (to be handled explicitly) 
     279      !!            :   wi      : now vertical velocity (for implicit treatment) 
     280      !! 
     281      !! Reference  : Leclair, M., and G. Madec, 2009, Ocean Modelling. 
     282      !!---------------------------------------------------------------------- 
     283      INTEGER, INTENT(in) ::   kt   ! time step 
     284      ! 
     285      INTEGER  ::   ji, jj, jk   ! dummy loop indices 
     286      REAL(wp)             ::   zCu, zcff, z1_e3w                     ! local scalars 
     287      REAL(wp) , PARAMETER ::   Cu_min = 0.15_wp                      ! local parameters 
     288      REAL(wp) , PARAMETER ::   Cu_max = 0.27                         ! local parameters 
     289      REAL(wp) , PARAMETER ::   Cu_cut = 2._wp*Cu_max - Cu_min        ! local parameters 
     290      REAL(wp) , PARAMETER ::   Fcu    = 4._wp*Cu_max*(Cu_max-Cu_min) ! local parameters 
     291      !!---------------------------------------------------------------------- 
     292      ! 
     293      IF( ln_timing )   CALL timing_start('wAimp') 
     294      ! 
     295      IF( kt == nit000 ) THEN 
     296         IF(lwp) WRITE(numout,*) 
     297         IF(lwp) WRITE(numout,*) 'wAimp : Courant number-based partitioning of now vertical velocity ' 
     298         IF(lwp) WRITE(numout,*) '~~~~~ ' 
     299         ! 
     300         Cu_adv(:,:,jpk) = 0._wp              ! bottom value : Cu_adv=0 (set once for all) 
     301      ENDIF 
     302      ! 
     303      DO jk = 1, jpkm1            ! calculate Courant numbers 
     304         DO jj = 2, jpjm1 
     305            DO ji = 2, fs_jpim1   ! vector opt. 
     306               z1_e3w = 1._wp / e3w_n(ji,jj,jk) 
     307               Cu_adv(ji,jj,jk) = r2dt * ( ( MAX( wn(ji,jj,jk) , 0._wp ) - MIN( wn(ji,jj,jk+1) , 0._wp ) )    & 
     308                  &                      + ( MAX( e2u(ji  ,jj)*e3uw_n(ji  ,jj,jk)*un(ji  ,jj,jk), 0._wp ) -   & 
     309                  &                          MIN( e2u(ji-1,jj)*e3uw_n(ji-1,jj,jk)*un(ji-1,jj,jk), 0._wp ) )   & 
     310                  &                        * r1_e1e2t(ji,jj)                                                  & 
     311                  &                      + ( MAX( e1v(ji,jj  )*e3vw_n(ji,jj  ,jk)*vn(ji,jj  ,jk), 0._wp ) -   & 
     312                  &                          MIN( e1v(ji,jj-1)*e3vw_n(ji,jj-1,jk)*vn(ji,jj-1,jk), 0._wp ) )   & 
     313                  &                        * r1_e1e2t(ji,jj)                                                  & 
     314                  &                      ) * z1_e3w 
     315            END DO 
     316         END DO 
     317      END DO 
     318      ! 
     319      CALL iom_put("Courant",Cu_adv) 
     320      ! 
     321      wi(:,:,:) = 0._wp                                 ! Includes top and bottom values set to zero 
     322      IF( MAXVAL( Cu_adv(:,:,:) ) > Cu_min ) THEN       ! Quick check if any breaches anywhere 
     323         DO jk = 1, jpkm1                               ! or scan Courant criterion and partition 
     324            DO jj = 2, jpjm1                            ! w where necessary 
     325               DO ji = 2, fs_jpim1   ! vector opt. 
     326                  ! 
     327                  zCu = MAX( Cu_adv(ji,jj,jk) , Cu_adv(ji,jj,jk+1) ) 
     328                  ! 
     329                  IF( zCu < Cu_min ) THEN               !<-- Fully explicit 
     330                     zcff = 0._wp 
     331                  ELSEIF( zCu < Cu_cut ) THEN           !<-- Mixed explicit 
     332                     zcff = ( zCu - Cu_min )**2 
     333                     zcff = zcff / ( Fcu + zcff ) 
     334                  ELSE                                  !<-- Mostly implicit 
     335                     zcff = ( zCu - Cu_max )/ zCu 
     336                  ENDIF 
     337                  zcff = MIN(1._wp, zcff) 
     338                  ! 
     339                  wi(ji,jj,jk) =           zcff   * wn(ji,jj,jk) 
     340                  wn(ji,jj,jk) = ( 1._wp - zcff ) * wn(ji,jj,jk) 
     341                  ! 
     342                  Cu_adv(ji,jj,jk) = zcff               ! Reuse array to output coefficient 
     343               END DO 
     344            END DO 
     345         END DO 
     346      ELSE 
     347         ! Fully explicit everywhere 
     348         Cu_adv = 0.0_wp                                ! Reuse array to output coefficient 
     349      ENDIF 
     350      CALL iom_put("wimp",wi)  
     351      CALL iom_put("wi_cff",Cu_adv) 
     352      CALL iom_put("wexp",wn) 
     353      ! 
     354      IF( ln_timing )   CALL timing_stop('wAimp') 
     355      ! 
     356   END SUBROUTINE wAimp 
    267357   !!====================================================================== 
    268358END MODULE sshwzv 
  • NEMO/branches/2018/dev_r10164_HPC09_ESIWACE_PREP_MERGE/src/OCE/IOM/iom.F90

    r10345 r10368  
    209209          zw_bnds(2,1:jpkm1  ) = gdepw_1d(jkmin:jpk) 
    210210          zw_bnds(2,jpk:     ) = gdepw_1d(jpk) + e3t_1d(jpk) 
    211           CALL iom_set_axis_attr( "deptht", bounds=zt_bnds ) 
    212           CALL iom_set_axis_attr( "depthu", bounds=zt_bnds ) 
    213           CALL iom_set_axis_attr( "depthv", bounds=zt_bnds ) 
    214           CALL iom_set_axis_attr( "depthw", bounds=zw_bnds ) 
     211          CALL iom_set_axis_attr( "deptht", bounds=zw_bnds ) 
     212          CALL iom_set_axis_attr( "depthu", bounds=zw_bnds ) 
     213          CALL iom_set_axis_attr( "depthv", bounds=zw_bnds ) 
     214          CALL iom_set_axis_attr( "depthw", bounds=zt_bnds ) 
    215215          ! 
    216216# if defined key_floats 
  • NEMO/branches/2018/dev_r10164_HPC09_ESIWACE_PREP_MERGE/src/OCE/LBC/lib_mpp.F90

    r10359 r10368  
    6464 
    6565   INTERFACE mpp_nfd 
    66       MODULE PROCEDURE   mpp_nfd_2d      , mpp_nfd_3d      , mpp_nfd_4d 
     66      MODULE PROCEDURE   mpp_nfd_2d    , mpp_nfd_3d    , mpp_nfd_4d 
    6767      MODULE PROCEDURE   mpp_nfd_2d_ptr, mpp_nfd_3d_ptr, mpp_nfd_4d_ptr 
    6868   END INTERFACE 
    6969 
    7070   ! Interface associated to the mpp_lnk_... routines is defined in lbclnk 
    71    PUBLIC   mpp_lnk_2d      , mpp_lnk_3d      , mpp_lnk_4d 
     71   PUBLIC   mpp_lnk_2d    , mpp_lnk_3d    , mpp_lnk_4d 
    7272   PUBLIC   mpp_lnk_2d_ptr, mpp_lnk_3d_ptr, mpp_lnk_4d_ptr 
    7373   ! 
     
    165165   LOGICAL, PUBLIC                               ::   l_full_nf_update = .FALSE.   !: logical for a full (2lines) update of bc at North fold report 
    166166   INTEGER,                    PARAMETER, PUBLIC ::   nbdelay = 1       !: number of delayed operations 
     167!$AGRIF_DO_NOT_TREAT 
    167168   CHARACTER(len=32), DIMENSION(nbdelay), PUBLIC ::   c_delaylist       !: name (used as id) of allreduce-delayed operations 
    168169   DATA c_delaylist(1) / 'advumx_delay' / 
     170!$AGRIF_END_DO_NOT_TREAT 
    169171   REAL(wp),          DIMENSION(nbdelay), PUBLIC ::   todelay           !: current value of the delayed operations 
    170172   INTEGER,           DIMENSION(nbdelay), PUBLIC ::   ndelayid = -1     !: mpi request id of the delayed operations 
     
    18091811      ! 
    18101812      iost=0 
    1811       IF( cdacce(1:6) == 'DIRECT' )  THEN 
     1813      IF( cdacce(1:6) == 'DIRECT' )  THEN         ! cdacce has always more than 6 characters 
    18121814         OPEN( UNIT=knum, FILE=clfile, FORM=cdform, ACCESS=cdacce, STATUS=cdstat, RECL=klengh         , ERR=100, IOSTAT=iost ) 
    1813       ELSE IF( cdstat(1:6) == 'APPEND' )  THEN 
     1815      ELSE IF( TRIM(cdstat) == 'APPEND' )  THEN   ! cdstat can have less than 6 characters 
    18141816         OPEN( UNIT=knum, FILE=clfile, FORM=cdform, ACCESS=cdacce, STATUS='UNKNOWN', POSITION='APPEND', ERR=100, IOSTAT=iost ) 
    18151817      ELSE 
  • NEMO/branches/2018/dev_r10164_HPC09_ESIWACE_PREP_MERGE/src/OCE/OBS/obsinter_h2d.h90

    r10068 r10368  
    11241124            WRITE(numout,*)'Grid lons    : ',plammm, plammp, plampm, plampp 
    11251125            WRITE(numout,*)'Current i,j  : ',ziguess, zjguess 
    1126             WRITE(numout,*)'Weights      : ',pbiwmm, pbiwmp, pbiwpm, pbiwpp 
    11271126            WRITE(numout,*)'jiter        = ',jiter 
    11281127            WRITE(numout,*)'zeps         = ',zeps 
  • NEMO/branches/2018/dev_r10164_HPC09_ESIWACE_PREP_MERGE/src/OCE/TRA/tradmp.F90

    r10297 r10368  
    123123            DO jj = 2, jpjm1 
    124124               DO ji = fs_2, fs_jpim1   ! vector opt. 
    125                   IF( avt(ji,jj,jk) <= 5.e-4_wp ) THEN 
     125                  IF( avt(ji,jj,jk) <= avt_c ) THEN 
    126126                     tsa(ji,jj,jk,jp_tem) = tsa(ji,jj,jk,jp_tem)   & 
    127127                        &                 + resto(ji,jj,jk) * ( zts_dta(ji,jj,jk,jp_tem) - tsb(ji,jj,jk,jp_tem) ) 
  • NEMO/branches/2018/dev_r10164_HPC09_ESIWACE_PREP_MERGE/src/OCE/TRA/trazdf.F90

    r10170 r10368  
    136136      ! 
    137137      INTEGER  ::  ji, jj, jk, jn   ! dummy loop indices 
    138       REAL(wp) ::  zrhs            ! local scalars 
     138      REAL(wp) ::  zrhs, zzwi, zzws ! local scalars 
    139139      REAL(wp), DIMENSION(jpi,jpj,jpk) ::  zwi, zwt, zwd, zws 
    140140      !!--------------------------------------------------------------------- 
     
    177177            ! 
    178178            ! Diagonal, lower (i), upper (s)  (including the bottom boundary condition since avt is masked) 
    179             DO jk = 1, jpkm1 
    180                DO jj = 2, jpjm1 
    181                   DO ji = fs_2, fs_jpim1   ! vector opt. 
    182 !!gm BUG  I think, use e3w_a instead of e3w_n, not sure of that 
    183                      zwi(ji,jj,jk) = - p2dt * zwt(ji,jj,jk  ) / e3w_n(ji,jj,jk  ) 
    184                      zws(ji,jj,jk) = - p2dt * zwt(ji,jj,jk+1) / e3w_n(ji,jj,jk+1) 
    185                      zwd(ji,jj,jk) = e3t_a(ji,jj,jk) - zwi(ji,jj,jk) - zws(ji,jj,jk) 
    186                  END DO 
    187                END DO 
    188             END DO 
     179            IF( ln_zad_Aimp ) THEN         ! Adaptive implicit vertical advection 
     180               DO jk = 1, jpkm1 
     181                  DO jj = 2, jpjm1 
     182                     DO ji = fs_2, fs_jpim1   ! vector opt. (ensure same order of calculation as below if wi=0.) 
     183                        zzwi = - p2dt * zwt(ji,jj,jk  ) / e3w_n(ji,jj,jk  ) 
     184                        zzws = - p2dt * zwt(ji,jj,jk+1) / e3w_n(ji,jj,jk+1) 
     185                        zwd(ji,jj,jk) = e3t_a(ji,jj,jk) - zzwi - zzws   & 
     186                           &                 + p2dt * ( MAX( wi(ji,jj,jk  ) , 0._wp ) - MIN( wi(ji,jj,jk+1) , 0._wp ) ) 
     187                        zwi(ji,jj,jk) = zzwi + p2dt *   MIN( wi(ji,jj,jk  ) , 0._wp ) 
     188                        zws(ji,jj,jk) = zzws - p2dt *   MAX( wi(ji,jj,jk+1) , 0._wp ) 
     189                    END DO 
     190                  END DO 
     191               END DO 
     192            ELSE 
     193               DO jk = 1, jpkm1 
     194                  DO jj = 2, jpjm1 
     195                     DO ji = fs_2, fs_jpim1   ! vector opt. 
     196                        zwi(ji,jj,jk) = - p2dt * zwt(ji,jj,jk  ) / e3w_n(ji,jj,jk) 
     197                        zws(ji,jj,jk) = - p2dt * zwt(ji,jj,jk+1) / e3w_n(ji,jj,jk+1) 
     198                        zwd(ji,jj,jk) = e3t_a(ji,jj,jk) - zwi(ji,jj,jk) - zws(ji,jj,jk) 
     199                    END DO 
     200                  END DO 
     201               END DO 
     202            ENDIF 
    189203            ! 
    190204            !! Matrix inversion from the first level 
  • NEMO/branches/2018/dev_r10164_HPC09_ESIWACE_PREP_MERGE/src/OCE/ZDF/zdf_oce.F90

    r10068 r10368  
    1818 
    1919   !                            !!* namelist namzdf: vertical physics * 
     20   !                             ! Adaptive-implicit vertical advection flag 
     21   LOGICAL , PUBLIC ::   ln_zad_Aimp !: adaptive (Courant number-based) implicit vertical advection 
    2022   !                             ! vertical closure scheme flags 
    2123   LOGICAL , PUBLIC ::   ln_zdfcst   !: constant coefficients 
  • NEMO/branches/2018/dev_r10164_HPC09_ESIWACE_PREP_MERGE/src/OCE/ZDF/zdfmxl.F90

    r10297 r10368  
    3333 
    3434   REAL(wp), PUBLIC ::   rho_c = 0.01_wp    !: density criterion for mixed layer depth 
    35    REAL(wp)        ::   avt_c = 5.e-4_wp   ! Kz criterion for the turbocline depth 
     35   REAL(wp), PUBLIC ::   avt_c = 5.e-4_wp   ! Kz criterion for the turbocline depth 
    3636 
    3737   !!---------------------------------------------------------------------- 
  • NEMO/branches/2018/dev_r10164_HPC09_ESIWACE_PREP_MERGE/src/OCE/ZDF/zdfosm.F90

    r10297 r10368  
    3232 
    3333   !!---------------------------------------------------------------------- 
    34    !!   'key_zdfosm'                                             OSMOSIS scheme 
     34   !!   'ln_zdfosm'                                             OSMOSIS scheme 
    3535   !!---------------------------------------------------------------------- 
    3636   !!   zdf_osm       : update momentum and tracer Kz from osm scheme 
  • NEMO/branches/2018/dev_r10164_HPC09_ESIWACE_PREP_MERGE/src/OCE/ZDF/zdfphy.F90

    r10170 r10368  
    8080         &             ln_zdfswm,                                    &     ! surface  wave-induced mixing 
    8181         &             ln_zdfiwm,                                    &     ! internal  -      -      - 
     82         &             ln_zad_Aimp,                                  &     ! apdative-implicit vertical advection 
    8283         &             rn_avm0, rn_avt0, nn_avb, nn_havtb                  ! coefficients 
    8384      !!---------------------------------------------------------------------- 
     
    101102      IF(lwp) THEN                      ! Parameter print 
    102103         WRITE(numout,*) '   Namelist namzdf : set vertical mixing mixing parameters' 
     104         WRITE(numout,*) '      adaptive-implicit vertical advection' 
     105         WRITE(numout,*) '         Courant number targeted application   ln_zad_Aimp = ', ln_zad_Aimp 
    103106         WRITE(numout,*) '      vertical closure scheme' 
    104107         WRITE(numout,*) '         constant vertical mixing coefficient    ln_zdfcst = ', ln_zdfcst 
     
    127130      ENDIF 
    128131 
     132      IF( ln_zad_Aimp ) THEN 
     133         IF( zdf_phy_alloc() /= 0 )   & 
     134        &       CALL ctl_stop( 'STOP', 'zdf_phy_init : unable to allocate adaptive-implicit z-advection arrays' ) 
     135         wi(:,:,:) = 0._wp 
     136      ENDIF 
    129137      !                          !==  Background eddy viscosity and diffusivity  ==! 
    130138      IF( nn_avb == 0 ) THEN             ! Define avmb, avtb from namelist parameter 
     
    316324      ! 
    317325   END SUBROUTINE zdf_phy 
     326   INTEGER FUNCTION zdf_phy_alloc() 
     327      !!---------------------------------------------------------------------- 
     328      !!                 ***  FUNCTION zdf_phy_alloc  *** 
     329      !!---------------------------------------------------------------------- 
     330     ! Allocate wi array (declared in oce.F90) for use with the adaptive-implicit vertical velocity option 
     331     ALLOCATE(     wi(jpi,jpj,jpk), Cu_adv(jpi,jpj,jpk),  STAT= zdf_phy_alloc ) 
     332     IF( zdf_phy_alloc /= 0 )   CALL ctl_warn('zdf_phy_alloc: failed to allocate ln_zad_Aimp=T required arrays') 
     333     IF( lk_mpp             )   CALL mpp_sum ( 'zdfphy', zdf_phy_alloc ) 
     334   END FUNCTION zdf_phy_alloc 
    318335 
    319336   !!====================================================================== 
  • NEMO/branches/2018/dev_r10164_HPC09_ESIWACE_PREP_MERGE/src/OCE/nemogcm.F90

    r10360 r10368  
    421421      IF( ln_traqsr    )   CALL tra_qsr_init      ! penetrative solar radiation qsr 
    422422                           CALL tra_bbc_init      ! bottom heat flux 
    423       IF( ln_trabbl    )   CALL tra_bbl_init      ! advective (and/or diffusive) bottom boundary layer scheme 
     423                           CALL tra_bbl_init      ! advective (and/or diffusive) bottom boundary layer scheme 
    424424                           CALL tra_dmp_init      ! internal tracer damping 
    425425                           CALL tra_adv_init      ! horizontal & vertical advection 
  • NEMO/branches/2018/dev_r10164_HPC09_ESIWACE_PREP_MERGE/src/OCE/oce.F90

    r10068 r10368  
    2222   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:)   ::   vb   ,  vn    , va     !: j-horizontal velocity        [m/s] 
    2323   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:)   ::           wn             !: vertical velocity            [m/s] 
     24   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:)   ::           wi             !: vertical vel. (adaptive-implicit) [m/s] 
    2425   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:)   ::           hdivn          !: horizontal divergence        [s-1] 
    2526   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:,:) ::   tsb  ,  tsn   , tsa    !: 4D T-S fields                  [Celsius,psu]  
     
    2930   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:) ::   rhd    !: in situ density anomalie rhd=(rho-rau0)/rau0  [no units] 
    3031   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:) ::   rhop   !: potential volumic mass                           [kg/m3] 
     32   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:) ::   Cu_adv                   !: vertical Courant number (adaptive-implicit) 
    3133 
    3234   !! free surface                                      !  before  ! now    ! after  ! 
  • NEMO/branches/2018/dev_r10164_HPC09_ESIWACE_PREP_MERGE/src/OCE/step.F90

    r10068 r10368  
    160160 
    161161                            CALL ssh_nxt       ( kstp )  ! after ssh (includes call to div_hor) 
    162       IF(.NOT.ln_linssh )   CALL dom_vvl_sf_nxt( kstp )  ! after vertical scale factors  
     162      IF( .NOT.ln_linssh )  CALL dom_vvl_sf_nxt( kstp )  ! after vertical scale factors  
    163163                            CALL wzv           ( kstp )  ! now cross-level velocity  
     164      IF( ln_zad_Aimp )     CALL wAimp         ( kstp )  ! Adaptive-implicit vertical advection partitioning 
    164165                            CALL eos    ( tsn, rhd, rhop, gdept_n(:,:,:) )  ! now in situ density for hpg computation 
    165166                             
     
    198199         IF(.NOT.ln_linssh) CALL dom_vvl_sf_nxt( kstp, kcall=2 )  ! after vertical scale factors (update depth average component) 
    199200                            CALL wzv        ( kstp )              ! now cross-level velocity  
     201         IF( ln_zad_Aimp )  CALL wAimp      ( kstp )  ! Adaptive-implicit vertical advection partitioning 
    200202      ENDIF 
    201203       
  • NEMO/branches/2018/dev_r10164_HPC09_ESIWACE_PREP_MERGE/src/OCE/stpctl.F90

    r10358 r10368  
    2424   USE lbclnk          ! ocean lateral boundary conditions (or mpp link) 
    2525   USE lib_mpp         ! distributed memory computing 
     26   USE zdf_oce ,  ONLY : ln_zad_Aimp       ! ocean vertical physics variables 
    2627   USE wet_dry,   ONLY : ll_wd, ssh_ref    ! reference depth for negative bathy 
    2728 
     
    3233   PUBLIC stp_ctl           ! routine called by step.F90 
    3334 
    34    INTEGER  ::   idrun, idtime, idssh, idu, ids1, ids2, istatus 
     35   INTEGER  ::   idrun, idtime, idssh, idu, ids1, ids2, idt1, idt2, idc1, idw1, istatus 
    3536   LOGICAL  ::   lsomeoce 
    3637   !!---------------------------------------------------------------------- 
     
    6566      INTEGER, DIMENSION(3)  ::   iu, is1, is2        ! min/max loc indices 
    6667      REAL(wp)               ::   zzz                 ! local real  
    67       REAL(wp), DIMENSION(5) ::   zmax 
     68      REAL(wp), DIMENSION(9) ::   zmax 
    6869      CHARACTER(len=20) :: clname 
    6970      !!---------------------------------------------------------------------- 
     
    8687            istatus = NF90_DEF_VAR( idrun,       's_min', NF90_DOUBLE, (/ idtime /), ids1  ) 
    8788            istatus = NF90_DEF_VAR( idrun,       's_max', NF90_DOUBLE, (/ idtime /), ids2  ) 
     89            istatus = NF90_DEF_VAR( idrun,       't_min', NF90_DOUBLE, (/ idtime /), idt1  ) 
     90            istatus = NF90_DEF_VAR( idrun,       't_max', NF90_DOUBLE, (/ idtime /), idt2  ) 
     91            IF( ln_zad_Aimp ) THEN 
     92               istatus = NF90_DEF_VAR( idrun,   'abs_wi_max', NF90_DOUBLE, (/ idtime /), idw1  ) 
     93               istatus = NF90_DEF_VAR( idrun,       'Cu_max', NF90_DOUBLE, (/ idtime /), idc1  ) 
     94            ENDIF 
    8895            istatus = NF90_ENDDEF(idrun) 
     96            zmax(8:9) = 0._wp    ! initialise to zero in case ln_zad_Aimp option is not in use 
    8997         ENDIF 
    9098      ENDIF 
     
    105113      zmax(3) = MAXVAL( -tsn(:,:,:,jp_sal) , mask = tmask(:,:,:) == 1._wp )   ! minus salinity max 
    106114      zmax(4) = MAXVAL(  tsn(:,:,:,jp_sal) , mask = tmask(:,:,:) == 1._wp )   !       salinity max 
    107       zmax(5) = REAL( nstop , wp )                                            ! stop indicator 
     115      zmax(5) = MAXVAL( -tsn(:,:,:,jp_tem) , mask = tmask(:,:,:) == 1._wp )   ! minus temperature max 
     116      zmax(6) = MAXVAL(  tsn(:,:,:,jp_tem) , mask = tmask(:,:,:) == 1._wp )   !       temperature max 
     117      zmax(7) = REAL( nstop , wp )                                            ! stop indicator 
     118      IF( ln_zad_Aimp ) THEN 
     119         zmax(8) = MAXVAL(  ABS( wi(:,:,:) ) , mask = wmask(:,:,:) == 1._wp ) ! implicit vertical vel. max 
     120         zmax(9) = MAXVAL(   Cu_adv(:,:,:)   , mask = tmask(:,:,:) == 1._wp ) !       cell Courant no. max 
     121      ENDIF 
    108122      ! 
    109123      IF( lk_mpp .AND. ln_ctl ) THEN 
    110124         CALL mpp_max( "stpctl", zmax )          ! max over the global domain 
    111          nstop = NINT( zmax(5) )                 ! nstop indicator sheared among all local domains 
     125         nstop = NINT( zmax(7) )                 ! nstop indicator sheared among all local domains 
    112126      ENDIF 
    113127      !                                   !==  run statistics  ==!   ("run.stat" files) 
     
    118132         istatus = NF90_PUT_VAR( idrun,  ids1, (/-zmax(3)/), (/kt/), (/1/) ) 
    119133         istatus = NF90_PUT_VAR( idrun,  ids2, (/ zmax(4)/), (/kt/), (/1/) ) 
     134         istatus = NF90_PUT_VAR( idrun,  idt1, (/-zmax(5)/), (/kt/), (/1/) ) 
     135         istatus = NF90_PUT_VAR( idrun,  idt2, (/ zmax(6)/), (/kt/), (/1/) ) 
     136         IF( ln_zad_Aimp ) THEN 
     137            istatus = NF90_PUT_VAR( idrun,  idw1, (/ zmax(8)/), (/kt/), (/1/) ) 
     138            istatus = NF90_PUT_VAR( idrun,  idc1, (/ zmax(9)/), (/kt/), (/1/) ) 
     139         ENDIF 
    120140         IF( MOD( kt , 100 ) == 0 ) istatus = NF90_SYNC(idrun) 
    121141         IF( kt == nitend         ) istatus = NF90_CLOSE(idrun) 
  • NEMO/branches/2018/dev_r10164_HPC09_ESIWACE_PREP_MERGE/src/TOP/PISCES/P4Z/p4zfechem.F90

    r10069 r10368  
    2727   LOGICAL          ::   ln_fechem    !: boolean for complex iron chemistry following Tagliabue and voelker 
    2828   LOGICAL          ::   ln_ligvar    !: boolean for variable ligand concentration following Tagliabue and voelker 
    29    LOGICAL          ::   ln_fecolloid !: boolean for variable colloidal fraction 
    3029   REAL(wp), PUBLIC ::   xlam1        !: scavenging rate of Iron  
    3130   REAL(wp), PUBLIC ::   xlamdust     !: scavenging rate of Iron by dust  
     
    6867      REAL(wp) ::   za, zb, zc, zkappa1, zkappa2, za0, za1, za2 
    6968      REAL(wp) ::   zxs, zfunc, zp, zq, zd, zr, zphi, zfff, zp3, zq2 
    70       REAL(wp) ::   ztfe, zoxy, zhplus 
     69      REAL(wp) ::   ztfe, zoxy, zhplus, zxlam 
    7170      REAL(wp) ::   zaggliga, zaggligb 
    7271      REAL(wp) ::   dissol, zligco 
     72      REAL(wp) :: zrfact2 
    7373      CHARACTER (len=25) :: charout 
    7474      REAL(wp), DIMENSION(jpi,jpj,jpk) ::   zTL1, zFe3, ztotlig, precip, zFeL1 
     75      REAL(wp), DIMENSION(jpi,jpj,jpk) ::   zcoll3d, zscav3d, zlcoll3d 
    7576      REAL(wp), ALLOCATABLE, DIMENSION(:,:,:) ::   zFeL2, zTL2, zFe2, zFeP 
    7677      REAL(wp), ALLOCATABLE, DIMENSION(:,:  ) ::   zstrn, zstrn2 
     
    129130         ! ------------------------------------------------------------ 
    130131         DO jn = 1, 2 
    131           DO jk = 1, jpkm1 
    132             DO jj = 1, jpj 
    133                DO ji = 1, jpi 
    134                   zlight = etot(ji,jj,jk) * zstrn(ji,jj) * REAL( 2-jn, wp ) 
    135                   zzstrn2 = zstrn2(ji,jj) * REAL( 2-jn, wp ) + (1. - zstrn2(ji,jj) ) * REAL( jn-1, wp ) 
    136                   ! Calculate ligand concentrations : assume 2/3rd of excess goes to 
    137                   ! strong ligands (L1) and 1/3rd to weak ligands (L2) 
    138                   ztligand       = ztotlig(ji,jj,jk) - ligand * 1E9 
    139                   zTL1(ji,jj,jk) =                0.000001 + 0.67 * ztligand 
    140                   zTL2(ji,jj,jk) = ligand * 1E9 - 0.000001 + 0.33 * ztligand 
    141                   ! ionic strength from Millero et al. 1987 
    142                   zph    = -LOG10( MAX( hi(ji,jj,jk), rtrn) ) 
    143                   zoxy   = trb(ji,jj,jk,jpoxy) 
    144                   ! Fe2+ oxydation rate from Santana-Casiano et al. (2005) 
    145                   zkox   = 35.407 - 6.7109 * zph + 0.5342 * zph * zph - 5362.6 / ( tempis(ji,jj,jk) + 273.15 )  & 
    146                     &    - 0.04406 * SQRT( salinprac(ji,jj,jk) ) - 0.002847 * salinprac(ji,jj,jk) 
    147                   zkox   = ( 10.** zkox ) * spd 
    148                   zkox   = zkox * MAX( 1.e-6, zoxy) / ( chemo2(ji,jj,jk) + rtrn ) 
    149                   ! PHOTOREDUCTION of complexed iron : Tagliabue and Arrigo (2006) 
    150                   zkph2 = MAX( 0., 15. * zlight / ( zlight + 2. ) ) * (1. - fr_i(ji,jj)) 
    151                   zkph1 = zkph2 / 5. 
    152                   ! pass the dfe concentration from PISCES 
    153                   ztfe = trb(ji,jj,jk,jpfer) * 1e9 
    154                   ! ---------------------------------------------------------- 
    155                   ! ANALYTICAL SOLUTION OF ROOTS OF THE FE3+ EQUATION 
    156                   ! As shown in Tagliabue and Voelker (2009), Fe3+ is the root of a 3rd order polynom.  
    157                   ! ---------------------------------------------------------- 
    158                   ! calculate some parameters 
    159                   za = 1 + ks / kpr 
    160                   zb = 1 + ( zkph1 + kth ) / ( zkox + rtrn ) 
    161                   zc = 1 + zkph2 / ( zkox + rtrn ) 
    162                   zkappa1 = ( kb1 + zkph1 + kth ) / kl1 
    163                   zkappa2 = ( kb2 + zkph2 ) / kl2 
    164                   za2 = zTL1(ji,jj,jk) * zb / za + zTL2(ji,jj,jk) * zc / za + zkappa1 + zkappa2 - ztfe / za 
    165                   za1 = zkappa2 * zTL1(ji,jj,jk) * zb / za + zkappa1 * zTL2(ji,jj,jk) * zc / za & 
     132            DO jk = 1, jpkm1 
     133               DO jj = 1, jpj 
     134                  DO ji = 1, jpi 
     135                     zlight = etot(ji,jj,jk) * zstrn(ji,jj) * REAL( 2-jn, wp ) 
     136                     zzstrn2 = zstrn2(ji,jj) * REAL( 2-jn, wp ) + (1. - zstrn2(ji,jj) ) * REAL( jn-1, wp ) 
     137                     ! Calculate ligand concentrations : assume 2/3rd of excess goes to 
     138                     ! strong ligands (L1) and 1/3rd to weak ligands (L2) 
     139                     ztligand       = ztotlig(ji,jj,jk) - ligand * 1E9 
     140                     zTL1(ji,jj,jk) =                0.000001 + 0.67 * ztligand 
     141                     zTL2(ji,jj,jk) = ligand * 1E9 - 0.000001 + 0.33 * ztligand 
     142                     ! ionic strength from Millero et al. 1987 
     143                     zph    = -LOG10( MAX( hi(ji,jj,jk), rtrn) ) 
     144                     zoxy   = trb(ji,jj,jk,jpoxy) 
     145                     ! Fe2+ oxydation rate from Santana-Casiano et al. (2005) 
     146                     zkox   = 35.407 - 6.7109 * zph + 0.5342 * zph * zph - 5362.6 / ( tempis(ji,jj,jk) + 273.15 )  & 
     147                     &        - 0.04406 * SQRT( salinprac(ji,jj,jk) ) - 0.002847 * salinprac(ji,jj,jk) 
     148                     zkox   = ( 10.** zkox ) * spd 
     149                     zkox   = zkox * MAX( 1.e-6, zoxy) / ( chemo2(ji,jj,jk) + rtrn ) 
     150                     ! PHOTOREDUCTION of complexed iron : Tagliabue and Arrigo (2006) 
     151                     zkph2 = MAX( 0., 15. * zlight / ( zlight + 2. ) ) * (1. - fr_i(ji,jj)) 
     152                     zkph1 = zkph2 / 5. 
     153                     ! pass the dfe concentration from PISCES 
     154                     ztfe = trb(ji,jj,jk,jpfer) * 1e9 
     155                     ! ---------------------------------------------------------- 
     156                     ! ANALYTICAL SOLUTION OF ROOTS OF THE FE3+ EQUATION 
     157                     ! As shown in Tagliabue and Voelker (2009), Fe3+ is the root of a 3rd order polynom.  
     158                     ! ---------------------------------------------------------- 
     159                     ! calculate some parameters 
     160                     za = 1.0 + ks / kpr 
     161                     zb = 1.0 + zkph2 / ( zkox ) 
     162                     zc = 1.0 + ( zkph1 + kth ) / ( zkox ) 
     163                     zkappa1 = ( kb1 + zkph1 + kth ) / kl1 
     164                     zkappa2 = ( kb2 + zkph2 ) / kl2 
     165                     za2 = zTL2(ji,jj,jk) * zb / za + zTL2(ji,jj,jk) * zc / za + zkappa1 + zkappa2 - ztfe / za 
     166                     za1 = zkappa1 * zTL2(ji,jj,jk) * zb / za + zkappa2 * zTL1(ji,jj,jk) * zc / za & 
    166167                      & + zkappa1 * zkappa2 - ( zkappa1 + zkappa2 ) * ztfe / za 
    167                   za0 = -zkappa1 * zkappa2 * ztfe / za 
    168                   zp  = za1 - za2 * za2 / 3. 
    169                   zq  = za2 * za2 * za2 * 2. / 27. - za2 * za1 / 3. + za0 
    170                   zp3 = zp / 3. 
    171                   zq2 = zq / 2. 
    172                   zd  = zp3 * zp3 * zp3 + zq2 * zq2 
    173                   zr  = zq / ABS( zq ) * SQRT( ABS( zp ) / 3. ) 
    174                   ! compute the roots 
    175                   IF( zp > 0.) THEN 
    176                      ! zphi = ASINH( zq / ( 2. * zr * zr * zr ) ) 
    177                      zphi =  zq / ( 2. * zr * zr * zr )  
    178                      zphi = LOG( zphi + SQRT( zphi * zphi + 1 ) )  ! asinh(x) = log(x + sqrt(x^2+1)) 
    179                      zxs  = -2. * zr * SINH( zphi / 3. ) - za1 / 3. 
    180                   ELSE 
    181                      IF( zd > 0. ) THEN 
    182                         zfff = MAX( 1., zq / ( 2. * zr * zr * zr ) ) 
    183                         ! zphi = ACOSH( zfff ) 
    184                         zphi = LOG( zfff + SQRT( zfff * zfff - 1 ) )  ! acosh(x) = log(x + sqrt(x^2-1)) 
    185                         zxs = -2. * zr * COSH( zphi / 3. ) - za1 / 3. 
     168                     za0 = -zkappa1 * zkappa2 * ztfe / za 
     169                     zp  = za1 - za2 * za2 / 3. 
     170                     zq  = za2 * za2 * za2 * 2. / 27. - za2 * za1 / 3. + za0 
     171                     zp3 = zp / 3. 
     172                     zq2 = zq / 2. 
     173                     zd  = zp3 * zp3 * zp3 + zq2 * zq2 
     174                     zr  = zq / ABS( zq ) * SQRT( ABS( zp ) / 3. ) 
     175                     ! compute the roots 
     176                     IF( zp > 0.) THEN 
     177                        ! zphi = ASINH( zq / ( 2. * zr * zr * zr ) ) 
     178                        zphi =  zq / ( 2. * zr * zr * zr )  
     179                        zphi = LOG( zphi + SQRT( zphi * zphi + 1 ) )  ! asinh(x) = log(x + sqrt(x^2+1)) 
     180                        zxs  = -2. * zr * SINH( zphi / 3. ) - za1 / 3. 
    186181                     ELSE 
    187                         zfff = MIN( 1., zq / ( 2. * zr * zr * zr ) ) 
    188                         zphi = ACOS( zfff ) 
    189                         DO jic = 1, 3 
    190                            zfunc = -2 * zr * COS( zphi / 3. + 2. * REAL( jic - 1, wp ) * rpi / 3. ) - za2 / 3. 
    191                            IF( zfunc > 0. .AND. zfunc <= ztfe)  zxs = zfunc 
    192                         END DO 
     182                        IF( zd > 0. ) THEN 
     183                           zfff = MAX( 1., zq / ( 2. * zr * zr * zr ) ) 
     184                           ! zphi = ACOSH( zfff ) 
     185                           zphi = LOG( zfff + SQRT( zfff * zfff - 1 ) )  ! acosh(x) = log(x + sqrt(x^2-1)) 
     186                           zxs = -2. * zr * COSH( zphi / 3. ) - za1 / 3. 
     187                        ELSE 
     188                           zfff = MIN( 1., zq / ( 2. * zr * zr * zr ) ) 
     189                           zphi = ACOS( zfff ) 
     190                           DO jic = 1, 3 
     191                              zfunc = -2 * zr * COS( zphi / 3. + 2. * REAL( jic - 1, wp ) * rpi / 3. ) - za2 / 3. 
     192                              IF( zfunc > 0. .AND. zfunc <= ztfe)  zxs = zfunc 
     193                           END DO 
     194                        ENDIF 
    193195                     ENDIF 
    194                   ENDIF 
    195                   ! solve for the other Fe species 
    196                   zzFe3 = MAX( 0., zxs ) 
    197                   zzFep = MAX( 0., ( ks * zzFe3 / kpr ) ) 
    198                   zkappa2 = ( kb2 + zkph2 ) / kl2 
    199                   zzFeL2 = MAX( 0., ( zzFe3 * zTL2(ji,jj,jk) ) / ( zkappa2 + zzFe3 ) ) 
    200                   zzFeL1 = MAX( 0., ( ztfe / zb - za / zb * zzFe3 - zc / zb * zzFeL2 ) ) 
    201                   zzFe2  = MAX( 0., ( ( zkph1 * zzFeL1 + zkph2 * zzFeL2 ) / zkox ) ) 
    202                   zFe3(ji,jj,jk)  = zFe3(ji,jj,jk)  + zzFe3 * zzstrn2 
    203                   zFe2(ji,jj,jk)  = zFe2(ji,jj,jk)  + zzFe2 * zzstrn2 
    204                   zFeL2(ji,jj,jk) = zFeL2(ji,jj,jk) + zzFeL2 * zzstrn2 
    205                   zFeL1(ji,jj,jk) = zFeL1(ji,jj,jk) + zzFeL1 * zzstrn2 
    206                   zFeP(ji,jj,jk)  = zFeP(ji,jj,jk)  + zzFeP * zzstrn2 
     196                     ! solve for the other Fe species 
     197                     zzFe3  = MAX( 0., zxs ) 
     198                     zzFep  = MAX( 0., ( ks * zzFe3 / kpr ) ) 
     199                     zzFeL1 = MAX( 0., ( zzFe3 * zTL1(ji,jj,jk) ) / ( zkappa1 + zzFe3 ) ) 
     200                     zzFeL2 = (ztfe - za * zzFe3 - zc * zzFeL1 ) / zb 
     201                     zzFe2  = MAX( 0., ( ( ( zkph1 + kth ) * zzFeL1 + zkph2 * zzFeL2 ) / zkox ) ) 
     202                     zzFep  = ztfe - zzFe3 - zzFe2 - zzFeL1 - zzFeL2 
     203                     zFe3(ji,jj,jk)  = zFe3(ji,jj,jk)  + zzFe3 * zzstrn2 
     204                     zFe2(ji,jj,jk)  = zFe2(ji,jj,jk)  + zzFe2 * zzstrn2 
     205                     zFeL2(ji,jj,jk) = zFeL2(ji,jj,jk) + zzFeL2 * zzstrn2 
     206                     zFeL1(ji,jj,jk) = zFeL1(ji,jj,jk) + zzFeL1 * zzstrn2 
     207                     zFeP(ji,jj,jk)  = zFeP(ji,jj,jk)  + zzFeP * zzstrn2 
     208                  END DO 
    207209               END DO 
    208210            END DO 
    209          END DO 
    210211         END DO 
    211212      ELSE 
     
    218219            DO jj = 1, jpj 
    219220               DO ji = 1, jpi 
    220                   zTL1(ji,jj,jk) = ztotlig(ji,jj,jk) 
    221                   zkeq           = fekeq(ji,jj,jk) 
    222                   zfesatur       = zTL1(ji,jj,jk) * 1E-9 
    223                   ztfe           = trb(ji,jj,jk,jpfer)  
     221                  zTL1(ji,jj,jk)  = ztotlig(ji,jj,jk) 
     222                  zkeq            = fekeq(ji,jj,jk) 
     223                  zfesatur        = zTL1(ji,jj,jk) * 1E-9 
     224                  ztfe            = trb(ji,jj,jk,jpfer)  
    224225                  ! Fe' is the root of a 2nd order polynom 
    225226                  zFe3 (ji,jj,jk) = ( -( 1. + zfesatur * zkeq - zkeq * ztfe )               & 
    226                      &             + SQRT( ( 1. + zfesatur * zkeq - zkeq * ztfe )**2       & 
    227                      &               + 4. * ztfe * zkeq) ) / ( 2. * zkeq ) 
     227                     &              + SQRT( ( 1. + zfesatur * zkeq - zkeq * ztfe )**2       & 
     228                     &              + 4. * ztfe * zkeq) ) / ( 2. * zkeq ) 
    228229                  zFe3 (ji,jj,jk) = zFe3(ji,jj,jk) * 1E9 
    229230                  zFeL1(ji,jj,jk) = MAX( 0., trb(ji,jj,jk,jpfer) * 1E9 - zFe3(ji,jj,jk) ) 
     
    242243               ! Scavenging onto dust is also included as evidenced from the DUNE experiments. 
    243244               ! -------------------------------------------------------------------------------------- 
     245               zhplus  = max( rtrn, hi(ji,jj,jk) ) 
     246               fe3sol  = fesol(ji,jj,jk,1) * ( zhplus**3 + fesol(ji,jj,jk,2) * zhplus**2  & 
     247               &         + fesol(ji,jj,jk,3) * zhplus + fesol(ji,jj,jk,4)     & 
     248               &         + fesol(ji,jj,jk,5) / zhplus ) 
    244249               IF( ln_fechem ) THEN 
    245250                  zfeequi = ( zFe3(ji,jj,jk) + zFe2(ji,jj,jk) + zFeP(ji,jj,jk) ) * 1E-9 
    246251                  zfecoll = ( 0.3 * zFeL1(ji,jj,jk) + 0.5 * zFeL2(ji,jj,jk) ) * 1E-9 
     252                  precip(ji,jj,jk) = 0.0 
    247253               ELSE 
    248254                  zfeequi = zFe3(ji,jj,jk) * 1E-9 
    249                   IF (ln_fecolloid) THEN 
    250                      zhplus   = max( rtrn, hi(ji,jj,jk) ) 
    251                      fe3sol  = fesol(ji,jj,jk,1) * ( zhplus**3 + fesol(ji,jj,jk,2) * zhplus**2  & 
    252                      &         + fesol(ji,jj,jk,3) * zhplus + fesol(ji,jj,jk,4)     & 
    253                      &         + fesol(ji,jj,jk,5) / zhplus ) 
    254                      zfecoll = max( ( 0.1 * zFeL1(ji,jj,jk) * 1E-9 ), ( zFeL1(ji,jj,jk) * 1E-9 -fe3sol ) ) 
    255                   ELSE 
    256                      zfecoll = 0.5 * zFeL1(ji,jj,jk) * 1E-9 
    257                      fe3sol  = 0. 
    258                   ENDIF 
     255                  zhplus  = max( rtrn, hi(ji,jj,jk) ) 
     256                  fe3sol  = fesol(ji,jj,jk,1) * ( zhplus**3 + fesol(ji,jj,jk,2) * zhplus**2  & 
     257                  &         + fesol(ji,jj,jk,3) * zhplus + fesol(ji,jj,jk,4)     & 
     258                  &         + fesol(ji,jj,jk,5) / zhplus ) 
     259                  zfecoll = 0.5 * zFeL1(ji,jj,jk) * 1E-9 
     260                  ! precipitation of Fe3+, creation of nanoparticles 
     261                  precip(ji,jj,jk) = MAX( 0., ( zFe3(ji,jj,jk) * 1E-9 - fe3sol ) ) * kfep * xstep 
    259262               ENDIF 
    260263               ! 
    261264               ztrc   = ( trb(ji,jj,jk,jppoc) + trb(ji,jj,jk,jpgoc) + trb(ji,jj,jk,jpcal) + trb(ji,jj,jk,jpgsi) ) * 1.e6  
    262                IF( ln_dust )  zdust  = dust(ji,jj) / ( wdust / rday ) * tmask(ji,jj,jk) ! dust in kg/m2/s 
    263                zlam1b = 3.e-5 + xlamdust * zdust + xlam1 * ztrc 
     265               IF( ln_dust )  zdust  = dust(ji,jj) / ( wdust / rday ) * tmask(ji,jj,jk) & 
     266               &  * EXP( -gdept_n(ji,jj,jk) / 540. ) 
     267               IF (ln_ligand) THEN 
     268                  zxlam  = xlam1 * MAX( 1.E-3, EXP(-2 * etot(ji,jj,jk) / 10. ) * (1. - EXP(-2 * trb(ji,jj,jk,jpoxy) / 100.E-6 ) )) 
     269               ELSE 
     270                  zxlam  = xlam1 * 1.0 
     271               ENDIF 
     272               zlam1b = 3.e-5 + xlamdust * zdust + zxlam * ztrc 
    264273               zscave = zfeequi * zlam1b * xstep 
    265274 
     
    267276               ! to later allocate scavenged iron to the different organic pools 
    268277               ! --------------------------------------------------------- 
    269                zdenom1 = xlam1 * trb(ji,jj,jk,jppoc) / zlam1b 
    270                zdenom2 = xlam1 * trb(ji,jj,jk,jpgoc) / zlam1b 
     278               zdenom1 = zxlam * trb(ji,jj,jk,jppoc) / zlam1b 
     279               zdenom2 = zxlam * trb(ji,jj,jk,jpgoc) / zlam1b 
    271280 
    272281               !  Increased scavenging for very high iron concentrations found near the coasts  
     
    276285               zlamfac = MIN( 1.  , zlamfac ) 
    277286               zdep    = MIN( 1., 1000. / gdept_n(ji,jj,jk) ) 
    278                zlam1b  = xlam1 * MAX( 0.e0, ( trb(ji,jj,jk,jpfer) * 1.e9 - ztotlig(ji,jj,jk) ) ) 
    279                zcoag   = zfeequi * zlam1b * xstep + 1E-4 * ( 1. - zlamfac ) * zdep * xstep * trb(ji,jj,jk,jpfer) 
     287               zcoag   = 1E-4 * ( 1. - zlamfac ) * zdep * xstep * trb(ji,jj,jk,jpfer) 
    280288 
    281289               !  Compute the coagulation of colloidal iron. This parameterization  
     
    283291               !  It requires certainly some more work as it is very poorly constrained. 
    284292               !  ---------------------------------------------------------------- 
    285                zlam1a  = ( 0.369  * 0.3 * trb(ji,jj,jk,jpdoc) + 102.4  * trb(ji,jj,jk,jppoc) ) * xdiss(ji,jj,jk)    & 
    286                    &   + ( 114.   * 0.3 * trb(ji,jj,jk,jpdoc) ) 
     293               zlam1a   = ( 0.369  * 0.3 * trb(ji,jj,jk,jpdoc) + 102.4  * trb(ji,jj,jk,jppoc) ) * xdiss(ji,jj,jk)    & 
     294                   &      + ( 114.   * 0.3 * trb(ji,jj,jk,jpdoc) ) 
    287295               zaggdfea = zlam1a * xstep * zfecoll 
    288296               ! 
    289                zlam1b = 3.53E3 *  trb(ji,jj,jk,jpgoc) * xdiss(ji,jj,jk) 
     297               zlam1b   = 3.53E3 * trb(ji,jj,jk,jpgoc) * xdiss(ji,jj,jk) 
    290298               zaggdfeb = zlam1b * xstep * zfecoll 
    291                ! 
    292                ! precipitation of Fe3+, creation of nanoparticles 
    293                precip(ji,jj,jk) = MAX( 0., ( zfeequi - fe3sol ) ) * kfep * xstep 
    294299               ! 
    295300               tra(ji,jj,jk,jpfer) = tra(ji,jj,jk,jpfer) - zscave - zaggdfea - zaggdfeb & 
     
    297302               tra(ji,jj,jk,jpsfe) = tra(ji,jj,jk,jpsfe) + zscave * zdenom1 + zaggdfea 
    298303               tra(ji,jj,jk,jpbfe) = tra(ji,jj,jk,jpbfe) + zscave * zdenom2 + zaggdfeb 
     304               zscav3d(ji,jj,jk)   = zscave 
     305               zcoll3d(ji,jj,jk)   = zaggdfea + zaggdfeb 
    299306               ! 
    300307            END DO 
     
    317324                  ! 
    318325                  zlam1b   = 3.53E3 *   trb(ji,jj,jk,jpgoc) * xdiss(ji,jj,jk) 
    319                   zligco   = MAX( ( 0.1 * trb(ji,jj,jk,jplgw) ), ( trb(ji,jj,jk,jplgw) - fe3sol ) ) 
     326                  zligco   = 0.5 * trn(ji,jj,jk,jplgw) 
    320327                  zaggliga = zlam1a * xstep * zligco 
    321328                  zaggligb = zlam1b * xstep * zligco 
    322329                  tra(ji,jj,jk,jpfep) = tra(ji,jj,jk,jpfep) + precip(ji,jj,jk) 
    323330                  tra(ji,jj,jk,jplgw) = tra(ji,jj,jk,jplgw) - zaggliga - zaggligb 
     331                  zlcoll3d(ji,jj,jk)  = zaggliga + zaggligb 
    324332               END DO 
    325333            END DO 
     
    328336         IF( .NOT.ln_fechem) THEN 
    329337            plig(:,:,:) =  MAX( 0., ( ( zFeL1(:,:,:) * 1E-9 ) / ( trb(:,:,:,jpfer) +rtrn ) ) ) 
    330             plig(:,:,:) =  MAX( 0. , plig(:,:,:) ) 
    331338         ENDIF 
    332339         ! 
     
    336343      IF( lk_iomput ) THEN 
    337344         IF( knt == nrdttrc ) THEN 
     345         zrfact2 = 1.e3 * rfact2r  ! conversion from mol/L/timestep into mol/m3/s 
    338346         IF( iom_use("Fe3")    )  CALL iom_put("Fe3"    , zFe3   (:,:,:)       * tmask(:,:,:) )   ! Fe3+ 
    339347         IF( iom_use("FeL1")   )  CALL iom_put("FeL1"   , zFeL1  (:,:,:)       * tmask(:,:,:) )   ! FeL1 
    340348         IF( iom_use("TL1")    )  CALL iom_put("TL1"    , zTL1   (:,:,:)       * tmask(:,:,:) )   ! TL1 
    341349         IF( iom_use("Totlig") )  CALL iom_put("Totlig" , ztotlig(:,:,:)       * tmask(:,:,:) )   ! TL 
    342          IF( iom_use("Biron")  )  CALL iom_put("Biron"  , biron  (:,:,:) * 1e9 * tmask(:,:,:) )   ! biron 
     350         IF( iom_use("Biron")  )  CALL iom_put("Biron"  , biron  (:,:,:)  * 1e9 * tmask(:,:,:) )   ! biron 
     351         IF( iom_use("FESCAV") )  CALL iom_put("FESCAV" , zscav3d(:,:,:)  * 1e9 * tmask(:,:,:) * zrfact2 ) 
     352         IF( iom_use("FECOLL") )  CALL iom_put("FECOLL" , zcoll3d(:,:,:)  * 1e9 * tmask(:,:,:) * zrfact2 ) 
     353         IF( iom_use("LGWCOLL"))  CALL iom_put("LGWCOLL", zlcoll3d(:,:,:) * 1e9 * tmask(:,:,:) * zrfact2 ) 
    343354         IF( ln_fechem ) THEN 
    344355            IF( iom_use("Fe2")  ) CALL iom_put("Fe2"    , zFe2   (:,:,:)       * tmask(:,:,:) )   ! Fe2+ 
     
    380391      INTEGER ::   ios   ! Local integer  
    381392      !! 
    382       NAMELIST/nampisfer/ ln_fechem, ln_ligvar, ln_fecolloid, xlam1, xlamdust, ligand, kfep  
     393      NAMELIST/nampisfer/ ln_fechem, ln_ligvar, xlam1, xlamdust, ligand, kfep  
    383394      !!---------------------------------------------------------------------- 
    384395      ! 
     
    401412         WRITE(numout,*) '      enable complex iron chemistry scheme      ln_fechem    =', ln_fechem 
    402413         WRITE(numout,*) '      variable concentration of ligand          ln_ligvar    =', ln_ligvar 
    403          WRITE(numout,*) '      Variable colloidal fraction of Fe3+       ln_fecolloid =', ln_fecolloid 
    404414         WRITE(numout,*) '      scavenging rate of Iron                   xlam1        =', xlam1 
    405415         WRITE(numout,*) '      scavenging rate of Iron by dust           xlamdust     =', xlamdust 
     
    407417         WRITE(numout,*) '      rate constant for nanoparticle formation  kfep         =', kfep 
    408418      ENDIF 
     419      !  
     420      IF (ln_ligand .AND. ln_fechem) CALL ctl_stop( 'STOP', 'p4z_fechem_init: ln_ligand and ln_fechem are incompatible') 
    409421      ! 
    410422      IF( ln_fechem ) THEN             ! set some constants used by the complexe chemistry scheme 
  • NEMO/branches/2018/dev_r10164_HPC09_ESIWACE_PREP_MERGE/src/TOP/PISCES/P4Z/p4zligand.F90

    r10069 r10368  
    1313   USE sms_pisces      ! PISCES Source Minus Sink variables 
    1414   USE prtctl_trc      ! print control for debugging 
     15   USE iom             !  I/O manager 
    1516 
    1617   IMPLICIT NONE 
     
    4344      INTEGER  ::   ji, jj, jk 
    4445      REAL(wp) ::   zlgwp, zlgwpr, zlgwr, zlablgw, zrfepa, zfepr 
     46      REAL(wp), DIMENSION(jpi,jpj,jpk) :: zligrem, zligpr, zrligprod 
     47      REAL(wp), ALLOCATABLE, DIMENSION(:,:,:) ::   zw3d 
    4548      CHARACTER (len=25) ::   charout 
    4649      !!--------------------------------------------------------------------- 
     
    5659               ! ------------------------------------------------------------------ 
    5760               ! production from remineralisation of organic matter 
    58                zlgwp  = orem(ji,jj,jk) * rlig 
     61               zlgwp = orem(ji,jj,jk) * rlig 
    5962               ! decay of weak ligand 
    6063               ! This is based on the idea that as LGW is lower 
    6164               ! there is a larger fraction of refractory OM 
    6265               zlgwr = max( rlgs , rlgw * exp( -2 * (trb(ji,jj,jk,jplgw)*1e9) ) ) ! years 
    63                zlgwr = 1. / zlgwr * tgfunc(ji,jj,jk) * ( xstep / nyear_len(1) ) * trb(ji,jj,jk,jplgw) 
     66               zlgwr = 1. / zlgwr * tgfunc(ji,jj,jk) * ( xstep / nyear_len(1) ) * blim(ji,jj,jk) * trb(ji,jj,jk,jplgw) 
    6467               ! photochem loss of weak ligand 
    6568               zlgwpr = prlgw * xstep * etot(ji,jj,jk) * trb(ji,jj,jk,jplgw) * (1. - fr_i(ji,jj)) 
    6669               tra(ji,jj,jk,jplgw) = tra(ji,jj,jk,jplgw) + zlgwp - zlgwr - zlgwpr 
     70               zligrem(ji,jj,jk)   = zlgwr 
     71               zligpr(ji,jj,jk)    = zlgwpr 
     72               zrligprod(ji,jj,jk) = zlgwp 
    6773               ! 
    6874               ! ---------------------------------------------------------- 
     
    7480               zrfepa = rfep * ( 1. - EXP( -1. * etot(ji,jj,jk) / 25. ) ) * (1.- fr_i(ji,jj)) 
    7581               zrfepa = MAX( (zrfepa / 10.0), zrfepa ) ! min of 10 days lifetime 
    76                zfepr = rfep * xstep * trb(ji,jj,jk,jpfep) 
     82               zfepr  = rfep * xstep * trb(ji,jj,jk,jpfep) 
    7783               tra(ji,jj,jk,jpfep) = tra(ji,jj,jk,jpfep) - zfepr 
    7884               tra(ji,jj,jk,jpfer) = tra(ji,jj,jk,jpfer) + zfepr 
     
    8187         END DO 
    8288      END DO 
     89      ! 
     90      !  Output of some diagnostics variables 
     91      !     --------------------------------- 
     92      IF( lk_iomput .AND. knt == nrdttrc ) THEN 
     93         ALLOCATE( zw3d(jpi,jpj,jpk) ) 
     94         IF( iom_use( "LIGREM" ) ) THEN 
     95            zw3d(:,:,:) = zligrem(:,:,:) * 1e9 * 1.e+3 * rfact2r * tmask(:,:,:) 
     96            CALL iom_put( "LIGREM", zw3d ) 
     97         ENDIF 
     98         IF( iom_use( "LIGPR" ) ) THEN 
     99            zw3d(:,:,:) = zligpr(:,:,:) * 1e9 * 1.e+3 * rfact2r * tmask(:,:,:)  
     100            CALL iom_put( "LIGPR", zw3d ) 
     101         ENDIF 
     102         IF( iom_use( "LPRODR" ) ) THEN 
     103            zw3d(:,:,:) = zrligprod(:,:,:) * 1e9 * 1.e+3 * rfact2r * tmask(:,:,:)  
     104            CALL iom_put( "LPRODR", zw3d ) 
     105         ENDIF 
     106         DEALLOCATE( zw3d ) 
     107      ENDIF 
    83108      ! 
    84109      IF(ln_ctl)   THEN  ! print mean trends (used for debugging) 
  • NEMO/branches/2018/dev_r10164_HPC09_ESIWACE_PREP_MERGE/src/TOP/PISCES/P4Z/p4zmeso.F90

    r10345 r10368  
    2525 
    2626   REAL(wp), PUBLIC ::  part2        !: part of calcite not dissolved in mesozoo guts 
    27    REAL(wp), PUBLIC ::  xprefc       !: mesozoo preference for POC  
    28    REAL(wp), PUBLIC ::  xprefp       !: mesozoo preference for nanophyto 
    29    REAL(wp), PUBLIC ::  xprefz       !: mesozoo preference for diatoms 
    30    REAL(wp), PUBLIC ::  xprefpoc     !: mesozoo preference for POC  
     27   REAL(wp), PUBLIC ::  xpref2d      !: mesozoo preference for diatoms 
     28   REAL(wp), PUBLIC ::  xpref2n      !: mesozoo preference for nanophyto 
     29   REAL(wp), PUBLIC ::  xpref2z      !: mesozoo preference for microzooplankton 
     30   REAL(wp), PUBLIC ::  xpref2c      !: mesozoo preference for POC  
    3131   REAL(wp), PUBLIC ::  xthresh2zoo  !: zoo feeding threshold for mesozooplankton  
    3232   REAL(wp), PUBLIC ::  xthresh2dia  !: diatoms feeding threshold for mesozooplankton  
     
    4040   REAL(wp), PUBLIC ::  unass2       !: Efficicency of mesozoo growth  
    4141   REAL(wp), PUBLIC ::  sigma2       !: Fraction of mesozoo excretion as DOM  
    42    REAL(wp), PUBLIC ::  epsher2      !: half sturation constant for grazing 2 
     42   REAL(wp), PUBLIC ::  epsher2      !: growth efficiency 
     43   REAL(wp), PUBLIC ::  epsher2min   !: minimum growth efficiency at high food for grazing 2 
    4344   REAL(wp), PUBLIC ::  grazflux     !: mesozoo flux feeding rate 
    4445 
     
    6364      REAL(wp) :: zcompadi, zcompaph, zcompapoc, zcompaz, zcompam 
    6465      REAL(wp) :: zgraze2 , zdenom, zdenom2 
    65       REAL(wp) :: zfact   , zfood, zfoodlim, zproport 
     66      REAL(wp) :: zfact   , zfood, zfoodlim, zproport, zbeta 
    6667      REAL(wp) :: zmortzgoc, zfrac, zfracfe, zratio, zratio2, zfracal, zgrazcal 
    67       REAL(wp) :: zepshert, zepsherv, zgrarsig, zgraztot, zgraztotn, zgraztotf 
    68       REAL(wp) :: zgrarem2, zgrafer2, zgrapoc2, zprcaca, zmortz2, zgrasrat, zgrasratn 
    69       REAL(wp) :: zrespz2, ztortz2, zgrazd, zgrazz, zgrazpof 
     68      REAL(wp) :: zepsherf, zepshert, zepsherv, zgrarsig, zgraztotc, zgraztotn, zgraztotf 
     69      REAL(wp) :: zgrarem2, zgrafer2, zgrapoc2, zprcaca, zmortz, zgrasrat, zgrasratn 
     70      REAL(wp) :: zrespz, ztortz, zgrazd, zgrazz, zgrazpof 
    7071      REAL(wp) :: zgrazn, zgrazpoc, zgraznf, zgrazf 
    7172      REAL(wp) :: zgrazfffp, zgrazfffg, zgrazffep, zgrazffeg 
    7273      CHARACTER (len=25) :: charout 
    73       REAL(wp), DIMENSION(jpi,jpj,jpk) :: zgrazing 
    74       REAL(wp), ALLOCATABLE, DIMENSION(:,:,:) ::   zw3d 
     74      REAL(wp), DIMENSION(jpi,jpj,jpk) :: zgrazing, zfezoo2 
     75      REAL(wp), ALLOCATABLE, DIMENSION(:,:,:) ::   zw3d, zz2ligprod 
    7576      !!--------------------------------------------------------------------- 
    7677      ! 
     
    7879      ! 
    7980      zgrazing(:,:,:) = 0._wp 
    80  
     81      zfezoo2 (:,:,:) = 0._wp 
     82      ! 
     83      IF (ln_ligand) THEN 
     84         ALLOCATE( zz2ligprod(jpi,jpj,jpk) ) 
     85         zz2ligprod(:,:,:) = 0._wp 
     86      ENDIF 
     87      ! 
    8188      DO jk = 1, jpkm1 
    8289         DO jj = 1, jpj 
     
    8794               !  Respiration rates of both zooplankton 
    8895               !  ------------------------------------- 
    89                zrespz2   = resrat2 * zfact * trb(ji,jj,jk,jpmes) / ( xkmort + trb(ji,jj,jk,jpmes) )  & 
    90                   &      + resrat2 * zfact * 3. * nitrfac(ji,jj,jk) 
     96               zrespz    = resrat2 * zfact * ( trb(ji,jj,jk,jpmes) / ( xkmort + trb(ji,jj,jk,jpmes) )  & 
     97               &           + 3. * nitrfac(ji,jj,jk) ) 
    9198 
    9299               !  Zooplankton mortality. A square function has been selected with 
    93100               !  no real reason except that it seems to be more stable and may mimic predation 
    94101               !  --------------------------------------------------------------- 
    95                ztortz2   = mzrat2 * 1.e6 * zfact * trb(ji,jj,jk,jpmes)  * (1. - nitrfac(ji,jj,jk) ) 
     102               ztortz    = mzrat2 * 1.e6 * zfact * trb(ji,jj,jk,jpmes)  * (1. - nitrfac(ji,jj,jk) ) 
    96103               ! 
    97104               zcompadi  = MAX( ( trb(ji,jj,jk,jpdia) - xthresh2dia ), 0.e0 ) 
    98105               zcompaz   = MAX( ( trb(ji,jj,jk,jpzoo) - xthresh2zoo ), 0.e0 ) 
     106               zcompapoc = MAX( ( trb(ji,jj,jk,jppoc) - xthresh2poc ), 0.e0 ) 
    99107               ! Size effect of nanophytoplankton on grazing : the smaller it is, the less prone 
    100108               ! it is to predation by mesozooplankton 
     
    102110               zcompaph  = MAX( ( trb(ji,jj,jk,jpphy) - xthresh2phy ), 0.e0 ) & 
    103111                  &      * MIN(1., MAX( 0., ( quotan(ji,jj,jk) - 0.2) / 0.3 ) ) 
    104                zcompapoc = MAX( ( trb(ji,jj,jk,jppoc) - xthresh2poc ), 0.e0 ) 
    105  
    106                zfood     = xprefc * zcompadi + xprefz * zcompaz + xprefp * zcompaph + xprefpoc * zcompapoc  
     112 
     113               !   Mesozooplankton grazing 
     114               !   ------------------------ 
     115               zfood     = xpref2d * zcompadi + xpref2z * zcompaz + xpref2n * zcompaph + xpref2c * zcompapoc  
    107116               zfoodlim  = MAX( 0., zfood - MIN( 0.5 * zfood, xthresh2 ) ) 
    108117               zdenom    = zfoodlim / ( xkgraz2 + zfoodlim ) 
     
    110119               zgraze2   = grazrat2 * xstep * tgfunc2(ji,jj,jk) * trb(ji,jj,jk,jpmes) * (1. - nitrfac(ji,jj,jk))  
    111120 
    112                zgrazd    = zgraze2  * xprefc   * zcompadi  * zdenom2  
    113                zgrazz    = zgraze2  * xprefz   * zcompaz   * zdenom2  
    114                zgrazn    = zgraze2  * xprefp   * zcompaph  * zdenom2  
    115                zgrazpoc  = zgraze2  * xprefpoc * zcompapoc * zdenom2  
     121               zgrazd    = zgraze2  * xpref2d  * zcompadi  * zdenom2  
     122               zgrazz    = zgraze2  * xpref2z  * zcompaz   * zdenom2  
     123               zgrazn    = zgraze2  * xpref2n  * zcompaph  * zdenom2  
     124               zgrazpoc  = zgraze2  * xpref2c * zcompapoc * zdenom2  
    116125 
    117126               zgraznf   = zgrazn   * trb(ji,jj,jk,jpnfe) / ( trb(ji,jj,jk,jpphy) + rtrn) 
     
    129138               &           * (1. - nitrfac(ji,jj,jk)) 
    130139               zgrazfffp = zgrazffep * trb(ji,jj,jk,jpsfe) / (trb(ji,jj,jk,jppoc) + rtrn) 
    131               ! 
    132               zgraztot = zgrazd + zgrazz + zgrazn + zgrazpoc + zgrazffep + zgrazffeg 
    133               ! Compute the proportion of filter feeders 
    134               zproport  = (zgrazffep + zgrazffeg)/(rtrn + zgraztot) 
    135               ! Compute fractionation of aggregates. It is assumed that  
    136               ! diatoms based aggregates are more prone to fractionation 
    137               ! since they are more porous (marine snow instead of fecal pellets) 
    138               zratio    = trb(ji,jj,jk,jpgsi) / ( trb(ji,jj,jk,jpgoc) + rtrn ) 
    139               zratio2   = zratio * zratio 
    140               zfrac     = zproport * grazflux  * xstep * wsbio4(ji,jj,jk)      & 
     140               ! 
     141               zgraztotc = zgrazd + zgrazz + zgrazn + zgrazpoc + zgrazffep + zgrazffeg 
     142               ! Compute the proportion of filter feeders 
     143               zproport  = (zgrazffep + zgrazffeg)/(rtrn + zgraztotc) 
     144               ! Compute fractionation of aggregates. It is assumed that  
     145               ! diatoms based aggregates are more prone to fractionation 
     146               ! since they are more porous (marine snow instead of fecal pellets) 
     147               zratio    = trb(ji,jj,jk,jpgsi) / ( trb(ji,jj,jk,jpgoc) + rtrn ) 
     148               zratio2   = zratio * zratio 
     149               zfrac     = zproport * grazflux  * xstep * wsbio4(ji,jj,jk)      & 
    141150               &          * trb(ji,jj,jk,jpgoc) * trb(ji,jj,jk,jpmes)          & 
    142151               &          * ( 0.2 + 3.8 * zratio2 / ( 1.**2 + zratio2 ) ) 
    143               zfracfe   = zfrac * trb(ji,jj,jk,jpbfe) / (trb(ji,jj,jk,jpgoc) + rtrn) 
    144  
    145               zgrazffep = zproport * zgrazffep 
    146               zgrazffeg = zproport * zgrazffeg 
    147               zgrazfffp = zproport * zgrazfffp 
    148               zgrazfffg = zproport * zgrazfffg 
    149               zgraztot = zgrazd + zgrazz + zgrazn + zgrazpoc + zgrazffep + zgrazffeg 
    150               zgraztotn = zgrazd * quotad(ji,jj,jk) + zgrazz + zgrazn * quotan(ji,jj,jk)   & 
    151               &   + zgrazpoc + zgrazffep + zgrazffeg 
    152               zgraztotf = zgrazf + zgraznf + zgrazz * ferat3 + zgrazpof + zgrazfffp + zgrazfffg 
    153  
    154               ! Total grazing ( grazing by microzoo is already computed in p4zmicro ) 
    155               zgrazing(ji,jj,jk) = zgraztot 
    156  
    157               !    Mesozooplankton efficiency 
    158               !    -------------------------- 
    159                zgrasrat  =  ( zgraztotf +rtrn )/ ( zgraztot + rtrn ) 
    160                zgrasratn =  ( zgraztotn +rtrn )/ ( zgraztot + rtrn ) 
     152               zfracfe   = zfrac * trb(ji,jj,jk,jpbfe) / (trb(ji,jj,jk,jpgoc) + rtrn) 
     153 
     154               zgrazffep = zproport * zgrazffep 
     155               zgrazffeg = zproport * zgrazffeg 
     156               zgrazfffp = zproport * zgrazfffp 
     157               zgrazfffg = zproport * zgrazfffg 
     158               zgraztotc = zgrazd + zgrazz + zgrazn + zgrazpoc + zgrazffep + zgrazffeg 
     159               zgraztotn = zgrazd * quotad(ji,jj,jk) + zgrazz + zgrazn * quotan(ji,jj,jk)   & 
     160               &   + zgrazpoc + zgrazffep + zgrazffeg 
     161               zgraztotf = zgrazf + zgraznf + zgrazz * ferat3 + zgrazpof + zgrazfffp + zgrazfffg 
     162 
     163               ! Total grazing ( grazing by microzoo is already computed in p4zmicro ) 
     164               zgrazing(ji,jj,jk) = zgraztotc 
     165 
     166               !    Mesozooplankton efficiency 
     167               !    -------------------------- 
     168               zgrasrat  =  ( zgraztotf + rtrn )/ ( zgraztotc + rtrn ) 
     169               zgrasratn =  ( zgraztotn + rtrn )/ ( zgraztotc + rtrn ) 
    161170               zepshert  = MIN( 1., zgrasratn, zgrasrat / ferat3) 
    162                zepsherv  = zepshert * MIN( epsher2, (1. - unass2) * zgrasrat / ferat3, (1. - unass2) * zgrasratn ) 
    163                zgrarem2  = zgraztot * ( 1. - zepsherv - unass2 ) & 
    164                 &       + ( 1. - epsher2 - unass2 ) / ( 1. - epsher2 ) * ztortz2 
    165                zgrafer2  = zgraztot * MAX( 0. , ( 1. - unass2 ) * zgrasrat - ferat3 * zepsherv )    & 
    166                 &       + ferat3 * ( ( 1. - epsher2 - unass2 ) /( 1. - epsher2 ) * ztortz2 ) 
    167                zgrapoc2  = zgraztot * unass2 
     171               zbeta     = MAX(0., (epsher2 - epsher2min) ) 
     172               zepsherf  = epsher2min + zbeta / ( 1.0 + 0.04E6 * 12. * zfood * zbeta )  
     173               zepsherv  = zepsherf * zepshert  
     174 
     175               zgrarem2  = zgraztotc * ( 1. - zepsherv - unass2 ) & 
     176               &         + ( 1. - epsher2 - unass2 ) / ( 1. - epsher2 ) * ztortz 
     177               zgrafer2  = zgraztotc * MAX( 0. , ( 1. - unass2 ) * zgrasrat - ferat3 * zepsherv )    & 
     178               &         + ferat3 * ( ( 1. - epsher2 - unass2 ) /( 1. - epsher2 ) * ztortz ) 
     179               zgrapoc2  = zgraztotc * unass2 
    168180 
    169181               !   Update the arrays TRA which contain the biological sources and sinks 
     
    173185               tra(ji,jj,jk,jpdoc) = tra(ji,jj,jk,jpdoc) + zgrarem2 - zgrarsig 
    174186               ! 
    175                IF( ln_ligand ) tra(ji,jj,jk,jplgw) = tra(ji,jj,jk,jplgw) + (zgrarem2 - zgrarsig) * ldocz 
     187               IF( ln_ligand ) THEN  
     188                  tra(ji,jj,jk,jplgw) = tra(ji,jj,jk,jplgw) + (zgrarem2 - zgrarsig) * ldocz 
     189                  zz2ligprod(ji,jj,jk) = (zgrarem2 - zgrarsig) * ldocz 
     190               ENDIF 
    176191               ! 
    177192               tra(ji,jj,jk,jpoxy) = tra(ji,jj,jk,jpoxy) - o2ut * zgrarsig 
    178193               tra(ji,jj,jk,jpfer) = tra(ji,jj,jk,jpfer) + zgrafer2 
     194               zfezoo2(ji,jj,jk)   = zgrafer2 
    179195               tra(ji,jj,jk,jpdic) = tra(ji,jj,jk,jpdic) + zgrarsig 
    180196               tra(ji,jj,jk,jptal) = tra(ji,jj,jk,jptal) + rno3 * zgrarsig               
    181197 
    182                zmortz2 = ztortz2 + zrespz2 
    183                zmortzgoc = unass2 / ( 1. - epsher2 ) * ztortz2 + zrespz2 
    184                tra(ji,jj,jk,jpmes) = tra(ji,jj,jk,jpmes) - zmortz2 + zepsherv * zgraztot  
     198               zmortz = ztortz + zrespz 
     199               zmortzgoc = unass2 / ( 1. - epsher2 ) * ztortz + zrespz 
     200               tra(ji,jj,jk,jpmes) = tra(ji,jj,jk,jpmes) - zmortz + zepsherv * zgraztotc  
    185201               tra(ji,jj,jk,jpdia) = tra(ji,jj,jk,jpdia) - zgrazd 
    186202               tra(ji,jj,jk,jpzoo) = tra(ji,jj,jk,jpzoo) - zgrazz 
     
    193209               tra(ji,jj,jk,jpdfe) = tra(ji,jj,jk,jpdfe) - zgrazf 
    194210 
    195               tra(ji,jj,jk,jppoc) = tra(ji,jj,jk,jppoc) - zgrazpoc - zgrazffep + zfrac 
    196               prodpoc(ji,jj,jk) = prodpoc(ji,jj,jk) + zfrac 
    197               conspoc(ji,jj,jk) = conspoc(ji,jj,jk) - zgrazpoc - zgrazffep 
    198               tra(ji,jj,jk,jpgoc) = tra(ji,jj,jk,jpgoc) + zmortzgoc - zgrazffeg + zgrapoc2 - zfrac 
    199               prodgoc(ji,jj,jk) = prodgoc(ji,jj,jk) + zmortzgoc + zgrapoc2 
    200               consgoc(ji,jj,jk) = consgoc(ji,jj,jk) - zgrazffeg - zfrac 
    201               tra(ji,jj,jk,jpsfe) = tra(ji,jj,jk,jpsfe) - zgrazpof - zgrazfffp + zfracfe 
    202               tra(ji,jj,jk,jpbfe) = tra(ji,jj,jk,jpbfe) + ferat3 * zmortzgoc - zgrazfffg     & 
     211               tra(ji,jj,jk,jppoc) = tra(ji,jj,jk,jppoc) - zgrazpoc - zgrazffep + zfrac 
     212               prodpoc(ji,jj,jk) = prodpoc(ji,jj,jk) + zfrac 
     213               conspoc(ji,jj,jk) = conspoc(ji,jj,jk) - zgrazpoc - zgrazffep 
     214               tra(ji,jj,jk,jpgoc) = tra(ji,jj,jk,jpgoc) + zmortzgoc - zgrazffeg + zgrapoc2 - zfrac 
     215               prodgoc(ji,jj,jk) = prodgoc(ji,jj,jk) + zmortzgoc + zgrapoc2 
     216               consgoc(ji,jj,jk) = consgoc(ji,jj,jk) - zgrazffeg - zfrac 
     217               tra(ji,jj,jk,jpsfe) = tra(ji,jj,jk,jpsfe) - zgrazpof - zgrazfffp + zfracfe 
     218               tra(ji,jj,jk,jpbfe) = tra(ji,jj,jk,jpbfe) + ferat3 * zmortzgoc - zgrazfffg     & 
    203219                 &                + zgraztotf * unass2 - zfracfe 
    204               zfracal = trb(ji,jj,jk,jpcal) / (trb(ji,jj,jk,jppoc) + trb(ji,jj,jk,jpgoc) + rtrn ) 
    205               zgrazcal = (zgrazffeg + zgrazpoc) * (1. - part2) * zfracal 
    206               ! calcite production 
    207               zprcaca = xfracal(ji,jj,jk) * zgrazn 
    208               prodcal(ji,jj,jk) = prodcal(ji,jj,jk) + zprcaca  ! prodcal=prodcal(nanophy)+prodcal(microzoo)+prodcal(mesozoo) 
    209               ! 
    210               zprcaca = part2 * zprcaca 
    211               tra(ji,jj,jk,jpdic) = tra(ji,jj,jk,jpdic) + zgrazcal - zprcaca 
    212               tra(ji,jj,jk,jptal) = tra(ji,jj,jk,jptal) - 2. * ( zgrazcal + zprcaca ) 
    213               tra(ji,jj,jk,jpcal) = tra(ji,jj,jk,jpcal) - zgrazcal + zprcaca 
     220               zfracal = trb(ji,jj,jk,jpcal) / (trb(ji,jj,jk,jppoc) + trb(ji,jj,jk,jpgoc) + rtrn ) 
     221               zgrazcal = (zgrazffeg + zgrazpoc) * (1. - part2) * zfracal 
     222               ! calcite production 
     223               zprcaca = xfracal(ji,jj,jk) * zgrazn 
     224               prodcal(ji,jj,jk) = prodcal(ji,jj,jk) + zprcaca  ! prodcal=prodcal(nanophy)+prodcal(microzoo)+prodcal(mesozoo) 
     225               ! 
     226               zprcaca = part2 * zprcaca 
     227               tra(ji,jj,jk,jpdic) = tra(ji,jj,jk,jpdic) + zgrazcal - zprcaca 
     228               tra(ji,jj,jk,jptal) = tra(ji,jj,jk,jptal) - 2. * ( zgrazcal + zprcaca ) 
     229               tra(ji,jj,jk,jpcal) = tra(ji,jj,jk,jpcal) - zgrazcal + zprcaca 
    214230            END DO 
    215231         END DO 
     
    226242            CALL iom_put( "PCAL", zw3d )   
    227243         ENDIF 
     244         IF( iom_use( "FEZOO2" ) ) THEN 
     245            zw3d(:,:,:) = zfezoo2(:,:,:) * 1e9 * 1.e+3 * rfact2r * tmask(:,:,:)   ! 
     246            CALL iom_put( "FEZOO2", zw3d ) 
     247         ENDIF 
     248         IF( iom_use( "LPRODZ2" ) .AND. ln_ligand )  THEN 
     249            zw3d(:,:,:) = zz2ligprod(:,:,:) * 1e9 * 1.e+3 * rfact2r * tmask(:,:,:) 
     250            CALL iom_put( "LPRODZ2"  , zw3d ) 
     251         ENDIF 
    228252         DEALLOCATE( zw3d ) 
    229253      ENDIF 
     
    253277      INTEGER ::   ios   ! Local integer 
    254278      ! 
    255       NAMELIST/namp4zmes/ part2, grazrat2, resrat2, mzrat2, xprefc, xprefp, xprefz,   & 
    256          &                xprefpoc, xthresh2dia, xthresh2phy, xthresh2zoo, xthresh2poc, & 
    257          &                xthresh2, xkgraz2, epsher2, sigma2, unass2, grazflux 
     279      NAMELIST/namp4zmes/ part2, grazrat2, resrat2, mzrat2, xpref2n, xpref2d, xpref2z,   & 
     280         &                xpref2c, xthresh2dia, xthresh2phy, xthresh2zoo, xthresh2poc, & 
     281         &                xthresh2, xkgraz2, epsher2, epsher2min, sigma2, unass2, grazflux 
    258282      !!---------------------------------------------------------------------- 
    259283      ! 
     
    275299         WRITE(numout,*) '   Namelist : namp4zmes' 
    276300         WRITE(numout,*) '      part of calcite not dissolved in mesozoo guts  part2        =', part2 
    277          WRITE(numout,*) '      mesozoo preference for phyto                   xprefc       =', xprefc 
    278          WRITE(numout,*) '      mesozoo preference for POC                     xprefp       =', xprefp 
    279          WRITE(numout,*) '      mesozoo preference for zoo                     xprefz       =', xprefz 
    280          WRITE(numout,*) '      mesozoo preference for poc                     xprefpoc     =', xprefpoc 
     301         WRITE(numout,*) '      mesozoo preference for phyto                   xpref2n      =', xpref2n 
     302         WRITE(numout,*) '      mesozoo preference for diatoms                 xpref2d      =', xpref2d 
     303         WRITE(numout,*) '      mesozoo preference for zoo                     xpref2z      =', xpref2z 
     304         WRITE(numout,*) '      mesozoo preference for poc                     xpref2c      =', xpref2c 
    281305         WRITE(numout,*) '      microzoo feeding threshold  for mesozoo        xthresh2zoo  =', xthresh2zoo 
    282306         WRITE(numout,*) '      diatoms feeding threshold  for mesozoo         xthresh2dia  =', xthresh2dia 
     
    289313         WRITE(numout,*) '      mesozoo flux feeding rate                      grazflux     =', grazflux 
    290314         WRITE(numout,*) '      non assimilated fraction of P by mesozoo       unass2       =', unass2 
    291          WRITE(numout,*) '      Efficicency of Mesozoo growth                  epsher2      =', epsher2 
     315         WRITE(numout,*) '      Efficiency of Mesozoo growth                   epsher2      =', epsher2 
     316         WRITE(numout,*) '      Minimum Efficiency of Mesozoo growth           epsher2min  =', epsher2min 
    292317         WRITE(numout,*) '      Fraction of mesozoo excretion as DOM           sigma2       =', sigma2 
    293318         WRITE(numout,*) '      half sturation constant for grazing 2          xkgraz2      =', xkgraz2 
  • NEMO/branches/2018/dev_r10164_HPC09_ESIWACE_PREP_MERGE/src/TOP/PISCES/P4Z/p4zmicro.F90

    r10345 r10368  
    2626 
    2727   REAL(wp), PUBLIC ::   part        !: part of calcite not dissolved in microzoo guts 
    28    REAL(wp), PUBLIC ::   xpref2c     !: microzoo preference for POC  
    29    REAL(wp), PUBLIC ::   xpref2p     !: microzoo preference for nanophyto 
    30    REAL(wp), PUBLIC ::   xpref2d     !: microzoo preference for diatoms 
     28   REAL(wp), PUBLIC ::   xprefc      !: microzoo preference for POC  
     29   REAL(wp), PUBLIC ::   xprefn      !: microzoo preference for nanophyto 
     30   REAL(wp), PUBLIC ::   xprefd      !: microzoo preference for diatoms 
    3131   REAL(wp), PUBLIC ::   xthreshdia  !: diatoms feeding threshold for microzooplankton  
    3232   REAL(wp), PUBLIC ::   xthreshphy  !: nanophyto threshold for microzooplankton  
     
    3636   REAL(wp), PUBLIC ::   mzrat       !: microzooplankton mortality rate  
    3737   REAL(wp), PUBLIC ::   grazrat     !: maximal microzoo grazing rate 
    38    REAL(wp), PUBLIC ::   xkgraz      !: non assimilated fraction of P by microzoo  
    39    REAL(wp), PUBLIC ::   unass       !: Efficicency of microzoo growth  
     38   REAL(wp), PUBLIC ::   xkgraz      !: Half-saturation constant of assimilation 
     39   REAL(wp), PUBLIC ::   unass       !: Non-assimilated part of food 
    4040   REAL(wp), PUBLIC ::   sigma1      !: Fraction of microzoo excretion as DOM  
    41    REAL(wp), PUBLIC ::   epsher      !: half sturation constant for grazing 1  
     41   REAL(wp), PUBLIC ::   epsher      !: growth efficiency for grazing 1  
     42   REAL(wp), PUBLIC ::   epshermin   !: minimum growth efficiency for grazing 1 
    4243 
    4344   !!---------------------------------------------------------------------- 
     
    6263      REAL(wp) :: zcompadi, zcompaz , zcompaph, zcompapoc 
    6364      REAL(wp) :: zgraze  , zdenom, zdenom2 
    64       REAL(wp) :: zfact   , zfood, zfoodlim 
    65       REAL(wp) :: zepshert, zepsherv, zgrarsig, zgraztot, zgraztotn, zgraztotf 
     65      REAL(wp) :: zfact   , zfood, zfoodlim, zbeta 
     66      REAL(wp) :: zepsherf, zepshert, zepsherv, zgrarsig, zgraztotc, zgraztotn, zgraztotf 
    6667      REAL(wp) :: zgrarem, zgrafer, zgrapoc, zprcaca, zmortz 
    6768      REAL(wp) :: zrespz, ztortz, zgrasrat, zgrasratn 
    6869      REAL(wp) :: zgrazp, zgrazm, zgrazsd 
    6970      REAL(wp) :: zgrazmf, zgrazsf, zgrazpf 
    70       REAL(wp), DIMENSION(jpi,jpj,jpk) :: zgrazing  
    71       REAL(wp), DIMENSION(:,:,:), ALLOCATABLE :: zw3d 
     71      REAL(wp), DIMENSION(jpi,jpj,jpk) :: zgrazing, zfezoo 
     72      REAL(wp), DIMENSION(:,:,:), ALLOCATABLE :: zw3d, zzligprod 
    7273      CHARACTER (len=25) :: charout 
    7374      !!--------------------------------------------------------------------- 
    7475      ! 
    7576      IF( ln_timing )   CALL timing_start('p4z_micro') 
     77      ! 
     78      IF (ln_ligand) THEN 
     79         ALLOCATE( zzligprod(jpi,jpj,jpk) ) 
     80         zzligprod(:,:,:) = 0._wp 
     81      ENDIF 
    7682      ! 
    7783      DO jk = 1, jpkm1 
     
    97103               !     Microzooplankton grazing 
    98104               !     ------------------------ 
    99                zfood     = xpref2p * zcompaph + xpref2c * zcompapoc + xpref2d * zcompadi 
     105               zfood     = xprefn * zcompaph + xprefc * zcompapoc + xprefd * zcompadi 
    100106               zfoodlim  = MAX( 0. , zfood - min(xthresh,0.5*zfood) ) 
    101107               zdenom    = zfoodlim / ( xkgraz + zfoodlim ) 
     
    103109               zgraze    = grazrat * xstep * tgfunc2(ji,jj,jk) * trb(ji,jj,jk,jpzoo) * (1. - nitrfac(ji,jj,jk)) 
    104110 
    105                zgrazp    = zgraze  * xpref2p * zcompaph  * zdenom2  
    106                zgrazm    = zgraze  * xpref2c * zcompapoc * zdenom2  
    107                zgrazsd   = zgraze  * xpref2d * zcompadi  * zdenom2  
     111               zgrazp    = zgraze  * xprefn * zcompaph  * zdenom2  
     112               zgrazm    = zgraze  * xprefc * zcompapoc * zdenom2  
     113               zgrazsd   = zgraze  * xprefd * zcompadi  * zdenom2  
    108114 
    109115               zgrazpf   = zgrazp  * trb(ji,jj,jk,jpnfe) / (trb(ji,jj,jk,jpphy) + rtrn) 
     
    111117               zgrazsf   = zgrazsd * trb(ji,jj,jk,jpdfe) / (trb(ji,jj,jk,jpdia) + rtrn) 
    112118               ! 
    113                zgraztot  = zgrazp  + zgrazm  + zgrazsd  
     119               zgraztotc = zgrazp  + zgrazm  + zgrazsd  
    114120               zgraztotf = zgrazpf + zgrazsf + zgrazmf  
    115121               zgraztotn = zgrazp * quotan(ji,jj,jk) + zgrazm + zgrazsd * quotad(ji,jj,jk) 
    116122 
    117123               ! Grazing by microzooplankton 
    118                zgrazing(ji,jj,jk) = zgraztot 
     124               zgrazing(ji,jj,jk) = zgraztotc 
    119125 
    120126               !    Various remineralization and excretion terms 
    121127               !    -------------------------------------------- 
    122                zgrasrat  = ( zgraztotf + rtrn ) / ( zgraztot + rtrn ) 
    123                zgrasratn = ( zgraztotn + rtrn ) / ( zgraztot + rtrn ) 
     128               zgrasrat  = ( zgraztotf + rtrn ) / ( zgraztotc + rtrn ) 
     129               zgrasratn = ( zgraztotn + rtrn ) / ( zgraztotc + rtrn ) 
    124130               zepshert  =  MIN( 1., zgrasratn, zgrasrat / ferat3) 
    125                zepsherv  = zepshert * MIN( epsher, (1. - unass) * zgrasrat / ferat3, (1. - unass) * zgrasratn ) 
    126                zgrafer   = zgraztot * MAX( 0. , ( 1. - unass ) * zgrasrat - ferat3 * zepsherv )  
    127                zgrarem   = zgraztot * ( 1. - zepsherv - unass ) 
    128                zgrapoc   = zgraztot * unass 
     131               zbeta     = MAX(0., (epsher - epshermin) ) 
     132               zepsherf  = epshermin + zbeta / ( 1.0 + 0.04E6 * 12. * zfood * zbeta ) 
     133               zepsherv  = zepsherf * zepshert  
     134 
     135               zgrafer   = zgraztotc * MAX( 0. , ( 1. - unass ) * zgrasrat - ferat3 * zepsherv )  
     136               zgrarem   = zgraztotc * ( 1. - zepsherv - unass ) 
     137               zgrapoc   = zgraztotc * unass 
    129138 
    130139               !  Update of the TRA arrays 
     
    135144               tra(ji,jj,jk,jpdoc) = tra(ji,jj,jk,jpdoc) + zgrarem - zgrarsig 
    136145               ! 
    137                IF( ln_ligand ) tra(ji,jj,jk,jplgw) = tra(ji,jj,jk,jplgw) + (zgrarem - zgrarsig) * ldocz 
     146               IF( ln_ligand ) THEN 
     147                  tra(ji,jj,jk,jplgw) = tra(ji,jj,jk,jplgw) + (zgrarem - zgrarsig) * ldocz 
     148                  zzligprod(ji,jj,jk) = (zgrarem - zgrarsig) * ldocz 
     149               ENDIF 
    138150               ! 
    139151               tra(ji,jj,jk,jpoxy) = tra(ji,jj,jk,jpoxy) - o2ut * zgrarsig 
    140152               tra(ji,jj,jk,jpfer) = tra(ji,jj,jk,jpfer) + zgrafer 
     153               zfezoo(ji,jj,jk)    = zgrafer 
    141154               tra(ji,jj,jk,jppoc) = tra(ji,jj,jk,jppoc) + zgrapoc 
    142                prodpoc(ji,jj,jk) = prodpoc(ji,jj,jk) + zgrapoc 
     155               prodpoc(ji,jj,jk)   = prodpoc(ji,jj,jk) + zgrapoc 
    143156               tra(ji,jj,jk,jpsfe) = tra(ji,jj,jk,jpsfe) + zgraztotf * unass 
    144157               tra(ji,jj,jk,jpdic) = tra(ji,jj,jk,jpdic) + zgrarsig 
     
    147160               !   -------------------------------------------------------------------- 
    148161               zmortz = ztortz + zrespz 
    149                tra(ji,jj,jk,jpzoo) = tra(ji,jj,jk,jpzoo) - zmortz + zepsherv * zgraztot  
     162               tra(ji,jj,jk,jpzoo) = tra(ji,jj,jk,jpzoo) - zmortz + zepsherv * zgraztotc  
    150163               tra(ji,jj,jk,jpphy) = tra(ji,jj,jk,jpphy) - zgrazp 
    151164               tra(ji,jj,jk,jpdia) = tra(ji,jj,jk,jpdia) - zgrazsd 
     
    180193              CALL iom_put( "GRAZ1", zw3d ) 
    181194           ENDIF 
     195           IF( iom_use( "FEZOO" ) ) THEN 
     196              zw3d(:,:,:) = zfezoo(:,:,:) * 1e9 * 1.e+3 * rfact2r * tmask(:,:,:)   ! 
     197              CALL iom_put( "FEZOO", zw3d ) 
     198           ENDIF 
     199           IF( iom_use( "LPRODZ" ) .AND. ln_ligand )  THEN 
     200              zw3d(:,:,:) = zzligprod(:,:,:) * 1e9 * 1.e+3 * rfact2r * tmask(:,:,:) 
     201              CALL iom_put( "LPRODZ"  , zw3d ) 
     202           ENDIF 
    182203           DEALLOCATE( zw3d ) 
    183204         ENDIF 
     
    209230      INTEGER ::   ios   ! Local integer 
    210231      ! 
    211       NAMELIST/namp4zzoo/ part, grazrat, resrat, mzrat, xpref2c, xpref2p, & 
    212          &                xpref2d,  xthreshdia,  xthreshphy,  xthreshpoc, & 
    213          &                xthresh, xkgraz, epsher, sigma1, unass 
     232      NAMELIST/namp4zzoo/ part, grazrat, resrat, mzrat, xprefn, xprefc, & 
     233         &                xprefd,  xthreshdia,  xthreshphy,  xthreshpoc, & 
     234         &                xthresh, xkgraz, epsher, epshermin, sigma1, unass 
    214235      !!---------------------------------------------------------------------- 
    215236      ! 
     
    231252         WRITE(numout,*) '   Namelist : namp4zzoo' 
    232253         WRITE(numout,*) '      part of calcite not dissolved in microzoo guts  part        =', part 
    233          WRITE(numout,*) '      microzoo preference for POC                     xpref2c     =', xpref2c 
    234          WRITE(numout,*) '      microzoo preference for nano                    xpref2p     =', xpref2p 
    235          WRITE(numout,*) '      microzoo preference for diatoms                 xpref2d     =', xpref2d 
     254         WRITE(numout,*) '      microzoo preference for POC                     xprefc      =', xprefc 
     255         WRITE(numout,*) '      microzoo preference for nano                    xprefn      =', xprefn 
     256         WRITE(numout,*) '      microzoo preference for diatoms                 xprefd      =', xprefd 
    236257         WRITE(numout,*) '      diatoms feeding threshold  for microzoo         xthreshdia  =', xthreshdia 
    237258         WRITE(numout,*) '      nanophyto feeding threshold for microzoo        xthreshphy  =', xthreshphy 
     
    243264         WRITE(numout,*) '      non assimilated fraction of P by microzoo       unass       =', unass 
    244265         WRITE(numout,*) '      Efficicency of microzoo growth                  epsher      =', epsher 
     266         WRITE(numout,*) '      Minimum efficicency of microzoo growth          epshermin   =', epshermin 
    245267         WRITE(numout,*) '      Fraction of microzoo excretion as DOM           sigma1      =', sigma1 
    246268         WRITE(numout,*) '      half sturation constant for grazing 1           xkgraz      =', xkgraz 
  • NEMO/branches/2018/dev_r10164_HPC09_ESIWACE_PREP_MERGE/src/TOP/PISCES/P4Z/p4zopt.F90

    r10069 r10368  
    109109         DO jk = 1, nksrp       
    110110            etot_ndcy(:,:,jk) =        ze1(:,:,jk) +        ze2(:,:,jk) +       ze3(:,:,jk) 
    111             enano    (:,:,jk) =  2.1 * ze1(:,:,jk) + 0.42 * ze2(:,:,jk) + 0.4 * ze3(:,:,jk) 
    112             ediat    (:,:,jk) =  1.6 * ze1(:,:,jk) + 0.69 * ze2(:,:,jk) + 0.7 * ze3(:,:,jk) 
     111            enano    (:,:,jk) =  1.85 * ze1(:,:,jk) + 0.69 * ze2(:,:,jk) + 0.46 * ze3(:,:,jk) 
     112            ediat    (:,:,jk) =  1.62 * ze1(:,:,jk) + 0.74 * ze2(:,:,jk) + 0.63 * ze3(:,:,jk) 
    113113         END DO 
    114114         IF( ln_p5z ) THEN 
    115115            DO jk = 1, nksrp       
    116               epico  (:,:,jk) =  2.1 * ze1(:,:,jk) + 0.42 * ze2(:,:,jk) + 0.4 * ze3(:,:,jk) 
     116              epico  (:,:,jk) =  1.94 * ze1(:,:,jk) + 0.66 * ze2(:,:,jk) + 0.4 * ze3(:,:,jk) 
    117117            END DO 
    118118         ENDIF 
     
    134134         DO jk = 1, nksrp       
    135135            etot (:,:,jk) =        ze1(:,:,jk) +        ze2(:,:,jk) +       ze3(:,:,jk) 
    136             enano(:,:,jk) =  2.1 * ze1(:,:,jk) + 0.42 * ze2(:,:,jk) + 0.4 * ze3(:,:,jk) 
    137             ediat(:,:,jk) =  1.6 * ze1(:,:,jk) + 0.69 * ze2(:,:,jk) + 0.7 * ze3(:,:,jk) 
     136            enano(:,:,jk) =  1.85 * ze1(:,:,jk) + 0.69 * ze2(:,:,jk) + 0.46 * ze3(:,:,jk) 
     137            ediat(:,:,jk) =  1.62 * ze1(:,:,jk) + 0.74 * ze2(:,:,jk) + 0.63 * ze3(:,:,jk) 
    138138         END DO 
    139139         IF( ln_p5z ) THEN 
    140140            DO jk = 1, nksrp       
    141               epico(:,:,jk) =  2.1 * ze1(:,:,jk) + 0.42 * ze2(:,:,jk) + 0.4 * ze3(:,:,jk) 
     141              epico(:,:,jk) =  1.94 * ze1(:,:,jk) + 0.66 * ze2(:,:,jk) + 0.4 * ze3(:,:,jk) 
    142142            END DO 
    143143         ENDIF 
     
    182182      zetmp1 (:,:)   = 0.e0 
    183183      zetmp2 (:,:)   = 0.e0 
    184       zetmp3 (:,:)   = 0.e0 
    185       zetmp4 (:,:)   = 0.e0 
    186184 
    187185      DO jk = 1, nksrp 
     
    191189                  zetmp1 (ji,jj) = zetmp1 (ji,jj) + etot     (ji,jj,jk) * e3t_n(ji,jj,jk) ! remineralisation 
    192190                  zetmp2 (ji,jj) = zetmp2 (ji,jj) + etot_ndcy(ji,jj,jk) * e3t_n(ji,jj,jk) ! production 
    193                   zetmp3 (ji,jj) = zetmp3 (ji,jj) + enano    (ji,jj,jk) * e3t_n(ji,jj,jk) ! production 
    194                   zetmp4 (ji,jj) = zetmp4 (ji,jj) + ediat    (ji,jj,jk) * e3t_n(ji,jj,jk) ! production 
    195191                  zdepmoy(ji,jj) = zdepmoy(ji,jj) +                       e3t_n(ji,jj,jk) 
    196192               ENDIF 
     
    209205                  emoy (ji,jj,jk) = zetmp1(ji,jj) * z1_dep 
    210206                  zpar (ji,jj,jk) = zetmp2(ji,jj) * z1_dep 
    211                   enano(ji,jj,jk) = zetmp3(ji,jj) * z1_dep 
    212                   ediat(ji,jj,jk) = zetmp4(ji,jj) * z1_dep 
     207               ENDIF 
     208            END DO 
     209         END DO 
     210      END DO 
     211      ! 
     212      zdepmoy(:,:)   = 0.e0 
     213      zetmp3 (:,:)   = 0.e0 
     214      zetmp4 (:,:)   = 0.e0 
     215      ! 
     216      DO jk = 1, nksrp 
     217         DO jj = 1, jpj 
     218            DO ji = 1, jpi 
     219               IF( gdepw_n(ji,jj,jk+1) <= MIN(hmld(ji,jj), heup_01(ji,jj)) ) THEN 
     220                  zetmp3 (ji,jj) = zetmp3 (ji,jj) + enano    (ji,jj,jk) * e3t_n(ji,jj,jk) ! production 
     221                  zetmp4 (ji,jj) = zetmp4 (ji,jj) + ediat    (ji,jj,jk) * e3t_n(ji,jj,jk) ! production 
     222                  zdepmoy(ji,jj) = zdepmoy(ji,jj) +                       e3t_n(ji,jj,jk) 
     223               ENDIF 
     224            END DO 
     225         END DO 
     226      END DO 
     227      enanom(:,:,:) = enano(:,:,:) 
     228      ediatm(:,:,:) = ediat(:,:,:) 
     229      ! 
     230      DO jk = 1, nksrp 
     231         DO jj = 1, jpj 
     232            DO ji = 1, jpi 
     233               IF( gdepw_n(ji,jj,jk+1) <= hmld(ji,jj) ) THEN 
     234                  z1_dep = 1. / ( zdepmoy(ji,jj) + rtrn ) 
     235                  enanom(ji,jj,jk) = zetmp3(ji,jj) * z1_dep 
     236                  ediatm(ji,jj,jk) = zetmp4(ji,jj) * z1_dep 
    213237               ENDIF 
    214238            END DO 
     
    221245            DO jj = 1, jpj 
    222246               DO ji = 1, jpi 
    223                   IF( gdepw_n(ji,jj,jk+1) <= hmld(ji,jj) ) THEN  
     247                  IF( gdepw_n(ji,jj,jk+1) <= MIN(hmld(ji,jj), heup_01(ji,jj)) ) THEN 
     248                     zetmp5(ji,jj)  = zetmp5 (ji,jj) + epico(ji,jj,jk) * e3t_n(ji,jj,jk) ! production 
     249                  ENDIF 
     250               END DO 
     251            END DO 
     252         END DO 
     253         ! 
     254         epicom(:,:,:) = epico(:,:,:) 
     255         ! 
     256         DO jk = 1, nksrp 
     257            DO jj = 1, jpj 
     258               DO ji = 1, jpi 
     259                  IF( gdepw_n(ji,jj,jk+1) <= hmld(ji,jj) ) THEN 
    224260                     z1_dep = 1. / ( zdepmoy(ji,jj) + rtrn ) 
    225                      zetmp5(ji,jj)  = zetmp5 (ji,jj) + epico(ji,jj,jk) * e3t_n(ji,jj,jk) ! production 
    226                      epico(ji,jj,jk) = zetmp5(ji,jj) * z1_dep 
     261                     epicom(ji,jj,jk) = zetmp5(ji,jj) * z1_dep 
    227262                  ENDIF 
    228263               END DO 
  • NEMO/branches/2018/dev_r10164_HPC09_ESIWACE_PREP_MERGE/src/TOP/PISCES/P4Z/p4zpoc.F90

    r10069 r10368  
    6262      CHARACTER (len=25) :: charout 
    6363      REAL(wp), DIMENSION(jpi,jpj  )   :: totprod, totthick, totcons  
    64       REAL(wp), DIMENSION(jpi,jpj,jpk)   :: zremipoc, zremigoc, zorem3, ztremint 
     64      REAL(wp), DIMENSION(jpi,jpj,jpk)   :: zremipoc, zremigoc, zorem3, ztremint, zfolimi 
    6565      REAL(wp), DIMENSION(jpi,jpj,jpk,jcpoc) :: alphag 
    6666      !!--------------------------------------------------------------------- 
     
    9090      orem  (:,:,:)   = 0. 
    9191      ztremint(:,:,:) = 0. 
     92      zfolimi (:,:,:) = 0. 
    9293 
    9394      DO jn = 1, jcpoc 
     
    209210                  tra(ji,jj,jk,jpdoc) = tra(ji,jj,jk,jpdoc) + zorem2 
    210211                  tra(ji,jj,jk,jpfer) = tra(ji,jj,jk,jpfer) + zofer2 
     212                  zfolimi(ji,jj,jk)   = zofer2 
    211213               END DO 
    212214            END DO 
     
    239241                  tra(ji,jj,jk,jpgop) = tra(ji,jj,jk,jpgop) - zopop2 * (1. + solgoc) 
    240242                  tra(ji,jj,jk,jpbfe) = tra(ji,jj,jk,jpbfe) - zofer2 * (1. + solgoc) 
     243                  zfolimi(ji,jj,jk)   = zofer2 
    241244               END DO 
    242245            END DO 
     
    414417                    tra(ji,jj,jk,jppoc) = tra(ji,jj,jk,jppoc) - zorem 
    415418                    tra(ji,jj,jk,jpsfe) = tra(ji,jj,jk,jpsfe) - zofer 
     419                    zfolimi(ji,jj,jk)   = zfolimi(ji,jj,jk) + zofer 
    416420                  ENDIF 
    417421               END DO 
     
    439443                tra(ji,jj,jk,jpdop) = tra(ji,jj,jk,jpdop) + zopop  
    440444                tra(ji,jj,jk,jpfer) = tra(ji,jj,jk,jpfer) + zofer  
     445                zfolimi(ji,jj,jk)   = zfolimi(ji,jj,jk) + zofer 
    441446             END DO 
    442447           END DO 
     
    449454          CALL iom_put( "REMINP" , zremipoc(:,:,:) * tmask(:,:,:) )  ! Remineralisation rate 
    450455          CALL iom_put( "REMING" , zremigoc(:,:,:) * tmask(:,:,:) )  ! Remineralisation rate 
     456          CALL iom_put( "REMINF" , zfolimi(:,:,:)  * tmask(:,:,:)  * 1.e+9 * zrfact2 )  ! Remineralisation rate 
    451457        ENDIF 
    452458     ENDIF 
  • NEMO/branches/2018/dev_r10164_HPC09_ESIWACE_PREP_MERGE/src/TOP/PISCES/P4Z/p4zprod.F90

    r10345 r10368  
    4040   REAL(wp), PUBLIC ::   grosip       !: 
    4141 
    42    REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:) ::   prmax    !: optimal production = f(temperature) 
    4342   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:) ::   quotan   !: proxy of N quota in Nanophyto 
    4443   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:) ::   quotad   !: proxy of N quota in diatomee 
     
    7877      REAL(wp), ALLOCATABLE, DIMENSION(:,:,:) :: zw3d 
    7978      REAL(wp), DIMENSION(jpi,jpj    ) :: zstrn, zmixnano, zmixdiat 
     79      REAL(wp), DIMENSION(jpi,jpj,jpk) :: zprmaxn,zprmaxd 
    8080      REAL(wp), DIMENSION(jpi,jpj,jpk) :: zpislopeadn, zpislopeadd, zysopt   
    8181      REAL(wp), DIMENSION(jpi,jpj,jpk) :: zprdia, zprbio, zprdch, zprnch    
     
    8383      REAL(wp), DIMENSION(jpi,jpj,jpk) :: zpronewn, zpronewd 
    8484      REAL(wp), DIMENSION(jpi,jpj,jpk) :: zmxl_fac, zmxl_chl 
     85      REAL(wp), DIMENSION(jpi,jpj,jpk) :: zpligprod1, zpligprod2 
    8586      !!--------------------------------------------------------------------- 
    8687      ! 
     
    9697 
    9798      ! Computation of the optimal production 
    98       prmax(:,:,:) = 0.8_wp * r1_rday * tgfunc(:,:,:)  
     99      zprmaxn(:,:,:) = 0.8_wp * r1_rday * tgfunc(:,:,:) 
     100      zprmaxd(:,:,:) = zprmaxn(:,:,:) 
    99101 
    100102      ! compute the day length depending on latitude and the day 
     
    128130      END DO 
    129131 
    130       zprbio(:,:,:) = prmax(:,:,:) * zmxl_fac(:,:,:) 
    131       zprdia(:,:,:) = zprbio(:,:,:) 
     132      zprbio(:,:,:) = zprmaxn(:,:,:) * zmxl_fac(:,:,:) 
     133      zprdia(:,:,:) = zprmaxd(:,:,:) * zmxl_fac(:,:,:) 
    132134 
    133135      ! Maximum light intensity 
     
    169171                      !  Computation of production function for Chlorophyll 
    170172                      !-------------------------------------------------- 
    171                       zpislopen = zpislopeadn(ji,jj,jk) / ( prmax(ji,jj,jk) * zmxl_chl(ji,jj,jk) * rday + rtrn ) 
    172                       zpisloped = zpislopeadd(ji,jj,jk) / ( prmax(ji,jj,jk) * zmxl_chl(ji,jj,jk) * rday + rtrn ) 
    173                       zprnch(ji,jj,jk) = prmax(ji,jj,jk) * ( 1.- EXP( -zpislopen * enano(ji,jj,jk) ) ) 
    174                       zprdch(ji,jj,jk) = prmax(ji,jj,jk) * ( 1.- EXP( -zpisloped * ediat(ji,jj,jk) ) ) 
     173                      zpislopen = zpislopeadn(ji,jj,jk) / ( zprmaxn(ji,jj,jk) * zmxl_chl(ji,jj,jk) * rday + rtrn ) 
     174                      zpisloped = zpislopeadd(ji,jj,jk) / ( zprmaxd(ji,jj,jk) * zmxl_chl(ji,jj,jk) * rday + rtrn ) 
     175                      zprnch(ji,jj,jk) = zprmaxn(ji,jj,jk) * ( 1.- EXP( -zpislopen * enanom(ji,jj,jk) ) ) 
     176                      zprdch(ji,jj,jk) = zprmaxd(ji,jj,jk) * ( 1.- EXP( -zpisloped * ediatm(ji,jj,jk) ) ) 
    175177                  ENDIF 
    176178               END DO 
     
    192194                      zpislopen = zpislopen * zmxl_fac(ji,jj,jk) / ( zmxl_chl(ji,jj,jk) + rtrn ) 
    193195                      zpisloped = zpisloped * zmxl_fac(ji,jj,jk) / ( zmxl_chl(ji,jj,jk) + rtrn ) 
    194                       zprnch(ji,jj,jk) = prmax(ji,jj,jk) * ( 1.- EXP( -zpislopen * enano(ji,jj,jk) ) ) 
    195                       zprdch(ji,jj,jk) = prmax(ji,jj,jk) * ( 1.- EXP( -zpisloped * ediat(ji,jj,jk) ) ) 
     196                      zprnch(ji,jj,jk) = zprmaxn(ji,jj,jk) * ( 1.- EXP( -zpislopen * enanom(ji,jj,jk) ) ) 
     197                      zprdch(ji,jj,jk) = zprmaxd(ji,jj,jk) * ( 1.- EXP( -zpisloped * ediatm(ji,jj,jk) ) ) 
    196198                  ENDIF 
    197199               END DO 
     
    206208            DO ji = 1, jpi 
    207209                zval = MIN( xnanopo4(ji,jj,jk), ( xnanonh4(ji,jj,jk) + xnanono3(ji,jj,jk) ) )   & 
    208                 &      * prmax(ji,jj,jk) / ( zprbio(ji,jj,jk) + rtrn ) 
     210                &      * zprmaxn(ji,jj,jk) / ( zprbio(ji,jj,jk) + rtrn ) 
    209211                quotan(ji,jj,jk) = MIN( 1., 0.2 + 0.8 * zval ) 
    210212                zval = MIN( xdiatpo4(ji,jj,jk), ( xdiatnh4(ji,jj,jk) + xdiatno3(ji,jj,jk) ) )   & 
    211                 &      * prmax(ji,jj,jk) / ( zprdia(ji,jj,jk) + rtrn ) 
     213                &      * zprmaxd(ji,jj,jk) / ( zprdia(ji,jj,jk) + rtrn ) 
    212214                quotad(ji,jj,jk) = MIN( 1., 0.2 + 0.8 * zval ) 
    213215            END DO 
     
    227229                   !    to mimic the very high ratios observed in the Southern Ocean (silpot2) 
    228230                  zlim  = trb(ji,jj,jk,jpsil) / ( trb(ji,jj,jk,jpsil) + xksi1 ) 
    229                   zsilim = MIN( zprdia(ji,jj,jk) / ( prmax(ji,jj,jk) + rtrn ), xlimsi(ji,jj,jk) ) 
     231                  zsilim = MIN( zprdia(ji,jj,jk) / ( zprmaxd(ji,jj,jk) + rtrn ), xlimsi(ji,jj,jk) ) 
    230232                  zsilfac = 4.4 * EXP( -4.23 * zsilim ) * MAX( 0.e0, MIN( 1., 2.2 * ( zlim - 0.5 ) )  ) + 1.e0 
    231233                  zsiborn = trb(ji,jj,jk,jpsil) * trb(ji,jj,jk,jpsil) * trb(ji,jj,jk,jpsil) 
     
    266268                  zratio = trb(ji,jj,jk,jpnfe) / ( trb(ji,jj,jk,jpphy) * fecnm + rtrn ) 
    267269                  zmax   = MAX( 0., ( 1. - zratio ) / ABS( 1.05 - zratio ) )  
    268                   zprofen(ji,jj,jk) = fecnm * prmax(ji,jj,jk) * ( 1.0 - fr_i(ji,jj) )  & 
     270                  zprofen(ji,jj,jk) = fecnm * zprmaxn(ji,jj,jk) * ( 1.0 - fr_i(ji,jj) )  & 
    269271                  &             * ( 4. - 4.5 * xlimnfe(ji,jj,jk) / ( xlimnfe(ji,jj,jk) + 0.5 ) )    & 
    270272                  &             * biron(ji,jj,jk) / ( biron(ji,jj,jk) + concnfe(ji,jj,jk) )  & 
     
    276278                  zratio = trb(ji,jj,jk,jpdfe) / ( trb(ji,jj,jk,jpdia) * fecdm + rtrn ) 
    277279                  zmax   = MAX( 0., ( 1. - zratio ) / ABS( 1.05 - zratio ) )  
    278                   zprofed(ji,jj,jk) = fecdm * prmax(ji,jj,jk) * ( 1.0 - fr_i(ji,jj) )  & 
     280                  zprofed(ji,jj,jk) = fecdm * zprmaxd(ji,jj,jk) * ( 1.0 - fr_i(ji,jj) )  & 
    279281                  &             * ( 4. - 4.5 * xlimdfe(ji,jj,jk) / ( xlimdfe(ji,jj,jk) + 0.5 ) )    & 
    280282                  &             * biron(ji,jj,jk) / ( biron(ji,jj,jk) + concdfe(ji,jj,jk) )  & 
     
    291293               IF( etot_ndcy(ji,jj,jk) > 1.E-3 ) THEN 
    292294                  !  production terms for nanophyto. ( chlorophyll ) 
    293                   znanotot = enano(ji,jj,jk) / ( zmxl_chl(ji,jj,jk) + rtrn ) 
     295                  znanotot = enanom(ji,jj,jk) / ( zmxl_chl(ji,jj,jk) + rtrn ) 
    294296                  zprod    = rday * zprorcan(ji,jj,jk) * zprnch(ji,jj,jk) * xlimphy(ji,jj,jk) 
    295297                  zprochln = chlcmin * 12. * zprorcan (ji,jj,jk) 
     
    298300                                        & (  zpislopeadn(ji,jj,jk) * znanotot +rtrn) 
    299301                  !  production terms for diatoms ( chlorophyll ) 
    300                   zdiattot = ediat(ji,jj,jk) / ( zmxl_chl(ji,jj,jk) + rtrn ) 
     302                  zdiattot = ediatm(ji,jj,jk) / ( zmxl_chl(ji,jj,jk) + rtrn ) 
    301303                  zprod    = rday * zprorcad(ji,jj,jk) * zprdch(ji,jj,jk) * xlimdia(ji,jj,jk) 
    302304                  zprochld = chlcmin * 12. * zprorcad(ji,jj,jk) 
     
    351353                    zfeup    = texcretn * zprofen(ji,jj,jk) + texcretd * zprofed(ji,jj,jk) 
    352354                    tra(ji,jj,jk,jplgw) = tra(ji,jj,jk,jplgw) + zdocprod * ldocp - zfeup * plig(ji,jj,jk) * lthet 
     355                    zpligprod1(ji,jj,jk) = zdocprod * ldocp 
     356                    zpligprod2(ji,jj,jk) = zfeup * plig(ji,jj,jk) * lthet 
    353357                 ENDIF 
    354358              END DO 
     
    392396              CALL iom_put( "PFeD"  , zw3d ) 
    393397          ENDIF 
     398          IF( iom_use( "LPRODP" ) )  THEN 
     399              zw3d(:,:,:) = zpligprod1(:,:,:) * 1e9 * zfact * tmask(:,:,:) 
     400              CALL iom_put( "LPRODP"  , zw3d ) 
     401          ENDIF 
     402          IF( iom_use( "LDETP" ) )  THEN 
     403              zw3d(:,:,:) = zpligprod2(:,:,:) * 1e9 * zfact * tmask(:,:,:) 
     404              CALL iom_put( "LDETP"  , zw3d ) 
     405          ENDIF 
    394406          IF( iom_use( "Mumax" ) )  THEN 
    395               zw3d(:,:,:) = prmax(:,:,:) * tmask(:,:,:)   ! Maximum growth rate 
     407              zw3d(:,:,:) = zprmaxn(:,:,:) * tmask(:,:,:)   ! Maximum growth rate 
    396408              CALL iom_put( "Mumax"  , zw3d ) 
    397409          ENDIF 
     
    404416          ENDIF 
    405417          IF( iom_use( "LNlight" ) .OR. iom_use( "LDlight" ) )  THEN 
    406               zw3d(:,:,:) = zprbio (:,:,:) / (prmax(:,:,:) + rtrn) * tmask(:,:,:) ! light limitation term 
     418              zw3d(:,:,:) = zprbio (:,:,:) / (zprmaxn(:,:,:) + rtrn) * tmask(:,:,:) ! light limitation term 
    407419              CALL iom_put( "LNlight"  , zw3d ) 
    408420              ! 
    409               zw3d(:,:,:) =  zprdia (:,:,:) / (prmax(:,:,:) + rtrn) * tmask(:,:,:)  ! light limitation term 
     421              zw3d(:,:,:) = zprdia (:,:,:) / (zprmaxd(:,:,:) + rtrn) * tmask(:,:,:)  ! light limitation term 
    410422              CALL iom_put( "LDlight"  , zw3d ) 
    411423          ENDIF 
     
    542554      !!                     ***  ROUTINE p4z_prod_alloc  *** 
    543555      !!---------------------------------------------------------------------- 
    544       ALLOCATE( prmax(jpi,jpj,jpk), quotan(jpi,jpj,jpk), quotad(jpi,jpj,jpk), STAT = p4z_prod_alloc ) 
     556      ALLOCATE( quotan(jpi,jpj,jpk), quotad(jpi,jpj,jpk), STAT = p4z_prod_alloc ) 
    545557      ! 
    546558      IF( p4z_prod_alloc /= 0 ) CALL ctl_warn('p4z_prod_alloc : failed to allocate arrays.') 
  • NEMO/branches/2018/dev_r10164_HPC09_ESIWACE_PREP_MERGE/src/TOP/PISCES/P4Z/p4zrem.F90

    r10345 r10368  
    6363      REAL(wp) ::   zsatur, zsatur2, znusil, znusil2, zdep, zdepmin, zfactdep 
    6464      REAL(wp) ::   zbactfer, zolimit, zonitr, zrfact2 
    65       REAL(wp) ::   zammonic, zoxyrem 
     65      REAL(wp) ::   zammonic, zoxyremc, zoxyremn, zoxyremp 
    6666      REAL(wp) ::   zosil, ztem, zdenitnh4, zolimic, zolimin, zolimip, zdenitrn, zdenitrp 
    6767      CHARACTER (len=25) :: charout 
    6868      REAL(wp), DIMENSION(jpi,jpj    ) :: ztempbac 
    69       REAL(wp), DIMENSION(jpi,jpj,jpk) :: zdepbac, zolimi, zdepprod, zfacsi, zfacsib 
     69      REAL(wp), DIMENSION(jpi,jpj,jpk) :: zdepbac, zolimi, zdepprod, zfacsi, zfacsib, zdepeff, zfebact 
    7070      REAL(wp), ALLOCATABLE, DIMENSION(:,:,:) :: zw3d 
    7171      !!--------------------------------------------------------------------- 
     
    7575      ! Initialisation of arrys 
    7676      zdepprod(:,:,:) = 1._wp 
     77      zdepeff (:,:,:) = 0.3_wp 
    7778      ztempbac(:,:)   = 0._wp 
    7879      zfacsib(:,:,:)  = xsilab / ( 1.0 - xsilab ) 
     80      zfebact(:,:,:)  = 0._wp 
    7981      zfacsi(:,:,:)   = xsilab 
    8082 
     
    9597                  zdepbac (ji,jj,jk) = zdepmin**0.683 * ztempbac(ji,jj) 
    9698                  zdepprod(ji,jj,jk) = zdepmin**0.273 
     99                  zdepeff (ji,jj,jk) = zdepeff(ji,jj,jk) * zdepmin**0.3 
    97100               ENDIF 
    98101            END DO 
     
    117120                  denitr(ji,jj,jk)  = zammonic * ( 1. - nitrfac2(ji,jj,jk) ) 
    118121                  denitr(ji,jj,jk)  = MIN( ( trb(ji,jj,jk,jpno3) - rtrn ) / rdenit, denitr(ji,jj,jk) ) 
    119                   zoxyrem           = zammonic - denitr(ji,jj,jk) 
     122                  zoxyremc          = zammonic - denitr(ji,jj,jk) 
    120123                  ! 
    121124                  zolimi (ji,jj,jk) = MAX( 0.e0, zolimi (ji,jj,jk) ) 
    122125                  denitr (ji,jj,jk) = MAX( 0.e0, denitr (ji,jj,jk) ) 
    123                   zoxyrem           = MAX( 0.e0, zoxyrem ) 
     126                  zoxyremc          = MAX( 0.e0, zoxyremc ) 
    124127 
    125128                  ! 
    126                   tra(ji,jj,jk,jppo4) = tra(ji,jj,jk,jppo4) + zolimi (ji,jj,jk) + denitr(ji,jj,jk) + zoxyrem 
    127                   tra(ji,jj,jk,jpnh4) = tra(ji,jj,jk,jpnh4) + zolimi (ji,jj,jk) + denitr(ji,jj,jk) + zoxyrem 
     129                  tra(ji,jj,jk,jppo4) = tra(ji,jj,jk,jppo4) + zolimi (ji,jj,jk) + denitr(ji,jj,jk) + zoxyremc 
     130                  tra(ji,jj,jk,jpnh4) = tra(ji,jj,jk,jpnh4) + zolimi (ji,jj,jk) + denitr(ji,jj,jk) + zoxyremc 
    128131                  tra(ji,jj,jk,jpno3) = tra(ji,jj,jk,jpno3) - denitr (ji,jj,jk) * rdenit 
    129                   tra(ji,jj,jk,jpdoc) = tra(ji,jj,jk,jpdoc) - zolimi (ji,jj,jk) - denitr(ji,jj,jk) - zoxyrem 
     132                  tra(ji,jj,jk,jpdoc) = tra(ji,jj,jk,jpdoc) - zolimi (ji,jj,jk) - denitr(ji,jj,jk) - zoxyremc 
    130133                  tra(ji,jj,jk,jpoxy) = tra(ji,jj,jk,jpoxy) - zolimi (ji,jj,jk) * o2ut 
    131                   tra(ji,jj,jk,jpdic) = tra(ji,jj,jk,jpdic) + zolimi (ji,jj,jk) + denitr(ji,jj,jk) + zoxyrem 
    132                   tra(ji,jj,jk,jptal) = tra(ji,jj,jk,jptal) + rno3 * ( zolimi(ji,jj,jk) + zoxyrem    & 
     134                  tra(ji,jj,jk,jpdic) = tra(ji,jj,jk,jpdic) + zolimi (ji,jj,jk) + denitr(ji,jj,jk) + zoxyremc 
     135                  tra(ji,jj,jk,jptal) = tra(ji,jj,jk,jptal) + rno3 * ( zolimi(ji,jj,jk) + zoxyremc    & 
    133136                  &                     + ( rdenit + 1.) * denitr(ji,jj,jk) ) 
    134137               END DO 
     
    159162                  ! Ammonification in suboxic waters with denitrification 
    160163                  ! ------------------------------------------------------- 
    161                   zolimit = zremikc * nitrfac(ji,jj,jk) * trb(ji,jj,jk,jpdoc) 
    162                   denitr(ji,jj,jk)  = MIN(  ( trb(ji,jj,jk,jpno3) - rtrn ) / rdenit, zolimit ) 
    163                   denitr(ji,jj,jk) = MAX( 0.e0, denitr(ji,jj,jk) ) 
     164                  zammonic = zremikc * nitrfac(ji,jj,jk) * trb(ji,jj,jk,jpdoc) 
     165                  denitr(ji,jj,jk)  = zammonic * ( 1. - nitrfac2(ji,jj,jk) ) 
     166                  denitr(ji,jj,jk)  = MAX(0., MIN(  ( trb(ji,jj,jk,jpno3) - rtrn ) / rdenit, denitr(ji,jj,jk) ) ) 
     167                  zoxyremc          = MAX(0., zammonic - denitr(ji,jj,jk)) 
    164168                  zdenitrn  = zremikn * denitr(ji,jj,jk) * trb(ji,jj,jk,jpdon) / ( trb(ji,jj,jk,jpdoc) + rtrn ) 
    165169                  zdenitrp  = zremikp * denitr(ji,jj,jk) * trb(ji,jj,jk,jpdop) / ( trb(ji,jj,jk,jpdoc) + rtrn ) 
    166  
    167                   tra(ji,jj,jk,jppo4) = tra(ji,jj,jk,jppo4) + zolimip + zdenitrp 
    168                   tra(ji,jj,jk,jpnh4) = tra(ji,jj,jk,jpnh4) + zolimin + zdenitrn 
     170                  zoxyremn  = zremikn * zoxyremc * trb(ji,jj,jk,jpdon) / ( trb(ji,jj,jk,jpdoc) + rtrn ) 
     171                  zoxyremp  = zremikp * zoxyremc * trb(ji,jj,jk,jpdop) / ( trb(ji,jj,jk,jpdoc) + rtrn ) 
     172 
     173                  tra(ji,jj,jk,jppo4) = tra(ji,jj,jk,jppo4) + zolimip + zdenitrp + zoxyremp 
     174                  tra(ji,jj,jk,jpnh4) = tra(ji,jj,jk,jpnh4) + zolimin + zdenitrn + zoxyremn 
    169175                  tra(ji,jj,jk,jpno3) = tra(ji,jj,jk,jpno3) - denitr(ji,jj,jk) * rdenit 
    170                   tra(ji,jj,jk,jpdoc) = tra(ji,jj,jk,jpdoc) - zolimic - denitr(ji,jj,jk) 
    171                   tra(ji,jj,jk,jpdon) = tra(ji,jj,jk,jpdon) - zolimin - zdenitrn 
    172                   tra(ji,jj,jk,jpdop) = tra(ji,jj,jk,jpdop) - zolimip - zdenitrp 
     176                  tra(ji,jj,jk,jpdoc) = tra(ji,jj,jk,jpdoc) - zolimic - denitr(ji,jj,jk) - zoxyremc 
     177                  tra(ji,jj,jk,jpdon) = tra(ji,jj,jk,jpdon) - zolimin - zdenitrn - zoxyremn 
     178                  tra(ji,jj,jk,jpdop) = tra(ji,jj,jk,jpdop) - zolimip - zdenitrp - zoxyremp 
    173179                  tra(ji,jj,jk,jpoxy) = tra(ji,jj,jk,jpoxy) - zolimic * o2ut 
    174                   tra(ji,jj,jk,jpdic) = tra(ji,jj,jk,jpdic) + zolimic + denitr(ji,jj,jk) 
    175                   tra(ji,jj,jk,jptal) = tra(ji,jj,jk,jptal) + rno3 * ( zolimin + ( rdenit + 1.) * zdenitrn ) 
     180                  tra(ji,jj,jk,jpdic) = tra(ji,jj,jk,jpdic) + zolimic + denitr(ji,jj,jk) + zoxyremc 
     181                  tra(ji,jj,jk,jptal) = tra(ji,jj,jk,jptal) + rno3 * ( zolimin + zoxyremn + ( rdenit + 1.) * zdenitrn ) 
    176182               END DO 
    177183            END DO 
     
    215221               ! studies (especially at Papa) have shown this uptake to be significant 
    216222               ! ---------------------------------------------------------- 
    217                zbactfer = feratb *  rfact2 * prmax(ji,jj,jk) * xlimbacl(ji,jj,jk)             & 
     223               zbactfer = feratb *  rfact2 * 0.6_wp / rday * tgfunc(ji,jj,jk) * xlimbacl(ji,jj,jk)     & 
    218224                  &              * trb(ji,jj,jk,jpfer) / ( xkferb + trb(ji,jj,jk,jpfer) )    & 
    219                   &              * zdepprod(ji,jj,jk) * zdepbac(ji,jj,jk) 
    220                tra(ji,jj,jk,jpfer) = tra(ji,jj,jk,jpfer) - zbactfer*0.16 
    221                tra(ji,jj,jk,jpsfe) = tra(ji,jj,jk,jpsfe) + zbactfer*0.12 
    222                tra(ji,jj,jk,jpbfe) = tra(ji,jj,jk,jpbfe) + zbactfer*0.04 
     225                  &              * zdepprod(ji,jj,jk) * zdepeff(ji,jj,jk) * zdepbac(ji,jj,jk) 
     226               tra(ji,jj,jk,jpfer) = tra(ji,jj,jk,jpfer) - zbactfer*0.33 
     227               tra(ji,jj,jk,jpsfe) = tra(ji,jj,jk,jpsfe) + zbactfer*0.25 
     228               tra(ji,jj,jk,jpbfe) = tra(ji,jj,jk,jpbfe) + zbactfer*0.08 
     229               zfebact(ji,jj,jk)   = zbactfer * 0.33 
     230               blim(ji,jj,jk)      = xlimbacl(ji,jj,jk)  * zdepbac(ji,jj,jk) / 1.e-6 * zdepprod(ji,jj,jk) 
    223231            END DO 
    224232         END DO 
     
    256264               tra(ji,jj,jk,jpgsi) = tra(ji,jj,jk,jpgsi) - zosil 
    257265               tra(ji,jj,jk,jpsil) = tra(ji,jj,jk,jpsil) + zosil 
    258                ! 
    259266            END DO 
    260267         END DO 
     
    268275 
    269276      IF( knt == nrdttrc ) THEN 
     277          zrfact2 = 1.e3 * rfact2r 
    270278          ALLOCATE( zw3d(jpi,jpj,jpk) ) 
    271279          zfact = 1.e+3 * rfact2r  !  conversion from mol/l/kt to  mol/m3/s 
     
    278286              zw3d(:,:,:) = denitr(:,:,:) * rdenit * rno3 * tmask(:,:,:) * zfact ! Denitrification 
    279287              CALL iom_put( "DENIT"  , zw3d ) 
     288          ENDIF 
     289          IF( iom_use( "BACT" ) )  THEN 
     290               zw3d(:,:,:) = zdepbac(:,:,:) * 1.E6 * tmask(:,:,:)  ! Bacterial biomass 
     291               CALL iom_put( "BACT", zw3d ) 
     292          ENDIF 
     293          IF( iom_use( "FEBACT" ) )  THEN 
     294               zw3d(:,:,:) = zfebact(:,:,:) * 1E9 * tmask(:,:,:) * zrfact2   ! Bacterial iron consumption 
     295               CALL iom_put( "FEBACT" , zw3d ) 
    280296          ENDIF 
    281297          ! 
  • NEMO/branches/2018/dev_r10164_HPC09_ESIWACE_PREP_MERGE/src/TOP/PISCES/P4Z/p4zsbc.F90

    r10314 r10368  
    3131   REAL(wp), PUBLIC ::   dustsolub    !: Solubility of the dust 
    3232   REAL(wp), PUBLIC ::   mfrac        !: Mineral Content of the dust 
     33   REAL(wp), PUBLIC ::   rdustfep     !: Fraction of dust that is dissolvable 
    3334   REAL(wp), PUBLIC ::   icefeinput   !: Iron concentration in sea ice 
    3435   REAL(wp), PUBLIC ::   wdust        !: Sinking speed of the dust  
     
    134135                  DO ji = 1, jpi 
    135136                     zcoef = ryyss * e1e2t(ji,jj) * h_rnf(ji,jj)  
    136                      rivalk(ji,jj) =   sf_river(jr_dic)%fnow(ji,jj,1)                                    & 
     137                     rivalk(ji,jj) =   sf_river(jr_dic)%fnow(ji,jj,1)  & 
    137138                        &              * 1.E3        / ( 12. * zcoef + rtrn ) 
    138                      rivdic(ji,jj) = ( sf_river(jr_dic)%fnow(ji,jj,1) + sf_river(jr_doc)%fnow(ji,jj,1) ) & 
     139                     rivdic(ji,jj) =   sf_river(jr_dic)%fnow(ji,jj,1) & 
    139140                        &              * 1.E3         / ( 12. * zcoef + rtrn ) 
    140                      rivdin(ji,jj) = ( sf_river(jr_din)%fnow(ji,jj,1) + sf_river(jr_don)%fnow(ji,jj,1) ) & 
     141                     rivdin(ji,jj) =   sf_river(jr_din)%fnow(ji,jj,1) & 
    141142                        &              * 1.E3 / rno3 / ( 14. * zcoef + rtrn ) 
    142                      rivdip(ji,jj) = ( sf_river(jr_dip)%fnow(ji,jj,1) + sf_river(jr_dop)%fnow(ji,jj,1) ) & 
     143                     rivdip(ji,jj) =   sf_river(jr_dip)%fnow(ji,jj,1) & 
    143144                        &              * 1.E3 / po4r / ( 31. * zcoef + rtrn ) 
    144                      rivdsi(ji,jj) =   sf_river(jr_dsi)%fnow(ji,jj,1)                                    & 
     145                     rivdsi(ji,jj) =   sf_river(jr_dsi)%fnow(ji,jj,1)  & 
    145146                        &              * 1.E3        / ( 28.1 * zcoef + rtrn ) 
     147                     rivdoc(ji,jj) =   sf_river(jr_doc)%fnow(ji,jj,1)  & 
     148                        &              * 1.E3        / ( 12. * zcoef + rtrn )  
    146149                  END DO 
    147150               END DO 
     
    158161                     rivdip(ji,jj) = ( sf_river(jr_dip)%fnow(ji,jj,1) ) & 
    159162                        &              * 1.E3 / po4r / ( 31. * zcoef + rtrn ) * tmask(ji,jj,1) 
    160                      rivdoc(ji,jj) = ( sf_river(jr_doc)%fnow(ji,jj,1) ) & 
    161                         &              * 1.E3 / ( 12. * zcoef + rtrn ) * tmask(ji,jj,1) 
    162163                     rivdon(ji,jj) = ( sf_river(jr_don)%fnow(ji,jj,1) ) & 
    163164                        &              * 1.E3 / rno3 / ( 14. * zcoef + rtrn ) * tmask(ji,jj,1) 
    164165                     rivdop(ji,jj) = ( sf_river(jr_dop)%fnow(ji,jj,1) ) & 
    165166                        &              * 1.E3 / po4r / ( 31. * zcoef + rtrn ) * tmask(ji,jj,1) 
     167                     rivdsi(ji,jj) =   sf_river(jr_dsi)%fnow(ji,jj,1)  & 
     168                        &              * 1.E3        / ( 28.1 * zcoef + rtrn ) 
     169                     rivdoc(ji,jj) =   sf_river(jr_doc)%fnow(ji,jj,1)  & 
     170                        &              * 1.E3        / ( 12. * zcoef + rtrn ) 
    166171                  END DO 
    167172               END DO 
     
    223228        &                ln_dust, ln_solub, ln_river, ln_ndepo, ln_ironsed, ln_ironice, ln_hydrofe,    & 
    224229        &                sedfeinput, distcoast, dustsolub, icefeinput, wdust, mfrac, nitrfix, diazolight, concfediaz, & 
    225         &                hratio, fep_rats, fep_rath, lgw_rath 
     230        &                hratio, fep_rats, fep_rath, rdustfep, lgw_rath 
    226231      !!---------------------------------------------------------------------- 
    227232      ! 
     
    262267            WRITE(numout,*) '      Fep/Fer ratio from sed sources            fep_rats   = ', fep_rats 
    263268            WRITE(numout,*) '      Fep/Fer ratio from sed hydro sources      fep_rath   = ', fep_rath 
     269            WRITE(numout,*) '      Fraction of dust that is dissolvable      rdustfep   = ', rdustfep 
    264270            WRITE(numout,*) '      Weak ligand ratio from sed hydro sources  lgw_rath   = ', lgw_rath 
    265271         ENDIF 
     
    343349         slf_river(jr_dsi) = sn_riverdsi   
    344350         ! 
    345          ALLOCATE( rivdic(jpi,jpj), rivalk(jpi,jpj), rivdin(jpi,jpj), rivdip(jpi,jpj), rivdsi(jpi,jpj) )  
    346          IF( ln_p5z )  ALLOCATE( rivdon(jpi,jpj), rivdop(jpi,jpj), rivdoc(jpi,jpj) ) 
     351         ALLOCATE( rivdic(jpi,jpj), rivalk(jpi,jpj), rivdin(jpi,jpj), rivdip(jpi,jpj), rivdsi(jpi,jpj), rivdoc(jpi,jpj) )  
     352         IF( ln_p5z )  ALLOCATE( rivdon(jpi,jpj), rivdop(jpi,jpj) ) 
    347353         ! 
    348354         ALLOCATE( sf_river(jpriv), rivinput(jpriv), STAT=ierr1 )    !* allocate and fill sf_river (forcing structure) with sn_river_ 
  • NEMO/branches/2018/dev_r10164_HPC09_ESIWACE_PREP_MERGE/src/TOP/PISCES/P4Z/p4zsed.F90

    r10345 r10368  
    6464      CHARACTER (len=25) :: charout 
    6565      REAL(wp), DIMENSION(jpi,jpj    ) :: zdenit2d, zbureff, zwork 
    66       REAL(wp), DIMENSION(jpi,jpj    ) :: zwsbio3, zwsbio4, zwscal 
     66      REAL(wp), DIMENSION(jpi,jpj    ) :: zwsbio3, zwsbio4 
    6767      REAL(wp), DIMENSION(jpi,jpj    ) :: zsedcal, zsedsi, zsedc 
    6868      REAL(wp), DIMENSION(jpi,jpj,jpk) :: zsoufer, zlight 
     
    8585      ! 
    8686      ! Allocate temporary workspace 
    87       IF( ln_p5z )    ALLOCATE( ztrpo4(jpi,jpj,jpk), ztrdop(jpi,jpj,jpk) ) 
     87      ALLOCATE( ztrpo4(jpi,jpj,jpk) ) 
     88      IF( ln_p5z )    ALLOCATE( ztrdop(jpi,jpj,jpk) ) 
    8889      IF( ln_ligand ) ALLOCATE( zwsfep(jpi,jpj) ) 
    89  
    9090 
    9191      zdenit2d(:,:) = 0.e0 
     
    9595      zsedcal (:,:) = 0.e0 
    9696      zsedc   (:,:) = 0.e0 
    97  
    9897 
    9998      ! Iron input/uptake due to sea ice : Crude parameterization based on Lancelot et al. 
     
    132131         ELSE 
    133132            zirondep(:,:,1) = dustsolub  * dust(:,:) * mfrac * rfact2 / e3t_n(:,:,1) / 55.85 + 3.e-10 * r1_ryyss  
     133         ENDIF 
     134         IF ( ln_ligand ) THEN 
     135            IF( ln_solub ) THEN 
     136               tra(:,:,1,jpfep) = tra(:,:,1,jpfep) + rdustfep * (1.0 - solub(:,:)) * dust(:,:) * mfrac * rfact2 / e3t_n(:,:,1) / 55.85 
     137            ELSE 
     138               tra(:,:,1,jpfep) = tra(:,:,1,jpfep) + rdustfep * (1.0 - dustsolub) * dust(:,:) * mfrac * rfact2 / e3t_n(:,:,1) / 55.85 
     139            ENDIF 
    134140         ENDIF 
    135141         zsidep(:,:)   = 8.8 * 0.075 * dust(:,:) * mfrac * rfact2 / e3t_n(:,:,1) / 28.1  
     
    173179                  tra(ji,jj,jk,jpdic) = tra(ji,jj,jk,jpdic) +  rivdic(ji,jj) * rfact2 
    174180                  tra(ji,jj,jk,jptal) = tra(ji,jj,jk,jptal) +  ( rivalk(ji,jj) - rno3 * rivdin(ji,jj) ) * rfact2 
     181                  tra(ji,jj,jk,jpdoc) = tra(ji,jj,jk,jpdoc) +  rivdoc(ji,jj) * rfact2 
    175182               ENDDO 
    176183            ENDDO 
    177184         ENDDO 
     185         IF (ln_ligand) THEN 
     186            DO jj = 1, jpj 
     187               DO ji = 1, jpi 
     188                  DO jk = 1, nk_rnf(ji,jj) 
     189                     tra(ji,jj,jk,jplgw) = tra(ji,jj,jk,jplgw) +  rivdic(ji,jj) * 5.e-5 * rfact2 
     190                  ENDDO 
     191               ENDDO 
     192            ENDDO 
     193         ENDIF 
    178194         IF( ln_p5z ) THEN 
    179195            DO jj = 1, jpj 
     
    182198                     tra(ji,jj,jk,jpdop) = tra(ji,jj,jk,jpdop) + rivdop(ji,jj) * rfact2 
    183199                     tra(ji,jj,jk,jpdon) = tra(ji,jj,jk,jpdon) + rivdon(ji,jj) * rfact2 
    184                      tra(ji,jj,jk,jpdoc) = tra(ji,jj,jk,jpdoc) + rivdoc(ji,jj) * rfact2 
    185200                  ENDDO 
    186201               ENDDO 
     
    216231            zdep = e3t_n(ji,jj,ikt) / xstep 
    217232            zwsbio4(ji,jj) = MIN( 0.99 * zdep, wsbio4(ji,jj,ikt) ) 
    218             zwscal (ji,jj) = MIN( 0.99 * zdep, wscal (ji,jj,ikt) ) 
    219233            zwsbio3(ji,jj) = MIN( 0.99 * zdep, wsbio3(ji,jj,ikt) ) 
    220234         END DO 
     
    278292            ikt  = mbkt(ji,jj) 
    279293            zdep = xstep / e3t_n(ji,jj,ikt)  
    280             zwsc = zwscal (ji,jj) * zdep 
     294            zwsc = zwsbio4(ji,jj) * zdep 
    281295            zsiloss = trb(ji,jj,ikt,jpgsi) * zwsc 
    282296            zcaloss = trb(ji,jj,ikt,jpcal) * zwsc 
     
    292306               ikt  = mbkt(ji,jj) 
    293307               zdep = xstep / e3t_n(ji,jj,ikt)  
    294                zwsc = zwscal (ji,jj) * zdep 
     308               zwsc = zwsbio4(ji,jj) * zdep 
    295309               zsiloss = trb(ji,jj,ikt,jpgsi) * zwsc 
    296310               zcaloss = trb(ji,jj,ikt,jpcal) * zwsc 
     
    393407               DO ji = 1, jpi 
    394408                  !                      ! Potential nitrogen fixation dependant on temperature and iron 
    395                   zlim = ( 1.- xnanono3(ji,jj,jk) - xnanonh4(ji,jj,jk) ) 
    396                   IF( zlim <= 0.2 )   zlim = 0.01 
     409                  ztemp = tsn(ji,jj,jk,jp_tem) 
     410                  zmudia = MAX( 0.,-0.001096*ztemp**2 + 0.057*ztemp -0.637 ) * 7.625 
     411                  !       Potential nitrogen fixation dependant on temperature and iron 
     412                  xdianh4 = trb(ji,jj,jk,jpnh4) / ( concnnh4 + trb(ji,jj,jk,jpnh4) ) 
     413                  xdiano3 = trb(ji,jj,jk,jpno3) / ( concnno3 + trb(ji,jj,jk,jpno3) ) * (1. - xdianh4) 
     414                  zlim = ( 1.- xdiano3 - xdianh4 ) 
     415                  IF( zlim <= 0.1 )   zlim = 0.01 
    397416                  zfact = zlim * rfact2 
    398  
    399                   ztrfer  = biron(ji,jj,jk)       / ( concfediaz + biron(ji,jj,jk)       ) 
    400                   ztrpo4s = trb  (ji,jj,jk,jppo4) / ( concnnh4   + trb  (ji,jj,jk,jppo4) )  
    401                   nitrpot(ji,jj,jk) =  MAX( 0.e0, ( 0.6 * tgfunc(ji,jj,jk) - 2.15 ) * r1_rday ) & 
    402                     &                *  zfact * MIN( ztrfer, ztrpo4s ) * zlight(ji,jj,jk) 
     417                  ztrfer = biron(ji,jj,jk) / ( concfediaz + biron(ji,jj,jk) ) 
     418                  ztrpo4(ji,jj,jk) = trb(ji,jj,jk,jppo4) / ( 1E-6 + trb(ji,jj,jk,jppo4) ) 
     419                  ztrdp = ztrpo4(ji,jj,jk) 
     420                  nitrpot(ji,jj,jk) =  zmudia * r1_rday * zfact * MIN( ztrfer, ztrdp ) * zlight(ji,jj,jk) 
    403421               END DO 
    404422            END DO 
     
    434452               DO ji = 1, jpi 
    435453                  zfact = nitrpot(ji,jj,jk) * nitrfix 
    436                   tra(ji,jj,jk,jpnh4) = tra(ji,jj,jk,jpnh4) +             zfact 
    437                   tra(ji,jj,jk,jptal) = tra(ji,jj,jk,jptal) + rno3      * zfact 
    438                   tra(ji,jj,jk,jpoxy) = tra(ji,jj,jk,jpoxy) + o2nit     * zfact  
     454                  tra(ji,jj,jk,jpnh4) = tra(ji,jj,jk,jpnh4) + zfact / 3.0 
     455                  tra(ji,jj,jk,jptal) = tra(ji,jj,jk,jptal) + rno3 * zfact / 3.0 
     456                  tra(ji,jj,jk,jppo4) = tra(ji,jj,jk,jppo4) - zfact * 2.0 / 3.0 
     457                  tra(ji,jj,jk,jpdoc) = tra(ji,jj,jk,jpdoc) + zfact * 1.0 / 3.0 
     458                  tra(ji,jj,jk,jppoc) = tra(ji,jj,jk,jppoc) + zfact * 1.0 / 3.0 * 2.0 / 3.0 
     459                  tra(ji,jj,jk,jpgoc) = tra(ji,jj,jk,jpgoc) + zfact * 1.0 / 3.0 * 1.0 / 3.0 
     460                  tra(ji,jj,jk,jpoxy) = tra(ji,jj,jk,jpoxy) + ( o2ut + o2nit ) * zfact * 2.0 / 3.0 + o2nit * zfact / 3.0 
     461                  tra(ji,jj,jk,jpfer) = tra(ji,jj,jk,jpfer) - 30E-6 * zfact * 1.0 / 3.0 
     462                  tra(ji,jj,jk,jpsfe) = tra(ji,jj,jk,jpsfe) + 30E-6 * zfact * 1.0 / 3.0 * 2.0 / 3.0 
     463                  tra(ji,jj,jk,jpbfe) = tra(ji,jj,jk,jpbfe) + 30E-6 * zfact * 1.0 / 3.0 * 1.0 / 3.0 
     464                  tra(ji,jj,jk,jpfer) = tra(ji,jj,jk,jpfer) + 0.002 * 4E-10 * zsoufer(ji,jj,jk) * rfact2 / rday 
    439465                  tra(ji,jj,jk,jppo4) = tra(ji,jj,jk,jppo4) + concdnh4 / ( concdnh4 + trb(ji,jj,jk,jppo4) ) & 
    440                   &                     * 0.002 * trb(ji,jj,jk,jpdoc) * xstep 
    441                   tra(ji,jj,jk,jpfer) = tra(ji,jj,jk,jpfer) + 0.002 * 4E-10 * zsoufer(ji,jj,jk) * xstep 
     466                  &                     * 0.001 * trb(ji,jj,jk,jpdoc) * xstep 
    442467              END DO 
    443468            END DO  
  • NEMO/branches/2018/dev_r10164_HPC09_ESIWACE_PREP_MERGE/src/TOP/PISCES/P4Z/p4zsink.F90

    r10314 r10368  
    143143      END DO 
    144144 
    145       wscal (:,:,:) = wsbio4(:,:,:) 
    146  
    147145      !  Initializa to zero all the sinking arrays  
    148146      !   ----------------------------------------- 
     
    165163        CALL p4z_sink2( wsbio4, sinkfer2, jpbfe, iiter2 ) 
    166164        CALL p4z_sink2( wsbio4, sinksil , jpgsi, iiter2 ) 
    167         CALL p4z_sink2( wscal , sinkcal , jpcal, iiter2 ) 
     165        CALL p4z_sink2( wsbio4, sinkcal , jpcal, iiter2 ) 
    168166      END DO 
    169167 
  • NEMO/branches/2018/dev_r10164_HPC09_ESIWACE_PREP_MERGE/src/TOP/PISCES/P4Z/p5zlim.F90

    r10345 r10368  
    1414   USE oce_trc         ! Shared ocean-passive tracers variables 
    1515   USE trc             ! Tracers defined 
     16   USE p4zlim 
    1617   USE sms_pisces      ! PISCES variables 
    1718   USE iom             !  I/O manager 
     
    2526 
    2627   !! * Shared module variables 
    27    REAL(wp), PUBLIC ::  concnno3    !:  NO3, PO4 half saturation    
    2828   REAL(wp), PUBLIC ::  concpno3    !:  NO3, PO4 half saturation    
    29    REAL(wp), PUBLIC ::  concdno3    !:  Phosphate half saturation for diatoms   
    30    REAL(wp), PUBLIC ::  concnnh4    !:  NH4 half saturation for phyto   
    3129   REAL(wp), PUBLIC ::  concpnh4    !:  NH4 half saturation for phyto   
    32    REAL(wp), PUBLIC ::  concdnh4    !:  NH4 half saturation for diatoms 
    3330   REAL(wp), PUBLIC ::  concnpo4    !:  NH4 half saturation for diatoms 
    3431   REAL(wp), PUBLIC ::  concppo4    !:  NH4 half saturation for diatoms 
    3532   REAL(wp), PUBLIC ::  concdpo4    !:  NH4 half saturation for diatoms 
    36    REAL(wp), PUBLIC ::  concnfer    !:  Iron half saturation for nanophyto  
    3733   REAL(wp), PUBLIC ::  concpfer    !:  Iron half saturation for nanophyto  
    38    REAL(wp), PUBLIC ::  concdfer    !:  Iron half saturation for diatoms   
    39    REAL(wp), PUBLIC ::  concbno3    !:  NO3 half saturation  for bacteria  
    40    REAL(wp), PUBLIC ::  concbnh4    !:  NH4 half saturation for bacteria 
    4134   REAL(wp), PUBLIC ::  concbpo4    !:  PO4 half saturation for bacteria 
    42    REAL(wp), PUBLIC ::  xsizedia    !:  Minimum size criteria for diatoms 
    4335   REAL(wp), PUBLIC ::  xsizepic    !:  Minimum size criteria for diatoms 
    44    REAL(wp), PUBLIC ::  xsizephy    !:  Minimum size criteria for nanophyto 
    45    REAL(wp), PUBLIC ::  xsizern     !:  Size ratio for nanophytoplankton 
    4636   REAL(wp), PUBLIC ::  xsizerp     !:  Size ratio for nanophytoplankton 
    47    REAL(wp), PUBLIC ::  xsizerd     !:  Size ratio for diatoms 
    48    REAL(wp), PUBLIC ::  xksi1       !:  half saturation constant for Si uptake  
    49    REAL(wp), PUBLIC ::  xksi2       !:  half saturation constant for Si/C  
    50    REAL(wp), PUBLIC ::  xkdoc       !:  2nd half-sat. of DOC remineralization   
    51    REAL(wp), PUBLIC ::  concbfe     !:  Fe half saturation for bacteria  
    52    REAL(wp), PUBLIC ::  oxymin      !:  half saturation constant for anoxia 
    5337   REAL(wp), PUBLIC ::  qfnopt      !:  optimal Fe quota for nanophyto 
    5438   REAL(wp), PUBLIC ::  qfpopt      !:  optimal Fe quota for nanophyto 
    5539   REAL(wp), PUBLIC ::  qfdopt      !:  optimal Fe quota for diatoms 
    56    REAL(wp), PUBLIC ::  caco3r      !:  mean rainratio  
    5740   REAL(wp), PUBLIC ::  qnnmin      !:  optimal Fe quota for diatoms 
    5841   REAL(wp), PUBLIC ::  qnnmax      !:  optimal Fe quota for diatoms 
     
    8972 
    9073   !!* Phytoplankton limitation terms 
    91    REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:)  ::   xnanono3   !: ??? 
    92    REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:)  ::   xdiatno3   !: ??? 
    93    REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:)  ::   xnanonh4   !: ??? 
    94    REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:)  ::   xdiatnh4   !: ??? 
    95    REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:)  ::   xnanopo4   !: ??? 
    96    REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:)  ::   xdiatpo4   !: ??? 
    97    REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:)  ::   xlimphy    !: ??? 
    98    REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:)  ::   xlimdia    !: ??? 
    99    REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:)  ::   xlimnfe    !: ??? 
    100    REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:)  ::   xlimdfe    !: ??? 
    101    REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:)  ::   xlimsi     !: ??? 
    102    REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:)  ::   xlimbac    !: ?? 
    103    REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:)  ::   xlimbacl   !: ?? 
    10474   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:)  ::   xpicono3   !: ??? 
    10575   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:)  ::   xpiconh4   !: ??? 
     
    574544         &      xpicopo4(jpi,jpj,jpk), xpicodop(jpi,jpj,jpk),       & 
    575545         &      xnanodop(jpi,jpj,jpk), xdiatdop(jpi,jpj,jpk),       & 
    576          &      xnanono3(jpi,jpj,jpk), xdiatno3(jpi,jpj,jpk),       & 
    577          &      xnanonh4(jpi,jpj,jpk), xdiatnh4(jpi,jpj,jpk),       & 
    578          &      xnanopo4(jpi,jpj,jpk), xdiatpo4(jpi,jpj,jpk),       & 
    579          &      xlimphy (jpi,jpj,jpk), xlimdia (jpi,jpj,jpk),       & 
    580          &      xlimnfe (jpi,jpj,jpk), xlimdfe (jpi,jpj,jpk),       & 
    581          &      xlimbac (jpi,jpj,jpk), xlimbacl(jpi,jpj,jpk),       & 
    582546         &      xnanofer(jpi,jpj,jpk), xdiatfer(jpi,jpj,jpk),       & 
    583547         &      xpicofer(jpi,jpj,jpk), xlimpfe (jpi,jpj,jpk),       & 
    584548         &      fvnuptk (jpi,jpj,jpk), fvduptk (jpi,jpj,jpk),       & 
    585          &      fvpuptk (jpi,jpj,jpk), xlimpic (jpi,jpj,jpk),       & 
    586          &      xlimsi  (jpi,jpj,jpk), STAT=ierr(1) ) 
     549         &      fvpuptk (jpi,jpj,jpk), xlimpic (jpi,jpj,jpk),    STAT=ierr(1) ) 
    587550         ! 
    588551      !*  Minimum/maximum quotas of phytoplankton 
  • NEMO/branches/2018/dev_r10164_HPC09_ESIWACE_PREP_MERGE/src/TOP/PISCES/P4Z/p5zmeso.F90

    r10070 r10368  
    2727   REAL(wp), PUBLIC ::  part2        !: part of calcite not dissolved in mesozoo guts 
    2828   REAL(wp), PUBLIC ::  xpref2c      !: mesozoo preference for POC  
    29    REAL(wp), PUBLIC ::  xpref2p      !: mesozoo preference for nanophyto 
     29   REAL(wp), PUBLIC ::  xpref2n      !: mesozoo preference for nanophyto 
    3030   REAL(wp), PUBLIC ::  xpref2z      !: mesozoo preference for zooplankton 
    3131   REAL(wp), PUBLIC ::  xpref2d      !: mesozoo preference for Diatoms  
     
    4040   REAL(wp), PUBLIC ::  mzrat2       !: microzooplankton mortality rate  
    4141   REAL(wp), PUBLIC ::  grazrat2     !: maximal mesozoo grazing rate 
    42    REAL(wp), PUBLIC ::  xkgraz2      !: non assimilated fraction of P by mesozoo 
    43    REAL(wp), PUBLIC ::  unass2c      !: Efficicency of mesozoo growth  
    44    REAL(wp), PUBLIC ::  unass2n      !: Efficicency of mesozoo growth  
    45    REAL(wp), PUBLIC ::  unass2p      !: Efficicency of mesozoo growth  
    46    REAL(wp), PUBLIC ::  epsher2      !: half sturation constant for grazing 2 
     42   REAL(wp), PUBLIC ::  xkgraz2      !: Half-saturation constant of assimilation 
     43   REAL(wp), PUBLIC ::  unass2c      !: Non-assimilated fraction of food 
     44   REAL(wp), PUBLIC ::  unass2n      !: Non-assimilated fraction of food 
     45   REAL(wp), PUBLIC ::  unass2p      !: Non-assimilated fraction of food 
     46   REAL(wp), PUBLIC ::  epsher2      !: Growth efficiency of mesozoo 
     47   REAL(wp), PUBLIC ::  epsher2min   !: Minimum growth efficiency of mesozoo 
    4748   REAL(wp), PUBLIC ::  ssigma2      !: Fraction excreted as semi-labile DOM 
    4849   REAL(wp), PUBLIC ::  srespir2     !: Active respiration 
     
    8586      CHARACTER (len=25) :: charout 
    8687      REAL(wp) :: zrfact2, zmetexcess 
    87       REAL(wp), DIMENSION(jpi,jpj,jpk) :: zgrazing 
    88       REAL(wp), ALLOCATABLE, DIMENSION(:,:,:) :: zw3d 
     88      REAL(wp), DIMENSION(jpi,jpj,jpk) :: zgrazing, zfezoo2 
     89      REAL(wp), ALLOCATABLE, DIMENSION(:,:,:) :: zw3d, zz2ligprod 
    8990 
    9091      !!--------------------------------------------------------------------- 
     
    9293      IF( ln_timing )   CALL timing_start('p5z_meso') 
    9394      ! 
     95 
    9496      zgrazing(:,:,:) = 0._wp 
     97      zfezoo2 (:,:,:) = 0._wp 
     98      ! 
     99      IF (ln_ligand) THEN 
     100         ALLOCATE( zz2ligprod(jpi,jpj,jpk) ) 
     101         zz2ligprod(:,:,:) = 0._wp 
     102      ENDIF 
    95103 
    96104      zmetexcess = 0.0 
     
    111119               !   no real reason except that it seems to be more stable and may mimic predation 
    112120               !   --------------------------------------------------------------- 
    113                ztortz   = mzrat2 * 1.e6 * zfact * trb(ji,jj,jk,jpmes) 
     121               ztortz   = mzrat2 * 1.e6 * zfact * trb(ji,jj,jk,jpmes) * (1. - nitrfac(ji,jj,jk)) 
    114122 
    115123               !   Computation of the abundance of the preys 
     
    124132               !   Mesozooplankton grazing 
    125133               !   ------------------------ 
    126                zfood     = xpref2d * zcompadi + xpref2z * zcompaz + xpref2p * zcompaph + xpref2c * zcompapoc   & 
     134               zfood     = xpref2d * zcompadi + xpref2z * zcompaz + xpref2n * zcompaph + xpref2c * zcompapoc   & 
    127135               &           + xpref2m * zcompames  
    128136               zfoodlim  = MAX( 0., zfood - MIN( 0.5 * zfood, xthresh2 ) ) 
    129137               zdenom    = zfoodlim / ( xkgraz2 + zfoodlim ) 
    130                zgraze2   = grazrat2 * xstep * tgfunc2(ji,jj,jk) * trb(ji,jj,jk,jpmes)  
     138               zgraze2   = grazrat2 * xstep * tgfunc2(ji,jj,jk) * trb(ji,jj,jk,jpmes) * (1. - nitrfac(ji,jj,jk))  
    131139 
    132140               !   An active switching parameterization is used here. 
     
    138146               !   most abundant species 
    139147               !   ------------------------------------------------------------   
    140                ztmp1 = xpref2p * zcompaph**1.5 
     148               ztmp1 = xpref2n * zcompaph**1.5 
    141149               ztmp2 = xpref2m * zcompames**1.5 
    142150               ztmp3 = xpref2c * zcompapoc**1.5 
     
    170178               !   ---------------------------------- 
    171179               zgrazffeg = grazflux  * xstep * wsbio4(ji,jj,jk)      & 
    172                &           * tgfunc2(ji,jj,jk) * trb(ji,jj,jk,jpgoc) * trb(ji,jj,jk,jpmes) 
     180               &           * tgfunc2(ji,jj,jk) * trb(ji,jj,jk,jpgoc) * trb(ji,jj,jk,jpmes)  & 
     181               &           * (1. - nitrfac(ji,jj,jk)) 
    173182               zgrazfffg = zgrazffeg * trb(ji,jj,jk,jpbfe) / (trb(ji,jj,jk,jpgoc) + rtrn) 
    174183               zgrazffng = zgrazffeg * trb(ji,jj,jk,jpgon) / (trb(ji,jj,jk,jpgoc) + rtrn) 
    175184               zgrazffpg = zgrazffeg * trb(ji,jj,jk,jpgop) / (trb(ji,jj,jk,jpgoc) + rtrn) 
    176185               zgrazffep = grazflux  * xstep *  wsbio3(ji,jj,jk)     & 
    177                &           * tgfunc2(ji,jj,jk) * trb(ji,jj,jk,jppoc) * trb(ji,jj,jk,jpmes) 
     186               &           * tgfunc2(ji,jj,jk) * trb(ji,jj,jk,jppoc) * trb(ji,jj,jk,jpmes)   & 
     187               &           * (1. - nitrfac(ji,jj,jk)) 
    178188               zgrazfffp = zgrazffep * trb(ji,jj,jk,jpsfe) / (trb(ji,jj,jk,jppoc) + rtrn) 
    179189               zgrazffnp = zgrazffep * trb(ji,jj,jk,jppon) / (trb(ji,jj,jk,jppoc) + rtrn) 
     
    226236               !   --------------------------------------------------- 
    227237               zepshert  = MIN( 1., zgrasratn/ no3rat3, zgrasratp/ po4rat3, zgrasratf / ferat3) 
    228                zbeta = 1./ (epsher2 - 0.2) 
    229                zepsherf = 0.2 + 1./ (zbeta + 0.04 * 12. * zfood *1E6 ) 
     238               zbeta     = MAX(0., (epsher2 - epsher2min) ) 
     239               zepsherf  = epsher2min + zbeta / ( 1.0 + 0.04E6 * 12. * zfood * zbeta ) 
    230240               zepsherv  = zepsherf * zepshert 
    231241 
     
    290300               tra(ji,jj,jk,jpdoc) = tra(ji,jj,jk,jpdoc) + zgradoc 
    291301               ! 
    292                IF( ln_ligand ) tra(ji,jj,jk,jplgw) = tra(ji,jj,jk,jplgw) + zgradoc * ldocz 
     302               IF( ln_ligand ) THEN 
     303                  tra(ji,jj,jk,jplgw)  = tra(ji,jj,jk,jplgw) + zgradoc * ldocz 
     304                  zz2ligprod(ji,jj,jk) = zgradoc * ldocz 
     305               ENDIF 
    293306               ! 
    294307               tra(ji,jj,jk,jpdon) = tra(ji,jj,jk,jpdon) + zgradon 
     
    296309               tra(ji,jj,jk,jpoxy) = tra(ji,jj,jk,jpoxy) - o2ut * zgrarem 
    297310               tra(ji,jj,jk,jpfer) = tra(ji,jj,jk,jpfer) + zgraref 
     311               zfezoo2(ji,jj,jk)   = zgraref 
    298312               tra(ji,jj,jk,jpdic) = tra(ji,jj,jk,jpdic) + zgrarem 
    299313               tra(ji,jj,jk,jptal) = tra(ji,jj,jk,jptal) + rno3 * zgraren 
     
    351365            CALL iom_put( "PCAL", zw3d ) 
    352366         ENDIF 
     367         IF( iom_use( "FEZOO2" ) ) THEN 
     368            zw3d(:,:,:) = zfezoo2(:,:,:) * 1e9 * 1.e+3 * rfact2r * tmask(:,:,:)   ! 
     369            CALL iom_put( "FEZOO2", zw3d ) 
     370         ENDIF 
     371         IF( iom_use( "LPRODZ2" ) .AND. ln_ligand )  THEN 
     372            zw3d(:,:,:) = zz2ligprod(:,:,:) * 1e9 * 1.e+3 * rfact2r * tmask(:,:,:) 
     373            CALL iom_put( "LPRODZ2"  , zw3d ) 
     374         ENDIF 
    353375         DEALLOCATE( zw3d ) 
    354376      ENDIF 
     
    379401      INTEGER :: ios                 ! Local integer output status for namelist read 
    380402      !! 
    381       NAMELIST/namp5zmes/part2, bmetexc2, grazrat2, resrat2, mzrat2, xpref2c, xpref2p, xpref2z, & 
     403      NAMELIST/namp5zmes/part2, bmetexc2, grazrat2, resrat2, mzrat2, xpref2c, xpref2n, xpref2z, & 
    382404         &                xpref2m, xpref2d, xthresh2dia, xthresh2phy, xthresh2zoo, xthresh2poc, & 
    383          &                xthresh2mes, xthresh2, xkgraz2, epsher2, ssigma2, unass2c, & 
     405         &                xthresh2mes, xthresh2, xkgraz2, epsher2, epsher2min, ssigma2, unass2c, & 
    384406         &                unass2n, unass2p, srespir2, grazflux 
    385407      !!---------------------------------------------------------------------- 
     
    399421         WRITE(numout,*) ' ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~' 
    400422         WRITE(numout,*) '    part of calcite not dissolved in mesozoo guts  part2       = ', part2 
    401          WRITE(numout,*) '    mesozoo preference for nano.                   xpref2p     = ', xpref2p 
     423         WRITE(numout,*) '    mesozoo preference for nano.                   xpref2n     = ', xpref2n 
    402424         WRITE(numout,*) '    mesozoo preference for diatoms                 xpref2d     = ', xpref2d 
    403425         WRITE(numout,*) '    mesozoo preference for zoo                     xpref2z     = ', xpref2z 
     
    415437         WRITE(numout,*) '    mesozoo flux feeding rate                      grazflux    = ', grazflux 
    416438         WRITE(numout,*) '    C egested fraction of food by mesozoo          unass2c     = ', unass2c 
    417          WRITE(numout,*) '    C egested fraction of food by mesozoo          unass2n     = ', unass2n 
    418          WRITE(numout,*) '    C egested fraction of food by mesozoo          unass2p     = ', unass2p 
     439         WRITE(numout,*) '    N egested fraction of food by mesozoo          unass2n     = ', unass2n 
     440         WRITE(numout,*) '    P egested fraction of food by mesozoo          unass2p     = ', unass2p 
    419441         WRITE(numout,*) '    Efficicency of Mesozoo growth                  epsher2     = ', epsher2 
     442         WRITE(numout,*) '    Minimum Efficiency of Mesozoo growth           epsher2min  =', epsher2min 
    420443         WRITE(numout,*) '    Fraction excreted as semi-labile DOM           ssigma2     = ', ssigma2 
    421444         WRITE(numout,*) '    Active respiration                             srespir2    = ', srespir2 
  • NEMO/branches/2018/dev_r10164_HPC09_ESIWACE_PREP_MERGE/src/TOP/PISCES/P4Z/p5zmicro.F90

    r10345 r10368  
    1515   USE trc             !  passive tracers common variables  
    1616   USE sms_pisces      !  PISCES Source Minus Sink variables 
     17   USE p4zlim 
    1718   USE p5zlim          !  Phytoplankton limitation terms 
    1819   USE iom             !  I/O manager 
     
    4142   REAL(wp), PUBLIC ::  mzrat       !: microzooplankton mortality rate  
    4243   REAL(wp), PUBLIC ::  grazrat     !: maximal microzoo grazing rate 
    43    REAL(wp), PUBLIC ::  xkgraz      !: non assimilated fraction of P by microzoo 
    44    REAL(wp), PUBLIC ::  unassc      !: Efficicency of microzoo growth  
    45    REAL(wp), PUBLIC ::  unassn      !: Efficicency of microzoo growth  
    46    REAL(wp), PUBLIC ::  unassp      !: Efficicency of microzoo growth  
    47    REAL(wp), PUBLIC ::  epsher      !: half sturation constant for grazing 1  
     44   REAL(wp), PUBLIC ::  xkgraz      !: Half-saturation constant of assimilation 
     45   REAL(wp), PUBLIC ::  unassc      !: Non-assimilated part of food 
     46   REAL(wp), PUBLIC ::  unassn      !: Non-assimilated part of food 
     47   REAL(wp), PUBLIC ::  unassp      !: Non-assimilated part of food 
     48   REAL(wp), PUBLIC ::  epsher      !: Growth efficiency for microzoo 
     49   REAL(wp), PUBLIC ::  epshermin   !: Minimum growth efficiency for microzoo 
    4850   REAL(wp), PUBLIC ::  srespir     !: half sturation constant for grazing 1  
    4951   REAL(wp), PUBLIC ::  ssigma      !: Fraction excreted as semi-labile DOM 
     
    8284      REAL(wp) :: zgrazdc, zgrazdn, zgrazdp, zgrazdf, zgraznf, zgrazz 
    8385      REAL(wp) :: zgrazpc, zgrazpn, zgrazpp, zgrazpf, zbeta, zrfact2, zmetexcess 
    84       REAL(wp), DIMENSION(jpi,jpj,jpk) :: zgrazing 
    85       REAL(wp), ALLOCATABLE, DIMENSION(:,:,:) :: zw3d 
     86      REAL(wp), DIMENSION(jpi,jpj,jpk) :: zgrazing, zfezoo 
     87      REAL(wp), ALLOCATABLE, DIMENSION(:,:,:) :: zw3d, zzligprod 
    8688      CHARACTER (len=25) :: charout 
    8789      !!--------------------------------------------------------------------- 
    8890      ! 
    8991      IF( ln_timing )   CALL timing_start('p5z_micro') 
     92      ! 
     93      IF (ln_ligand) THEN 
     94         ALLOCATE( zzligprod(jpi,jpj,jpk) ) 
     95         zzligprod(:,:,:) = 0._wp 
     96      ENDIF 
    9097      ! 
    9198      zmetexcess = 0.0 
     
    106113               !   no real reason except that it seems to be more stable and may mimic predation. 
    107114               !   ------------------------------------------------------------------------------ 
    108                ztortz = mzrat * 1.e6 * zfact * trb(ji,jj,jk,jpzoo) 
     115               ztortz = mzrat * 1.e6 * zfact * trb(ji,jj,jk,jpzoo) * (1. - nitrfac(ji,jj,jk)) 
    109116 
    110117               !   Computation of the abundance of the preys 
     
    123130               zfoodlim  = MAX( 0. , zfood - min(xthresh,0.5*zfood) ) 
    124131               zdenom    = zfoodlim / ( xkgraz + zfoodlim ) 
    125                zgraze    = grazrat * xstep * tgfunc2(ji,jj,jk) * trb(ji,jj,jk,jpzoo)  
     132               zgraze    = grazrat * xstep * tgfunc2(ji,jj,jk) * trb(ji,jj,jk,jpzoo) * (1. - nitrfac(ji,jj,jk))  
    126133 
    127134               !   An active switching parameterization is used here. 
     
    183190               !   --------------------------------------------------- 
    184191               zepshert  = MIN( 1., zgrasratn/ no3rat3, zgrasratp/ po4rat3, zgrasratf / ferat3) 
    185                zbeta = 1./ (epsher - 0.2) 
    186                zepsherf = 0.2 + 1./ (zbeta + 0.04 * 12. * zfood * 1E6 ) 
     192               zbeta     = MAX( 0., (epsher - epshermin) ) 
     193               zepsherf  = epshermin + zbeta / ( 1.0 + 0.04E6 * 12. * zfood * zbeta ) 
    187194               zepsherv  = zepsherf * zepshert 
    188195 
     
    244251               tra(ji,jj,jk,jpdoc) = tra(ji,jj,jk,jpdoc) + zgradoc 
    245252               ! 
    246                IF( ln_ligand ) tra(ji,jj,jk,jplgw) = tra(ji,jj,jk,jplgw) + zgradoc * ldocz 
     253               IF( ln_ligand ) THEN  
     254                  tra(ji,jj,jk,jplgw) = tra(ji,jj,jk,jplgw) + zgradoc * ldocz 
     255                  zzligprod(ji,jj,jk) = zgradoc * ldocz 
     256               ENDIF 
    247257               ! 
    248258               tra(ji,jj,jk,jpdon) = tra(ji,jj,jk,jpdon) + zgradon 
     
    250260               tra(ji,jj,jk,jpoxy) = tra(ji,jj,jk,jpoxy) - o2ut * zgrarem  
    251261               tra(ji,jj,jk,jpfer) = tra(ji,jj,jk,jpfer) + zgraref 
     262               zfezoo(ji,jj,jk)    = zgraref 
    252263               tra(ji,jj,jk,jpzoo) = tra(ji,jj,jk,jpzoo) + zepsherv * zgraztotc - zrespirc - ztortz - zgrazz 
    253264               tra(ji,jj,jk,jpphy) = tra(ji,jj,jk,jpphy) - zgraznc 
     
    288299      END DO 
    289300      ! 
    290       IF( lk_iomput .AND. knt == nrdttrc ) THEN 
    291          ALLOCATE( zw3d(jpi,jpj,jpk) ) 
    292          IF( iom_use( "GRAZ1" ) ) THEN 
    293             zw3d(:,:,:) = zgrazing(:,:,:) * 1.e+3 * rfact2r * tmask(:,:,:)  !  Total grazing of phyto by zooplankton 
    294             CALL iom_put( "GRAZ1", zw3d ) 
     301      IF( lk_iomput ) THEN 
     302         IF( knt == nrdttrc ) THEN 
     303            ALLOCATE( zw3d(jpi,jpj,jpk) ) 
     304            IF( iom_use( "GRAZ1" ) ) THEN 
     305               zw3d(:,:,:) = zgrazing(:,:,:) * 1.e+3 * rfact2r * tmask(:,:,:)  !  Total grazing of phyto by zooplankton 
     306               CALL iom_put( "GRAZ1", zw3d ) 
     307            ENDIF 
     308            IF( iom_use( "FEZOO" ) ) THEN 
     309               zw3d(:,:,:) = zfezoo(:,:,:) * 1e9 * 1.e+3 * rfact2r * tmask(:,:,:)   ! 
     310               CALL iom_put( "FEZOO", zw3d ) 
     311            ENDIF 
     312            IF( iom_use( "LPRODZ" ) .AND. ln_ligand )  THEN 
     313               zw3d(:,:,:) = zzligprod(:,:,:) * 1e9 * 1.e+3 * rfact2r * tmask(:,:,:) 
     314               CALL iom_put( "LPRODZ"  , zw3d ) 
     315            ENDIF 
     316            DEALLOCATE( zw3d ) 
    295317         ENDIF 
    296          DEALLOCATE( zw3d ) 
    297318      ENDIF 
    298319      ! 
     
    325346         &                xprefp, xprefd, xprefz, xthreshdia, xthreshphy, & 
    326347         &                xthreshpic, xthreshpoc, xthreshzoo, xthresh, xkgraz, & 
    327          &                epsher, ssigma, srespir, unassc, unassn, unassp 
     348         &                epsher, epshermin, ssigma, srespir, unassc, unassn, unassp 
    328349      !!---------------------------------------------------------------------- 
    329350      ! 
     
    360381         WRITE(numout,*) '    P egested fraction of fodd by microzoo          unassp      =', unassp 
    361382         WRITE(numout,*) '    Efficicency of microzoo growth                  epsher      =', epsher 
     383         WRITE(numout,*) '    Minimum Efficiency of Microzoo growth           epshermin   =', epshermin 
    362384         WRITE(numout,*) '    Fraction excreted as semi-labile DOM            ssigma      =', ssigma 
    363385         WRITE(numout,*) '    Active respiration                              srespir     =', srespir 
  • NEMO/branches/2018/dev_r10164_HPC09_ESIWACE_PREP_MERGE/src/TOP/PISCES/P4Z/p5zmort.F90

    r10345 r10368  
    1414   USE trc             !  passive tracers common variables  
    1515   USE sms_pisces      !  PISCES Source Minus Sink variables 
     16   USE p4zlim 
    1617   USE p5zlim          !  Phytoplankton limitation terms 
    1718   USE prtctl_trc      !  print control for debugging 
  • NEMO/branches/2018/dev_r10164_HPC09_ESIWACE_PREP_MERGE/src/TOP/PISCES/P4Z/p5zprod.F90

    r10345 r10368  
    1616   USE trc             !  passive tracers common variables  
    1717   USE sms_pisces      !  PISCES Source Minus Sink variables 
     18   USE p4zlim 
    1819   USE p5zlim          !  Co-limitations of differents nutrients 
    1920   USE prtctl_trc      !  print control for debugging 
     
    4243   REAL(wp), PUBLIC ::  grosip          !: 
    4344 
    44    REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:) ::   prmaxn   !: optimal production = f(temperature) 
    45    REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:) ::   prmaxp   !: optimal production = f(temperature) 
    46    REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:) ::   prmaxd   !: optimal production = f(temperature) 
    4745   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:)   ::   zdaylen 
    4846    
     
    8381      REAL(wp), DIMENSION(jpi,jpj    ) :: zmixnano, zmixpico, zmixdiat, zstrn 
    8482      REAL(wp), DIMENSION(jpi,jpj,jpk) :: zpislopeadn, zpislopeadp, zpislopeadd 
     83      REAL(wp), DIMENSION(jpi,jpj,jpk) :: zprnut, zprmaxp, zprmaxn, zprmaxd 
    8584      REAL(wp), DIMENSION(jpi,jpj,jpk) :: zprbio, zprpic, zprdia, zysopt 
    8685      REAL(wp), DIMENSION(jpi,jpj,jpk) :: zprchln, zprchlp, zprchld 
     
    9190      REAL(wp), DIMENSION(jpi,jpj,jpk) :: zpropo4n, zpropo4p, zpropo4d 
    9291      REAL(wp), DIMENSION(jpi,jpj,jpk) :: zprodopn, zprodopp, zprodopd 
    93       REAL(wp), DIMENSION(jpi,jpj,jpk) :: zrespn, zrespp, zrespd, zprnut 
     92      REAL(wp), DIMENSION(jpi,jpj,jpk) :: zrespn, zrespp, zrespd 
    9493      REAL(wp), DIMENSION(jpi,jpj,jpk) :: zcroissn, zcroissp, zcroissd 
    9594      REAL(wp), DIMENSION(jpi,jpj,jpk) :: zmxl_fac, zmxl_chl 
     95      REAL(wp), DIMENSION(jpi,jpj,jpk) :: zpligprod1, zpligprod2 
    9696      REAL(wp), ALLOCATABLE, DIMENSION(:,:,:) :: zw3d 
    9797      REAL(wp), ALLOCATABLE, DIMENSION(:,:  ) :: zw2d 
     
    106106      zpropo4n(:,:,:) = 0._wp ; zpropo4p(:,:,:) = 0._wp ; zpropo4d(:,:,:) = 0._wp 
    107107      zprdia  (:,:,:) = 0._wp ; zprpic  (:,:,:) = 0._wp ; zprbio  (:,:,:) = 0._wp 
     108      zprodopn(:,:,:) = 0._wp ; zprodopp(:,:,:) = 0._wp ; zprodopd(:,:,:) = 0._wp 
    108109      zysopt  (:,:,:) = 0._wp 
    109110      zrespn  (:,:,:) = 0._wp ; zrespp  (:,:,:) = 0._wp ; zrespd  (:,:,:) = 0._wp  
    110111 
    111112      ! Computation of the optimal production 
    112       prmaxn(:,:,:) = ( 0.65_wp * (1. + zpsino3 * qnpmax ) ) * r1_rday * tgfunc(:,:,:) 
    113       prmaxp(:,:,:) = 0.5 / 0.65 * prmaxn(:,:,:)  
    114       prmaxd(:,:,:) = prmaxn(:,:,:)  
    115       zprnut(:,:,:) = 0.65_wp * r1_rday * tgfunc(:,:,:) 
     113      zprnut (:,:,:) = 0.65_wp * r1_rday * tgfunc(:,:,:) 
     114      zprmaxn(:,:,:) = ( 0.65_wp * (1. + zpsino3 * qnpmax ) ) * r1_rday * tgfunc(:,:,:) 
     115      zprmaxp(:,:,:) = 0.5 / 0.65 * zprmaxn(:,:,:)  
     116      zprmaxd(:,:,:) = zprmaxn(:,:,:)  
    116117 
    117118      ! compute the day length depending on latitude and the day 
     
    145146      END DO 
    146147 
    147       zprbio(:,:,:) = prmaxn(:,:,:) * zmxl_fac(:,:,:) 
    148       zprdia(:,:,:) = prmaxd(:,:,:) * zmxl_fac(:,:,:) 
    149       zprpic(:,:,:) = prmaxp(:,:,:) * zmxl_fac(:,:,:) 
     148      zprbio(:,:,:) = zprmaxn(:,:,:) * zmxl_fac(:,:,:) 
     149      zprdia(:,:,:) = zprmaxd(:,:,:) * zmxl_fac(:,:,:) 
     150      zprpic(:,:,:) = zprmaxp(:,:,:) * zmxl_fac(:,:,:) 
    150151 
    151152 
     
    184185                  zpisloped = zpisloped * zmxl_fac(ji,jj,jk) / ( zmxl_chl(ji,jj,jk) + rtrn ) 
    185186                  zpislopep = zpislopep * zmxl_fac(ji,jj,jk) / ( zmxl_chl(ji,jj,jk) + rtrn ) 
    186                   zprchln(ji,jj,jk) = prmaxn(ji,jj,jk) * ( 1.- EXP( -zpislopen * enano(ji,jj,jk) )  ) 
    187                   zprchlp(ji,jj,jk) = prmaxp(ji,jj,jk) * ( 1.- EXP( -zpislopep * epico(ji,jj,jk) )  ) 
    188                   zprchld(ji,jj,jk) = prmaxd(ji,jj,jk) * ( 1.- EXP( -zpisloped * ediat(ji,jj,jk) )  ) 
     187                  zprchln(ji,jj,jk) = zprmaxn(ji,jj,jk) * ( 1.- EXP( -zpislopen * enanom(ji,jj,jk) )  ) 
     188                  zprchlp(ji,jj,jk) = zprmaxp(ji,jj,jk) * ( 1.- EXP( -zpislopep * epicom(ji,jj,jk) )  ) 
     189                  zprchld(ji,jj,jk) = zprmaxd(ji,jj,jk) * ( 1.- EXP( -zpisloped * ediatm(ji,jj,jk) )  ) 
    189190               ENDIF 
    190191            END DO 
     
    203204                  !    to mimic the very high ratios observed in the Southern Ocean (silpot2) 
    204205                  zlim  = trb(ji,jj,jk,jpsil) / ( trb(ji,jj,jk,jpsil) + xksi1 ) 
    205                   zsilim = MIN( zprdia(ji,jj,jk) / ( prmaxd(ji,jj,jk) + rtrn ), xlimsi(ji,jj,jk) ) 
     206                  zsilim = MIN( zprdia(ji,jj,jk) / ( zprmaxd(ji,jj,jk) + rtrn ), xlimsi(ji,jj,jk) ) 
    206207                  zsilfac = 3.4 * EXP( -4.23 * zsilim ) * MAX( 0.e0, MIN( 1., 2.2 * ( zlim - 0.5 ) )  ) + 1.e0 
    207208                  zsiborn = trb(ji,jj,jk,jpsil) * trb(ji,jj,jk,jpsil) * trb(ji,jj,jk,jpsil) 
     
    347348               IF( etot_ndcy(ji,jj,jk) > 1.E-3 ) THEN 
    348349                     !  production terms for nanophyto. ( chlorophyll ) 
    349                   znanotot = enano(ji,jj,jk) / ( zmxl_chl(ji,jj,jk) + rtrn ) 
     350                  znanotot = enanom(ji,jj,jk) / ( zmxl_chl(ji,jj,jk) + rtrn ) 
    350351                  zprod = rday * (zpronewn(ji,jj,jk) + zproregn(ji,jj,jk)) * zprchln(ji,jj,jk) * xlimphy(ji,jj,jk) 
    351352                  thetannm_n   = MIN ( thetannm, ( thetannm / (1. - 1.14 / 43.4 *tsn(ji,jj,jk,jp_tem)))   & 
     
    354355                  zprochln = MAX(zprochln, chlcmin * 12. * zprorcan (ji,jj,jk) ) 
    355356                     !  production terms for picophyto. ( chlorophyll ) 
    356                   zpicotot = epico(ji,jj,jk) / ( zmxl_chl(ji,jj,jk) + rtrn ) 
     357                  zpicotot = epicom(ji,jj,jk) / ( zmxl_chl(ji,jj,jk) + rtrn ) 
    357358                  zprod = rday * (zpronewp(ji,jj,jk) + zproregp(ji,jj,jk)) * zprchlp(ji,jj,jk) * xlimpic(ji,jj,jk) 
    358359                  thetanpm_n   = MIN ( thetanpm, ( thetanpm / (1. - 1.14 / 43.4 *tsn(ji,jj,jk,jp_tem)))   & 
     
    361362                  zprochlp = MAX(zprochlp, chlcmin * 12. * zprorcap(ji,jj,jk) ) 
    362363                  !  production terms for diatomees ( chlorophyll ) 
    363                   zdiattot = ediat(ji,jj,jk) / ( zmxl_chl(ji,jj,jk) + rtrn ) 
     364                  zdiattot = ediatm(ji,jj,jk) / ( zmxl_chl(ji,jj,jk) + rtrn ) 
    364365                  zprod = rday * (zpronewd(ji,jj,jk) + zproregd(ji,jj,jk)) * zprchld(ji,jj,jk) * xlimdia(ji,jj,jk) 
    365366                  thetandm_n   = MIN ( thetandm, ( thetandm / (1. - 1.14 / 43.4 *tsn(ji,jj,jk,jp_tem)))   & 
     
    449450                 zfeup    = texcretn * zprofen(ji,jj,jk) + texcretd * zprofed(ji,jj,jk) + texcretp * zprofep(ji,jj,jk) 
    450451                 tra(ji,jj,jk,jplgw) = tra(ji,jj,jk,jplgw) + zdocprod * ldocp - zfeup * plig(ji,jj,jk) * lthet 
     452                 zpligprod1(ji,jj,jk) = zdocprod * ldocp 
     453                 zpligprod2(ji,jj,jk) = zfeup * plig(ji,jj,jk) * lthet 
    451454              END DO 
    452455           END DO 
     
    500503              CALL iom_put( "PFeD"  , zw3d ) 
    501504          ENDIF 
     505          IF( iom_use( "LPRODP" ) )  THEN 
     506              zw3d(:,:,:) = zpligprod1(:,:,:) * 1e9 * zfact * tmask(:,:,:) 
     507              CALL iom_put( "LPRODP"  , zw3d ) 
     508          ENDIF 
     509          IF( iom_use( "LDETP" ) )  THEN 
     510              zw3d(:,:,:) = zpligprod2(:,:,:) * 1e9 * zfact * tmask(:,:,:) 
     511              CALL iom_put( "LDETP"  , zw3d ) 
     512          ENDIF 
    502513          IF( iom_use( "Mumax" ) )  THEN 
    503               zw3d(:,:,:) = prmaxn(:,:,:) * tmask(:,:,:)   ! Maximum growth rate 
     514              zw3d(:,:,:) = zprmaxn(:,:,:) * tmask(:,:,:)   ! Maximum growth rate 
    504515              CALL iom_put( "Mumax"  , zw3d ) 
    505516          ENDIF 
     
    515526          ENDIF 
    516527          IF( iom_use( "LNlight" ) .OR. iom_use( "LDlight" ) .OR. iom_use( "LPlight" ) )  THEN 
    517               zw3d(:,:,:) = zprbio (:,:,:) / (prmaxn(:,:,:) + rtrn) * tmask(:,:,:) ! light limitation term 
     528              zw3d(:,:,:) = zprbio (:,:,:) / (zprmaxn(:,:,:) + rtrn) * tmask(:,:,:) ! light limitation term 
    518529              CALL iom_put( "LNlight"  , zw3d ) 
    519530              ! 
    520               zw3d(:,:,:) = zprpic (:,:,:) / (prmaxp(:,:,:) + rtrn) * tmask(:,:,:) ! light limitation term 
     531              zw3d(:,:,:) = zprpic (:,:,:) / (zprmaxp(:,:,:) + rtrn) * tmask(:,:,:) ! light limitation term 
    521532              CALL iom_put( "LPlight"  , zw3d ) 
    522533              ! 
    523               zw3d(:,:,:) =  zprdia (:,:,:) / (prmaxd(:,:,:) + rtrn) * tmask(:,:,:)  ! light limitation term 
     534              zw3d(:,:,:) =  zprdia (:,:,:) / (zprmaxd(:,:,:) + rtrn) * tmask(:,:,:)  ! light limitation term 
    524535              CALL iom_put( "LDlight"  , zw3d ) 
    525536          ENDIF 
     
    611622      !!                     ***  ROUTINE p5z_prod_alloc  *** 
    612623      !!---------------------------------------------------------------------- 
    613       ALLOCATE( prmaxn(jpi,jpj,jpk), prmaxp(jpi,jpj,jpk), prmaxd(jpi,jpj,jpk),  & 
    614       &          zdaylen(jpi,jpj), STAT = p5z_prod_alloc ) 
     624      ALLOCATE( zdaylen(jpi,jpj), STAT = p5z_prod_alloc ) 
    615625      ! 
    616626      IF( p5z_prod_alloc /= 0 ) CALL ctl_warn('p5z_prod_alloc : failed to allocate arrays.') 
  • NEMO/branches/2018/dev_r10164_HPC09_ESIWACE_PREP_MERGE/src/TOP/PISCES/SED/oce_sed.F90

    r10345 r10368  
    1111   USE par_trc 
    1212 
    13    USE dom_oce , ONLY :   nidom     =>   nidom          !: 
    1413   USE dom_oce , ONLY :   glamt     =>   glamt          !: longitude of t-point (degre) 
    1514   USE dom_oce , ONLY :   gphit     =>   gphit          !: latitude  of t-point (degre) 
     
    2120   USE dom_oce , ONLY :   rdt       =>   rdt            !: time step for the dynamics 
    2221   USE dom_oce , ONLY :   nyear     =>   nyear          !: Current year 
    23    USE dom_oce , ONLY :   nmonth    =>   nmonth         !: Current month 
    24    USE dom_oce , ONLY :   nday      =>   nday           !: Current day 
    2522   USE dom_oce , ONLY :   ndastp    =>   ndastp         !: time step date in year/month/day aammjj 
    26    USE dom_oce , ONLY :   nday_year =>   nday_year      !: curent day counted from jan 1st of the current year 
    2723   USE dom_oce , ONLY :   adatrj    =>   adatrj         !: number of elapsed days since the begining of the run 
    2824   USE trc     , ONLY :  nittrc000  =>   nittrc000 
     
    3531   USE sms_pisces, ONLY : wsbio4    =>   wsbio4          !: sinking flux for POC 
    3632   USE sms_pisces, ONLY : wsbio3    =>   wsbio3          !: sinking flux for GOC 
    37    USE sms_pisces, ONLY : wscal     =>   wscal           !: sinking flux for calcite 
    3833   USE sms_pisces, ONLY : wsbio2    =>   wsbio2           !: sinking flux for calcite 
    3934   USE sms_pisces, ONLY : wsbio     =>   wsbio           !: sinking flux for calcite 
     35   USE sms_pisces, ONLY : ln_p5z    =>   ln_p5z          !: PISCES-QUOTA flag 
    4036   USE p4zche, ONLY     : akb3      =>   akb3            !: Chemical constants   
    4137   USE sms_pisces, ONLY : ak13      =>   ak13            !: Chemical constants   
  • NEMO/branches/2018/dev_r10164_HPC09_ESIWACE_PREP_MERGE/src/TOP/PISCES/SED/sedchem.F90

    r10345 r10368  
    678678         &         - 0.07711*zsal + 0.0041249*zsal15 
    679679 
     680         ! CONVERT FROM DIFFERENT PH SCALES 
     681         total2free  = 1.0/(1.0 + zst/zcks) 
     682         free2SWS    = 1. + zst/zcks + zft/(zckf*total2free) 
     683         total2SWS   = total2free * free2SWS 
     684         SWS2total   = 1.0 / total2SWS 
     685 
     686 
    680687         ! K1, K2 OF CARBONIC ACID, KB OF BORIC ACID, KW (H2O) (LIT.?) 
    681688         zak1    = 10**(zck1) * total2SWS 
     
    743750         aksis(ji) = zaksi * EXP( zbuf1 * zcpexp + zbuf2 * zcpexp2 ) 
    744751 
    745          ! CONVERT FROM DIFFERENT PH SCALES 
    746          total2free  = 1.0/(1.0 + zst/aks3s(ji)) 
    747          free2SWS    = 1. + zst/aks3s(ji) + zft/akf3s(ji) 
    748          total2SWS   = total2free * free2SWS 
    749          SWS2total   = 1.0 / total2SWS 
    750  
    751752         ! Convert to total scale 
    752753         ak1s(ji)  = ak1s(ji)  * SWS2total 
  • NEMO/branches/2018/dev_r10164_HPC09_ESIWACE_PREP_MERGE/src/TOP/PISCES/SED/seddsr.F90

    r10345 r10368  
    591591            zgamma = pwcp(ji,jk,jwalk) - 2.0 * zsedtra(4) 
    592592            zdelta = pwcp(ji,jk,jwpo4) - redfep * zsedtra(4) 
    593             zsedtra(4) = ( zsedtra(4) * zalpha ) / ( 0.25 * zsedtra(4) * ( exp( reac_fe2 * zalpha * dtsed2 / 2. ) - 1.0 )  & 
    594             &            + zalpha * exp( reac_fe2 * zalpha * dtsed2 / 2. ) + rtrn ) 
     593            IF ( zalpha == 0. ) THEN 
     594               zsedtra(4) = zsedtra(4) / ( 1.0 + zsedtra(4) * reac_fe2 * dtsed2 / 2.0 ) 
     595            ELSE 
     596               zsedtra(4) = ( zsedtra(4) * zalpha ) / ( 0.25 * zsedtra(4) * ( exp( reac_fe2 * zalpha * dtsed2 / 2. ) - 1.0 )  & 
     597               &            + zalpha * exp( reac_fe2 * zalpha * dtsed2 / 2. ) ) 
     598            ENDIF 
    595599            zsedtra(1) = zalpha + 0.25 * zsedtra(4)  
    596600            zsedtra(5) = zbeta  - zsedtra(4) 
     
    601605            zbeta  = pwcp(ji,jk,jwso4) + zsedtra(2) 
    602606            zgamma = pwcp(ji,jk,jwalk) - 2.0 * zsedtra(2) 
    603             zsedtra(2) = ( zsedtra(2) * zalpha ) / ( 2.0 * zsedtra(2) * ( exp( reac_h2s * zalpha * dtsed2 / 2. ) - 1.0 )  & 
    604             &            + zalpha * exp( reac_h2s * zalpha * dtsed2 / 2. ) + rtrn ) 
     607            IF ( zalpha == 0. ) THEN 
     608               zsedtra(2) = zsedtra(2) / ( 1.0 + zsedtra(2) * reac_h2s * dtsed2 / 2.0 ) 
     609            ELSE 
     610               zsedtra(2) = ( zsedtra(2) * zalpha ) / ( 2.0 * zsedtra(2) * ( exp( reac_h2s * zalpha * dtsed2 / 2. ) - 1.0 )  & 
     611               &            + zalpha * exp( reac_h2s * zalpha * dtsed2 / 2. ) ) 
     612            ENDIF 
    605613            zsedtra(1) = zalpha + 2.0 * zsedtra(2) 
    606614            pwcp(ji,jk,jwalk) = zgamma + 2.0 * zsedtra(2) 
     
    609617            zalpha = zsedtra(1) - 2.0 * zsedtra(3) 
    610618            zgamma = pwcp(ji,jk,jwalk) - 2.0 * zsedtra(3) 
    611             zsedtra(3) = ( zsedtra(3) * zalpha ) / ( 2.0 * zsedtra(3) * ( exp( reac_nh4 * zadsnh4 * zalpha * dtsed2 / 2. ) - 1.0 )  & 
    612             &            + zalpha * exp( reac_nh4 * zadsnh4 * zalpha * dtsed2 /2. ) + rtrn ) 
     619            IF ( zalpha == 0. ) THEN 
     620               zsedtra(3) = zsedtra(3) / ( 1.0 + zsedtra(3) * reac_nh4 * zadsnh4 * dtsed2 / 2.0 ) 
     621            ELSE 
     622               zsedtra(3) = ( zsedtra(3) * zalpha ) / ( 2.0 * zsedtra(3) * ( exp( reac_nh4 * zadsnh4 * zalpha * dtsed2 / 2. ) - 1.0 )  & 
     623               &            + zalpha * exp( reac_nh4 * zadsnh4 * zalpha * dtsed2 /2. ) ) 
     624            ENDIF 
    613625            zsedtra(1) = zalpha + 2.0 * zsedtra(3) 
    614626            pwcp(ji,jk,jwalk) = zgamma + 2.0 * zsedtra(3) 
     
    617629            zbeta  = zsedtra(4) + zsedtra(6) 
    618630            zgamma = pwcp(ji,jk,jwso4) + zsedtra(6) 
    619             zsedtra(6) = ( zsedtra(6) * zalpha ) / ( 2.0 * zsedtra(6) * ( exp( reac_feso * zalpha * dtsed2 / 2. ) - 1.0 )  & 
    620             &            + zalpha * exp( reac_feso * zalpha * dtsed2 /2. ) + rtrn ) 
     631            IF ( zalpha == 0. ) THEN 
     632               zsedtra(6) = zsedtra(6) / ( 1.0 + zsedtra(6) * reac_feso * dtsed2 / 2.0 ) 
     633            ELSE 
     634               zsedtra(6) = ( zsedtra(6) * zalpha ) / ( 2.0 * zsedtra(6) * ( exp( reac_feso * zalpha * dtsed2 / 2. ) - 1.0 )  & 
     635               &            + zalpha * exp( reac_feso * zalpha * dtsed2 /2. ) ) 
     636            ENDIF 
    621637            zsedtra(1) = zalpha + 2.0 * zsedtra(6) 
    622638            zsedtra(4) = zbeta  - zsedtra(6) 
     
    626642            zbeta  = zsedtra(4) + zsedtra(6) 
    627643            zgamma = pwcp(ji,jk,jwalk) - 2.0 * zsedtra(4) 
    628             zsedtra(4) = ( zsedtra(4) * zalpha ) / ( zsedtra(4) * ( exp( reac_fes * zalpha * dtsed2 / 2. ) - 1.0 )  & 
    629             &            + zalpha * exp( reac_fes * zalpha * dtsed2 /2. ) + rtrn ) 
     644            IF ( zalpha == 0. ) THEN 
     645               zsedtra(4) = zsedtra(4) / ( 1.0 + zsedtra(4) * reac_fes * dtsed2 / 2.0 ) 
     646            ELSE 
     647               zsedtra(4) = ( zsedtra(4) * zalpha ) / ( zsedtra(4) * ( exp( reac_fes * zalpha * dtsed2 / 2. ) - 1.0 )  & 
     648               &            + zalpha * exp( reac_fes * zalpha * dtsed2 /2. ) ) 
     649            ENDIF 
    630650            zsedtra(2) = zalpha + zsedtra(4) 
    631651            zsedtra(6) = zbeta  - zsedtra(4) 
     
    637657            zdelta = pwcp(ji,jk,jwso4) + zsedtra(2) 
    638658            zepsi  = pwcp(ji,jk,jwpo4) + redfep * zsedtra(5) 
    639             zsedtra(2) = ( zsedtra(2) * zalpha ) / ( 2.0 * zsedtra(2) * ( exp( reac_feh2s * zalpha * dtsed2 ) - 1.0 )  & 
    640             &            + zalpha * exp( reac_feh2s * zalpha * dtsed2 ) + rtrn ) 
     659            IF ( zalpha == 0. ) THEN 
     660               zsedtra(2) = zsedtra(2) / ( 1.0 + zsedtra(2) * reac_feh2s * dtsed2 ) 
     661            ELSE 
     662               zsedtra(2) = ( zsedtra(2) * zalpha ) / ( 2.0 * zsedtra(2) * ( exp( reac_feh2s * zalpha * dtsed2 ) - 1.0 )  & 
     663               &            + zalpha * exp( reac_feh2s * zalpha * dtsed2 ) ) 
     664            ENDIF 
    641665            zsedtra(5) = zalpha + 2.0 * zsedtra(2) 
    642666            zsedtra(4) = zbeta  - zsedtra(5) 
     
    648672            zbeta  = zsedtra(4) + zsedtra(6) 
    649673            zgamma = pwcp(ji,jk,jwalk) - 2.0 * zsedtra(4) 
    650             zsedtra(4) = ( zsedtra(4) * zalpha ) / ( zsedtra(4) * ( exp( reac_fes * zalpha * dtsed2 / 2. ) - 1.0 )  & 
    651             &            + zalpha * exp( reac_fes * zalpha * dtsed2 /2. ) + rtrn ) 
     674            IF ( zalpha == 0. ) THEN 
     675               zsedtra(4) = zsedtra(4) / ( 1.0 + zsedtra(4) * reac_fes * dtsed2 / 2.0 ) 
     676            ELSE 
     677               zsedtra(4) = ( zsedtra(4) * zalpha ) / ( zsedtra(4) * ( exp( reac_fes * zalpha * dtsed2 / 2. ) - 1.0 )  & 
     678               &            + zalpha * exp( reac_fes * zalpha * dtsed2 /2. ) ) 
     679            ENDIF 
    652680            zsedtra(2) = zalpha + zsedtra(4) 
    653681            zsedtra(6) = zbeta  - zsedtra(4) 
     
    657685            zbeta  = zsedtra(4) + zsedtra(6) 
    658686            zgamma = pwcp(ji,jk,jwso4) + zsedtra(6) 
    659             zsedtra(6) = ( zsedtra(6) * zalpha ) / ( 2.0 * zsedtra(6) * ( exp( reac_feso * zalpha * dtsed2 / 2. ) - 1.0 )  & 
    660             &            + zalpha * exp( reac_feso * zalpha * dtsed2 /2. ) + rtrn ) 
     687            IF (zalpha == 0.) THEN 
     688               zsedtra(6) = zsedtra(6) / ( 1.0 + zsedtra(6) * reac_feso * dtsed2 / 2. ) 
     689            ELSE 
     690               zsedtra(6) = ( zsedtra(6) * zalpha ) / ( 2.0 * zsedtra(6) * ( exp( reac_feso * zalpha * dtsed2 / 2. ) - 1.0 )  & 
     691               &            + zalpha * exp( reac_feso * zalpha * dtsed2 /2. ) ) 
     692            ENDIF 
    661693            zsedtra(1) = zalpha + 2.0 * zsedtra(6) 
    662694            zsedtra(4) = zbeta  - zsedtra(6) 
     
    665697            zalpha = zsedtra(1) - 2.0 * zsedtra(3) 
    666698            zgamma = pwcp(ji,jk,jwalk) - 2.0 * zsedtra(3) 
    667             zsedtra(3) = ( zsedtra(3) * zalpha ) / ( 2.0 * zsedtra(3) * ( exp( reac_nh4 * zadsnh4 * zalpha * dtsed2 / 2. ) - 1.0 )  & 
    668             &            + zalpha * exp( reac_nh4 * zadsnh4 * zalpha * dtsed2 /2. ) + rtrn ) 
     699            IF (zalpha == 0.) THEN 
     700               zsedtra(3) = zsedtra(3) / ( 1.0 + zsedtra(3) * reac_nh4 * zadsnh4 * dtsed2 / 2.0)  
     701            ELSE 
     702               zsedtra(3) = ( zsedtra(3) * zalpha ) / ( 2.0 * zsedtra(3) * ( exp( reac_nh4 * zadsnh4 * zalpha * dtsed2 / 2. ) - 1.0 )  & 
     703               &            + zalpha * exp( reac_nh4 * zadsnh4 * zalpha * dtsed2 /2. ) ) 
     704            ENDIF 
    669705            zsedtra(1) = zalpha + 2.0 * zsedtra(3) 
    670706            pwcp(ji,jk,jwalk) = zgamma + 2.0 * zsedtra(3) 
     
    673709            zbeta  = pwcp(ji,jk,jwso4) + zsedtra(2) 
    674710            zgamma = pwcp(ji,jk,jwalk) - 2.0 * zsedtra(2) 
    675             zsedtra(2) = ( zsedtra(2) * zalpha ) / ( 2.0 * zsedtra(2) * ( exp( reac_h2s * zalpha * dtsed2 / 2. ) - 1.0 )  & 
    676             &            + zalpha * exp( reac_h2s * zalpha * dtsed2 / 2. ) + rtrn ) 
     711            IF ( zalpha == 0. ) THEN 
     712               zsedtra(2) = zsedtra(2) / ( 1.0 + zsedtra(2) * reac_h2s * dtsed2 / 2.0 ) 
     713            ELSE 
     714               zsedtra(2) = ( zsedtra(2) * zalpha ) / ( 2.0 * zsedtra(2) * ( exp( reac_h2s * zalpha * dtsed2 / 2. ) - 1.0 )  & 
     715               &            + zalpha * exp( reac_h2s * zalpha * dtsed2 / 2. ) ) 
     716            ENDIF 
    677717            zsedtra(1) = zalpha + 2.0 * zsedtra(2) 
    678718            pwcp(ji,jk,jwso4) = zbeta - zsedtra(2) 
     
    683723            zgamma = pwcp(ji,jk,jwalk) - 2.0 * zsedtra(4) 
    684724            zdelta = pwcp(ji,jk,jwpo4) - redfep * zsedtra(4) 
    685             zsedtra(4) = ( zsedtra(4) * zalpha ) / ( 0.25 * zsedtra(4) * ( exp( reac_fe2 * zalpha * dtsed2 / 2. ) - 1.0 )  & 
    686             &            + zalpha * exp( reac_fe2 * zalpha * dtsed2 / 2. ) + rtrn ) 
     725            IF ( zalpha == 0. ) THEN 
     726               zsedtra(4) = zsedtra(4) / ( 1.0 + zsedtra(4) * reac_fe2 * dtsed2 / 2.0 ) 
     727            ELSE 
     728               zsedtra(4) = ( zsedtra(4) * zalpha ) / ( 0.25 * zsedtra(4) * ( exp( reac_fe2 * zalpha * dtsed2 / 2. ) - 1.0 )  & 
     729               &            + zalpha * exp( reac_fe2 * zalpha * dtsed2 / 2. ) ) 
     730            ENDIF 
    687731            zsedtra(1) = zalpha + 0.25 * zsedtra(4) 
    688732            zsedtra(5) = zbeta  - zsedtra(4) 
  • NEMO/branches/2018/dev_r10164_HPC09_ESIWACE_PREP_MERGE/src/TOP/PISCES/SED/seddta.F90

    r10345 r10368  
    4949 
    5050      REAL(wp), DIMENSION(jpoce) :: zdtap, zdtag 
    51       REAL(wp), DIMENSION(jpi,jpj) :: zwsbio4, zwsbio3, zwscal 
     51      REAL(wp), DIMENSION(jpi,jpj) :: zwsbio4, zwsbio3 
    5252      REAL(wp) :: zf0, zf1, zf2, zkapp, zratio, zdep 
    5353 
     
    9797               zwsbio4(ji,jj) = wsbio2 / rday 
    9898               zwsbio3(ji,jj) = wsbio  / rday 
    99                zwscal (ji,jj) = wsbio2 / rday 
    10099            END DO 
    101100         END DO 
     
    106105               zdep = e3t_n(ji,jj,ikt) / r2dttrc 
    107106               zwsbio4(ji,jj) = MIN( 0.99 * zdep, wsbio4(ji,jj,ikt) / rday ) 
    108                zwscal (ji,jj) = MIN( 0.99 * zdep, wscal (ji,jj,ikt) / rday ) 
    109107               zwsbio3(ji,jj) = MIN( 0.99 * zdep, wsbio3(ji,jj,ikt) / rday ) 
    110108            END DO 
     
    130128               trc_data(ji,jj,12 ) = MIN(trb(ji,jj,ikt,jppoc), 1E-4) * zwsbio3(ji,jj) * 1E3 
    131129               trc_data(ji,jj,13 ) = MIN(trb(ji,jj,ikt,jpgoc), 1E-4) * zwsbio4(ji,jj) * 1E3 
    132                trc_data(ji,jj,14)  = MIN(trb(ji,jj,ikt,jpcal), 1E-4) * zwscal (ji,jj) * 1E3 
     130               trc_data(ji,jj,14)  = MIN(trb(ji,jj,ikt,jpcal), 1E-4) * zwsbio4(ji,jj) * 1E3 
    133131               trc_data(ji,jj,15)  = tsn(ji,jj,ikt,jp_tem) 
    134132               trc_data(ji,jj,16)  = tsn(ji,jj,ikt,jp_sal) 
  • NEMO/branches/2018/dev_r10164_HPC09_ESIWACE_PREP_MERGE/src/TOP/PISCES/SED/sedini.F90

    r10345 r10368  
    472472      ENDIF 
    473473 
     474      IF ( ln_p5z .AND. ln_sed_2way ) CALL ctl_stop( '2 ways coupling with sediment cannot be activated with PISCES-QUOTA' ) 
     475 
    474476      REWIND( numnamsed_ref )              ! Namelist nam_geom in reference namelist : Pisces variables 
    475477      READ  ( numnamsed_ref, nam_geom, IOSTAT = ios, ERR = 903) 
  • NEMO/branches/2018/dev_r10164_HPC09_ESIWACE_PREP_MERGE/src/TOP/PISCES/sms_pisces.F90

    r10345 r10368  
    7676   REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:,:) ::  etot_ndcy      !: PAR over 24h in case of diurnal cycle 
    7777   REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:,:) ::  enano, ediat   !: PAR for phyto, nano and diat  
     78   REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:,:) ::  enanom, ediatm !: PAR for phyto, nano and diat  
    7879   REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:,:) ::  epico          !: PAR for pico 
     80   REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:,:) ::  epicom         !: PAR for pico 
    7981   REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:,:) ::  emoy           !: averaged PAR in the mixed layer 
    8082   REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:)   ::  heup_01 !: Absolute euphotic layer depth 
     
    8991   REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:,:) ::   wsbio3   !: POC sinking speed  
    9092   REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:,:) ::   wsbio4   !: GOC sinking speed 
    91    REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:,:) ::   wscal    !: Calcite and BSi sinking speeds 
    9293   REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:,:) ::   wsfep 
    9394 
     
    105106   REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:,:) ::   prodgoc    !: Calcite production 
    106107   REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:,:) ::   consgoc    !: Calcite production 
    107  
     108   REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:,:) ::   blim       !: bacterial production factor 
    108109   REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:,:) ::   sizen      !: size of diatoms  
    109110   REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:,:) ::   sizep      !: size of diatoms  
     
    147148         !*  Biological fluxes for light  
    148149         ALLOCATE(  enano(jpi,jpj,jpk)    , ediat(jpi,jpj,jpk) ,   & 
     150           &        enanom(jpi,jpj,jpk)   , ediatm(jpi,jpj,jpk),   & 
    149151           &        etot_ndcy(jpi,jpj,jpk), emoy(jpi,jpj,jpk)  ,  STAT=ierr(2) )  
    150152 
     
    157159            &      prodcal(jpi,jpj,jpk) , xdiss   (jpi,jpj,jpk),    & 
    158160            &      prodpoc(jpi,jpj,jpk) , conspoc(jpi,jpj,jpk) ,    & 
    159             &      prodgoc(jpi,jpj,jpk) , consgoc(jpi,jpj,jpk) ,  STAT=ierr(4) ) 
     161            &      prodgoc(jpi,jpj,jpk) , consgoc(jpi,jpj,jpk) ,    & 
     162            &      blim   (jpi,jpj,jpk) ,                         STAT=ierr(4) ) 
    160163 
    161164         !* Variable for chemistry of the CO2 cycle 
     
    170173         !* Sinkong speed 
    171174         ALLOCATE( wsbio3 (jpi,jpj,jpk) , wsbio4 (jpi,jpj,jpk),     & 
    172             &      wscal(jpi,jpj,jpk)                         ,   STAT=ierr(7) )    
     175            &                             STAT=ierr(7) )    
    173176         !  
    174177         IF( ln_ligand ) THEN 
     
    180183      IF( ln_p5z ) THEN 
    181184         !        
    182          ALLOCATE( epico(jpi,jpj,jpk)                         ,   STAT=ierr(9) )  
     185         ALLOCATE( epico(jpi,jpj,jpk)   , epicom(jpi,jpj,jpk) ,   STAT=ierr(9) )  
    183186 
    184187         !*  Size of phytoplankton cells 
  • NEMO/branches/2018/dev_r10164_HPC09_ESIWACE_PREP_MERGE/src/TOP/PISCES/trcini_pisces.F90

    r10345 r10368  
    109109      ierr = ierr +  p4z_flx_alloc() 
    110110      ierr = ierr +  p4z_sed_alloc() 
    111       ierr = ierr +  p4z_rem_alloc() 
     111      ierr = ierr +  p4z_lim_alloc() 
    112112      IF( ln_p4z ) THEN 
    113          ierr = ierr +  p4z_lim_alloc() 
    114113         ierr = ierr +  p4z_prod_alloc() 
    115114      ELSE 
     
    117116         ierr = ierr +  p5z_prod_alloc() 
    118117      ENDIF 
     118      ierr = ierr +  p4z_rem_alloc() 
    119119      ! 
    120120      IF( lk_mpp    )   CALL mpp_sum( 'trcini_pisces', ierr ) 
     
    165165        IF( cltra == 'PICN'     )   jpnpi = jn      !: Picophytoplankton N biomass 
    166166        IF( cltra == 'PICP'     )   jpppi = jn      !: Picophytoplankton P biomass 
     167        IF( cltra == 'PCHL'     )   jppch = jn      !: Diatoms Chlorophyll Concentration 
    167168        IF( cltra == 'PFe'      )   jppfe = jn      !: Picophytoplankton Fe biomass 
    168169        IF( cltra == 'LGW'      )   jplgw = jn      !: Weak ligands 
     
    238239         xksi(:,:)    = 2.e-6 
    239240         xksimax(:,:) = xksi(:,:) 
     241         IF( ln_p5z ) THEN 
     242            sized(:,:,:) = 1.0 
     243            sizen(:,:,:) = 1.0 
     244            sized(:,:,:) = 1.0 
     245         ENDIF 
    240246      END IF 
    241247 
  • NEMO/branches/2018/dev_r10164_HPC09_ESIWACE_PREP_MERGE/src/TOP/TRP/trcdmp.F90

    r10345 r10368  
    122122                     DO jj = 2, jpjm1 
    123123                        DO ji = fs_2, fs_jpim1   ! vector opt. 
    124                            IF( avs(ji,jj,jk) <= 5.e-4_wp )  THEN  
     124                           IF( avt(ji,jj,jk) <= avt_c )  THEN  
    125125                              tra(ji,jj,jk,jn) = tra(ji,jj,jk,jn) + restotr(ji,jj,jk) * ( ztrcdta(ji,jj,jk) - trb(ji,jj,jk,jn) ) 
    126126                           ENDIF 
  • NEMO/branches/2018/dev_r10164_HPC09_ESIWACE_PREP_MERGE/src/TOP/oce_trc.F90

    r10097 r10368  
    8787 
    8888   !* vertical diffusion * 
    89    USE zdf_oce , ONLY :   avs        =>   avs         !: vert. diffusivity coef. for salinity (w-point) 
     89   USE zdf_oce , ONLY :   avs        =>   avs         !: vert. diffusivity coef. for salinity    (w-point) 
     90   USE zdf_oce , ONLY :   avt        =>   avt         !: vert. diffusivity coef. for temperature (w-point) 
    9091 
    9192   !* mixing & mixed layer depth * 
     
    9495   USE zdfmxl , ONLY :   hmlp        =>   hmlp        !: mixed layer depth  (rho=rho0+zdcrit) (m) 
    9596   USE zdfmxl , ONLY :   hmlpt       =>   hmlpt       !: mixed layer depth at t-points (m) 
     97   USE zdfmxl , ONLY :   avt_c       =>   avt_c       !: Kz criterion for the turbocline depth 
    9698 
    9799END MODULE oce_trc 
Note: See TracChangeset for help on using the changeset viewer.