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 2457 – NEMO

Changeset 2457


Ignore:
Timestamp:
2010-12-07T10:51:47+01:00 (13 years ago)
Author:
cetlod
Message:

Improve TOP & OFF components in v3.3beta, see ticket #774

Location:
branches/nemo_v3_3_beta/NEMOGCM/NEMO
Files:
18 edited

Legend:

Unmodified
Added
Removed

  • branches/nemo_v3_3_beta/NEMOGCM/NEMO/OFF_SRC/domrea.F90

    r2444 r2457  
    1313   USE dom_oce         ! ocean space and time domain 
    1414   USE dommsk          ! domain: masks 
     15   USE lbclnk          !  
    1516   USE in_out_manager  ! I/O manager 
    1617 
    … …  
    6970      INTEGER  ::   ik, inum0 , inum1 , inum2 , inum3 , inum4   ! local integers 
    7071      REAL(wp) ::   zrefdep         ! local real 
    71       REAL(wp), DIMENSION(jpi,jpj) ::   zprt   ! 2D workspace 
     72      REAL(wp), DIMENSION(jpi,jpj) ::   zmbk   ! 2D workspace 
    7273      !!---------------------------------------------------------------------- 
    7374 
    … …  
    7677      IF(lwp) WRITE(numout,*) '~~~~~~~' 
    7778 
    78       zprt(:,:) = 0._wp 
     79      zmbk(:,:) = 0._wp 
    7980 
    8081      SELECT CASE (nmsh) 
    … …  
    156157         CALL iom_get( inum3, jpdom_data, 'ff', ff ) 
    157158 
    158          CALL iom_get( inum4, jpdom_data, 'mbathy', zprt ) 
     159         CALL iom_get( inum4, jpdom_data, 'mbathy', zmbk ) 
    159160      
    160          DO jj = 1, jpj 
    161             DO ji = 1, jpi 
    162                mbathy(ji,jj) = MAX( zprt(ji,jj) * tmask(ji,jj,1), 1._wp ) + 1 
    163             ENDDO 
    164          ENDDO 
    165  
     161         ! 
     162         mbkt(:,:) = MAX( zmbk(:,:) * tmask(:,:,1), 1._wp )    ! bottom k-index of T-level (=1 over land) 
     163         !   
     164         DO jj = 1, jpjm1                      ! bottom k-index of u- (v-) level 
     165            DO ji = 1, jpim1 
     166               mbku(ji,jj) = MIN(  mbkt(ji+1,jj  ) , mbkt(ji,jj)  ) 
     167               mbkv(ji,jj) = MIN(  mbkt(ji  ,jj+1) , mbkt(ji,jj)  ) 
     168            END DO 
     169         END DO 
     170         ! converte into REAL to use lbc_lnk ; impose a min value of 1 as a zero can be set in lbclnk  
     171         zmbk(:,:) = REAL( mbku(:,:), wp )   ;   CALL lbc_lnk(zmbk,'U',1.)   ;   mbku  (:,:) = MAX( INT( zmbk(:,:) ), 1 ) 
     172         zmbk(:,:) = REAL( mbkv(:,:), wp )   ;   CALL lbc_lnk(zmbk,'V',1.)   ;   mbkv  (:,:) = MAX( INT( zmbk(:,:) ), 1 ) 
     173         !                                                     ! bottom k-index of W-level = mbkt+1 
     174         ! 
    166175         IF( ln_sco ) THEN                                         ! s-coordinate 
    167176            CALL iom_get( inum4, jpdom_data, 'hbatt', hbatt ) 
    … …  
    217226              DO jj = 1, jpj 
    218227                DO ji = 1, jpi 
    219                   ik = mbathy(ji,jj) - 1 
     228                  ik = mbkt(ji,jj) 
    220229                  ! ocean point only  
    221230                  IF( ik > 0 ) THEN 

  • branches/nemo_v3_3_beta/NEMOGCM/NEMO/OFF_SRC/opa.F90

    r2444 r2457  
    2525   USE traqsr          ! solar radiation penetration    (tra_qsr_init routine) 
    2626   USE trabbl          ! bottom boundary layer          (tra_bbl_init routine) 
    27    USE zpshde          ! partial step: hor. derivative  (zps_hde_init routine) 
    2827   USE zdfini          ! vertical physics: initialization 
    2928   USE phycst          ! physical constant                  (par_cst routine) 
    … …  
    180179                            CALL     dom_cfg    ! Domain configuration 
    181180                            CALL     dom_init   ! Domain 
    182  
    183       IF( ln_zps        )   CALL zps_hde_init   ! Partial steps:  horizontal derivative 
    184181                            CALL  istate_init   ! ocean initial state (Dynamics and tracers) 
    185182 

  • branches/nemo_v3_3_beta/NEMOGCM/NEMO/TOP_SRC/LOBSTER/trcexp.F90

    r2300 r2457  
    5555      INTEGER, INTENT( in ) ::   kt      ! ocean time-step index       
    5656      !! 
    57       INTEGER  ::   ji, jj, jk, jl, ikbot 
     57      INTEGER  ::   ji, jj, jk, jl, ikt 
    5858      REAL(wp) ::   zgeolpoc, zfact, zwork, ze3t, zsedpocd 
    5959      REAL(wp), DIMENSION(:,:,:), ALLOCATABLE ::   ztrbio 
    … …  
    9595      DO jj = 2, jpjm1 
    9696         DO ji = fs_2, fs_jpim1 
    97             ikbot = mbathy(ji,jj) - 1 
    98             tra(ji,jj,ikbot,jp_lob_no3) = tra(ji,jj,ikbot,jp_lob_no3) + sedlam * sedpocn(ji,jj) / fse3t(ji,jj,ikbot)  
     97            ikt = mbkt(ji,jj)  
     98            tra(ji,jj,ikt,jp_lob_no3) = tra(ji,jj,ikt,jp_lob_no3) + sedlam * sedpocn(ji,jj) / fse3t(ji,jj,ikt)  
    9999            ! Deposition of organic matter in the sediment 
    100             zwork = vsed * trn(ji,jj,ikbot,jp_lob_det) 
     100            zwork = vsed * trn(ji,jj,ikt,jp_lob_det) 
    101101            sedpoca(ji,jj) = ( zwork + dminl(ji,jj) * fbod(ji,jj)   & 
    102102               &           - sedlam * sedpocn(ji,jj) - sedlostpoc * sedpocn(ji,jj) ) * rdt 

  • branches/nemo_v3_3_beta/NEMOGCM/NEMO/TOP_SRC/LOBSTER/trcini_lobster.F90

    r2287 r2457  
    2222   USE lbclnk  
    2323   USE lib_mpp  
     24   USE lib_fortran  
    2425 
    2526   IMPLICIT NONE 
    … …  
    137138      ! Coastal surface 
    138139      ! --------------- 
    139       areacot = 0.e0 
    140       DO ji = 1, jpi 
    141          DO jj = 1, jpj 
    142             areacot = areacot + e1t(ji,jj) * e2t(ji,jj) * cmask(ji,jj) 
    143          END DO 
    144       END DO 
    145       ! 
    146       IF( lk_mpp ) CALL mpp_sum( areacot )   ! sum over the global domain 
     140      areacot = glob_sum( e1t(:,:) * e2t(:,:) * cmask(:,:) ) 
    147141 
    148142      ! Initialization of tracer concentration in case of  no restart  

  • branches/nemo_v3_3_beta/NEMOGCM/NEMO/TOP_SRC/PISCES/p4zflx.F90

    r2287 r2457  
    2828#endif 
    2929   USE lib_mpp 
     30   USE lib_fortran 
    3031 
    3132   IMPLICIT NONE 
    … …  
    3536   PUBLIC   p4z_flx_init   
    3637 
    37    REAL(wp) :: &  ! pre-industrial atmospheric [co2] (ppm)   
    38       atcox  = 0.20946 ,    &  !: 
    39       atcco2 = 278.            !: 
    40  
    41    REAL(wp) :: & 
    42       xconv  = 0.01/3600      !: coefficients for conversion  
    43  
    44    INTEGER  ::  nspyr         !: number of timestep per year 
    45  
    46 #if defined key_cpl_carbon_cycle 
    47    REAL(wp), PUBLIC, DIMENSION(jpi,jpj) ::  & 
    48       oce_co2            !: ocean carbon flux 
    49    REAL(wp) :: & 
    50       t_atm_co2_flx,  &  !: Total atmospheric carbon flux per year 
    51       t_oce_co2_flx      !: Total ocean carbon flux per year 
    52 #endif 
     38   REAL(wp), PUBLIC, DIMENSION(jpi,jpj) ::  oce_co2            !: ocean carbon flux  
     39   REAL(wp), PUBLIC, DIMENSION(jpi,jpj) ::  satmco2            !: atmospheric pco2 
     40   REAL(wp)                             ::  t_oce_co2_flx      !: Total ocean carbon flux  
     41   REAL(wp)                             ::  t_atm_co2_flx      !: global mean of atmospheric pco2 
     42   REAL(wp)                             ::  area               !: ocean surface 
     43   REAL(wp)                             ::  atcco2 = 278.      !: pre-industrial atmospheric [co2] (ppm)     
     44   REAL(wp)                             ::  atcox  = 0.20946   !: 
     45   REAL(wp)                             ::  xconv  = 0.01/3600 !: coefficients for conversion  
    5346 
    5447   !!* Substitution 
    … …  
    7770      REAL(wp), DIMENSION(jpi,jpj) ::   zkgco2, zkgo2, zh2co3 
    7871#if defined key_diatrc && defined key_iomput 
    79       REAL(wp), DIMENSION(jpi,jpj) ::  zcflx, zoflx, zkg, zdpco2, zdpo2 
     72      REAL(wp), DIMENSION(jpi,jpj) ::  zoflx, zkg, zdpco2, zdpo2 
    8073#endif 
    8174      CHARACTER (len=25) :: charout 
    … …  
    8679      !     SURFACE LAYER); THE RESULT OF THIS CALCULATION 
    8780      !     IS USED TO COMPUTE AIR-SEA FLUX OF CO2 
     81 
     82#if defined key_cpl_carbon_cycle 
     83      satmco2(:,:) = atm_co2(:,:) 
     84#endif 
    8885 
    8986      DO jrorr = 1, 10 
    … …  
    150147         DO ji = 1, jpi 
    151148            ! Compute CO2 flux for the sea and air 
    152 #if ! defined key_cpl_carbon_cycle 
    153             zfld = atcco2 * tmask(ji,jj,1) * chemc(ji,jj,1) * zkgco2(ji,jj) 
     149            zfld = satmco2(ji,jj) * tmask(ji,jj,1) * chemc(ji,jj,1) * zkgco2(ji,jj) 
    154150            zflu = zh2co3(ji,jj) * tmask(ji,jj,1) * zkgco2(ji,jj) 
    155 #else 
    156             zfld = atm_co2(ji,jj) * tmask(ji,jj,1) * chemc(ji,jj,1) * zkgco2(ji,jj) 
    157             zflu = zh2co3(ji,jj) * tmask(ji,jj,1) * zkgco2(ji,jj) 
    158             ! compute flux of carbon 
    159151            oce_co2(ji,jj) = ( zfld - zflu ) * rfact & 
    160152               &             * e1t(ji,jj) * e2t(ji,jj) * tmask(ji,jj,1) * 1000. 
    161 #endif 
     153            ! compute the trend 
    162154            tra(ji,jj,1,jpdic) = tra(ji,jj,1,jpdic) + ( zfld - zflu ) / fse3t(ji,jj,1) 
    163155 
    … …  
    170162            ! Save diagnostics 
    171163#  if ! defined key_iomput 
    172             trc2d(ji,jj,jp_pcs0_2d    ) = ( zfld - zflu )     * 1000. * tmask(ji,jj,1) 
     164            zfact = 1. / ( e1t(ji,jj) * e2t(ji,jj) ) / rfact 
     165            trc2d(ji,jj,jp_pcs0_2d    ) = oce_co2(ji,jj) * zfact 
    173166            trc2d(ji,jj,jp_pcs0_2d + 1) = ( zfld16 - zflu16 ) * 1000. * tmask(ji,jj,1) 
    174167            trc2d(ji,jj,jp_pcs0_2d + 2) = zkgco2(ji,jj) * tmask(ji,jj,1) 
    175             trc2d(ji,jj,jp_pcs0_2d + 3) = ( atcco2 - zh2co3(ji,jj) / ( chemc(ji,jj,1) + rtrn ) ) & 
     168            trc2d(ji,jj,jp_pcs0_2d + 3) = ( satmco2(ji,jj) - zh2co3(ji,jj) / ( chemc(ji,jj,1) + rtrn ) ) & 
    176169               &                            * tmask(ji,jj,1) 
    177170#  else 
    178             zcflx(ji,jj) = ( zfld - zflu ) * 1000.  * tmask(ji,jj,1) 
    179171            zoflx(ji,jj) = ( zfld16 - zflu16 ) * 1000. * tmask(ji,jj,1) 
    180172            zkg  (ji,jj) = zkgco2(ji,jj) * tmask(ji,jj,1) 
    181             zdpco2(ji,jj) = ( atcco2 - zh2co3(ji,jj)      / ( chemc(ji,jj,1) + rtrn ) ) & 
    182               &             * tmask(ji,jj,1) 
    183             zdpo2 (ji,jj) = ( atcox  - trn(ji,jj,1,jpoxy) / ( chemc(ji,jj,2) + rtrn ) ) & 
    184               &             * tmask(ji,jj,1) 
     173            zdpco2(ji,jj) = ( satmco2(ji,jj) - zh2co3(ji,jj) / ( chemc(ji,jj,1) + rtrn ) ) * tmask(ji,jj,1) 
     174            zdpo2 (ji,jj) = ( atcox  - trn(ji,jj,1,jpoxy) / ( chemc(ji,jj,2) + rtrn ) ) * tmask(ji,jj,1) 
    185175#  endif 
    186176#endif 
    … …  
    188178      END DO 
    189179 
    190 #if defined key_cpl_carbon_cycle 
    191       ! Total Flux of Carbon 
    192       DO jj = 1, jpj  
    193         DO ji = 1, jpi 
    194            t_atm_co2_flx = t_atm_co2_flx + atm_co2(ji,jj) * tmask_i(ji,jj) 
    195            t_oce_co2_flx = t_oce_co2_flx + oce_co2(ji,jj) * tmask_i(ji,jj) 
    196         END DO 
    197       END DO 
    198  
    199       IF( MOD( kt, nspyr ) == 0 ) THEN 
    200         IF( lk_mpp ) THEN 
    201           CALL mpp_sum( t_atm_co2_flx )   ! sum over the global domain 
    202           CALL mpp_sum( t_oce_co2_flx )   ! sum over the global domain 
    203         ENDIF 
    204         ! Conversion in GtC/yr ; negative for outgoing from ocean 
    205         t_oce_co2_flx = (-1.) * t_oce_co2_flx  * 12. / 1.e15 
    206         ! 
    207         WRITE(numout,*) ' Atmospheric pCO2    :' 
    208         WRITE(numout,*) '-------------------- : ',kt,'  ',t_atm_co2_flx 
    209         WRITE(numout,*) '(ppm)' 
    210         WRITE(numout,*) 'Total Flux of Carbon out of the ocean :' 
    211         WRITE(numout,*) '-------------------- : ',t_oce_co2_flx 
    212         WRITE(numout,*) '(GtC/yr)' 
    213         t_atm_co2_flx = 0. 
    214         t_oce_co2_flx = 0. 
    215 # if defined key_iomput 
    216         CALL iom_put( "tatpco2" , t_atm_co2_flx  ) 
    217         CALL iom_put( "tco2flx" , t_oce_co2_flx  ) 
    218 #endif 
     180      t_oce_co2_flx = t_oce_co2_flx + glob_sum( oce_co2(:,:) )                     ! Cumulative Total Flux of Carbon 
     181      IF( kt == nitend ) THEN 
     182         t_atm_co2_flx = glob_sum( satmco2(:,:) * e1t(:,:) * e2t(:,:) )            ! Total atmospheric pCO2 
     183         ! 
     184         t_oce_co2_flx = (-1.) * t_oce_co2_flx  * 12. / 1.e15                      ! Conversion in PgC ; negative for out of the ocean 
     185         t_atm_co2_flx = t_atm_co2_flx  / area                                     ! global mean of atmospheric pCO2 
     186         ! 
     187         IF( lwp) THEN 
     188            WRITE(numout,*) 
     189            WRITE(numout,*) ' Global mean of atmospheric pCO2 (ppm) at it= ', kt, ' date= ', ndastp 
     190            WRITE(numout,*) '------------------------------------------------------- :  ',t_atm_co2_flx 
     191            WRITE(numout,*) 
     192            WRITE(numout,*) ' Cumulative total Flux of Carbon out of the ocean (PgC) :' 
     193            WRITE(numout,*) '-------------------------------------------------------  ',t_oce_co2_flx 
     194         ENDIF 
     195         ! 
    219196      ENDIF 
    220 #endif 
    221197 
    222198      IF(ln_ctl)   THEN  ! print mean trends (used for debugging) 
    … …  
    227203 
    228204# if defined key_diatrc && defined key_iomput 
    229       CALL iom_put( "Cflx" , zcflx  ) 
     205      CALL iom_put( "Cflx" , oce_co2(:,:) / ( e1t(:,:) * e2t(:,:) ) / rfact  ) 
    230206      CALL iom_put( "Oflx" , zoflx  ) 
    231207      CALL iom_put( "Kg"   , zkg    ) 
    … …  
    261237      ENDIF 
    262238 
    263       ! number of time step per year   
    264       nspyr = INT( nyear_len(1) * rday / rdt ) 
    265  
    266 #if defined key_cpl_carbon_cycle 
     239      ! interior global domain surface 
     240      area = glob_sum( e1t(:,:) * e2t(:,:) )   
     241 
    267242      ! Initialization of Flux of Carbon 
    268       oce_co2(:,:) = 0. 
    269       t_atm_co2_flx = 0. 
    270       t_oce_co2_flx = 0. 
    271 #endif 
     243      oce_co2(:,:)  = 0._wp 
     244      t_atm_co2_flx = 0._wp 
     245      ! Initialisation of atmospheric pco2 
     246      satmco2(:,:)  = atcco2 
     247      t_oce_co2_flx = 0._wp 
    272248 
    273249   END SUBROUTINE p4z_flx_init 

  • branches/nemo_v3_3_beta/NEMOGCM/NEMO/TOP_SRC/PISCES/p4zprod.F90

    r2287 r2457  
    2323 
    2424   USE lib_mpp 
     25   USE lib_fortran 
    2526 
    2627   IMPLICIT NONE 
    … …  
    4950      texcret2                   ,  &  !: 1 - excret2         
    5051      tpp                              !: Total primary production 
    51  
    52    INTEGER  ::  nspyr                  !: number of timesteps per year 
    5352 
    5453   !!* Substitution 
    … …  
    326325 
    327326     ! Total primary production per year 
    328      DO jk = 1, jpkm1 
    329         DO jj = 1, jpj 
    330           DO ji = 1, jpi 
    331              zvol = cvol(ji,jj,jk) 
     327 
    332328#if defined key_degrad 
    333              zvol = zvol * facvol(ji,jj,jk) 
     329     tpp = tpp + glob_sum( ( zprorca(:,:,:) + zprorcad(:,:,:) ) * cvol(:,:,:) * facvol(:,:,:) ) 
     330#else 
     331     tpp = tpp + glob_sum( ( zprorca(:,:,:) + zprorcad(:,:,:) ) * cvol(:,:,:) ) 
    334332#endif 
    335              tpp  = tpp + ( zprorca(ji,jj,jk) + zprorcad(ji,jj,jk) ) & 
    336                           * zvol * tmask(ji,jj,jk) * tmask_i(ji,jj) 
    337           END DO 
    338         END DO 
    339       END DO 
    340  
    341  
    342       IF( MOD( kt, nspyr ) == 0 .AND. jnt == nrdttrc ) THEN 
    343         IF( lk_mpp ) CALL mpp_sum( tpp ) 
    344         WRITE(numout,*) 'Total PP :' 
     333 
     334     IF( kt == nitend .AND. jnt == nrdttrc ) THEN 
     335        WRITE(numout,*) 'Total PP (Gtc) :' 
    345336        WRITE(numout,*) '-------------------- : ',tpp * 12. / 1.E12 
    346         WRITE(numout,*) '(GtC/yr)' 
    347         tpp = 0. 
     337        WRITE(numout,*)  
    348338      ENDIF 
    349339 
    … …  
    418408      ENDIF 
    419409 
    420       ! number of timesteps per year 
    421       nspyr  = INT( nyear_len(1) * rday / rdt ) 
    422  
    423410      texcret   = 1.0 - excret 
    424411      texcret2  = 1.0 - excret2 

  • branches/nemo_v3_3_beta/NEMOGCM/NEMO/TOP_SRC/PISCES/p4zrem.F90

    r2287 r2457  
    6565      REAL(wp) ::   zremip, zremik , zlam1b 
    6666      REAL(wp) ::   zkeq  , zfeequi, zsiremin 
    67       REAL(wp) ::   zsatur, zsatur1, zsatur2, zsatur22, znusil 
    68       REAL(wp) ::   ztem1, ztem2 
     67      REAL(wp) ::   zsatur, zsatur2, znusil 
    6968      REAL(wp) ::   zbactfer, zorem, zorem2, zofer 
    7069      REAL(wp) ::   zosil, zdenom1, zscave, zaggdfe 

  • branches/nemo_v3_3_beta/NEMOGCM/NEMO/TOP_SRC/PISCES/p4zsed.F90

    r2403 r2457  
    1919   USE sms_pisces 
    2020   USE lib_mpp 
     21   USE lib_fortran 
    2122   USE prtctl_trc 
    2223   USE p4zbio 
    … …  
    4849 
    4950   !! * Module variables 
    50    INTEGER ::                   & 
    51      ryyss,                     &  !: number of seconds per year 
    52      rmtss                         !: number of seconds per month 
    53  
     51   REAL(wp) :: ryyss               !: number of seconds per year  
     52   REAL(wp) :: ryyss1              !: inverse of ryyss 
     53   REAL(wp) :: rmtss               !: number of seconds per month 
     54   REAL(wp) :: rday1               !: inverse of rday 
     55 
     56   INTEGER , PARAMETER :: & 
     57        jpmth = 12, jpyr = 1 
    5458   INTEGER ::                   & 
    5559      numdust,                  &  !: logical unit for surface fluxes data 
    5660      nflx1 , nflx2,            &  !: first and second record used 
    5761      nflx11, nflx12      ! ??? 
    58    REAL(wp), DIMENSION(jpi,jpj,2) ::    &  !: 
    59      dustmo                                !: 2 consecutive set of dust fields  
    60    REAL(wp), DIMENSION(jpi,jpj)   ::    & 
    61      rivinp, cotdep, nitdep, dust 
    62    REAL(wp), DIMENSION(jpi,jpj,jpk)  ::   & 
    63      ironsed 
     62   REAL(wp), DIMENSION(jpi,jpj,jpmth) ::  dustmo    !: set of dust fields 
     63   REAL(wp), DIMENSION(jpi,jpj)      ::  rivinp, cotdep, nitdep, dust  
     64   REAL(wp), DIMENSION(jpi,jpj)      ::  e1e2t 
     65   REAL(wp), DIMENSION(jpi,jpj,jpk)  ::  ironsed  
    6466   REAL(wp) :: sumdepsi, rivalkinput, rivpo4input, nitdepinput 
    6567 
    … …  
    7476CONTAINS 
    7577 
    76    SUBROUTINE p4z_sed(kt, jnt) 
     78   SUBROUTINE p4z_sed( kt, jnt ) 
    7779      !!--------------------------------------------------------------------- 
    7880      !!                     ***  ROUTINE p4z_sed  *** 
    … …  
    8587      !!--------------------------------------------------------------------- 
    8688      INTEGER, INTENT(in) ::   kt, jnt ! ocean time step 
    87       INTEGER  ::   ji, jj, jk 
    88       INTEGER  ::   ikt 
     89      INTEGER  ::   ji, jj, jk, ikt 
    8990#if ! defined key_sed 
    9091      REAL(wp) ::   zsumsedsi, zsumsedpo4, zsumsedcal 
     92      REAL(wp) ::   zrivalk, zrivsil, zrivpo4 
    9193#endif 
    92       REAL(wp) ::   zconctmp , zdenitot  , znitrpottot 
    93       REAL(wp) ::   zlim, zconctmp2, zfact, zrivalk 
     94      REAL(wp) ::   zdenitot, znitrpottot, zlim, zfact 
     95      REAL(wp) ::   zwsbio3, zwsbio4, zwscal 
    9496      REAL(wp), DIMENSION(jpi,jpj)     ::   zsidep 
     97      REAL(wp), DIMENSION(jpi,jpj)     ::   zwork, zwork1 
    9598      REAL(wp), DIMENSION(jpi,jpj,jpk) ::   znitrpot, zirondep 
    96 #if defined key_diatrc  
    97       REAL(wp) :: zrfact2 
    98 # if defined key_iomput 
    99      REAL(wp), DIMENSION(jpi,jpj)    ::    zw2d  
    100 # endif 
    101 #endif 
    10299      CHARACTER (len=25) :: charout 
    103100      !!--------------------------------------------------------------------- 
    104101 
    105       IF( ( jnt == 1 ) .AND. ( ln_dustfer ) )  CALL p4z_sbc( kt ) 
    106  
    107       zirondep(:,:,:) = 0.e0          ! Initialisation of variables used to compute deposition 
    108       zsidep  (:,:)   = 0.e0 
     102      IF( jnt == 1  .AND.  ln_dustfer  )  CALL p4z_sbc( kt ) 
    109103 
    110104      ! Iron and Si deposition at the surface 
    … …  
    113107      DO jj = 1, jpj 
    114108         DO ji = 1, jpi 
    115             zirondep(ji,jj,1) = ( dustsolub * dust(ji,jj) / ( 55.85 * rmtss ) + 3.e-10 / ryyss )   & 
     109            zirondep(ji,jj,1) = ( dustsolub * dust(ji,jj) / ( 55.85 * rmtss ) + 3.e-10 * ryyss1 )   & 
    116110               &             * rfact2 / fse3t(ji,jj,1) 
    117111            zsidep  (ji,jj)   = 8.8 * 0.075 * dust(ji,jj) * rfact2 / ( fse3t(ji,jj,1) * 28.1 * rmtss ) 
    … …  
    147141 
    148142#if ! defined key_sed 
    149       ! Initialisation of variables used to compute Sinking Speed 
    150       zsumsedsi  = 0.e0 
    151       zsumsedpo4 = 0.e0 
    152       zsumsedcal = 0.e0 
    153  
    154143      ! Loss of biogenic silicon, Caco3 organic carbon in the sediments.  
    155144      ! First, the total loss is computed. 
    … …  
    158147      DO jj = 1, jpj 
    159148         DO ji = 1, jpi 
    160             ikt = MAX( mbathy(ji,jj)-1, 1 ) 
    161             zfact = e1t(ji,jj) * e2t(ji,jj) / rday * tmask_i(ji,jj) 
     149            ikt = mbkt(ji,jj)  
    162150# if defined key_kriest 
    163             zsumsedsi  = zsumsedsi  + zfact * trn(ji,jj,ikt,jpdsi) * wscal (ji,jj,ikt) 
    164             zsumsedpo4 = zsumsedpo4 + zfact * trn(ji,jj,ikt,jppoc) * wsbio3(ji,jj,ikt) 
     151            zwork (ji,jj) = trn(ji,jj,ikt,jpdsi) * wscal (ji,jj,ikt) 
     152            zwork1(ji,jj) = trn(ji,jj,ikt,jppoc) * wsbio3(ji,jj,ikt) 
    165153# else 
    166             zsumsedsi  = zsumsedsi  + zfact *  trn(ji,jj,ikt,jpdsi) * wsbio4(ji,jj,ikt) 
    167             zsumsedpo4 = zsumsedpo4 + zfact *( trn(ji,jj,ikt,jpgoc) * wsbio4(ji,jj,ikt)   & 
    168                &       + trn(ji,jj,ikt,jppoc) * wsbio3(ji,jj,ikt) ) 
     154            zwork (ji,jj) = trn(ji,jj,ikt,jpdsi) * wsbio4(ji,jj,ikt) 
     155            zwork1(ji,jj) = trn(ji,jj,ikt,jpgoc) * wsbio4(ji,jj,ikt) + trn(ji,jj,ikt,jppoc) * wsbio3(ji,jj,ikt)  
    169156# endif 
    170             zsumsedcal = zsumsedcal + zfact *  trn(ji,jj,ikt,jpcal) * wscal (ji,jj,ikt) * 2.e0 
    171          END DO 
    172       END DO 
    173  
    174       IF( lk_mpp ) THEN 
    175          CALL mpp_sum( zsumsedsi  )   ! sums over the global domain 
    176          CALL mpp_sum( zsumsedcal )   ! sums over the global domain 
    177          CALL mpp_sum( zsumsedpo4 )   ! sums over the global domain 
    178       ENDIF 
    179  
     157         END DO 
     158      END DO 
     159      zsumsedsi  = glob_sum( zwork (:,:) * e1e2t(:,:) ) * rday1 
     160      zsumsedpo4 = glob_sum( zwork1(:,:) * e1e2t(:,:) ) * rday1 
     161      DO jj = 1, jpj 
     162         DO ji = 1, jpi 
     163            ikt = mbkt(ji,jj)  
     164            zwork (ji,jj) = trn(ji,jj,ikt,jpcal) * wscal (ji,jj,ikt) 
     165         END DO 
     166      END DO 
     167      zsumsedcal = glob_sum( zwork (:,:) * e1e2t(:,:) ) * 2.0 * rday1 
    180168#endif 
    181169 
    … …  
    188176      DO jj = 1, jpj 
    189177         DO ji = 1, jpi 
    190             ikt = MAX( mbathy(ji,jj) - 1, 1 ) 
    191 # if ! defined key_kriest 
    192             zconctmp = trn(ji,jj,ikt,jpdsi) * xstep / fse3t(ji,jj,ikt) * wscal (ji,jj,ikt)  
     178            ikt = mbkt(ji,jj) 
     179            zfact = xstep / fse3t(ji,jj,ikt) 
     180            zwsbio3 = 1._wp - zfact * wsbio3(ji,jj,ikt) 
     181            zwsbio4 = 1._wp - zfact * wsbio4(ji,jj,ikt) 
     182            zwscal  = 1._wp - zfact * wscal (ji,jj,ikt) 
     183            ! 
     184# if defined key_kriest 
     185            trn(ji,jj,ikt,jpdsi) = trn(ji,jj,ikt,jpdsi) * zwsbio4 
     186            trn(ji,jj,ikt,jpnum) = trn(ji,jj,ikt,jpnum) * zwsbio4 
     187            trn(ji,jj,ikt,jppoc) = trn(ji,jj,ikt,jppoc) * zwsbio3 
     188            trn(ji,jj,ikt,jpsfe) = trn(ji,jj,ikt,jpsfe) * zwsbio3 
    193189# else 
    194             zconctmp = trn(ji,jj,ikt,jpdsi) * xstep / fse3t(ji,jj,ikt) * wsbio4(ji,jj,ikt) 
     190            trn(ji,jj,ikt,jpdsi) = trn(ji,jj,ikt,jpdsi) * zwscal  
     191            trn(ji,jj,ikt,jpgoc) = trn(ji,jj,ikt,jpgoc) * zwsbio4 
     192            trn(ji,jj,ikt,jppoc) = trn(ji,jj,ikt,jppoc) * zwsbio3 
     193            trn(ji,jj,ikt,jpbfe) = trn(ji,jj,ikt,jpbfe) * zwsbio4 
     194            trn(ji,jj,ikt,jpsfe) = trn(ji,jj,ikt,jpsfe) * zwsbio3 
    195195# endif 
    196             trn(ji,jj,ikt,jpdsi) = trn(ji,jj,ikt,jpdsi) - zconctmp 
     196            trn(ji,jj,ikt,jpcal) = trn(ji,jj,ikt,jpcal) * zwscal 
     197         END DO 
     198      END DO 
    197199 
    198200#if ! defined key_sed 
    199             zrivalk = ( 1.- ( sumdepsi + rivalkinput / ryyss / 6. ) / zsumsedsi ) 
    200             trn(ji,jj,ikt,jpsil) = trn(ji,jj,ikt,jpsil) + zconctmp  * zrivalk  
    201 #endif 
    202          END DO 
    203       END DO 
    204  
     201      zrivsil =  1._wp - ( sumdepsi + rivalkinput * ryyss1 / 6. ) / zsumsedsi  
     202      zrivalk =  1._wp - ( rivalkinput * ryyss1 ) / zsumsedcal  
     203      zrivpo4 =  1._wp - ( rivpo4input * ryyss1 ) / zsumsedpo4  
    205204      DO jj = 1, jpj 
    206205         DO ji = 1, jpi 
    207             ikt = MAX( mbathy(ji,jj) - 1, 1 ) 
    208             zconctmp = trn(ji,jj,ikt,jpcal) * wscal(ji,jj,ikt) * xstep / fse3t(ji,jj,ikt) 
    209             trn(ji,jj,ikt,jpcal) = trn(ji,jj,ikt,jpcal) - zconctmp 
    210 #if ! defined key_sed 
    211             zrivalk = ( 1.- ( rivalkinput / ryyss ) / zsumsedcal ) 
    212             trn(ji,jj,ikt,jptal) = trn(ji,jj,ikt,jptal) + zconctmp * zrivalk * 2.0 
    213             trn(ji,jj,ikt,jpdic) = trn(ji,jj,ikt,jpdic) + zconctmp * zrivalk  
    214 #endif 
    215          END DO 
    216       END DO 
    217  
    218       DO jj = 1, jpj 
    219          DO ji = 1, jpi 
    220             ikt = MAX( mbathy(ji,jj) - 1, 1 ) 
     206            ikt = mbkt(ji,jj) 
    221207            zfact = xstep / fse3t(ji,jj,ikt) 
    222 # if ! defined key_kriest 
    223             zconctmp  = trn(ji,jj,ikt,jpgoc) 
    224             zconctmp2 = trn(ji,jj,ikt,jppoc) 
    225             trn(ji,jj,ikt,jpgoc) = trn(ji,jj,ikt,jpgoc) - zconctmp  * wsbio4(ji,jj,ikt) * zfact 
    226             trn(ji,jj,ikt,jppoc) = trn(ji,jj,ikt,jppoc) - zconctmp2 * wsbio3(ji,jj,ikt) * zfact 
    227 #if ! defined key_sed 
    228             trn(ji,jj,ikt,jpdoc) = trn(ji,jj,ikt,jpdoc)    & 
    229             &      + ( zconctmp  * wsbio4(ji,jj,ikt) + zconctmp2 * wsbio3(ji,jj,ikt) ) * zfact   & 
    230             &      * ( 1.- rivpo4input / (ryyss * zsumsedpo4 ) ) 
    231 #endif 
    232             trn(ji,jj,ikt,jpbfe) = trn(ji,jj,ikt,jpbfe) - trn(ji,jj,ikt,jpbfe) * wsbio4(ji,jj,ikt) * zfact 
    233             trn(ji,jj,ikt,jpsfe) = trn(ji,jj,ikt,jpsfe) - trn(ji,jj,ikt,jpsfe) * wsbio3(ji,jj,ikt) * zfact 
    234  
     208            zwsbio3 = zfact * wsbio3(ji,jj,ikt) 
     209            zwsbio4 = zfact * wsbio4(ji,jj,ikt) 
     210            zwscal  = zfact * wscal (ji,jj,ikt) 
     211            trn(ji,jj,ikt,jptal) =  trn(ji,jj,ikt,jptal) + trn(ji,jj,ikt,jpcal) * zwscal  * zrivalk * 2.0 
     212            trn(ji,jj,ikt,jpdic) =  trn(ji,jj,ikt,jpdic) + trn(ji,jj,ikt,jpcal) * zwscal  * zrivalk 
     213# if defined key_kriest 
     214            trn(ji,jj,ikt,jpsil) =  trn(ji,jj,ikt,jpsil) + trn(ji,jj,ikt,jpdsi) * zwsbio4 * zrivsil  
     215            trn(ji,jj,ikt,jpdoc) =  trn(ji,jj,ikt,jpdoc) + trn(ji,jj,ikt,jppoc) * zwsbio3 * zrivpo4  
    235216# else 
    236             zconctmp  = trn(ji,jj,ikt,jpnum) 
    237             zconctmp2 = trn(ji,jj,ikt,jppoc) 
    238             trn(ji,jj,ikt,jpnum) = trn(ji,jj,ikt,jpnum) - zconctmp  * wsbio4(ji,jj,ikt) * zfact  
    239             trn(ji,jj,ikt,jppoc) = trn(ji,jj,ikt,jppoc) - zconctmp2 * wsbio3(ji,jj,ikt) * zfact  
    240 #if ! defined key_sed 
    241             trn(ji,jj,ikt,jpdoc) = trn(ji,jj,ikt,jpdoc) + ( zconctmp2 * wsbio3(ji,jj,ikt) )   
    242             &                     * zfact * ( 1.- rivpo4input / ( ryyss * zsumsedpo4 ) ) 
    243 #endif 
    244             trn(ji,jj,ikt,jpsfe) = trn(ji,jj,ikt,jpsfe) - trn(ji,jj,ikt,jpsfe) * wsbio3(ji,jj,ikt) * zfact  
     217            trn(ji,jj,ikt,jpsil) =  trn(ji,jj,ikt,jpsil) + trn(ji,jj,ikt,jpdsi) * zwscal  * zrivsil  
     218            trn(ji,jj,ikt,jpdoc) =  trn(ji,jj,ikt,jpdoc)   & 
     219            &                     + ( trn(ji,jj,ikt,jppoc) * zwsbio3 + trn(ji,jj,ikt,jpgoc) * zwsbio4 ) * zrivpo4 
    245220# endif 
    246221         END DO 
    247222      END DO 
     223# endif 
    248224 
    249225      ! Nitrogen fixation (simple parameterization). The total gain 
    … …  
    252228      ! ------------------------------------------------------------- 
    253229 
    254       zdenitot = 0.e0 
    255       DO jk = 1, jpkm1 
    256          DO jj = 1,jpj 
    257             DO ji = 1,jpi 
    258                zdenitot = zdenitot + denitr(ji,jj,jk) * rdenit * cvol(ji,jj,jk) * xnegtr(ji,jj,jk) 
    259             END DO 
    260          END DO 
    261       END DO 
    262  
    263       IF( lk_mpp )   CALL mpp_sum( zdenitot )      ! sum over the global domain 
     230      zdenitot = glob_sum( denitr(:,:,:)  * cvol(:,:,:) * xnegtr(:,:,:) ) * rdenit 
    264231 
    265232      ! Potential nitrogen fixation dependant on temperature and iron 
    … …  
    274241               zlim = ( 1.- xnanono3(ji,jj,jk) - xnanonh4(ji,jj,jk) ) 
    275242               IF( zlim <= 0.2 )   zlim = 0.01 
    276                znitrpot(ji,jj,jk) = MAX( 0.e0, ( 0.6 * tgfunc(ji,jj,jk) - 2.15 ) / rday )   & 
     243               znitrpot(ji,jj,jk) = MAX( 0.e0, ( 0.6 * tgfunc(ji,jj,jk) - 2.15 ) * rday1 )   & 
    277244# if defined key_degrad 
    278245               &                  * facvol(ji,jj,jk)   & 
    … …  
    284251      END DO 
    285252 
    286       znitrpottot = 0.e0 
    287       DO jk = 1, jpkm1 
    288          DO jj = 1, jpj 
    289             DO ji = 1, jpi 
    290                znitrpottot = znitrpottot + znitrpot(ji,jj,jk) * cvol(ji,jj,jk) 
    291             END DO 
    292          END DO 
    293       END DO 
    294  
    295       IF( lk_mpp )   CALL mpp_sum( znitrpottot )  ! sum over the global domain 
     253      znitrpottot = glob_sum( znitrpot(:,:,:) * cvol(:,:,:) ) 
    296254 
    297255      ! Nitrogen change due to nitrogen fixation 
    … …  
    301259         DO jj = 1, jpj 
    302260            DO ji = 1, jpi 
    303 # if ! defined key_c1d && ( defined key_orca_r4 || defined key_orca_r2 || defined key_orca_r05 || defined key_orca_r025 ) 
    304 !!             zfact = znitrpot(ji,jj,jk) * zdenitot / znitrpottot 
    305261               zfact = znitrpot(ji,jj,jk) * 1.e-7 
    306 # else 
    307                zfact = znitrpot(ji,jj,jk) * 1.e-7 
    308 # endif 
    309262               trn(ji,jj,jk,jpnh4) = trn(ji,jj,jk,jpnh4) + zfact 
    310263               trn(ji,jj,jk,jpoxy) = trn(ji,jj,jk,jpoxy) + zfact   * o2nit 
    … …  
    315268 
    316269#if defined key_diatrc 
    317       zrfact2 = 1.e+3 * rfact2r 
     270      zfact = 1.e+3 * rfact2r 
    318271#  if  ! defined key_iomput 
    319       trc2d(:,:,jp_pcs0_2d + 11) = zirondep(:,:,1)         * zrfact2 * fse3t(:,:,1) * tmask(:,:,1) 
    320       trc2d(:,:,jp_pcs0_2d + 12) = znitrpot(:,:,1) * 1.e-7 * zrfact2 * fse3t(:,:,1) * tmask(:,:,1) 
    321 # else 
    322       ! surface downward net flux of iron 
    323       zw2d(:,:)   =  ( zirondep(:,:,1) + ironsed(:,:,1) * rfact2 ) * zrfact2 * fse3t(:,:,1) * tmask(:,:,1)  
    324       IF( jnt == nrdttrc ) CALL iom_put( "Irondep", zw2d ) 
    325       ! nitrogen fixation at surface 
    326       zw2d(:,:)   =  znitrpot(:,:,1) * 1.e-7 * zrfact2  * fse3t(:,:,1) * tmask(:,:,1) 
    327       IF( jnt == nrdttrc ) CALL iom_put( "Nfix" , zw2d ) 
    328 # endif 
    329 # endif 
     272      trc2d(:,:,jp_pcs0_2d + 11) = zirondep(:,:,1)         * zfact * fse3t(:,:,1) * tmask(:,:,1) 
     273      trc2d(:,:,jp_pcs0_2d + 12) = znitrpot(:,:,1) * 1.e-7 * zfact * fse3t(:,:,1) * tmask(:,:,1) 
     274#  else 
     275      zwork (:,:)  =  ( zirondep(:,:,1) + ironsed(:,:,1) * rfact2 ) * zfact * fse3t(:,:,1) * tmask(:,:,1)  
     276      zwork1(:,:)  =  znitrpot(:,:,1) * 1.e-7                       * zfact * fse3t(:,:,1) * tmask(:,:,1) 
     277      IF( jnt == nrdttrc ) THEN 
     278         CALL iom_put( "Irondep", zwork  )  ! surface downward net flux of iron 
     279         CALL iom_put( "Nfix"   , zwork1 )  ! nitrogen fixation at surface 
     280      ENDIF 
     281#  endif 
     282#endif 
    330283      ! 
    331284       IF(ln_ctl)   THEN  ! print mean trends (used for debugging) 
    … …  
    337290   END SUBROUTINE p4z_sed 
    338291 
    339    SUBROUTINE p4z_sbc(kt) 
     292   SUBROUTINE p4z_sbc( kt ) 
    340293 
    341294      !!---------------------------------------------------------------------- 
    … …  
    354307 
    355308      !! * Local declarations 
    356       INTEGER ::   & 
    357          imois, imois2,       &  ! temporary integers 
    358          i15  , iman             !    "          " 
    359       REAL(wp) ::   & 
    360          zxy                     !    "         " 
    361  
     309      INTEGER :: imois, i15, iman  
     310      REAL(wp) :: zxy 
    362311 
    363312      !!--------------------------------------------------------------------- 
    … …  
    370319      imois = nmonth + i15 - 1 
    371320      IF( imois == 0 ) imois = iman 
    372       imois2 = nmonth 
    373  
    374       ! 1. first call kt=nit000 
    375       ! ----------------------- 
    376  
    377       IF( kt == nit000 ) THEN 
    378          ! initializations 
    379          nflx1  = 0 
    380          nflx11 = 0 
    381          ! open the file 
    382          IF(lwp) THEN 
    383             WRITE(numout,*) ' ' 
    384             WRITE(numout,*) ' **** Routine p4z_sbc' 
    385          ENDIF 
    386          CALL iom_open ( 'dust.orca.nc', numdust ) 
    387       ENDIF 
    388  
    389  
    390      ! Read monthly file 
    391       ! ---------------- 
    392  
     321 
     322      ! Calendar computation 
    393323      IF( kt == nit000 .OR. imois /= nflx1 ) THEN 
    394324 
    395          ! Calendar computation 
     325         IF( kt == nit000 )  nflx1  = 0 
    396326 
    397327         ! nflx1 number of the first file record used in the simulation 
    … …  
    399329 
    400330         nflx1 = imois 
    401          nflx2 = nflx1+1 
     331         nflx2 = nflx1 + 1 
    402332         nflx1 = MOD( nflx1, iman ) 
    403333         nflx2 = MOD( nflx2, iman ) 
    404334         IF( nflx1 == 0 )   nflx1 = iman 
    405335         IF( nflx2 == 0 )   nflx2 = iman 
    406          IF(lwp) WRITE(numout,*) 'first record file used nflx1 ',nflx1 
    407          IF(lwp) WRITE(numout,*) 'last  record file used nflx2 ',nflx2 
    408  
    409          ! Read monthly fluxes data 
    410  
    411          ! humidity 
    412          CALL iom_get ( numdust, jpdom_data, 'dust', dustmo(:,:,1), nflx1 ) 
    413          CALL iom_get ( numdust, jpdom_data, 'dust', dustmo(:,:,2), nflx2 ) 
     336         IF(lwp) WRITE(numout,*)  
     337         IF(lwp) WRITE(numout,*) ' p4z_sbc : first record file used nflx1 ',nflx1 
     338         IF(lwp) WRITE(numout,*) ' p4z_sbc : last  record file used nflx2 ',nflx2 
    414339 
    415340      ENDIF 
    416341 
    417      ! 3. at every time step interpolation of fluxes 
     342      ! 3. at every time step interpolation of fluxes 
    418343      ! --------------------------------------------- 
    419344 
    420345      zxy = FLOAT( nday + 15 - 30 * i15 ) / 30 
    421       dust(:,:) = ( (1.-zxy) * dustmo(:,:,1) + zxy * dustmo(:,:,2) ) 
    422  
    423       IF( kt == nitend ) CALL iom_close (numdust) 
     346      dust(:,:) = ( (1.-zxy) * dustmo(:,:,nflx1) + zxy * dustmo(:,:,nflx2) ) 
    424347 
    425348   END SUBROUTINE p4z_sbc 
    … …  
    440363      !!---------------------------------------------------------------------- 
    441364 
    442       INTEGER ::   ji, jj, jk, jm 
    443       INTEGER , PARAMETER ::   jpmois = 12, jpan = 1 
     365      INTEGER :: ji, jj, jk, jm 
    444366      INTEGER :: numriv, numbath, numdep 
    445367 
    … …  
    449371      REAL(wp) , DIMENSION (jpi,jpj)     ::   riverdoc, river, ndepo 
    450372      REAL(wp) , DIMENSION (jpi,jpj,jpk) ::   cmask 
    451       REAL(wp) , DIMENSION(jpi,jpj,12)    ::   zdustmo 
    452373 
    453374      NAMELIST/nampissed/ ln_dustfer, ln_river, ln_ndepo, ln_sedinput, sedfeinput, dustsolub 
    … …  
    475396         IF(lwp) WRITE(numout,*) '    ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~ ' 
    476397         CALL iom_open ( 'dust.orca.nc', numdust ) 
    477          DO jm = 1, jpmois 
    478             CALL iom_get( numdust, jpdom_data, 'dust', zdustmo(:,:,jm), jm ) 
     398         DO jm = 1, jpmth 
     399            CALL iom_get( numdust, jpdom_data, 'dust', dustmo(:,:,jm), jm ) 
    479400         END DO 
    480401         CALL iom_close( numdust ) 
    481402      ELSE 
    482          zdustmo(:,:,:) = 0.e0 
     403         dustmo(:,:,:) = 0.e0 
    483404         dust(:,:) = 0.0 
    484405      ENDIF 
    … …  
    490411         IF(lwp) WRITE(numout,*) '    ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~' 
    491412         CALL iom_open ( 'river.orca.nc', numriv ) 
    492          CALL iom_get  ( numriv, jpdom_data, 'riverdic', river   (:,:), jpan ) 
    493          CALL iom_get  ( numriv, jpdom_data, 'riverdoc', riverdoc(:,:), jpan ) 
     413         CALL iom_get  ( numriv, jpdom_data, 'riverdic', river   (:,:), jpyr ) 
     414         CALL iom_get  ( numriv, jpdom_data, 'riverdoc', riverdoc(:,:), jpyr ) 
    494415         CALL iom_close( numriv ) 
    495416      ELSE 
    … …  
    504425         IF(lwp) WRITE(numout,*) '    ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~' 
    505426         CALL iom_open ( 'ndeposition.orca.nc', numdep ) 
    506          CALL iom_get  ( numdep, jpdom_data, 'ndep', ndepo(:,:), jpan ) 
     427         CALL iom_get  ( numdep, jpdom_data, 'ndep', ndepo(:,:), jpyr ) 
    507428         CALL iom_close( numdep ) 
    508429      ELSE 
    … …  
    517438         IF(lwp) WRITE(numout,*) '       from bathy.orca.nc file ' 
    518439         CALL iom_open ( 'bathy.orca.nc', numbath ) 
    519          CALL iom_get  ( numbath, jpdom_data, 'bathy', cmask(:,:,:), jpan ) 
     440         CALL iom_get  ( numbath, jpdom_data, 'bathy', cmask(:,:,:), jpyr ) 
    520441         CALL iom_close( numbath ) 
    521442         ! 
    … …  
    547468 
    548469 
    549       ! Number of seconds per year and per month 
    550       ryyss = nyear_len(1) * rday 
    551       rmtss = ryyss / raamo 
     470      !                                    ! Number of seconds per year and per month 
     471      ryyss  = nyear_len(1) * rday 
     472      rmtss  = ryyss / raamo 
     473      rday1  = 1. / rday 
     474      ryyss1 = 1. / ryyss 
     475      !                                    ! ocean surface cell 
     476      e1e2t(:,:) = e1t(:,:) * e2t(:,:) 
    552477 
    553478      ! total atmospheric supply of Si 
    554479      ! ------------------------------ 
    555480      sumdepsi = 0.e0 
    556       DO jm = 1, jpmois 
    557          DO jj = 2, jpjm1 
    558             DO ji = fs_2, fs_jpim1 
    559                sumdepsi = sumdepsi + zdustmo(ji,jj,jm) / (12.*rmtss) * 8.8        & 
    560                   &     * 0.075/28.1 * e1t(ji,jj) * e2t(ji,jj) * tmask(ji,jj,1) * tmask_i(ji,jj) 
    561             END DO 
    562          END DO 
    563       END DO 
    564       IF( lk_mpp )  CALL mpp_sum( sumdepsi )  ! sum over the global domain 
     481      DO jm = 1, jpmth 
     482         zcoef = 1. / ( 12. * rmtss ) * 8.8 * 0.075 / 28.1         
     483         sumdepsi = sumdepsi + glob_sum( dustmo(:,:,jm) * e1e2t(:,:) ) * zcoef 
     484      ENDDO 
    565485 
    566486      ! N/P and Si releases due to coastal rivers 
    … …  
    568488      DO jj = 1, jpj 
    569489         DO ji = 1, jpi 
    570             zcoef = ryyss * e1t(ji,jj) * e2t(ji,jj) * fse3t(ji,jj,1) * tmask(ji,jj,1) * tmask_i(ji,jj) 
     490            zcoef = ryyss * e1e2t(ji,jj)  * fse3t(ji,jj,1) * tmask(ji,jj,1)  
    571491            cotdep(ji,jj) =  river(ji,jj)                  *1E9 / ( 12. * zcoef + rtrn ) 
    572492            rivinp(ji,jj) = (river(ji,jj)+riverdoc(ji,jj)) *1E9 / ( 31.6* zcoef + rtrn ) 
    … …  
    577497      CALL lbc_lnk( cotdep , 'T', 1. )  ;  CALL lbc_lnk( rivinp , 'T', 1. )  ;  CALL lbc_lnk( nitdep , 'T', 1. ) 
    578498 
    579       rivpo4input = 0.e0 
    580       rivalkinput = 0.e0 
    581       nitdepinput = 0.e0 
    582       DO jj = 2 , jpjm1 
    583          DO ji = fs_2, fs_jpim1 
    584             zcoef = cvol(ji,jj,1) * ryyss 
    585             rivpo4input = rivpo4input + rivinp(ji,jj) * zcoef 
    586             rivalkinput = rivalkinput + cotdep(ji,jj) * zcoef 
    587             nitdepinput = nitdepinput + nitdep(ji,jj) * zcoef 
    588          END DO 
    589      END DO 
    590       IF( lk_mpp ) THEN 
    591          CALL mpp_sum( rivpo4input )  ! sum over the global domain 
    592          CALL mpp_sum( rivalkinput )  ! sum over the global domain 
    593          CALL mpp_sum( nitdepinput )  ! sum over the global domain 
    594       ENDIF 
     499      rivpo4input = glob_sum( rivinp(:,:) * cvol(:,:,1) ) * ryyss 
     500      rivalkinput = glob_sum( cotdep(:,:) * cvol(:,:,1) ) * ryyss 
     501      nitdepinput = glob_sum( nitdep(:,:) * cvol(:,:,1) ) * ryyss 
    595502 
    596503 

  • branches/nemo_v3_3_beta/NEMOGCM/NEMO/TOP_SRC/PISCES/trcrst_pisces.F90

    r2431 r2457  
    2222   USE trcdta 
    2323   USE lib_mpp 
     24   USE lib_fortran 
    2425 
    2526   IMPLICIT NONE 
    … …  
    123124         ! set total alkalinity, phosphate, nitrate & silicate 
    124125 
    125          zalksum = 0.e0 
    126          zpo4sum = 0.e0 
    127          zno3sum = 0.e0 
    128          zsilsum = 0.e0 
    129          DO jk = 1, jpk 
    130             DO jj = 1, jpj 
    131                DO ji = 1, jpi 
    132                   zvol = cvol(ji,jj,jk) 
    133 #  if defined key_degrad 
    134                   zvol = zvol * facvol(ji,jj,jk) 
    135 #  endif 
    136                   zalksum = zalksum + trn(ji,jj,jk,jptal) * zvol 
    137                   zpo4sum = zpo4sum + trn(ji,jj,jk,jppo4) * zvol 
    138                   zno3sum = zno3sum + trn(ji,jj,jk,jpno3) * zvol 
    139                   zsilsum = zsilsum + trn(ji,jj,jk,jpsil) * zvol 
    140                END DO 
    141             END DO 
    142          END DO 
    143          IF( lk_mpp )   CALL mpp_sum( zalksum )     ! sum over the global domain 
    144          IF( lk_mpp )   CALL mpp_sum( zpo4sum )     ! sum over the global domain 
    145          IF( lk_mpp )   CALL mpp_sum( zno3sum )     ! sum over the global domain 
    146          IF( lk_mpp )   CALL mpp_sum( zsilsum )     ! sum over the global domain 
    147126         zarea   = 1. / areatot * 1.e6 
    148          zalksum = zalksum * zarea 
    149          zpo4sum = zpo4sum * zarea / 122. 
    150          zno3sum = zno3sum * zarea / 7.6 
    151          zsilsum = zsilsum * zarea 
     127# if defined key_degrad 
     128         zalksum = glob_sum( trn(:,:,:,jptal) * cvol(:,:,:) * facvol(:,:,:) ) * zarea 
     129         zpo4sum = glob_sum( trn(:,:,:,jppo4) * cvol(:,:,:) * facvol(:,:,:) ) * zarea / 122. 
     130         zno3sum = glob_sum( trn(:,:,:,jpno3) * cvol(:,:,:) * facvol(:,:,:) ) * zarea / 7.6 
     131         zsilsum = glob_sum( trn(:,:,:,jpsil) * cvol(:,:,:) * facvol(:,:,:) ) * zarea 
     132# else 
     133         zalksum = glob_sum( trn(:,:,:,jptal) * cvol(:,:,:)  ) * zarea 
     134         zpo4sum = glob_sum( trn(:,:,:,jppo4) * cvol(:,:,:)  ) * zarea / 122. 
     135         zno3sum = glob_sum( trn(:,:,:,jpno3) * cvol(:,:,:)  ) * zarea / 7.6 
     136         zsilsum = glob_sum( trn(:,:,:,jpsil) * cvol(:,:,:)  ) * zarea 
     137# endif 
    152138 
    153139         IF(lwp) WRITE(numout,*) '       TALK mean : ', zalksum 

  • branches/nemo_v3_3_beta/NEMOGCM/NEMO/TOP_SRC/SED/sed.F90

    r2281 r2457  
    1717      gphit    =>   gphit  ,  & !: latitude  of t-point (degre) 
    1818      e3t_0    =>   e3t_0  ,  & !: reference depth of t-points (m) 
    19       mbathy   =>   mbathy ,  & !: bathymetry 
     19      mbkt     =>   mbkt   ,  & !: vertical index of the bottom last T- ocean level 
    2020      tmask    =>   tmask  ,  & !: land/ocean mask at t-points 
    2121      rdt      =>   rdt         !: time step for the dynamics 

  • branches/nemo_v3_3_beta/NEMOGCM/NEMO/TOP_SRC/SED/sedchem.F90

    r2281 r2457  
    216216      DO jj = 1,jpj 
    217217         DO ji = 1, jpi 
    218             ikt = MAX( mbathy(ji,jj)-1, 1 ) 
     218            ikt = mbkt(ji,jj)  
    219219            IF ( tmask(ji,jj,ikt) == 1 ) THEN 
    220220               zchem_data(ji,jj,1) = ak13 (ji,jj,ikt) 

  • branches/nemo_v3_3_beta/NEMOGCM/NEMO/TOP_SRC/SED/seddta.F90

    r2281 r2457  
    118118         DO jj = 1,jpj 
    119119            DO ji = 1, jpi 
    120                ikt = MAX( mbathy(ji,jj)-1, 1 ) 
     120               ikt = mbkt(ji,jj) 
    121121               IF ( tmask(ji,jj,ikt) == 1 ) THEN 
    122122                  trc_data(ji,jj,1)  = trn  (ji,jj,ikt,jptal) 

  • branches/nemo_v3_3_beta/NEMOGCM/NEMO/TOP_SRC/SED/sedini.F90

    r2281 r2457  
    135135      DO jj = 1, jpj 
    136136         DO ji = 1, jpi 
    137             ikt = MAX( mbathy(ji,jj) - 1, 1 ) 
     137            ikt = mbkt(ji,jj)  
    138138            IF( tmask(ji,jj,ikt) == 1 ) epkbot(ji,jj) = e3t_0(ikt) 
    139139         ENDDO 

  • branches/nemo_v3_3_beta/NEMOGCM/NEMO/TOP_SRC/SED/sedsfc.F90

    r2281 r2457  
    5252      DO jj = 1,jpj 
    5353         DO ji = 1, jpi 
    54             ikt = MAX( mbathy(ji,jj)-1, 1 ) 
     54            ikt = mbkt(ji,jj) 
    5555            IF ( tmask(ji,jj,ikt) == 1 ) THEN 
    5656               trn(ji,jj,ikt,jptal) = trc_data(ji,jj,1) 

  • branches/nemo_v3_3_beta/NEMOGCM/NEMO/TOP_SRC/oce_trc.F90

    r2431 r2457  
    135135# endif 
    136136   !* masks, bathymetry * 
    137    USE dom_oce , ONLY :   mbathy     =>   mbathy     !: number of ocean level (=0,  & 1, ... , jpk-1)  
     137   USE dom_oce , ONLY :   mbkt       =>   mbkt       !: vertical index of the bottom last T- ocean level 
     138   USE dom_oce , ONLY :   mbku       =>   mbku       !: vertical index of the bottom last U- ocean level 
     139   USE dom_oce , ONLY :   mbkv       =>   mbkv       !: vertical index of the bottom last V- ocean level 
    138140   USE dom_oce , ONLY :   tmask_i    =>   tmask_i    !: Interior mask at t-points 
    139141   USE dom_oce , ONLY :   tmask      =>   tmask      !: land/ocean mask at t-points 

  • branches/nemo_v3_3_beta/NEMOGCM/NEMO/TOP_SRC/trcdta.F90

    r2287 r2457  
    138138                     DO jj = 1, jpj                ! interpolation of temperature at the last level 
    139139                        DO ji = 1, jpi 
    140                            ik = mbathy(ji,jj) - 1 
     140                           ik = mbkt(ji,jj) 
    141141                           IF( ik > 2 ) THEN 
    142142                              zl = ( gdept_0(ik) - fsdept_0(ji,jj,ik) ) / ( gdept_0(ik) - gdept_0(ik-1) ) 

  • branches/nemo_v3_3_beta/NEMOGCM/NEMO/TOP_SRC/trcini.F90

    r2431 r2457  
    3333   USE prtctl_trc      ! Print control passive tracers (prt_ctl_trc_init routine) 
    3434   USE lib_mpp         ! distributed memory computing library 
     35   USE lib_fortran     !  
    3536    
    3637   IMPLICIT NONE 
    … …  
    6768      !                 ! masked grid volume 
    6869      DO jk = 1, jpk 
    69          cvol(:,:,jk) = e1t(:,:) * e2t(:,:) * fse3t(:,:,jk) * tmask(:,:,jk) * tmask_i(:,:) 
     70         cvol(:,:,jk) = e1t(:,:) * e2t(:,:) * fse3t(:,:,jk) * tmask(:,:,jk)  
    7071      END DO 
    7172 
    7273      ! total volume of the ocean 
    7374#if ! defined key_degrad 
    74       areatot = SUM( cvol(:,:,:) ) 
     75      areatot = glob_sum( cvol(:,:,:) ) 
    7576#else 
    76       areatot = SUM( cvol(:,:,:) * facvol(:,:,:) )  ! degrad option: reduction by facvol 
     77      areatot = glob_sum( cvol(:,:,:) * facvol(:,:,:) )  ! degrad option: reduction by facvol 
    7778#endif 
    78       IF( lk_mpp )   CALL mpp_sum( areatot )        ! sum over the global domain   
    7979 
    8080                                  CALL trc_nam      ! read passive tracers namelists 
    … …  
    146146      DO jn = 1, jptra 
    147147#if ! defined key_degrad 
    148          trai = trai + SUM( trn(:,:,:,jn) * cvol(:,:,:) ) 
     148         trai = trai + glob_sum( trn(:,:,:,jn) * cvol(:,:,:) ) 
    149149#else 
    150          trai = trai + SUM( trn(:,:,:,jn) * cvol(:,:,:) * facvol(:,:,:) ) ! degrad option: reduction by facvol 
     150         trai = trai + glob_sum( trn(:,:,:,jn) * cvol(:,:,:) * facvol(:,:,:) ) ! degrad option: reduction by facvol 
    151151#endif 
    152152      END DO       
    153       IF( lk_mpp )   CALL mpp_sum( trai )     ! sum over the global domain   
    154  
    155153 
    156154      !                 ! control print 

  • branches/nemo_v3_3_beta/NEMOGCM/NEMO/TOP_SRC/trcrst.F90

    r2287 r2457  
    2727   USE trcnam_trp 
    2828   USE lib_mpp 
     29   USE lib_fortran 
    2930   USE iom 
    3031   USE trcrst_cfc      ! CFC       
    … …  
    317318      zdiag_tot = 0.e0 
    318319      DO jn = 1, jptra 
    319          zdiag_var    = 0.e0 
    320          zdiag_varmin = 0.e0 
    321          zdiag_varmax = 0.e0 
    322          DO jk = 1, jpk 
    323             DO jj = 1, jpj 
    324                DO ji = 1, jpi 
    325                   zvol = cvol(ji,jj,jk) 
    326320#  if defined key_degrad 
    327                   zvol = zvol * facvol(ji,jj,jk) 
     321         zdiag_var = glob_sum( trn(:,:,:,jn) * cvol(:,:,:) * facvol(:,:,:) ) 
     322#  else 
     323         zdiag_var = glob_sum( trn(:,:,:,jn) * cvol(:,:,:)  ) 
    328324#  endif 
    329                   zdiag_var = zdiag_var + trn(ji,jj,jk,jn) * zvol 
    330                END DO 
    331             END DO 
    332          END DO 
    333           
    334325         zdiag_varmin = MINVAL( trn(:,:,:,jn), mask= ((tmask*SPREAD(tmask_i,DIM=3,NCOPIES=jpk).NE.0.)) ) 
    335326         zdiag_varmax = MAXVAL( trn(:,:,:,jn), mask= ((tmask*SPREAD(tmask_i,DIM=3,NCOPIES=jpk).NE.0.)) ) 
    … …  
    337328            CALL mpp_min( zdiag_varmin )      ! min over the global domain 
    338329            CALL mpp_max( zdiag_varmax )      ! max over the global domain 
    339             CALL mpp_sum( zdiag_var    )      ! sum over the global domain 
    340330         END IF 
    341331         zdiag_tot = zdiag_tot + zdiag_var 
Note: See TracChangeset for help on using the changeset viewer.