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 8868 for branches/2017/dev_METO_2017/NEMOGCM/NEMO – NEMO

Ignore:
Timestamp:
2017-12-01T09:43:23+01:00 (6 years ago)
Author:
timgraham
Message:

Merged dev_r8789_sbc into branch

Location:
branches/2017/dev_METO_2017/NEMOGCM/NEMO
Files:
42 edited

Legend:

Unmodified
Added
Removed
  • branches/2017/dev_METO_2017/NEMOGCM/NEMO/LIM_SRC_3/limthd.F90

    r7813 r8868  
    205205      IF( .NOT. ln_limdO )  qlead(:,:) = 0._wp 
    206206      ! In case we bypass growing/melting from top and bottom: we suppose ice is impermeable => ocean is isolated from atmosphere 
    207       IF( .NOT. ln_limdH )  hfx_in(:,:) = pfrld(:,:) * ( qns_oce(:,:) + qsr_oce(:,:) ) + qemp_oce(:,:) 
    208       IF( .NOT. ln_limdH )  fhtur (:,:) = 0._wp  ;  fhld  (:,:) = 0._wp 
     207      IF( .NOT. ln_limdH ) THEN 
     208         hfx_in(:,:) = pfrld(:,:) * ( qns_oce(:,:) + qsr_oce(:,:) ) + qemp_oce(:,:) 
     209         fhtur (:,:) = 0._wp 
     210         fhld  (:,:) = 0._wp 
     211      ENDIF 
    209212 
    210213      ! --------------------------------------------------------------------- 
  • branches/2017/dev_METO_2017/NEMOGCM/NEMO/OFF_SRC/nemogcm.F90

    r7761 r8868  
    548548 
    549549      ! Find the factors of n. 
    550       IF( kn == 1 )   GOTO 20 
    551  
    552       ! nu holds the unfactorised part of the number. 
    553       ! knfax holds the number of factors found. 
    554       ! l points to the allowed factor list. 
    555       ! ifac holds the current factor. 
    556  
    557       inu   = kn 
    558       knfax = 0 
    559  
    560       DO jl = ntest, 1, -1 
    561          ! 
    562          ifac = ilfax(jl) 
    563          IF( ifac > inu )   CYCLE 
    564  
    565          ! Test whether the factor will divide. 
    566  
    567          IF( MOD(inu,ifac) == 0 ) THEN 
     550      IF( kn .NE. 1 ) THEN 
     551 
     552         ! nu holds the unfactorised part of the number. 
     553         ! knfax holds the number of factors found. 
     554         ! l points to the allowed factor list. 
     555         ! ifac holds the current factor. 
     556    
     557         inu   = kn 
     558         knfax = 0 
     559    
     560         DO jl = ntest, 1, -1 
    568561            ! 
    569             knfax = knfax + 1            ! Add the factor to the list 
    570             IF( knfax > kmaxfax ) THEN 
    571                kerr = 6 
    572                write (*,*) 'FACTOR: insufficient space in factor array ', knfax 
    573                return 
     562            ifac = ilfax(jl) 
     563            IF( ifac > inu )   CYCLE 
     564    
     565            ! Test whether the factor will divide. 
     566    
     567            IF( MOD(inu,ifac) == 0 ) THEN 
     568               ! 
     569               knfax = knfax + 1            ! Add the factor to the list 
     570               IF( knfax > kmaxfax ) THEN 
     571                  kerr = 6 
     572                  write (*,*) 'FACTOR: insufficient space in factor array ', knfax 
     573                  return 
     574               ENDIF 
     575               kfax(knfax) = ifac 
     576               ! Store the other factor that goes with this one 
     577               knfax = knfax + 1 
     578               kfax(knfax) = inu / ifac 
     579               !WRITE (*,*) 'ARPDBG, factors ',knfax-1,' & ',knfax,' are ', kfax(knfax-1),' and ',kfax(knfax) 
    574580            ENDIF 
    575             kfax(knfax) = ifac 
    576             ! Store the other factor that goes with this one 
    577             knfax = knfax + 1 
    578             kfax(knfax) = inu / ifac 
    579             !WRITE (*,*) 'ARPDBG, factors ',knfax-1,' & ',knfax,' are ', kfax(knfax-1),' and ',kfax(knfax) 
    580          ENDIF 
    581          ! 
    582       END DO 
    583  
    584    20 CONTINUE      ! Label 20 is the exit point from the factor search loop. 
     581            ! 
     582         END DO 
     583    
     584      ENDIF 
    585585      ! 
    586586   END SUBROUTINE factorise 
  • branches/2017/dev_METO_2017/NEMOGCM/NEMO/OPA_SRC/ASM/asminc.F90

    r8030 r8868  
    2020   !!   dyn_asm_inc    : Apply the dynamic (u and v) increments 
    2121   !!   ssh_asm_inc    : Apply the SSH increment 
     22   !!   ssh_asm_div    : Apply divergence associated with SSH increment 
    2223   !!   seaice_asm_inc : Apply the seaice increment 
    2324   !!---------------------------------------------------------------------- 
     
    4849   PUBLIC   dyn_asm_inc    !: Apply the dynamic (u and v) increments 
    4950   PUBLIC   ssh_asm_inc    !: Apply the SSH increment 
     51   PUBLIC   ssh_asm_div    !: Apply the SSH divergence 
    5052   PUBLIC   seaice_asm_inc !: Apply the seaice increment 
    5153 
     
    768770            ENDIF 
    769771            ! 
     772#if defined key_asminc 
     773         ELSE IF( kt == nitiaufin_r+1 ) THEN 
     774            ! 
     775            ssh_iau(:,:) = 0._wp 
     776            ! 
     777#endif 
    770778         ENDIF 
    771779         !                          !----------------------------------------- 
     
    792800   END SUBROUTINE ssh_asm_inc 
    793801 
     802   SUBROUTINE ssh_asm_div( kt, phdivn ) 
     803      !!---------------------------------------------------------------------- 
     804      !!                  ***  ROUTINE ssh_asm_div  *** 
     805      !! 
     806      !! ** Purpose :   ssh increment with z* is incorporated via a correction of the local divergence           
     807      !!                across all the water column 
     808      !! 
     809      !! ** Method  : 
     810      !!                CAUTION : sshiau is positive (inflow) decreasing the 
     811      !!                          divergence and expressed in m/s 
     812      !! 
     813      !! ** Action  :   phdivn   decreased by the ssh increment 
     814      !!---------------------------------------------------------------------- 
     815      INTEGER, INTENT(IN) :: kt                               ! ocean time-step index 
     816      REAL(wp), DIMENSION(:,:,:), INTENT(inout) ::   phdivn   ! horizontal divergence 
     817      !! 
     818      INTEGER  ::   jk                                        ! dummy loop index 
     819      REAL(wp), DIMENSION(:,:)  , POINTER       ::   ztim     ! local array 
     820      !!---------------------------------------------------------------------- 
     821      !  
     822#if defined key_asminc 
     823      CALL ssh_asm_inc( kt ) !==   (calculate increments) 
     824      ! 
     825      IF( ln_linssh ) THEN  
     826         phdivn(:,:,1) = phdivn(:,:,1) - ssh_iau(:,:) / e3t_n(:,:,1) * tmask(:,:,1) 
     827      ELSE  
     828         CALL wrk_alloc( jpi,jpj, ztim) 
     829         ztim(:,:) = ssh_iau(:,:) / ( ht_n(:,:) + 1.0 - ssmask(:,:) ) 
     830         DO jk = 1, jpkm1                                  
     831            phdivn(:,:,jk) = phdivn(:,:,jk) - ztim(:,:) * tmask(:,:,jk)  
     832         END DO 
     833         ! 
     834         CALL wrk_dealloc( jpi,jpj, ztim) 
     835      ENDIF 
     836#endif 
     837      ! 
     838   END SUBROUTINE ssh_asm_div 
    794839 
    795840   SUBROUTINE seaice_asm_inc( kt, kindic ) 
  • branches/2017/dev_METO_2017/NEMOGCM/NEMO/OPA_SRC/BDY/bdydta.F90

    r7861 r8868  
    266266                        IF( ln_full_vel_array(ib_bdy) ) THEN 
    267267                           CALL fld_read( kt=kt, kn_fsbc=1, sd=bf(jstart:jend), map=nbmap_ptr(jstart:jend),  & 
    268                                      & kit=jit, kt_offset=time_offset , jpk_bdy=nb_jpk_bdy, fvl=ln_full_vel_array(ib_bdy)  ) 
     268                                     & kit=jit, kt_offset=time_offset , jpk_bdy=nb_jpk_bdy,   & 
     269                                     & fvl=ln_full_vel_array(ib_bdy)  ) 
    269270                        ELSE 
    270271                           CALL fld_read( kt=kt, kn_fsbc=1, sd=bf(jstart:jend), map=nbmap_ptr(jstart:jend),  & 
     
    335336                     jend = jstart + dta%nread(1) - 1 
    336337                     CALL fld_read( kt=kt, kn_fsbc=1, sd=bf(jstart:jend), & 
    337                                   & map=nbmap_ptr(jstart:jend), kt_offset=time_offset, jpk_bdy=nb_jpk_bdy, fvl=ln_full_vel_array(ib_bdy) ) 
     338                                  & map=nbmap_ptr(jstart:jend), kt_offset=time_offset, jpk_bdy=nb_jpk_bdy,   & 
     339                                  & fvl=ln_full_vel_array(ib_bdy) ) 
    338340                  ENDIF 
    339341                  ! If full velocities in boundary data then split into barotropic and baroclinic data 
  • branches/2017/dev_METO_2017/NEMOGCM/NEMO/OPA_SRC/DIA/dia25h.F90

    r7753 r8868  
    316316#endif 
    317317            cnt_25h = 1 
    318             IF (lwp)  WRITE(numout,*) 'dia_wri_tide : After 25hr mean write, reset sum to current value and cnt_25h to one for overlapping average',cnt_25h 
     318            IF (lwp)  WRITE(numout,*) 'dia_wri_tide :   & 
     319        &    After 25hr mean write, reset sum to current value and cnt_25h to one for overlapping average',cnt_25h 
    319320 
    320321      ENDIF !  cnt_25h .EQ. 25 .AND.  MOD( kt, i_steps * 24) == 0 .AND. kt .NE. nn_it000 
  • branches/2017/dev_METO_2017/NEMOGCM/NEMO/OPA_SRC/DIA/diacfl.F90

    r7753 r8868  
    151151            WRITE(numout,*) 'dia_cfl     : Maximum Courant number information for the run:' 
    152152            WRITE(numout,*) '~~~~~~~~~~~~' 
    153             WRITE(numout,FMT='(12x,a12,7x,f6.4,5x,a16,i4,1x,i4,1x,i4,a1)') 'Run Max Cu', cu_max, 'at (i, j, k) = (', cu_loc(1), cu_loc(2), cu_loc(3), ')' 
     153            WRITE(numout,FMT='(12x,a12,7x,f6.4,5x,a16,i4,1x,i4,1x,i4,a1)') 'Run Max Cu', cu_max, 'at (i, j, k) =   & 
     154                          &   (', cu_loc(1), cu_loc(2), cu_loc(3), ')' 
    154155            WRITE(numout,FMT='(12x,a8,11x,f7.1)') ' => dt/C', dt*(1.0/cu_max) 
    155             WRITE(numout,FMT='(12x,a12,7x,f6.4,5x,a16,i4,1x,i4,1x,i4,a1)') 'Run Max Cv', cv_max, 'at (i, j, k) = (', cv_loc(1), cv_loc(2), cv_loc(3), ')' 
     156            WRITE(numout,FMT='(12x,a12,7x,f6.4,5x,a16,i4,1x,i4,1x,i4,a1)') 'Run Max Cv', cv_max, 'at (i, j, k) =   & 
     157                          &   (', cv_loc(1), cv_loc(2), cv_loc(3), ')' 
    156158            WRITE(numout,FMT='(12x,a8,11x,f7.1)') ' => dt/C', dt*(1.0/cv_max) 
    157             WRITE(numout,FMT='(12x,a12,7x,f6.4,5x,a16,i4,1x,i4,1x,i4,a1)') 'Run Max Cw', cw_max, 'at (i, j, k) = (', cw_loc(1), cw_loc(2), cw_loc(3), ')' 
     159            WRITE(numout,FMT='(12x,a12,7x,f6.4,5x,a16,i4,1x,i4,1x,i4,a1)') 'Run Max Cw', cw_max, 'at (i, j, k) =   & 
     160                          &   (', cw_loc(1), cw_loc(2), cw_loc(3), ')' 
    158161            WRITE(numout,FMT='(12x,a8,11x,f7.1)') ' => dt/C', dt*(1.0/cw_max) 
    159162 
  • branches/2017/dev_METO_2017/NEMOGCM/NEMO/OPA_SRC/DIA/diawri.F90

    r7753 r8868  
    282282      ! 
    283283      IF ( iom_use("eken") ) THEN 
    284          rke(:,:,jk) = 0._wp                               !      kinetic energy  
     284         rke(:,:,jpk) = 0._wp                               !      kinetic energy  
    285285         DO jk = 1, jpkm1 
    286286            DO jj = 2, jpjm1 
    287287               DO ji = fs_2, fs_jpim1   ! vector opt. 
    288                   zztmp   = 1._wp / ( e1e2t(ji,jj) * e3t_n(ji,jj,jk) ) 
    289                   zztmpx  = 0.5 * (  un(ji-1,jj,jk) * un(ji-1,jj,jk) * e2u(ji-1,jj) * e3u_n(ji-1,jj,jk)    & 
    290                      &             + un(ji  ,jj,jk) * un(ji  ,jj,jk) * e2u(ji  ,jj) * e3u_n(ji  ,jj,jk) )  & 
     288                  zztmp   = r1_e1e2t(ji,jj) / e3t_n(ji,jj,jk) 
     289                  zztmpx  = 0.5 * (  un(ji-1,jj,jk) * un(ji-1,jj,jk) * e1e2u(ji-1,jj) * e3u_n(ji-1,jj,jk)    & 
     290                     &             + un(ji  ,jj,jk) * un(ji  ,jj,jk) * e1e2u(ji  ,jj) * e3u_n(ji  ,jj,jk) )  & 
    291291                     &          *  zztmp  
    292292                  ! 
    293                   zztmpy  = 0.5 * (  vn(ji,jj-1,jk) * vn(ji,jj-1,jk) * e1v(ji,jj-1) * e3v_n(ji,jj-1,jk)    & 
    294                      &             + vn(ji,jj  ,jk) * vn(ji,jj  ,jk) * e1v(ji,jj  ) * e3v_n(ji,jj  ,jk) )  & 
     293                  zztmpy  = 0.5 * (  vn(ji,jj-1,jk) * vn(ji,jj-1,jk) * e1e2v(ji,jj-1) * e3v_n(ji,jj-1,jk)    & 
     294                     &             + vn(ji,jj  ,jk) * vn(ji,jj  ,jk) * e1e2v(ji,jj  ) * e3v_n(ji,jj  ,jk) )  & 
    295295                     &          *  zztmp  
    296296                  ! 
  • branches/2017/dev_METO_2017/NEMOGCM/NEMO/OPA_SRC/DOM/daymod.F90

    r7646 r8868  
    134134 
    135135      ! control print 
    136       IF(lwp) WRITE(numout,'(a,i6,a,i2,a,i2,a,i8,a,i8,a,i8,a,i8)')' =======>> 1/2 time step before the start of the run DATE Y/M/D = ',   & 
     136      IF(lwp) WRITE(numout,'(a,i6,a,i2,a,i2,a,i8,a,i8,a,i8,a,i8)')   & 
     137           &                   ' =======>> 1/2 time step before the start of the run DATE Y/M/D = ',   & 
    137138           &                   nyear, '/', nmonth, '/', nday, '  nsec_day:', nsec_day, '  nsec_week:', nsec_week, '  & 
    138139           &                   nsec_month:', nsec_month , '  nsec_year:' , nsec_year 
  • branches/2017/dev_METO_2017/NEMOGCM/NEMO/OPA_SRC/DOM/iscplrst.F90

    r7646 r8868  
    219219               DO ji=1,jpi 
    220220                  IF (tmask(ji,jj,1) == 0._wp .OR. ptmask_b(ji,jj,1) == 0._wp) THEN 
    221                      e3t_n(ji,jj,jk) = e3t_0(ji,jj,jk) * ( 1._wp + sshn(ji,jj) / ( ht_0(ji,jj) + 1._wp - ssmask(ji,jj) ) * tmask(ji,jj,jk) ) 
     221                     e3t_n(ji,jj,jk) = e3t_0(ji,jj,jk) * ( 1._wp + sshn(ji,jj) /   & 
     222                     &   ( ht_0(ji,jj) + 1._wp - ssmask(ji,jj) ) * tmask(ji,jj,jk) ) 
    222223                  ENDIF 
    223224               END DO 
     
    390391            DO jj = 1,jpj 
    391392               DO ji = 1,jpi 
    392                   IF (zwmaskn(ji,jj,jk) * zwmaskb(ji,jj,jk) == 1._wp .AND. (tmask(ji,jj,1)==0._wp .OR. ptmask_b(ji,jj,1)==0._wp) ) THEN 
     393                  IF (zwmaskn(ji,jj,jk) * zwmaskb(ji,jj,jk) == 1._wp .AND.   & 
     394                  &      (tmask(ji,jj,1)==0._wp .OR. ptmask_b(ji,jj,1)==0._wp) ) THEN 
    393395                     !compute weight 
    394396                     zdzp1 = MAX(0._wp,gdepw_n(ji,jj,jk+1) - pdepw_b(ji,jj,jk+1)) 
  • branches/2017/dev_METO_2017/NEMOGCM/NEMO/OPA_SRC/DYN/divhor.F90

    r7753 r8868  
    2525   USE iscplhsb        ! ice sheet / ocean coupling 
    2626   USE iscplini        ! ice sheet / ocean coupling 
     27#if defined key_asminc    
     28   USE asminc          ! Assimilation increment 
     29#endif 
    2730   ! 
    2831   USE in_out_manager  ! I/O manager 
     
    9295      IF( ln_rnf )   CALL sbc_rnf_div( hdivn )      !==  runoffs    ==!   (update hdivn field) 
    9396      ! 
     97#if defined key_asminc  
     98      IF( ln_sshinc .AND. ln_asmiau )   CALL ssh_asm_div( kt, hdivn )   !==  SSH assimilation  ==!   (update hdivn field) 
     99      !  
     100#endif 
    94101      IF( ln_isf )   CALL sbc_isf_div( hdivn )      !==  ice shelf  ==!   (update hdivn field) 
    95102      ! 
  • branches/2017/dev_METO_2017/NEMOGCM/NEMO/OPA_SRC/DYN/sshwzv.F90

    r7753 r8868  
    2727   USE agrif_opa_interp 
    2828#endif 
    29 #if defined key_asminc    
    30    USE   asminc       ! Assimilation increment 
    31 #endif 
    3229   ! 
    3330   USE in_out_manager ! I/O manager 
     
    121118         ENDIF 
    122119      ENDIF 
    123  
    124 #if defined key_asminc 
    125       IF( lk_asminc .AND. ln_sshinc .AND. ln_asmiau ) THEN     ! Include the IAU weighted SSH increment 
    126          CALL ssh_asm_inc( kt ) 
    127          ssha(:,:) = ssha(:,:) + z2dt * ssh_iau(:,:) 
    128       ENDIF 
    129 #endif 
    130120      !                                           !------------------------------! 
    131121      !                                           !           outputs            ! 
  • branches/2017/dev_METO_2017/NEMOGCM/NEMO/OPA_SRC/ICB/icbtrj.F90

    r8080 r8868  
    8080 
    8181      ! define trajectory output name 
    82       IF( lk_mpp ) THEN   ;   WRITE(cl_filename,'("trajectory_icebergs_",A,"-",A,"_",I4.4,".nc")') TRIM(ADJUSTL(cldate_ini)), TRIM(ADJUSTL(cldate_end)), narea-1 
    83       ELSE                ;   WRITE(cl_filename,'("trajectory_icebergs_",A,"-",A         ,".nc")') TRIM(ADJUSTL(cldate_ini)), TRIM(ADJUSTL(cldate_end)) 
     82      IF( lk_mpp ) THEN   ;   WRITE(cl_filename,'("trajectory_icebergs_",A,"-",A,"_",I4.4,".nc")')   & 
     83                          &   TRIM(ADJUSTL(cldate_ini)), TRIM(ADJUSTL(cldate_end)), narea-1 
     84      ELSE                ;   WRITE(cl_filename,'("trajectory_icebergs_",A,"-",A         ,".nc")')   & 
     85                          &   TRIM(ADJUSTL(cldate_ini)), TRIM(ADJUSTL(cldate_end)) 
    8486      ENDIF 
    8587      IF ( lwp .AND. nn_verbose_level >= 0) WRITE(numout,'(2a)') 'icebergs, icb_trj_init: creating ',TRIM(cl_filename) 
  • branches/2017/dev_METO_2017/NEMOGCM/NEMO/OPA_SRC/IOM/iom.F90

    r7768 r8868  
    774774         istart(idmspc+1) = itime 
    775775 
    776          IF( PRESENT(kstart) .AND. .NOT. ll_depth_spec ) THEN ; istart(1:idmspc) = kstart(1:idmspc) ; icnt(1:idmspc) = kcount(1:idmspc) 
     776         IF( PRESENT(kstart) .AND. .NOT. ll_depth_spec ) THEN  
     777            istart(1:idmspc) = kstart(1:idmspc) ; icnt(1:idmspc) = kcount(1:idmspc) 
    777778         ELSE 
    778             IF(           idom == jpdom_unknown ) THEN                                                ; icnt(1:idmspc) = idimsz(1:idmspc) 
     779            IF(idom == jpdom_unknown ) THEN 
     780               icnt(1:idmspc) = idimsz(1:idmspc) 
    779781            ELSE  
    780782               IF( .NOT. PRESENT(pv_r1d) ) THEN   !   not a 1D array 
     
    15401542      ALLOCATE( zlon(ni*nj) )       ;       zlon(:) = 0. 
    15411543 
    1542       CALL dom_ngb( 180., 90., ix, iy, 'T' ) !  i-line that passes near the North Pole : Reference latitude (used in plots) 
     1544      CALL dom_ngb( -168.53, 65.03, ix, iy, 'T' ) !  i-line that passes through Bering Strait: Reference latitude (used in plots) 
    15431545      CALL iom_set_domain_attr("gznl", ni_glo=jpiglo, nj_glo=jpjglo, ibegin=nimpp+nldi-2, jbegin=njmpp+nldj-2, ni=ni, nj=nj) 
    15441546      CALL iom_set_domain_attr("gznl", data_dim=2, data_ibegin = 1-nldi, data_ni = jpi, data_jbegin = 1-nldj, data_nj = jpj) 
     
    15951597      ! frequency of the call of iom_put (attribut: freq_op) 
    15961598      f_op%timestep = 1        ;  f_of%timestep = 0  ; CALL iom_set_field_attr('field_definition', freq_op=f_op, freq_offset=f_of) 
     1599      f_op%timestep = 2        ;  f_of%timestep = 0  ; CALL iom_set_field_attr('trendT_even'     , freq_op=f_op, freq_offset=f_of) 
     1600      f_op%timestep = 2        ;  f_of%timestep = -1 ; CALL iom_set_field_attr('trendT_odd'      , freq_op=f_op, freq_offset=f_of) 
    15971601      f_op%timestep = nn_fsbc  ;  f_of%timestep = 0  ; CALL iom_set_field_attr('SBC'             , freq_op=f_op, freq_offset=f_of) 
    15981602      f_op%timestep = nn_fsbc  ;  f_of%timestep = 0  ; CALL iom_set_field_attr('SBC_scalar'      , freq_op=f_op, freq_offset=f_of) 
  • branches/2017/dev_METO_2017/NEMOGCM/NEMO/OPA_SRC/SBC/sbcapr.F90

    r6140 r8868  
    2323   PRIVATE 
    2424 
    25    PUBLIC   sbc_apr    ! routine called in sbcmod 
     25   PUBLIC   sbc_apr       ! routine called in sbcmod 
     26   PUBLIC   sbc_apr_init  ! routine called in sbcmod 
    2627    
    2728   !                                !!* namsbc_apr namelist (Atmospheric PRessure) * 
     
    4647CONTAINS 
    4748 
     49   SUBROUTINE sbc_apr_init 
     50      !!--------------------------------------------------------------------- 
     51      !!                     ***  ROUTINE sbc_apr  *** 
     52      !! 
     53      !! ** Purpose :   read atmospheric pressure fields in netcdf files. 
     54      !! 
     55      !! ** Method  : - Read namelist namsbc_apr 
     56      !!              - Read Patm fields in netcdf files  
     57      !!              - Compute reference atmospheric pressure 
     58      !!              - Compute inverse barometer ssh 
     59      !! ** action  :   apr      : atmospheric pressure at kt 
     60      !!                ssh_ib   : inverse barometer ssh at kt 
     61      !!--------------------------------------------------------------------- 
     62      INTEGER            ::   ierror  ! local integer  
     63      INTEGER            ::   ios     ! Local integer output status for namelist read 
     64      !! 
     65      CHARACTER(len=100) ::  cn_dir   ! Root directory for location of ssr files 
     66      TYPE(FLD_N)        ::  sn_apr   ! informations about the fields to be read 
     67      !! 
     68      NAMELIST/namsbc_apr/ cn_dir, sn_apr, ln_ref_apr, rn_pref, ln_apr_obc 
     69      !!---------------------------------------------------------------------- 
     70      REWIND( numnam_ref )              ! Namelist namsbc_apr in reference namelist : File for atmospheric pressure forcing 
     71      READ  ( numnam_ref, namsbc_apr, IOSTAT = ios, ERR = 901) 
     72901   IF( ios /= 0 ) CALL ctl_nam ( ios , 'namsbc_apr in reference namelist', lwp ) 
     73 
     74      REWIND( numnam_cfg )              ! Namelist namsbc_apr in configuration namelist : File for atmospheric pressure forcing 
     75      READ  ( numnam_cfg, namsbc_apr, IOSTAT = ios, ERR = 902 ) 
     76902   IF( ios /= 0 ) CALL ctl_nam ( ios , 'namsbc_apr in configuration namelist', lwp ) 
     77      IF(lwm) WRITE ( numond, namsbc_apr ) 
     78      ! 
     79      ALLOCATE( sf_apr(1), STAT=ierror )           !* allocate and fill sf_sst (forcing structure) with sn_sst 
     80      IF( ierror > 0 )   CALL ctl_stop( 'STOP', 'sbc_apr: unable to allocate sf_apr structure' ) 
     81      ! 
     82      CALL fld_fill( sf_apr, (/ sn_apr /), cn_dir, 'sbc_apr', 'Atmospheric pressure ', 'namsbc_apr' ) 
     83                                ALLOCATE( sf_apr(1)%fnow(jpi,jpj,1)   ) 
     84      IF( sn_apr%ln_tint )   ALLOCATE( sf_apr(1)%fdta(jpi,jpj,1,2) ) 
     85                             ALLOCATE( ssh_ib(jpi,jpj) , ssh_ibb(jpi,jpj) ) 
     86                             ALLOCATE( apr (jpi,jpj) ) 
     87      ! 
     88      IF( lwp )THEN                                 !* control print 
     89         WRITE(numout,*) 
     90         WRITE(numout,*) '   Namelist namsbc_apr : Atmospheric PRessure as extrenal forcing' 
     91         WRITE(numout,*) '      ref. pressure: global mean Patm (T) or a constant (F)  ln_ref_apr = ', ln_ref_apr 
     92      ENDIF 
     93      ! 
     94      IF( ln_ref_apr ) THEN                        !* Compute whole inner domain mean masked ocean surface 
     95         tarea = glob_sum( e1e2t(:,:) ) 
     96         IF(lwp) WRITE(numout,*) '         Variable ref. Patm computed over a ocean surface of ', tarea*1e-6, 'km2' 
     97      ELSE 
     98         IF(lwp) WRITE(numout,*) '         Reference Patm used : ', rn_pref, ' N/m2' 
     99      ENDIF 
     100      ! 
     101      r1_grau = 1.e0 / (grav * rau0)               !* constant for optimization 
     102      ! 
     103      !                                            !* control check 
     104      IF ( ln_apr_obc  ) THEN 
     105         IF(lwp) WRITE(numout,*) '         Inverse barometer added to OBC ssh data' 
     106      ENDIF 
     107!jc: stop below should rather be a warning  
     108      IF( ln_apr_obc .AND. .NOT.ln_apr_dyn   )   & 
     109            CALL ctl_warn( 'sbc_apr: use inverse barometer ssh at open boundary ONLY requires ln_apr_dyn=T' ) 
     110      ! 
     111   END SUBROUTINE sbc_apr_init 
     112 
    48113   SUBROUTINE sbc_apr( kt ) 
    49114      !!--------------------------------------------------------------------- 
     
    61126      INTEGER, INTENT(in)::   kt   ! ocean time step 
    62127      ! 
    63       INTEGER            ::   ierror  ! local integer  
    64       INTEGER            ::   ios     ! Local integer output status for namelist read 
    65       !! 
    66       CHARACTER(len=100) ::  cn_dir   ! Root directory for location of ssr files 
    67       TYPE(FLD_N)        ::  sn_apr   ! informations about the fields to be read 
    68       !! 
    69       NAMELIST/namsbc_apr/ cn_dir, sn_apr, ln_ref_apr, rn_pref, ln_apr_obc 
    70128      !!---------------------------------------------------------------------- 
    71       ! 
    72       !                                         ! -------------------- ! 
    73       IF( kt == nit000 ) THEN                   ! First call kt=nit000 ! 
    74          !                                      ! -------------------- ! 
    75          REWIND( numnam_ref )              ! Namelist namsbc_apr in reference namelist : File for atmospheric pressure forcing 
    76          READ  ( numnam_ref, namsbc_apr, IOSTAT = ios, ERR = 901) 
    77 901      IF( ios /= 0 ) CALL ctl_nam ( ios , 'namsbc_apr in reference namelist', lwp ) 
    78  
    79          REWIND( numnam_cfg )              ! Namelist namsbc_apr in configuration namelist : File for atmospheric pressure forcing 
    80          READ  ( numnam_cfg, namsbc_apr, IOSTAT = ios, ERR = 902 ) 
    81 902      IF( ios /= 0 ) CALL ctl_nam ( ios , 'namsbc_apr in configuration namelist', lwp ) 
    82          IF(lwm) WRITE ( numond, namsbc_apr ) 
    83          ! 
    84          ALLOCATE( sf_apr(1), STAT=ierror )           !* allocate and fill sf_sst (forcing structure) with sn_sst 
    85          IF( ierror > 0 )   CALL ctl_stop( 'STOP', 'sbc_apr: unable to allocate sf_apr structure' ) 
    86          ! 
    87          CALL fld_fill( sf_apr, (/ sn_apr /), cn_dir, 'sbc_apr', 'Atmospheric pressure ', 'namsbc_apr' ) 
    88                                 ALLOCATE( sf_apr(1)%fnow(jpi,jpj,1)   ) 
    89          IF( sn_apr%ln_tint )   ALLOCATE( sf_apr(1)%fdta(jpi,jpj,1,2) ) 
    90                                 ALLOCATE( ssh_ib(jpi,jpj) , ssh_ibb(jpi,jpj) ) 
    91                                 ALLOCATE( apr (jpi,jpj) ) 
    92          ! 
    93          IF(lwp) THEN                                 !* control print 
    94             WRITE(numout,*) 
    95             WRITE(numout,*) '   Namelist namsbc_apr : Atmospheric PRessure as extrenal forcing' 
    96             WRITE(numout,*) '      ref. pressure: global mean Patm (T) or a constant (F)  ln_ref_apr = ', ln_ref_apr 
    97          ENDIF 
    98          ! 
    99          IF( ln_ref_apr ) THEN                        !* Compute whole inner domain mean masked ocean surface 
    100             tarea = glob_sum( e1e2t(:,:) ) 
    101             IF(lwp) WRITE(numout,*) '         Variable ref. Patm computed over a ocean surface of ', tarea*1e-6, 'km2' 
    102          ELSE 
    103             IF(lwp) WRITE(numout,*) '         Reference Patm used : ', rn_pref, ' N/m2' 
    104          ENDIF 
    105          ! 
    106          r1_grau = 1.e0 / (grav * rau0)               !* constant for optimization 
    107          ! 
    108          !                                            !* control check 
    109          IF ( ln_apr_obc  ) THEN 
    110             IF(lwp) WRITE(numout,*) '         Inverse barometer added to OBC ssh data' 
    111          ENDIF 
    112 !jc: stop below should rather be a warning  
    113          IF( ln_apr_obc .AND. .NOT.ln_apr_dyn   )   & 
    114             CALL ctl_stop( 'sbc_apr: use inverse barometer ssh at open boundary ONLY requires ln_apr_dyn=T' ) 
    115       ENDIF 
    116129 
    117130      !                                         ! ========================== ! 
  • branches/2017/dev_METO_2017/NEMOGCM/NEMO/OPA_SRC/SBC/sbcblk.F90

    r7753 r8868  
    186186      !                             !** initialization of the chosen bulk formulae (+ check) 
    187187      !                                   !* select the bulk chosen in the namelist and check the choice 
    188       ;                                                        ioptio = 0 
     188                                                               ioptio = 0 
    189189      IF( ln_NCAR      ) THEN   ;   nblk =  np_NCAR        ;   ioptio = ioptio + 1   ;   ENDIF 
    190190      IF( ln_COARE_3p0 ) THEN   ;   nblk =  np_COARE_3p0   ;   ioptio = ioptio + 1   ;   ENDIF 
     
    219219         ALLOCATE( sf(ifpr)%fnow(jpi,jpj,1) ) 
    220220         IF( slf_i(ifpr)%ln_tint )   ALLOCATE( sf(ifpr)%fdta(jpi,jpj,1,2) ) 
     221 
     222           IF( slf_i(ifpr)%nfreqh .GT. 0._wp .AND. MOD( 3600._wp * slf_i(ifpr)%nfreqh , REAL(nn_fsbc, wp) * rdt) .NE. 0._wp  )   & 
     223            &  CALL ctl_warn( 'sbcmod time step rdt * nn_fsbc is NOT a submultiple of atmospheric forcing frequency' ) 
     224 
    221225      END DO 
    222226      !                                      !- fill the bulk structure with namelist informations 
  • branches/2017/dev_METO_2017/NEMOGCM/NEMO/OPA_SRC/SBC/sbccpl.F90

    r7968 r8868  
    18411841#endif 
    18421842      ! outputs 
    1843       IF( srcv(jpr_cal)%laction )    CALL iom_put('hflx_cal_cea' , - frcv(jpr_cal)%z3(:,:,1) * lfus                                  ) ! latent heat from calving 
    1844       IF( srcv(jpr_icb)%laction )    CALL iom_put('hflx_icb_cea' , - frcv(jpr_icb)%z3(:,:,1) * lfus                                  ) ! latent heat from icebergs melting 
    1845       IF( iom_use('hflx_snow_cea') ) CALL iom_put('hflx_snow_cea',  sprecip(:,:) * ( zcptsnw(:,:) - Lfus )                           ) ! heat flux from snow (cell average) 
    1846       IF( iom_use('hflx_rain_cea') ) CALL iom_put('hflx_rain_cea',( tprecip(:,:) - sprecip(:,:) ) * zcptrain(:,:)                    ) ! heat flux from rain (cell average) 
    1847       IF( iom_use('hflx_evap_cea') ) CALL iom_put('hflx_evap_cea',(frcv(jpr_tevp)%z3(:,:,1) - frcv(jpr_ievp)%z3(:,:,1) * zicefr(:,:) & ! heat flux from from evap (cell average) 
     1843      IF( srcv(jpr_cal)%laction )    CALL iom_put('hflx_cal_cea' ,   & 
     1844                                             &   - frcv(jpr_cal)%z3(:,:,1) * lfus) ! latent heat from calving 
     1845      IF( srcv(jpr_icb)%laction )    CALL iom_put('hflx_icb_cea' , - frcv(jpr_icb)%z3(:,:,1) * lfus) ! latent heat from icebergs melting 
     1846      IF( iom_use('hflx_snow_cea') ) CALL iom_put('hflx_snow_cea',  sprecip(:,:) * ( zcptsnw(:,:) - Lfus )) ! heat flux from snow (cell average) 
     1847      IF( iom_use('hflx_rain_cea') ) CALL iom_put('hflx_rain_cea',( tprecip(:,:) - sprecip(:,:) ) * zcptrain(:,:)) ! heat flux from rain (cell average) 
     1848      IF( iom_use('hflx_evap_cea') ) CALL iom_put('hflx_evap_cea',(frcv(jpr_tevp)%z3(:,:,1)   & 
     1849                                             &    - frcv(jpr_ievp)%z3(:,:,1) * zicefr(:,:) & ! heat flux from from evap (cell average) 
    18481850         &                                                        ) * zcptn(:,:) * tmask(:,:,1) ) 
    1849       IF( iom_use('hflx_snow_ao_cea') ) CALL iom_put('hflx_snow_ao_cea',sprecip(:,:) * (zcptsnw(:,:) - Lfus) * (1._wp - zsnw(:,:))   ) ! heat flux from snow (over ocean) 
    1850       IF( iom_use('hflx_snow_ai_cea') ) CALL iom_put('hflx_snow_ai_cea',sprecip(:,:) * (zcptsnw(:,:) - Lfus) *          zsnw(:,:)    ) ! heat flux from snow (over ice) 
     1851      IF( iom_use('hflx_snow_ao_cea') ) CALL iom_put('hflx_snow_ao_cea',sprecip(:,:)   & 
     1852                                             & * (zcptsnw(:,:) - Lfus) * (1._wp - zsnw(:,:))   ) ! heat flux from snow (over ocean) 
     1853      IF( iom_use('hflx_snow_ai_cea') ) CALL iom_put('hflx_snow_ai_cea',sprecip(:,:) * (zcptsnw(:,:) - Lfus)   & 
     1854                                             & *          zsnw(:,:)    ) ! heat flux from snow (over ice) 
    18511855      ! note: hflx for runoff and iceshelf are done in sbcrnf and sbcisf resp. 
    18521856      ! 
  • branches/2017/dev_METO_2017/NEMOGCM/NEMO/OPA_SRC/SBC/sbcisf.F90

    r7968 r8868  
    192192                     zqlatisf3d(ji,jj,jk) = zqlatisf3d(ji,jj,jk) + qisf     (ji,jj) * r1_hisf_tbl(ji,jj) * e3t_n(ji,jj,jk) 
    193193                  END DO 
    194                   zfwfisf3d (ji,jj,jk) = zfwfisf3d (ji,jj,jk) + fwfisf   (ji,jj) * r1_hisf_tbl(ji,jj) * ralpha(ji,jj) * e3t_n(ji,jj,jk) 
    195                   zqhcisf3d (ji,jj,jk) = zqhcisf3d (ji,jj,jk) + zqhcisf2d(ji,jj) * r1_hisf_tbl(ji,jj) * ralpha(ji,jj) * e3t_n(ji,jj,jk) 
    196                   zqlatisf3d(ji,jj,jk) = zqlatisf3d(ji,jj,jk) + qisf     (ji,jj) * r1_hisf_tbl(ji,jj) * ralpha(ji,jj) * e3t_n(ji,jj,jk) 
     194                  zfwfisf3d (ji,jj,jk) = zfwfisf3d (ji,jj,jk) + fwfisf   (ji,jj) * r1_hisf_tbl(ji,jj)   &  
     195                     &                                                                   * ralpha(ji,jj) * e3t_n(ji,jj,jk) 
     196                  zqhcisf3d (ji,jj,jk) = zqhcisf3d (ji,jj,jk) + zqhcisf2d(ji,jj) * r1_hisf_tbl(ji,jj)   &  
     197                     &                                                                   * ralpha(ji,jj) * e3t_n(ji,jj,jk) 
     198                  zqlatisf3d(ji,jj,jk) = zqlatisf3d(ji,jj,jk) + qisf     (ji,jj) * r1_hisf_tbl(ji,jj)   &   
     199                     &                                                                   * ralpha(ji,jj) * e3t_n(ji,jj,jk) 
    197200               END DO 
    198201            END DO 
  • branches/2017/dev_METO_2017/NEMOGCM/NEMO/OPA_SRC/SBC/sbcmod.F90

    r7822 r8868  
    4141   USE sbcssr         ! surface boundary condition: sea surface restoring 
    4242   USE sbcrnf         ! surface boundary condition: runoffs 
     43   USE sbcapr         ! surface boundary condition: atmo pressure  
    4344   USE sbcisf         ! surface boundary condition: ice shelf 
    4445   USE sbcfwb         ! surface boundary condition: freshwater budget 
     
    332333                          CALL sbc_rnf_init            ! Runof initialization 
    333334      ! 
     335      IF( ln_apr_dyn )    CALL sbc_apr_init            ! Atmo Pressure Forcing initialization 
     336      ! 
    334337      IF( nn_ice == 3 )   CALL sbc_lim_init            ! LIM3 initialization 
    335338      ! 
  • branches/2017/dev_METO_2017/NEMOGCM/NEMO/OPA_SRC/TRA/trabbl.F90

    r7753 r8868  
    545545      CALL wrk_dealloc( jpi, jpj, zmbk ) 
    546546 
    547       !                                 !* sign of grad(H) at u- and v-points 
    548       mgrhu(jpi,:) = 0   ;   mgrhu(:,jpj) = 0   ;   mgrhv(jpi,:) = 0   ;   mgrhv(:,jpj) = 0 
     547                                        !* sign of grad(H) at u- and v-points; zero if grad(H) = 0 
     548      mgrhu(:,:) = 0   ;   mgrhv(:,:) = 0 
    549549      DO jj = 1, jpjm1 
    550550         DO ji = 1, jpim1 
    551             mgrhu(ji,jj) = INT(  SIGN( 1.e0, gdept_0(ji+1,jj,mbkt(ji+1,jj)) - gdept_0(ji,jj,mbkt(ji,jj)) )  ) 
    552             mgrhv(ji,jj) = INT(  SIGN( 1.e0, gdept_0(ji,jj+1,mbkt(ji,jj+1)) - gdept_0(ji,jj,mbkt(ji,jj)) )  ) 
     551            IF( gdept_0(ji+1,jj,mbkt(ji+1,jj)) - gdept_0(ji,jj,mbkt(ji,jj)) /= 0._wp ) THEN 
     552               mgrhu(ji,jj) = INT(  SIGN( 1.e0, gdept_0(ji+1,jj,mbkt(ji+1,jj)) - gdept_0(ji,jj,mbkt(ji,jj)) )  ) 
     553            ENDIF 
     554            ! 
     555            IF( gdept_0(ji,jj+1,mbkt(ji,jj+1)) - gdept_0(ji,jj,mbkt(ji,jj)) /= 0._wp ) THEN 
     556               mgrhv(ji,jj) = INT(  SIGN( 1.e0, gdept_0(ji,jj+1,mbkt(ji,jj+1)) - gdept_0(ji,jj,mbkt(ji,jj)) )  ) 
     557            ENDIF 
    553558         END DO 
    554559      END DO 
  • branches/2017/dev_METO_2017/NEMOGCM/NEMO/OPA_SRC/TRA/tranxt.F90

    r7753 r8868  
    121121      IF( l_trdtra )   THEN                     
    122122         CALL wrk_alloc( jpi, jpj, jpk, ztrdt, ztrds ) 
    123          ztrdt(:,:,jk) = 0._wp 
    124          ztrds(:,:,jk) = 0._wp 
     123         ztrdt(:,:,jpk) = 0._wp 
     124         ztrds(:,:,jpk) = 0._wp 
    125125         IF( ln_traldf_iso ) THEN              ! diagnose the "pure" Kz diffusive trend  
    126126            CALL trd_tra( kt, 'TRA', jp_tem, jptra_zdfp, ztrdt ) 
     
    128128         ENDIF 
    129129         ! total trend for the non-time-filtered variables.  
    130             zfact = 1.0 / rdt 
     130         zfact = 1.0 / rdt 
     131         ! G Nurser 23 Mar 2017. Recalculate trend as Delta(e3t*T)/e3tn; e3tn cancel from tsn terms 
    131132         DO jk = 1, jpkm1 
    132             ztrdt(:,:,jk) = ( tsa(:,:,jk,jp_tem) - tsn(:,:,jk,jp_tem) ) * zfact  
    133             ztrds(:,:,jk) = ( tsa(:,:,jk,jp_sal) - tsn(:,:,jk,jp_sal) ) * zfact  
     133            ztrdt(:,:,jk) = ( tsa(:,:,jk,jp_tem)*e3t_a(:,:,jk) / e3t_n(:,:,jk) - tsn(:,:,jk,jp_tem)) * zfact 
     134            ztrds(:,:,jk) = ( tsa(:,:,jk,jp_sal)*e3t_a(:,:,jk) / e3t_n(:,:,jk) - tsn(:,:,jk,jp_sal)) * zfact 
    134135         END DO 
    135136         CALL trd_tra( kt, 'TRA', jp_tem, jptra_tot, ztrdt ) 
    136137         CALL trd_tra( kt, 'TRA', jp_sal, jptra_tot, ztrds ) 
    137          ! Store now fields before applying the Asselin filter  
    138          ! in order to calculate Asselin filter trend later. 
    139          ztrdt(:,:,:) = tsn(:,:,:,jp_tem)  
    140          ztrds(:,:,:) = tsn(:,:,:,jp_sal) 
     138         IF( ln_linssh ) THEN  
     139            ! Store now fields before applying the Asselin filter  
     140            ! in order to calculate Asselin filter trend later. 
     141            ztrdt(:,:,:) = tsn(:,:,:,jp_tem)  
     142            ztrds(:,:,:) = tsn(:,:,:,jp_sal) 
     143         ENDIF 
    141144      ENDIF 
    142145 
     
    147150            END DO 
    148151         END DO 
     152         IF (l_trdtra .AND. .NOT. ln_linssh) THEN  ! Zero Asselin filter contribution must be explicitly written out since for vvl 
     153                                                   ! Asselin filter is output by tra_nxt_vvl that is not called on this time step 
     154            ztrdt(:,:,:) = 0._wp 
     155            ztrds(:,:,:) = 0._wp 
     156            CALL trd_tra( kt, 'TRA', jp_tem, jptra_atf, ztrdt ) 
     157            CALL trd_tra( kt, 'TRA', jp_sal, jptra_atf, ztrds ) 
     158         END IF 
    149159         ! 
    150160      ELSE                                            ! Leap-Frog + Asselin filter time stepping 
     
    162172      ENDIF      
    163173      ! 
    164       IF( l_trdtra ) THEN      ! trend of the Asselin filter (tb filtered - tb)/dt      
     174      IF( l_trdtra .AND. ln_linssh ) THEN      ! trend of the Asselin filter (tb filtered - tb)/dt      
     175         zfact = 1._wp / r2dt              
    165176         DO jk = 1, jpkm1 
    166             zfact = 1._wp / r2dt              
    167177            ztrdt(:,:,jk) = ( tsb(:,:,jk,jp_tem) - ztrdt(:,:,jk) ) * zfact 
    168178            ztrds(:,:,jk) = ( tsb(:,:,jk,jp_sal) - ztrds(:,:,jk) ) * zfact 
     
    170180         CALL trd_tra( kt, 'TRA', jp_tem, jptra_atf, ztrdt ) 
    171181         CALL trd_tra( kt, 'TRA', jp_sal, jptra_atf, ztrds ) 
    172          CALL wrk_dealloc( jpi, jpj, jpk, ztrdt, ztrds ) 
    173182      END IF 
     183      IF( l_trdtra ) CALL wrk_dealloc( jpi, jpj, jpk, ztrdt, ztrds ) 
    174184      ! 
    175185      !                        ! control print 
     
    259269      LOGICAL  ::   ll_traqsr, ll_rnf, ll_isf   ! local logical 
    260270      INTEGER  ::   ji, jj, jk, jn              ! dummy loop indices 
    261       REAL(wp) ::   zfact1, ztc_a , ztc_n , ztc_b , ztc_f , ztc_d    ! local scalar 
     271      REAL(wp) ::   zfact, zfact1, ztc_a , ztc_n , ztc_b , ztc_f , ztc_d    ! local scalar 
    262272      REAL(wp) ::   zfact2, ze3t_b, ze3t_n, ze3t_a, ze3t_f, ze3t_d   !   -      - 
     273      REAL(wp), POINTER, DIMENSION(:,:,:,:) :: ztrd_atf 
    263274      !!---------------------------------------------------------------------- 
    264275      ! 
     
    279290      ENDIF 
    280291      ! 
     292      IF( ( l_trdtra .and. cdtype == 'TRA' ) .OR. ( l_trdtrc .and. cdtype == 'TRC' ) )   THEN 
     293         CALL wrk_alloc( jpi, jpj, jpk, kjpt, ztrd_atf ) 
     294         ztrd_atf(:,:,:,:) = 0.0_wp 
     295      ENDIF 
     296      zfact = 1._wp / r2dt 
    281297      DO jn = 1, kjpt       
    282298         DO jk = 1, jpkm1 
     
    331347                  ptn(ji,jj,jk,jn) = pta(ji,jj,jk,jn)     ! ptn <-- pta 
    332348                  ! 
     349                  IF( ( l_trdtra .and. cdtype == 'TRA' ) .OR. ( l_trdtrc .and. cdtype == 'TRC' ) ) THEN 
     350                     ztrd_atf(ji,jj,jk,jn) = (ztc_f - ztc_n) * zfact/ze3t_n 
     351                  ENDIF 
     352                  ! 
    333353               END DO 
    334354            END DO 
     
    337357      END DO 
    338358      ! 
     359      IF( l_trdtra .and. cdtype == 'TRA' ) THEN  
     360         CALL trd_tra( kt, cdtype, jp_tem, jptra_atf, ztrd_atf(:,:,:,jp_tem) ) 
     361         CALL trd_tra( kt, cdtype, jp_sal, jptra_atf, ztrd_atf(:,:,:,jp_sal) ) 
     362         CALL wrk_dealloc( jpi, jpj, jpk, kjpt, ztrd_atf ) 
     363      ENDIF 
     364      IF( l_trdtrc .and. cdtype == 'TRC' ) THEN 
     365         DO jn = 1, kjpt 
     366            CALL trd_tra( kt, cdtype, jn, jptra_atf, ztrd_atf(:,:,:,jn) ) 
     367         END DO 
     368         CALL wrk_dealloc( jpi, jpj, jpk, kjpt, ztrd_atf ) 
     369      ENDIF 
     370      ! 
    339371   END SUBROUTINE tra_nxt_vvl 
    340372 
  • branches/2017/dev_METO_2017/NEMOGCM/NEMO/OPA_SRC/TRA/trasbc.F90

    r7788 r8868  
    2727   USE trd_oce        ! trends: ocean variables 
    2828   USE trdtra         ! trends manager: tracers  
     29#if defined key_asminc    
     30   USE asminc         ! Assimilation increment 
     31#endif 
    2932   ! 
    3033   USE in_out_manager ! I/O manager 
     
    7275      INTEGER, INTENT(in) ::   kt   ! ocean time-step index 
    7376      ! 
    74       INTEGER  ::   ji, jj, jk, jn        ! dummy loop indices   
    75       INTEGER  ::   ikt, ikb              ! local integers 
    76       REAL(wp) ::   zfact, z1_e3t, zdep   ! local scalar 
     77      INTEGER  ::   ji, jj, jk, jn              ! dummy loop indices   
     78      INTEGER  ::   ikt, ikb                    ! local integers 
     79      REAL(wp) ::   zfact, z1_e3t, zdep, ztim   ! local scalar 
    7780      REAL(wp), POINTER, DIMENSION(:,:,:) ::  ztrdt, ztrds 
    7881      !!---------------------------------------------------------------------- 
     
    208211      IF( iom_use('rnf_x_sss') )   CALL iom_put( "rnf_x_sss", rnf*tsn(:,:,1,jp_sal) )   ! runoff term on sss 
    209212 
     213#if defined key_asminc 
     214      ! 
     215      !---------------------------------------- 
     216      !        Assmilation effects 
     217      !---------------------------------------- 
     218      ! 
     219      IF( ln_sshinc ) THEN         ! input of heat and salt due to assimilation 
     220          ! 
     221         IF( ln_linssh ) THEN  
     222            DO jj = 2, jpj  
     223               DO ji = fs_2, fs_jpim1 
     224                  ztim = ssh_iau(ji,jj) / e3t_n(ji,jj,1) 
     225                  tsa(ji,jj,1,jp_tem) = tsa(ji,jj,1,jp_tem) + tsn(ji,jj,1,jp_tem) * ztim 
     226                  tsa(ji,jj,1,jp_sal) = tsa(ji,jj,1,jp_sal) + tsn(ji,jj,1,jp_sal) * ztim 
     227               END DO 
     228            END DO 
     229         ELSE 
     230            DO jj = 2, jpj  
     231               DO ji = fs_2, fs_jpim1 
     232                  ztim = ssh_iau(ji,jj) / ( ht_n(ji,jj) + 1. - ssmask(ji, jj) ) 
     233                  tsa(ji,jj,:,jp_tem) = tsa(ji,jj,:,jp_tem) + tsn(ji,jj,:,jp_tem) * ztim 
     234                  tsa(ji,jj,:,jp_sal) = tsa(ji,jj,:,jp_sal) + tsn(ji,jj,:,jp_sal) * ztim 
     235               END DO   
     236            END DO   
     237         ENDIF 
     238         ! 
     239      ENDIF 
     240      ! 
     241#endif 
     242 
    210243      ! 
    211244      !---------------------------------------- 
  • branches/2017/dev_METO_2017/NEMOGCM/NEMO/OPA_SRC/TRA/trazdf.F90

    r7753 r8868  
    8989      IF( l_trdtra )   THEN                      ! save the vertical diffusive trends for further diagnostics 
    9090         DO jk = 1, jpkm1 
    91             ztrdt(:,:,jk) = ( ( tsa(:,:,jk,jp_tem) - tsb(:,:,jk,jp_tem) ) / r2dt ) - ztrdt(:,:,jk) 
    92             ztrds(:,:,jk) = ( ( tsa(:,:,jk,jp_sal) - tsb(:,:,jk,jp_sal) ) / r2dt ) - ztrds(:,:,jk) 
     91            ztrdt(:,:,jk) = ( ( tsa(:,:,jk,jp_tem)*e3t_a(:,:,jk) - tsb(:,:,jk,jp_tem)*e3t_b(:,:,jk) ) & 
     92                 & / (e3t_n(:,:,jk)*r2dt) ) - ztrdt(:,:,jk) 
     93            ztrds(:,:,jk) = ( ( tsa(:,:,jk,jp_sal)*e3t_a(:,:,jk) - tsb(:,:,jk,jp_sal)*e3t_b(:,:,jk) ) & 
     94                 & / (e3t_n(:,:,jk)*r2dt) ) - ztrds(:,:,jk) 
    9395         END DO 
    9496!!gm this should be moved in trdtra.F90 and done on all trends 
  • branches/2017/dev_METO_2017/NEMOGCM/NEMO/OPA_SRC/TRD/trdtra.F90

    r7646 r8868  
    104104                                 ztrds(:,:,:) = 0._wp 
    105105                                 CALL trd_tra_mng( trdt, ztrds, ktrd, kt ) 
     106         CASE( jptra_evd )   ;   avt_evd(:,:,:) = ptrd(:,:,:) * tmask(:,:,:) 
    106107         CASE DEFAULT                 ! other trends: masked trends 
    107108            trdt(:,:,:) = ptrd(:,:,:) * tmask(:,:,:)              ! mask & store 
     
    311312!!gm Rq: mask the trends already masked in trd_tra, but lbc_lnk should probably be added 
    312313      ! 
     314      ! Trends evaluated every time step that could go to the standard T file and can be output every ts into a 1ts file if 1ts output is selected 
    313315      SELECT CASE( ktrd ) 
    314       CASE( jptra_xad  )   ;   CALL iom_put( "ttrd_xad" , ptrdx )        ! x- horizontal advection 
    315                                CALL iom_put( "strd_xad" , ptrdy ) 
    316       CASE( jptra_yad  )   ;   CALL iom_put( "ttrd_yad" , ptrdx )        ! y- horizontal advection 
    317                                CALL iom_put( "strd_yad" , ptrdy ) 
    318       CASE( jptra_zad  )   ;   CALL iom_put( "ttrd_zad" , ptrdx )        ! z- vertical   advection 
    319                                CALL iom_put( "strd_zad" , ptrdy ) 
    320                                IF( ln_linssh ) THEN                   ! cst volume : adv flux through z=0 surface 
    321                                   CALL wrk_alloc( jpi, jpj, z2dx, z2dy ) 
    322                                   z2dx(:,:) = wn(:,:,1) * tsn(:,:,1,jp_tem) / e3t_n(:,:,1) 
    323                                   z2dy(:,:) = wn(:,:,1) * tsn(:,:,1,jp_sal) / e3t_n(:,:,1) 
    324                                   CALL iom_put( "ttrd_sad", z2dx ) 
    325                                   CALL iom_put( "strd_sad", z2dy ) 
    326                                   CALL wrk_dealloc( jpi, jpj, z2dx, z2dy ) 
    327                                ENDIF 
    328       CASE( jptra_totad  ) ;   CALL iom_put( "ttrd_totad" , ptrdx )      ! total   advection 
    329                                CALL iom_put( "strd_totad" , ptrdy ) 
    330       CASE( jptra_ldf  )   ;   CALL iom_put( "ttrd_ldf" , ptrdx )        ! lateral diffusion 
    331                                CALL iom_put( "strd_ldf" , ptrdy ) 
    332       CASE( jptra_zdf  )   ;   CALL iom_put( "ttrd_zdf" , ptrdx )        ! vertical diffusion (including Kz contribution) 
    333                                CALL iom_put( "strd_zdf" , ptrdy ) 
    334       CASE( jptra_zdfp )   ;   CALL iom_put( "ttrd_zdfp", ptrdx )        ! PURE vertical diffusion (no isoneutral contribution) 
    335                                CALL iom_put( "strd_zdfp", ptrdy ) 
    336       CASE( jptra_evd )    ;   CALL iom_put( "ttrd_evd", ptrdx )         ! EVD trend (convection) 
    337                                CALL iom_put( "strd_evd", ptrdy ) 
    338       CASE( jptra_dmp  )   ;   CALL iom_put( "ttrd_dmp" , ptrdx )        ! internal restoring (damping) 
    339                                CALL iom_put( "strd_dmp" , ptrdy ) 
    340       CASE( jptra_bbl  )   ;   CALL iom_put( "ttrd_bbl" , ptrdx )        ! bottom boundary layer 
    341                                CALL iom_put( "strd_bbl" , ptrdy ) 
    342       CASE( jptra_npc  )   ;   CALL iom_put( "ttrd_npc" , ptrdx )        ! static instability mixing 
    343                                CALL iom_put( "strd_npc" , ptrdy ) 
    344       CASE( jptra_nsr  )   ;   CALL iom_put( "ttrd_qns" , ptrdx(:,:,1) )        ! surface forcing + runoff (ln_rnf=T) 
    345                                CALL iom_put( "strd_cdt" , ptrdy(:,:,1) )        ! output as 2D surface fields 
    346       CASE( jptra_qsr  )   ;   CALL iom_put( "ttrd_qsr" , ptrdx )        ! penetrative solar radiat. (only on temperature) 
    347       CASE( jptra_bbc  )   ;   CALL iom_put( "ttrd_bbc" , ptrdx )        ! geothermal heating   (only on temperature) 
    348       CASE( jptra_atf  )   ;   CALL iom_put( "ttrd_atf" , ptrdx )        ! asselin time Filter 
    349                                CALL iom_put( "strd_atf" , ptrdy ) 
    350       CASE( jptra_tot  )   ;   CALL iom_put( "ttrd_tot" , ptrdx )        ! model total trend 
    351                                CALL iom_put( "strd_tot" , ptrdy ) 
     316      ! This total trend is done every time step 
     317      CASE( jptra_tot  )   ;   CALL iom_put( "ttrd_tot" , ptrdx )           ! model total trend 
     318         CALL iom_put( "strd_tot" , ptrdy ) 
    352319      END SELECT 
     320 
     321      ! These trends are done every second time step. When 1ts output is selected must go different (2ts) file from standard T-file 
     322      IF( MOD( kt, 2 ) == 0 ) THEN 
     323         SELECT CASE( ktrd ) 
     324         CASE( jptra_xad  )   ;   CALL iom_put( "ttrd_xad" , ptrdx )        ! x- horizontal advection 
     325            CALL iom_put( "strd_xad" , ptrdy ) 
     326         CASE( jptra_yad  )   ;   CALL iom_put( "ttrd_yad" , ptrdx )        ! y- horizontal advection 
     327            CALL iom_put( "strd_yad" , ptrdy ) 
     328         CASE( jptra_zad  )   ;   CALL iom_put( "ttrd_zad" , ptrdx )        ! z- vertical   advection 
     329            CALL iom_put( "strd_zad" , ptrdy ) 
     330            IF( ln_linssh ) THEN                   ! cst volume : adv flux through z=0 surface 
     331               CALL wrk_alloc( jpi, jpj, z2dx, z2dy ) 
     332               z2dx(:,:) = wn(:,:,1) * tsn(:,:,1,jp_tem) / e3t_n(:,:,1) 
     333               z2dy(:,:) = wn(:,:,1) * tsn(:,:,1,jp_sal) / e3t_n(:,:,1) 
     334               CALL iom_put( "ttrd_sad", z2dx ) 
     335               CALL iom_put( "strd_sad", z2dy ) 
     336               CALL wrk_dealloc( jpi, jpj, z2dx, z2dy ) 
     337            ENDIF 
     338         CASE( jptra_totad  ) ;   CALL iom_put( "ttrd_totad" , ptrdx )      ! total   advection 
     339            CALL iom_put( "strd_totad" , ptrdy ) 
     340         CASE( jptra_ldf  )   ;   CALL iom_put( "ttrd_ldf" , ptrdx )        ! lateral diffusion 
     341            CALL iom_put( "strd_ldf" , ptrdy ) 
     342         CASE( jptra_zdf  )   ;   CALL iom_put( "ttrd_zdf" , ptrdx )        ! vertical diffusion (including Kz contribution) 
     343            CALL iom_put( "strd_zdf" , ptrdy ) 
     344         CASE( jptra_zdfp )   ;   CALL iom_put( "ttrd_zdfp", ptrdx )        ! PURE vertical diffusion (no isoneutral contribution) 
     345            CALL iom_put( "strd_zdfp", ptrdy ) 
     346         CASE( jptra_evd )    ;   CALL iom_put( "ttrd_evd", ptrdx )         ! EVD trend (convection) 
     347            CALL iom_put( "strd_evd", ptrdy ) 
     348         CASE( jptra_dmp  )   ;   CALL iom_put( "ttrd_dmp" , ptrdx )        ! internal restoring (damping) 
     349            CALL iom_put( "strd_dmp" , ptrdy ) 
     350         CASE( jptra_bbl  )   ;   CALL iom_put( "ttrd_bbl" , ptrdx )        ! bottom boundary layer 
     351            CALL iom_put( "strd_bbl" , ptrdy ) 
     352         CASE( jptra_npc  )   ;   CALL iom_put( "ttrd_npc" , ptrdx )        ! static instability mixing 
     353            CALL iom_put( "strd_npc" , ptrdy ) 
     354         CASE( jptra_bbc  )   ;   CALL iom_put( "ttrd_bbc" , ptrdx )        ! geothermal heating   (only on temperature) 
     355         CASE( jptra_nsr  )   ;   CALL iom_put( "ttrd_qns" , ptrdx(:,:,1) ) ! surface forcing + runoff (ln_rnf=T) 
     356            CALL iom_put( "strd_cdt" , ptrdy(:,:,1) )        ! output as 2D surface fields 
     357         CASE( jptra_qsr  )   ;   CALL iom_put( "ttrd_qsr" , ptrdx )        ! penetrative solar radiat. (only on temperature) 
     358         END SELECT 
     359         ! the Asselin filter trend  is also every other time step but needs to be lagged one time step 
     360         ! Even when 1ts output is selected can go to the same (2ts) file as the trends plotted every even time step. 
     361      ELSE IF( MOD( kt, 2 ) == 1 ) THEN 
     362         SELECT CASE( ktrd ) 
     363         CASE( jptra_atf  )   ;   CALL iom_put( "ttrd_atf" , ptrdx )        ! asselin time Filter 
     364            CALL iom_put( "strd_atf" , ptrdy ) 
     365         END SELECT 
     366      END IF 
    353367      ! 
    354368   END SUBROUTINE trd_tra_iom 
  • branches/2017/dev_METO_2017/NEMOGCM/NEMO/OPA_SRC/ZDF/zdfbfr.F90

    r7753 r8868  
    384384            END IF 
    385385         ENDIF 
     386         ! 
     387         bfrua(:,:) = - bfrcoef2d(:,:) 
     388         bfrva(:,:) = - bfrcoef2d(:,:) 
     389         ! 
    386390         ! 
    387391      CASE DEFAULT 
  • branches/2017/dev_METO_2017/NEMOGCM/NEMO/OPA_SRC/ZDF/zdftmx.F90

    r7779 r8868  
    225225 
    226226      !                             ! compute the form function using N2 at each time step 
     227      zdn2dz     (:,:,jpk) = 0.e0 
    227228      zempba_3d_1(:,:,jpk) = 0.e0 
    228229      zempba_3d_2(:,:,jpk) = 0.e0 
  • branches/2017/dev_METO_2017/NEMOGCM/NEMO/OPA_SRC/nemogcm.F90

    r7761 r8868  
    206206      ! 
    207207#if defined key_agrif 
    208       IF( .NOT. Agrif_Root() ) THEN 
    209                          CALL Agrif_ParentGrid_To_ChildGrid() 
    210          IF( ln_diaobs ) CALL dia_obs_wri 
    211          IF( nn_timing == 1 )   CALL timing_finalize 
    212                                 CALL Agrif_ChildGrid_To_ParentGrid() 
    213       ENDIF 
     208      CALL Agrif_ParentGrid_To_ChildGrid() 
     209      IF( ln_diaobs ) CALL dia_obs_wri 
     210      IF( nn_timing == 1 )   CALL timing_finalize 
     211      CALL Agrif_ChildGrid_To_ParentGrid() 
    214212#endif 
    215213      IF( nn_timing == 1 )   CALL timing_finalize 
     
    452450      !                                      ! external forcing  
    453451!!gm to be added : creation and call of sbc_apr_init 
     452!==> cbr: sbc_apr_init in sbcmod as sbc_rnf_init 
    454453                            CALL    tide_init   ! tidal harmonics 
    455454                            CALL     sbc_init   ! surface boundary conditions (including sea-ice) 
     
    751750      ! 
    752751      ! Find the factors of n. 
    753       IF( kn == 1 )   GOTO 20 
    754  
    755       ! nu holds the unfactorised part of the number. 
    756       ! knfax holds the number of factors found. 
    757       ! l points to the allowed factor list. 
    758       ! ifac holds the current factor. 
    759       ! 
    760       inu   = kn 
    761       knfax = 0 
    762       ! 
    763       DO jl = ntest, 1, -1 
     752      IF( kn .NE. 1 ) THEN 
     753 
     754         ! nu holds the unfactorised part of the number. 
     755         ! knfax holds the number of factors found. 
     756         ! l points to the allowed factor list. 
     757         ! ifac holds the current factor. 
    764758         ! 
    765          ifac = ilfax(jl) 
    766          IF( ifac > inu )   CYCLE 
    767  
    768          ! Test whether the factor will divide. 
    769  
    770          IF( MOD(inu,ifac) == 0 ) THEN 
     759         inu   = kn 
     760         knfax = 0 
     761         ! 
     762         DO jl = ntest, 1, -1 
    771763            ! 
    772             knfax = knfax + 1            ! Add the factor to the list 
    773             IF( knfax > kmaxfax ) THEN 
    774                kerr = 6 
    775                write (*,*) 'FACTOR: insufficient space in factor array ', knfax 
    776                return 
     764            ifac = ilfax(jl) 
     765            IF( ifac > inu )   CYCLE 
     766    
     767            ! Test whether the factor will divide. 
     768    
     769            IF( MOD(inu,ifac) == 0 ) THEN 
     770               ! 
     771               knfax = knfax + 1            ! Add the factor to the list 
     772               IF( knfax > kmaxfax ) THEN 
     773                  kerr = 6 
     774                  write (*,*) 'FACTOR: insufficient space in factor array ', knfax 
     775                  return 
     776               ENDIF 
     777               kfax(knfax) = ifac 
     778               ! Store the other factor that goes with this one 
     779               knfax = knfax + 1 
     780               kfax(knfax) = inu / ifac 
     781               !WRITE (*,*) 'ARPDBG, factors ',knfax-1,' & ',knfax,' are ', kfax(knfax-1),' and ',kfax(knfax) 
    777782            ENDIF 
    778             kfax(knfax) = ifac 
    779             ! Store the other factor that goes with this one 
    780             knfax = knfax + 1 
    781             kfax(knfax) = inu / ifac 
    782             !WRITE (*,*) 'ARPDBG, factors ',knfax-1,' & ',knfax,' are ', kfax(knfax-1),' and ',kfax(knfax) 
    783          ENDIF 
     783            ! 
     784         END DO 
    784785         ! 
    785       END DO 
    786       ! 
    787    20 CONTINUE      ! Label 20 is the exit point from the factor search loop. 
     786      ENDIF 
    788787      ! 
    789788   END SUBROUTINE factorise 
  • branches/2017/dev_METO_2017/NEMOGCM/NEMO/SAO_SRC/nemogcm.F90

    r7646 r8868  
    499499      ! 
    500500      ! Find the factors of n. 
    501       IF( kn == 1 )   GOTO 20 
    502  
    503       ! nu holds the unfactorised part of the number. 
    504       ! knfax holds the number of factors found. 
    505       ! l points to the allowed factor list. 
    506       ! ifac holds the current factor. 
    507       ! 
    508       inu   = kn 
    509       knfax = 0 
    510       ! 
    511       DO jl = ntest, 1, -1 
    512          ! 
    513          ifac = ilfax(jl) 
    514          IF( ifac > inu )   CYCLE 
    515  
    516          ! Test whether the factor will divide. 
    517  
    518          IF( MOD(inu,ifac) == 0 ) THEN 
     501      IF( kn .NE. 1 ) THEN 
     502 
     503         ! nu holds the unfactorised part of the number. 
     504         ! knfax holds the number of factors found. 
     505         ! l points to the allowed factor list. 
     506         ! ifac holds the current factor. 
     507         ! 
     508         inu   = kn 
     509         knfax = 0 
     510         ! 
     511         DO jl = ntest, 1, -1 
    519512            ! 
    520             knfax = knfax + 1            ! Add the factor to the list 
    521             IF( knfax > kmaxfax ) THEN 
    522                kerr = 6 
    523                write (*,*) 'FACTOR: insufficient space in factor array ', knfax 
    524                return 
     513            ifac = ilfax(jl) 
     514            IF( ifac > inu )   CYCLE 
     515    
     516            ! Test whether the factor will divide. 
     517    
     518            IF( MOD(inu,ifac) == 0 ) THEN 
     519               ! 
     520               knfax = knfax + 1            ! Add the factor to the list 
     521               IF( knfax > kmaxfax ) THEN 
     522                  kerr = 6 
     523                  write (*,*) 'FACTOR: insufficient space in factor array ', knfax 
     524                  return 
     525               ENDIF 
     526               kfax(knfax) = ifac 
     527               ! Store the other factor that goes with this one 
     528               knfax = knfax + 1 
     529               kfax(knfax) = inu / ifac 
     530               !WRITE (*,*) 'ARPDBG, factors ',knfax-1,' & ',knfax,' are ', kfax(knfax-1),' and ',kfax(knfax) 
    525531            ENDIF 
    526             kfax(knfax) = ifac 
    527             ! Store the other factor that goes with this one 
    528             knfax = knfax + 1 
    529             kfax(knfax) = inu / ifac 
    530             !WRITE (*,*) 'ARPDBG, factors ',knfax-1,' & ',knfax,' are ', kfax(knfax-1),' and ',kfax(knfax) 
    531          ENDIF 
    532          ! 
    533       END DO 
    534       ! 
    535    20 CONTINUE      ! Label 20 is the exit point from the factor search loop. 
     532            ! 
     533         END DO 
     534         ! 
     535      ENDIF 
    536536      ! 
    537537   END SUBROUTINE factorise 
  • branches/2017/dev_METO_2017/NEMOGCM/NEMO/SAS_SRC/daymod.F90

    r7761 r8868  
    22   !!====================================================================== 
    33   !!                       ***  MODULE  daymod  *** 
    4    !! Ocean        :  calendar  
     4   !! Ocean :   management of the model calendar 
    55   !!===================================================================== 
    66   !! History :  OPA  ! 1994-09  (M. Pontaud M. Imbard)  Original code 
    77   !!                 ! 1997-03  (O. Marti) 
    8    !!                 ! 1997-05  (G. Madec)  
     8   !!                 ! 1997-05  (G. Madec) 
    99   !!                 ! 1997-08  (M. Imbard) 
    1010   !!   NEMO     1.0  ! 2003-09  (G. Madec)  F90 + nyear, nmonth, nday 
    1111   !!                 ! 2004-01  (A.M. Treguier) new calculation based on adatrj 
    1212   !!                 ! 2006-08  (G. Madec)  surface module major update 
    13    !!----------------------------------------------------------------------       
     13   !!                 ! 2015-11  (D. Lea) Allow non-zero initial time of day 
     14   !!---------------------------------------------------------------------- 
    1415 
    1516   !!---------------------------------------------------------------------- 
    1617   !!   day        : calendar 
    17    !!   
    18    !!           ------------------------------- 
    19    !!           ----------- WARNING ----------- 
    20    !! 
    21    !!   we suppose that the time step is deviding the number of second of in a day 
    22    !!             ---> MOD( rday, rdt ) == 0 
    23    !! 
    24    !!           ----------- WARNING ----------- 
    25    !!           ------------------------------- 
    26    !!   
    27    !!---------------------------------------------------------------------- 
    28    USE dom_oce         ! ocean space and time domain 
    29    USE phycst          ! physical constants 
    30    USE in_out_manager  ! I/O manager 
    31    USE iom             !  
    32    USE ioipsl, ONLY :   ymds2ju   ! for calendar 
    33    USE prtctl          ! Print control 
    34    USE restart         !  
    35    USE timing          ! Timing 
     18   !!---------------------------------------------------------------------- 
     19   !!                    ----------- WARNING ----------- 
     20   !!                    ------------------------------- 
     21   !!   sbcmod assume that the time step is dividing the number of second of  
     22   !!   in a day, i.e. ===> MOD( rday, rdt ) == 0  
     23   !!   except when user defined forcing is used (see sbcmod.F90) 
     24   !!---------------------------------------------------------------------- 
     25   USE dom_oce        ! ocean space and time domain 
     26   USE phycst         ! physical constants 
     27   USE ioipsl  , ONLY :   ymds2ju      ! for calendar 
     28   USE trc_oce , ONLY :   l_offline   ! offline flag 
     29   ! 
     30   USE in_out_manager ! I/O manager 
     31   USE prtctl         ! Print control 
     32   USE iom            ! 
     33   USE timing         ! Timing 
     34   USE restart        ! restart 
    3635 
    3736   IMPLICIT NONE 
     
    4039   PUBLIC   day        ! called by step.F90 
    4140   PUBLIC   day_init   ! called by istate.F90 
    42  
    43    INTEGER ::   nsecd, nsecd05, ndt, ndt05 
    44  
    45    !!---------------------------------------------------------------------- 
    46    !! NEMO/OPA 3.3 , NEMO Consortium (2010) 
     41   PUBLIC   day_mth    ! Needed by TAM 
     42 
     43   INTEGER, PUBLIC ::   nsecd, nsecd05, ndt, ndt05   !: (PUBLIC for TAM) 
     44 
     45   !!---------------------------------------------------------------------- 
     46   !! NEMO/OPA 4.0 , NEMO Consortium (2016) 
    4747   !! $Id$ 
    4848   !! Software governed by the CeCILL licence     (NEMOGCM/NEMO_CeCILL.txt) 
     
    5353      !!---------------------------------------------------------------------- 
    5454      !!                   ***  ROUTINE day_init  *** 
    55       !!  
    56       !! ** Purpose :   Initialization of the calendar values to their values 1 time step before nit000  
     55      !! 
     56      !! ** Purpose :   Initialization of the calendar values to their values 1 time step before nit000 
    5757      !!                because day will be called at the beginning of step 
    5858      !! 
     
    6767      !!              - nmonth_len, nyear_len, nmonth_half, nmonth_end through day_mth 
    6868      !!---------------------------------------------------------------------- 
    69       INTEGER  ::   inbday, idweek 
    70       REAL(wp) ::   zjul 
     69      INTEGER  ::   inbday, idweek   ! local integers 
     70      REAL(wp) ::   zjul             ! local scalar 
    7171      !!---------------------------------------------------------------------- 
    7272      ! 
     
    7676            &           'You must do a restart at higher frequency (or remove this stop and recompile the code in I8)' ) 
    7777      ENDIF 
    78       ! all calendar staff is based on the fact that MOD( rday, rdt ) == 0 
    79       IF( MOD( rday , rdt ) /= 0. )   CALL ctl_stop( 'the time step must devide the number of second of in a day' ) 
    80       IF( MOD( rday , 2.  ) /= 0. )   CALL ctl_stop( 'the number of second of in a day must be an even number'    ) 
    81       IF( MOD( rdt  , 2.  ) /= 0. )   CALL ctl_stop( 'the time step (in second) must be an even number'           ) 
    82       nsecd   = NINT(rday       ) 
    83       nsecd05 = NINT(0.5 * rday ) 
    84       ndt     = NINT(      rdt  ) 
    85       ndt05   = NINT(0.5 * rdt  ) 
    86  
    87       ! ==> clem: here we read the ocean restart for the date (only if it exists) 
    88       !           It is not clean and another solution should be found 
    89       CALL day_rst( nit000, 'READ' ) 
    90       ! ==> 
    91  
    92       ! set the calendar from ndastp (read in restart file and namelist) 
    93  
     78      nsecd   = NINT( rday       ) 
     79      nsecd05 = NINT( 0.5 * rday ) 
     80      ndt     = NINT(       rdt  ) 
     81      ndt05   = NINT( 0.5 * rdt  ) 
     82 
     83      IF( .NOT. l_offline )   CALL day_rst( nit000, 'READ' ) 
     84 
     85      ! set the calandar from ndastp (read in restart file and namelist) 
    9486      nyear   =   ndastp / 10000 
    9587      nmonth  = ( ndastp - (nyear * 10000) ) / 100 
    96       nday    =   ndastp - (nyear * 10000) - ( nmonth * 100 )  
    97  
    98       CALL ymds2ju( nyear, nmonth, nday, 0.0, fjulday )  ! we assume that we start run at 00:00 
     88      nday    =   ndastp - (nyear * 10000) - ( nmonth * 100 ) 
     89 
     90      nhour   =   nn_time0 / 100 
     91      nminute = ( nn_time0 - nhour * 100 ) 
     92 
     93      CALL ymds2ju( nyear, nmonth, nday, nhour*3600._wp+nminute*60._wp, fjulday )   
    9994      IF( ABS(fjulday - REAL(NINT(fjulday),wp)) < 0.1 / rday )   fjulday = REAL(NINT(fjulday),wp)   ! avoid truncation error 
    100       fjulday = fjulday + 1.                             ! move back to the day at nit000 (and not at nit000 - 1) 
     95      IF( nn_time0*3600 - ndt05 .lt. 0 ) fjulday = fjulday + 1.                    ! move back to the day at nit000 (and not at nit000 - 1) 
    10196 
    10297      nsec1jan000 = 0 
    10398      CALL day_mth 
    104        
     99 
    105100      IF ( nday == 0 ) THEN     !   for ex if ndastp = ndate0 - 1 
    106          nmonth = nmonth - 1   
     101         nmonth = nmonth - 1 
    107102         nday = nmonth_len(nmonth) 
    108103      ENDIF 
     
    113108         IF( nleapy == 1 )   CALL day_mth 
    114109      ENDIF 
    115        
     110 
    116111      ! day since january 1st 
    117112      nday_year = nday + SUM( nmonth_len(1:nmonth - 1) ) 
    118113 
    119       !compute number of days between last monday and today       
     114      !compute number of days between last monday and today 
    120115      CALL ymds2ju( 1900, 01, 01, 0.0, zjul )  ! compute julian day value of 01.01.1900 (our reference that was a Monday) 
    121       inbday = NINT(fjulday - zjul)            ! compute nb day between  01.01.1900 and current day   
    122       idweek = MOD(inbday, 7)                  ! compute nb day between last monday and current day   
     116      inbday = FLOOR(fjulday - zjul)            ! compute nb day between  01.01.1900 and start of current day 
     117      idweek = MOD(inbday, 7)                  ! compute nb day between last monday and current day 
     118      IF (idweek .lt. 0) idweek=idweek+7       ! Avoid negative values for dates before 01.01.1900 
    123119 
    124120      ! number of seconds since the beginning of current year/month/week/day at the middle of the time-step 
    125       nsec_year  = nday_year * nsecd - ndt05   ! 1 time step before the middle of the first time step 
    126       nsec_month = nday      * nsecd - ndt05   ! because day will be called at the beginning of step 
    127       nsec_week  = idweek    * nsecd - ndt05 
    128       nsec_day   =             nsecd - ndt05 
     121      IF (nhour*3600+nminute*60-ndt05 .gt. 0) THEN 
     122         ! 1 timestep before current middle of first time step is still the same day 
     123         nsec_year  = (nday_year-1) * nsecd + nhour*3600+nminute*60 - ndt05  
     124         nsec_month = (nday-1)      * nsecd + nhour*3600+nminute*60 - ndt05     
     125      ELSE 
     126         ! 1 time step before the middle of the first time step is the previous day  
     127         nsec_year  = nday_year * nsecd + nhour*3600+nminute*60 - ndt05  
     128         nsec_month = nday      * nsecd + nhour*3600+nminute*60 - ndt05    
     129      ENDIF 
     130      nsec_week  = idweek    * nsecd + nhour*3600+nminute*60 - ndt05 
     131      nsec_day   =             nhour*3600+nminute*60 - ndt05  
     132      IF( nsec_day .lt. 0 ) nsec_day = nsec_day + nsecd 
     133      IF( nsec_week .lt. 0 ) nsec_week = nsec_week + nsecd*7 
    129134 
    130135      ! control print 
    131       IF(lwp) WRITE(numout,'(a,i6,a,i2,a,i2,a,i8,a,i8)')' =======>> 1/2 time step before the start of the run DATE Y/M/D = ',   & 
    132            &                   nyear, '/', nmonth, '/', nday, '  nsec_day:', nsec_day, '  nsec_week:', nsec_week 
     136      IF(lwp) WRITE(numout,'(a,i6,a,i2,a,i2,a,i8,a,i8,a,i8,a,i8)')   & 
     137           &                   ' =======>> 1/2 time step before the start of the run DATE Y/M/D = ',   & 
     138           &                   nyear, '/', nmonth, '/', nday, '  nsec_day:', nsec_day, '  nsec_week:', nsec_week, '  & 
     139           &                   nsec_month:', nsec_month , '  nsec_year:' , nsec_year 
    133140 
    134141      ! Up to now, calendar parameters are related to the end of previous run (nit000-1) 
     
    142149      !!---------------------------------------------------------------------- 
    143150      !!                   ***  ROUTINE day_init  *** 
    144       !!  
     151      !! 
    145152      !! ** Purpose :   calendar values related to the months 
    146153      !! 
     
    154161 
    155162      ! length of the month of the current year (from nleapy, read in namelist) 
    156       IF ( nleapy < 2 ) THEN  
     163      IF ( nleapy < 2 ) THEN 
    157164         nmonth_len(:) = (/ 31, 31, 28, 31, 30, 31, 30, 31, 31, 30, 31, 30, 31, 31 /) 
    158165         nyear_len(:) = 365 
     
    177184      ! time since Jan 1st   0     1     2    ...    11    12    13 
    178185      !          ---------*--|--*--|--*--| ... |--*--|--*--|--*--|-------------------------------------- 
    179       !                 <---> <---> <--->  ...  <---> <---> <--->         
     186      !                 <---> <---> <--->  ...  <---> <---> <---> 
    180187      ! month number      0     1     2    ...    11    12    13 
    181188      ! 
     
    190197         nmonth_end(jm) = nmonth_end(jm-1) + nsecd * nmonth_len(jm) 
    191198      END DO 
    192       !            
    193    END SUBROUTINE  
     199      ! 
     200   END SUBROUTINE 
    194201 
    195202 
     
    197204      !!---------------------------------------------------------------------- 
    198205      !!                      ***  ROUTINE day  *** 
    199       !!  
     206      !! 
    200207      !! ** Purpose :   Compute the date with a day iteration IF necessary. 
    201208      !! 
     
    209216      !!              - adatrj    : date in days since the beginning of the run 
    210217      !!              - nsec_year : current time of the year (in second since 00h, jan 1st) 
    211       !!----------------------------------------------------------------------       
     218      !!---------------------------------------------------------------------- 
    212219      INTEGER, INTENT(in) ::   kt        ! ocean time-step indices 
    213220      ! 
     
    220227      zprec = 0.1 / rday 
    221228      !                                                 ! New time-step 
    222       nsec_year  = nsec_year  + ndt  
    223       nsec_month = nsec_month + ndt                  
     229      nsec_year  = nsec_year  + ndt 
     230      nsec_month = nsec_month + ndt 
    224231      nsec_week  = nsec_week  + ndt 
    225       nsec_day   = nsec_day   + ndt                 
     232      nsec_day   = nsec_day   + ndt 
    226233      adatrj  = adatrj  + rdt / rday 
    227234      fjulday = fjulday + rdt / rday 
    228235      IF( ABS(fjulday - REAL(NINT(fjulday),wp)) < zprec )   fjulday = REAL(NINT(fjulday),wp)   ! avoid truncation error 
    229236      IF( ABS(adatrj  - REAL(NINT(adatrj ),wp)) < zprec )   adatrj  = REAL(NINT(adatrj ),wp)   ! avoid truncation error 
    230        
     237 
    231238      IF( nsec_day > nsecd ) THEN                       ! New day 
    232239         ! 
     
    261268 
    262269      IF( nsec_week > 7*nsecd )   nsec_week = ndt05     ! New week 
    263        
     270 
    264271      IF(ln_ctl) THEN 
    265272         WRITE(charout,FMT="('kt =', I4,'  d/m/y =',I2,I2,I4)") kt, nday, nmonth, nyear 
     
    267274      ENDIF 
    268275 
    269       ! since we no longer call rst_opn, need to define nitrst here, used by ice restart routine 
    270       IF( kt == nit000 )  THEN 
    271          nitrst = nitend 
    272          lrst_oce = .FALSE.  ! init restart ocean (done in rst_opn when not SAS) 
    273       ENDIF 
    274  
    275       IF( MOD( kt - 1, nstock ) == 0 ) THEN 
    276          ! we use kt - 1 and not kt - nit000 to keep the same periodicity from the beginning of the experiment 
    277          nitrst = kt + nstock - 1                  ! define the next value of nitrst for restart writing 
    278          IF( nitrst > nitend )   nitrst = nitend   ! make sure we write a restart at the end of the run 
    279       ENDIF 
    280  
     276      IF( .NOT. l_offline ) CALL rst_opn( kt )               ! Open the restart file if needed and control lrst_oce 
     277      IF( lrst_oce         ) CALL day_rst( kt, 'WRITE' )      ! write day restart information 
     278      ! 
    281279      IF( nn_timing == 1 )  CALL timing_stop('day') 
    282280      ! 
     
    312310      CHARACTER(len=*), INTENT(in) ::   cdrw       ! "READ"/"WRITE" flag 
    313311      ! 
    314       REAL(wp) ::   zkt, zndastp 
     312      REAL(wp) ::   zkt, zndastp, zdayfrac, ksecs, ktime 
     313      INTEGER  ::   ihour, iminute 
    315314      !!---------------------------------------------------------------------- 
    316315 
     
    337336            ! define ndastp and adatrj 
    338337            IF ( nrstdt == 2 ) THEN 
    339                ! read the parameters correspondting to nit000 - 1 (last time step of previous run) 
     338               ! read the parameters corresponding to nit000 - 1 (last time step of previous run) 
    340339               CALL iom_get( numror, 'ndastp', zndastp ) 
    341340               ndastp = NINT( zndastp ) 
    342341               CALL iom_get( numror, 'adatrj', adatrj  ) 
     342          CALL iom_get( numror, 'ntime', ktime ) 
     343          nn_time0=INT(ktime) 
     344               ! calculate start time in hours and minutes 
     345          zdayfrac=adatrj-INT(adatrj) 
     346          ksecs = NINT(zdayfrac*86400)        ! Nearest second to catch rounding errors in adatrj          
     347          ihour = INT(ksecs/3600) 
     348          iminute = ksecs/60-ihour*60 
     349            
     350               ! Add to nn_time0 
     351               nhour   =   nn_time0 / 100 
     352               nminute = ( nn_time0 - nhour * 100 ) 
     353          nminute=nminute+iminute 
     354           
     355          IF( nminute >= 60 ) THEN 
     356             nminute=nminute-60 
     357        nhour=nhour+1 
     358          ENDIF 
     359          nhour=nhour+ihour 
     360          IF( nhour >= 24 ) THEN 
     361        nhour=nhour-24 
     362             adatrj=adatrj+1 
     363          ENDIF           
     364          nn_time0 = nhour * 100 + nminute 
     365          adatrj = INT(adatrj)                    ! adatrj set to integer as nn_time0 updated           
    343366            ELSE 
    344                ! parameters correspondting to nit000 - 1 (as we start the step loop with a call to day) 
    345                ndastp = ndate0 - 1     ! ndate0 read in the namelist in dom_nam, we assume that we start run at 00:00 
     367               ! parameters corresponding to nit000 - 1 (as we start the step loop with a call to day) 
     368               ndastp = ndate0        ! ndate0 read in the namelist in dom_nam 
     369               nhour   =   nn_time0 / 100 
     370               nminute = ( nn_time0 - nhour * 100 ) 
     371               IF( nhour*3600+nminute*60-ndt05 .lt. 0 )  ndastp=ndastp-1      ! Start hour is specified in the namelist (default 0) 
    346372               adatrj = ( REAL( nit000-1, wp ) * rdt ) / rday 
    347373               ! note this is wrong if time step has changed during run 
    348374            ENDIF 
    349375         ELSE 
    350             ! parameters correspondting to nit000 - 1 (as we start the step loop with a call to day) 
    351             ndastp = ndate0 - 1        ! ndate0 read in the namelist in dom_nam, we assume that we start run at 00:00 
     376            ! parameters corresponding to nit000 - 1 (as we start the step loop with a call to day) 
     377            ndastp = ndate0           ! ndate0 read in the namelist in dom_nam 
     378            nhour   =   nn_time0 / 100 
     379       nminute = ( nn_time0 - nhour * 100 ) 
     380            IF( nhour*3600+nminute*60-ndt05 .lt. 0 )  ndastp=ndastp-1      ! Start hour is specified in the namelist (default 0) 
    352381            adatrj = ( REAL( nit000-1, wp ) * rdt ) / rday 
    353382         ENDIF 
     
    358387            WRITE(numout,*) '   date ndastp                                      : ', ndastp 
    359388            WRITE(numout,*) '   number of elapsed days since the begining of run : ', adatrj 
     389       WRITE(numout,*) '   nn_time0                                         : ',nn_time0 
    360390            WRITE(numout,*) 
    361391         ENDIF 
     
    373403         CALL iom_rstput( kt, nitrst, numrow, 'adatrj' , adatrj            )   ! number of elapsed days since 
    374404         !                                                                     ! the begining of the run [s] 
     405    CALL iom_rstput( kt, nitrst, numrow, 'ntime'  , REAL( nn_time0, wp) ) ! time 
    375406      ENDIF 
    376407      ! 
    377408   END SUBROUTINE day_rst 
     409 
    378410   !!====================================================================== 
    379411END MODULE daymod 
  • branches/2017/dev_METO_2017/NEMOGCM/NEMO/SAS_SRC/nemogcm.F90

    r7761 r8868  
    2525   USE usrdef_nam     ! user defined configuration 
    2626   USE daymod         ! calendar 
     27   USE restart        ! open  restart file 
    2728   USE step           ! NEMO time-stepping                 (stp     routine) 
    2829   USE cpl_oasis3     ! 
     
    364365      IF( ln_ctl      )     CALL prt_ctl_init   ! Print control 
    365366                            CALL day_init   ! model calendar (using both namelist and restart infos) 
     367      IF( ln_rstart )       CALL rst_read_open 
    366368 
    367369                            CALL sbc_init   ! Forcings : surface module  
     
    596598      ! 
    597599      ! Find the factors of n. 
    598       IF( kn == 1 )   GOTO 20 
    599  
    600       ! nu holds the unfactorised part of the number. 
    601       ! knfax holds the number of factors found. 
    602       ! l points to the allowed factor list. 
    603       ! ifac holds the current factor. 
    604       ! 
    605       inu   = kn 
    606       knfax = 0 
    607       ! 
    608       DO jl = ntest, 1, -1 
    609          ! 
    610          ifac = ilfax(jl) 
    611          IF( ifac > inu )   CYCLE 
    612  
    613          ! Test whether the factor will divide. 
    614  
    615          IF( MOD(inu,ifac) == 0 ) THEN 
     600      IF( kn .NE. 1 ) THEN 
     601 
     602         ! nu holds the unfactorised part of the number. 
     603         ! knfax holds the number of factors found. 
     604         ! l points to the allowed factor list. 
     605         ! ifac holds the current factor. 
     606         ! 
     607         inu   = kn 
     608         knfax = 0 
     609         ! 
     610         DO jl = ntest, 1, -1 
    616611            ! 
    617             knfax = knfax + 1            ! Add the factor to the list 
    618             IF( knfax > kmaxfax ) THEN 
    619                kerr = 6 
    620                write (*,*) 'FACTOR: insufficient space in factor array ', knfax 
    621                return 
     612            ifac = ilfax(jl) 
     613            IF( ifac > inu )   CYCLE 
     614    
     615            ! Test whether the factor will divide. 
     616    
     617            IF( MOD(inu,ifac) == 0 ) THEN 
     618               ! 
     619               knfax = knfax + 1            ! Add the factor to the list 
     620               IF( knfax > kmaxfax ) THEN 
     621                  kerr = 6 
     622                  write (*,*) 'FACTOR: insufficient space in factor array ', knfax 
     623                  return 
     624               ENDIF 
     625               kfax(knfax) = ifac 
     626               ! Store the other factor that goes with this one 
     627               knfax = knfax + 1 
     628               kfax(knfax) = inu / ifac 
     629               !WRITE (*,*) 'ARPDBG, factors ',knfax-1,' & ',knfax,' are ', kfax(knfax-1),' and ',kfax(knfax) 
    622630            ENDIF 
    623             kfax(knfax) = ifac 
    624             ! Store the other factor that goes with this one 
    625             knfax = knfax + 1 
    626             kfax(knfax) = inu / ifac 
    627             !WRITE (*,*) 'ARPDBG, factors ',knfax-1,' & ',knfax,' are ', kfax(knfax-1),' and ',kfax(knfax) 
    628          ENDIF 
    629          ! 
    630       END DO 
    631       ! 
    632    20 CONTINUE      ! Label 20 is the exit point from the factor search loop. 
     631            ! 
     632         END DO 
     633         ! 
     634      ENDIF 
    633635      ! 
    634636   END SUBROUTINE factorise 
  • branches/2017/dev_METO_2017/NEMOGCM/NEMO/TOP_SRC/CFC/trcsms_cfc.F90

    r7646 r8868  
    180180      ! 
    181181      IF( lk_iomput ) THEN 
     182         jl = 0 
    182183         DO jn = jp_cfc0, jp_cfc1 
    183             CALL iom_put( 'qtr_'//ctrcnm(jn) , qtr_cfc (:,:,jn) ) 
    184             CALL iom_put( 'qint_'//ctrcnm(jn), qint_cfc(:,:,jn) ) 
     184            jl = jl + 1 
     185            CALL iom_put( 'qtr_'//TRIM(ctrcnm(jn)) , qtr_cfc (:,:,jl) ) 
     186            CALL iom_put( 'qint_'//TRIM(ctrcnm(jn)), qint_cfc(:,:,jl) ) 
    185187         ENDDO 
    186188      END IF 
  • branches/2017/dev_METO_2017/NEMOGCM/NEMO/TOP_SRC/PISCES/P4Z/p4zlim.F90

    r7753 r8868  
    207207                  &                                / ( oxymin + trb(ji,jj,jk,jpoxy) )  ) 
    208208               nitrfac(ji,jj,jk) = MIN( 1., nitrfac(ji,jj,jk) ) 
     209               ! 
     210               ! denitrification factor computed from NO3 levels 
     211               nitrfac2(ji,jj,jk) = MAX( 0.e0,       ( 1.E-6 - trb(ji,jj,jk,jpno3) )  & 
     212                  &                                / ( 1.E-6 + trb(ji,jj,jk,jpno3) ) ) 
     213               nitrfac2(ji,jj,jk) = MIN( 1., nitrfac2(ji,jj,jk) ) 
    209214            END DO 
    210215         END DO 
  • branches/2017/dev_METO_2017/NEMOGCM/NEMO/TOP_SRC/PISCES/P4Z/p4zmeso.F90

    r7753 r8868  
    9595               !  no real reason except that it seems to be more stable and may mimic predation 
    9696               !  --------------------------------------------------------------- 
    97                ztortz2   = mzrat2 * 1.e6 * zfact * trb(ji,jj,jk,jpmes) 
     97               ztortz2   = mzrat2 * 1.e6 * zfact * trb(ji,jj,jk,jpmes)  * (1. - nitrfac(ji,jj,jk) ) 
    9898               ! 
    9999               zcompadi  = MAX( ( trb(ji,jj,jk,jpdia) - xthresh2dia ), 0.e0 ) 
     
    125125               !  ---------------------------------- 
    126126               zgrazffeg = grazflux  * xstep * wsbio4(ji,jj,jk)      & 
    127                &           * tgfunc2(ji,jj,jk) * trb(ji,jj,jk,jpgoc) * trb(ji,jj,jk,jpmes) 
     127               &           * tgfunc2(ji,jj,jk) * trb(ji,jj,jk,jpgoc) * trb(ji,jj,jk,jpmes) & 
     128               &           * (1. - nitrfac(ji,jj,jk)) 
    128129               zgrazfffg = zgrazffeg * trb(ji,jj,jk,jpbfe) / (trb(ji,jj,jk,jpgoc) + rtrn) 
    129130               zgrazffep = grazflux  * xstep *  wsbio3(ji,jj,jk)     & 
    130                &           * tgfunc2(ji,jj,jk) * trb(ji,jj,jk,jppoc) * trb(ji,jj,jk,jpmes) 
     131               &           * tgfunc2(ji,jj,jk) * trb(ji,jj,jk,jppoc) * trb(ji,jj,jk,jpmes) & 
     132               &           * (1. - nitrfac(ji,jj,jk)) 
    131133               zgrazfffp = zgrazffep * trb(ji,jj,jk,jpsfe) / (trb(ji,jj,jk,jppoc) + rtrn) 
    132134              ! 
  • branches/2017/dev_METO_2017/NEMOGCM/NEMO/TOP_SRC/PISCES/P4Z/p4zmicro.F90

    r7753 r8868  
    9393               !  no real reason except that it seems to be more stable and may mimic predation. 
    9494               !  --------------------------------------------------------------- 
    95                ztortz = mzrat * 1.e6 * zfact * trb(ji,jj,jk,jpzoo) 
     95               ztortz = mzrat * 1.e6 * zfact * trb(ji,jj,jk,jpzoo) * (1. - nitrfac(ji,jj,jk)) 
    9696 
    9797               zcompadi  = MIN( MAX( ( trb(ji,jj,jk,jpdia) - xthreshdia ), 0.e0 ), xsizedia ) 
     
    105105               zdenom    = zfoodlim / ( xkgraz + zfoodlim ) 
    106106               zdenom2   = zdenom / ( zfood + rtrn ) 
    107                zgraze    = grazrat * xstep * tgfunc2(ji,jj,jk) * trb(ji,jj,jk,jpzoo)  
     107               zgraze    = grazrat * xstep * tgfunc2(ji,jj,jk) * trb(ji,jj,jk,jpzoo) * (1. - nitrfac(ji,jj,jk)) 
    108108 
    109109               zgrazp    = zgraze  * xpref2p * zcompaph  * zdenom2  
  • branches/2017/dev_METO_2017/NEMOGCM/NEMO/TOP_SRC/PISCES/P4Z/p4zrem.F90

    r7753 r8868  
    6565      REAL(wp) ::   zsatur, zsatur2, znusil, znusil2, zdep, zdepmin, zfactdep 
    6666      REAL(wp) ::   zbactfer, zolimit, zonitr, zrfact2 
     67      REAL(wp) ::   zammonic, zoxyrem 
    6768      REAL(wp) ::   zosil, ztem, zdenitnh4, zolimic, zolimin, zolimip, zdenitrn, zdenitrp 
    6869      CHARACTER (len=25) :: charout 
     
    118119                  ! Ammonification in suboxic waters with denitrification 
    119120                  ! ------------------------------------------------------- 
    120                   denitr(ji,jj,jk)  = MIN(  ( trb(ji,jj,jk,jpno3) - rtrn ) / rdenit,   & 
    121                      &                     zremik * nitrfac(ji,jj,jk) * trb(ji,jj,jk,jpdoc)  ) 
     121                  zammonic = zremik * nitrfac(ji,jj,jk) * trb(ji,jj,jk,jpdoc) 
     122                  denitr(ji,jj,jk)  = zammonic * ( 1. - nitrfac2(ji,jj,jk) ) 
     123                  zoxyrem           = zammonic *        nitrfac2(ji,jj,jk) 
    122124                  ! 
    123125                  zolimi (ji,jj,jk) = MAX( 0.e0, zolimi (ji,jj,jk) ) 
    124126                  denitr (ji,jj,jk) = MAX( 0.e0, denitr (ji,jj,jk) ) 
     127                  zoxyrem           = MAX( 0.e0, zoxyrem  ) 
     128 
    125129                  ! 
    126                   tra(ji,jj,jk,jppo4) = tra(ji,jj,jk,jppo4) + zolimi (ji,jj,jk) + denitr(ji,jj,jk) 
    127                   tra(ji,jj,jk,jpnh4) = tra(ji,jj,jk,jpnh4) + zolimi (ji,jj,jk) + denitr(ji,jj,jk) 
     130                  tra(ji,jj,jk,jppo4) = tra(ji,jj,jk,jppo4) + zolimi (ji,jj,jk) + denitr(ji,jj,jk) + zoxyrem 
     131                  tra(ji,jj,jk,jpnh4) = tra(ji,jj,jk,jpnh4) + zolimi (ji,jj,jk) + denitr(ji,jj,jk) + zoxyrem 
    128132                  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) 
     133                  tra(ji,jj,jk,jpdoc) = tra(ji,jj,jk,jpdoc) - zolimi (ji,jj,jk) - denitr(ji,jj,jk) - zoxyrem 
    130134                  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) 
    132                   tra(ji,jj,jk,jptal) = tra(ji,jj,jk,jptal) + rno3 * ( zolimi(ji,jj,jk)    & 
     135                  tra(ji,jj,jk,jpdic) = tra(ji,jj,jk,jpdic) + zolimi (ji,jj,jk) + denitr(ji,jj,jk) + zoxyrem 
     136                  tra(ji,jj,jk,jptal) = tra(ji,jj,jk,jptal) + rno3 * ( zolimi(ji,jj,jk) + zoxyrem    & 
    133137                  &                     + ( rdenit + 1.) * denitr(ji,jj,jk) ) 
    134138               END DO 
  • branches/2017/dev_METO_2017/NEMOGCM/NEMO/TOP_SRC/PISCES/P4Z/p4zsed.F90

    r7753 r8868  
    5454      REAL(wp) ::  zwflux, zfminus, zfplus 
    5555      REAL(wp) ::  zlim, zfact, zfactcal 
    56       REAL(wp) ::  zo2, zno3, zflx, zpdenit, z1pdenit, zdenitt, zolimit 
     56      REAL(wp) ::  zo2, zno3, zflx, zpdenit, z1pdenit, zolimit 
    5757      REAL(wp) ::  zsiloss, zcaloss, zws3, zws4, zwsc, zdep 
    5858      REAL(wp) ::  zwstpoc, zwstpon, zwstpop 
     
    319319               tra(ji,jj,ikt,jptal) =  tra(ji,jj,ikt,jptal) + zcaloss * zrivalk * 2.0 
    320320               tra(ji,jj,ikt,jpdic) =  tra(ji,jj,ikt,jpdic) + zcaloss * zrivalk 
    321                zsedcal(ji,jj) = (1.0 - zrivalk) * zcaloss / zdep 
    322                zsedsi (ji,jj) = (1.0 - zrivsil) * zsiloss / zdep 
     321               zsedcal(ji,jj) = (1.0 - zrivalk) * zcaloss * e3t_n(ji,jj,ikt)  
     322               zsedsi (ji,jj) = (1.0 - zrivsil) * zsiloss * e3t_n(ji,jj,ikt)  
    323323            END DO 
    324324         END DO 
     
    365365 
    366366      IF( .NOT.lk_sed ) THEN 
    367          ! The 0.5 factor in zpdenit and zdenitt is to avoid negative NO3 concentration after both denitrification 
    368          ! in the sediments and just above the sediments. Not very clever, but simpliest option. 
     367         ! The 0.5 factor in zpdenit is to avoid negative NO3 concentration after 
     368         ! denitrification in the sediments. Not very clever, but simpliest option. 
    369369         DO jj = 1, jpj 
    370370            DO ji = 1, jpi 
     
    378378               z1pdenit = zwstpoc * zrivno3 - zpdenit 
    379379               zolimit = MIN( ( trb(ji,jj,ikt,jpoxy) - rtrn ) / o2ut, z1pdenit * ( 1.- nitrfac(ji,jj,ikt) ) ) 
    380                zdenitt = MIN(  0.5 * ( trb(ji,jj,ikt,jpno3) - rtrn ) / rdenit, z1pdenit * nitrfac(ji,jj,ikt) ) 
    381                tra(ji,jj,ikt,jpdoc) = tra(ji,jj,ikt,jpdoc) + z1pdenit - zolimit - zdenitt 
    382                tra(ji,jj,ikt,jppo4) = tra(ji,jj,ikt,jppo4) + zpdenit + zolimit + zdenitt 
    383                tra(ji,jj,ikt,jpnh4) = tra(ji,jj,ikt,jpnh4) + zpdenit + zolimit + zdenitt 
    384                tra(ji,jj,ikt,jpno3) = tra(ji,jj,ikt,jpno3) - rdenit * (zpdenit + zdenitt) 
     380               tra(ji,jj,ikt,jpdoc) = tra(ji,jj,ikt,jpdoc) + z1pdenit - zolimit 
     381               tra(ji,jj,ikt,jppo4) = tra(ji,jj,ikt,jppo4) + zpdenit + zolimit 
     382               tra(ji,jj,ikt,jpnh4) = tra(ji,jj,ikt,jpnh4) + zpdenit + zolimit 
     383               tra(ji,jj,ikt,jpno3) = tra(ji,jj,ikt,jpno3) - rdenit * zpdenit 
    385384               tra(ji,jj,ikt,jpoxy) = tra(ji,jj,ikt,jpoxy) - zolimit * o2ut 
    386                tra(ji,jj,ikt,jptal) = tra(ji,jj,ikt,jptal) + rno3 * (zolimit + (1.+rdenit) * (zpdenit + zdenitt) ) 
    387                tra(ji,jj,ikt,jpdic) = tra(ji,jj,ikt,jpdic) + zpdenit + zolimit + zdenitt 
     385               tra(ji,jj,ikt,jptal) = tra(ji,jj,ikt,jptal) + rno3 * (zolimit + (1.+rdenit) * zpdenit ) 
     386               tra(ji,jj,ikt,jpdic) = tra(ji,jj,ikt,jpdic) + zpdenit + zolimit  
    388387               sdenit(ji,jj) = rdenit * zpdenit * e3t_n(ji,jj,ikt) 
    389                zsedc(ji,jj)   = (1. - zrivno3) * zwstpoc / zdep 
     388               zsedc(ji,jj)   = (1. - zrivno3) * zwstpoc * e3t_n(ji,jj,ikt) 
    390389               IF( ln_p5z ) THEN 
    391390                  zwstpop              = trb(ji,jj,ikt,jpgop) * zws4 + trb(ji,jj,ikt,jppop) * zws3 
    392391                  zwstpon              = trb(ji,jj,ikt,jpgon) * zws4 + trb(ji,jj,ikt,jppon) * zws3 
    393                   tra(ji,jj,ikt,jpdon) = tra(ji,jj,ikt,jpdon) + (z1pdenit - zolimit - zdenitt) * zwstpon / (zwstpoc + rtrn) 
    394                   tra(ji,jj,ikt,jpdop) = tra(ji,jj,ikt,jpdop) + (z1pdenit - zolimit - zdenitt) * zwstpop / (zwstpoc + rtrn) 
     392                  tra(ji,jj,ikt,jpdon) = tra(ji,jj,ikt,jpdon) + ( z1pdenit - zolimit ) * zwstpon / (zwstpoc + rtrn) 
     393                  tra(ji,jj,ikt,jpdop) = tra(ji,jj,ikt,jpdop) + ( z1pdenit - zolimit ) * zwstpop / (zwstpoc + rtrn) 
    395394               ENDIF 
    396395            END DO 
     
    494493      IF( lk_iomput ) THEN 
    495494         IF( knt == nrdttrc ) THEN 
    496             zfact = 1.e+3 * rfact2r * rno3  !  conversion from molC/l/kt  to molN/m3/s 
    497             IF( iom_use("Nfix"   ) ) CALL iom_put( "Nfix", nitrpot(:,:,:) * nitrfix * zfact * tmask(:,:,:) )  ! nitrogen fixation  
     495            zfact = 1.e+3 * rfact2r !  conversion from molC/l/kt  to molN/m3/s 
     496            IF( iom_use("Nfix"   ) ) CALL iom_put( "Nfix", nitrpot(:,:,:) * nitrfix * rno3 * zfact * tmask(:,:,:) )  ! nitrogen fixation  
    498497            IF( iom_use("INTNFIX") ) THEN   ! nitrogen fixation rate in ocean ( vertically integrated ) 
    499498               zwork1(:,:) = 0. 
    500499               DO jk = 1, jpkm1 
    501                  zwork1(:,:) = zwork1(:,:) + nitrpot(:,:,jk) * nitrfix * zfact * e3t_n(:,:,jk) * tmask(:,:,jk) 
     500                 zwork1(:,:) = zwork1(:,:) + nitrpot(:,:,jk) * nitrfix * rno3 * zfact * e3t_n(:,:,jk) * tmask(:,:,jk) 
    502501               ENDDO 
    503502               CALL iom_put( "INTNFIX" , zwork1 )  
    504503            ENDIF 
    505             IF( iom_use("SedCal" ) ) CALL iom_put( "SedCal", zsedcal(:,:) * 1.e+3 ) 
    506             IF( iom_use("SedSi" ) )  CALL iom_put( "SedSi",  zsedsi (:,:) * 1.e+3 ) 
    507             IF( iom_use("SedC" ) )   CALL iom_put( "SedC",   zsedc  (:,:) * 1.e+3 ) 
    508             IF( iom_use("Sdenit" ) ) CALL iom_put( "Sdenit", sdenit (:,:) * 1.e+3 * rno3 ) 
     504            IF( iom_use("SedCal" ) ) CALL iom_put( "SedCal", zsedcal(:,:) * zfact ) 
     505            IF( iom_use("SedSi" ) )  CALL iom_put( "SedSi",  zsedsi (:,:) * zfact ) 
     506            IF( iom_use("SedC" ) )   CALL iom_put( "SedC",   zsedc  (:,:) * zfact ) 
     507            IF( iom_use("Sdenit" ) ) CALL iom_put( "Sdenit", sdenit (:,:) * zfact * rno3 ) 
    509508         ENDIF 
    510509      ENDIF 
  • branches/2017/dev_METO_2017/NEMOGCM/NEMO/TOP_SRC/PISCES/P4Z/p4zsms.F90

    r7753 r8868  
    431431 
    432432      IF( kt == nittrc000 ) THEN  
     433         xfact1 = rfact2r * 12. / 1.e15 * ryyss    ! conversion molC/kt --> PgC/yr 
     434         xfact2 = 1.e+3 * rno3 * 14. / 1.e12 * ryyss   ! conversion molC/l/s ----> TgN/m3/yr 
     435         xfact3 = 1.e+3 * rfact2r * rno3   ! conversion molC/l/kt ----> molN/m3/s 
    433436         IF( ln_check_mass .AND. lwp) THEN      !   Open budget file of NO3, ALK, Si, Fer 
    434437            CALL ctl_opn( numco2, 'carbon.budget'  , 'REPLACE', 'FORMATTED', 'SEQUENTIAL', -1, 6, .FALSE., narea ) 
    435438            CALL ctl_opn( numnut, 'nutrient.budget', 'REPLACE', 'FORMATTED', 'SEQUENTIAL', -1, 6, .FALSE., narea ) 
    436439            CALL ctl_opn( numnit, 'nitrogen.budget', 'REPLACE', 'FORMATTED', 'SEQUENTIAL', -1, 6, .FALSE., narea ) 
    437             xfact1 = rfact2r * 12. / 1.e15 * ryyss    ! conversion molC/kt --> PgC/yr 
    438             xfact2 = 1.e+3 * rno3 * 14. / 1.e12 * ryyss   ! conversion molC/l/s ----> TgN/m3/yr 
    439             xfact3 = 1.e+3 * rfact2r * rno3   ! conversion molC/l/kt ----> molN/m3/s 
    440440            cltxt='time-step   Alkalinity        Nitrate        Phosphorus         Silicate           Iron' 
    441441            IF( lwp ) WRITE(numnut,*)  TRIM(cltxt) 
     
    517517      IF( iom_use( "tnfix" ) .OR.  ( ln_check_mass .AND. kt == nitend )  ) THEN 
    518518         znitrpottot  = glob_sum ( nitrpot(:,:,:) * nitrfix * cvol(:,:,:) ) 
    519          CALL iom_put( "tnfix"  , znitrpottot * 1.e+3 * rno3 )  ! Global  nitrogen fixation molC/l  to molN/m3  
     519         CALL iom_put( "tnfix"  , znitrpottot * xfact3 )  ! Global  nitrogen fixation molC/l  to molN/m3  
    520520      ENDIF 
    521521      ! 
    522522      IF( iom_use( "tdenit" ) .OR.  ( ln_check_mass .AND. kt == nitend )  ) THEN 
    523          zrdenittot   = glob_sum ( denitr(:,:,:) * rdenit * xnegtr(:,:,:) * cvol(:,:,:) ) 
    524          CALL iom_put( "tdenit"  , zrdenittot * 1.e+3 * rno3 )  ! Total denitrification molC/l to molN/m3  
    525       ENDIF 
    526       ! 
    527       IF( iom_use( "Sdenit" ) .OR.  ( ln_check_mass .AND. kt == nitend )  ) THEN 
    528          zsdenittot   = glob_sum ( sdenit(:,:) * e1e2t(:,:) ) 
    529          CALL iom_put( "Sdenit", sdenit(:,:) * xfact3 * tmask(:,:,1) )  ! Nitrate reduction in the sediments 
    530       ENDIF 
    531  
     523         zrdenittot = glob_sum ( denitr(:,:,:) * rdenit * xnegtr(:,:,:) * cvol(:,:,:) ) 
     524         zsdenittot = glob_sum ( sdenit(:,:) * e1e2t(:,:) * tmask(:,:,1) ) 
     525         CALL iom_put( "tdenit" , ( zrdenittot + zsdenittot ) * xfact3 )  ! Total denitrification molC/l to molN/m3  
     526      ENDIF 
     527      ! 
    532528      IF( ln_check_mass .AND. kt == nitend ) THEN   ! Compute the budget of NO3, ALK, Si, Fer 
    533529         t_atm_co2_flx  = t_atm_co2_flx / glob_sum( e1e2t(:,:) ) 
  • branches/2017/dev_METO_2017/NEMOGCM/NEMO/TOP_SRC/PISCES/SED/sedco3.F90

    r5215 r8868  
    7575 
    7676      DO jk = 1, jpksed 
    77 10001    CONTINUE 
    78          IF( itime <= 2 ) THEN 
     77         DO WHILE( itime <= 2 ) 
    7978            lconv  = .FALSE. 
    8079            IF( itime > 0 ) THEN   
     
    154153!                     WRITE(numsed,*) '    with re-initialization of initial PH field '        
    155154                     itime = 2 
    156                      GOTO 10001 
    157155                  ELSE 
    158156!                     WRITE(numsed,*) ' convergence after iter =', jiter, ' iterations ;  res =',zresm  
     
    165163!                     &               '  after iter =', jiter, ' iterations ;  res =',zresm   
    166164!                  WRITE(numsed,*) ' ' 
    167                   itime = 0 
     165                  itime = 3 
    168166               ENDIF 
    169167            ELSE 
     
    172170               IF ( itime == 1 ) THEN 
    173171                  WRITE(numsed,*) ' try one more time with more iterations and higher relax. value' 
    174                   GOTO 10001 
    175172               ELSE IF ( itime == 2 ) THEN 
    176173                  WRITE(numsed,*) ' try one more time for with more iterations, higher relax. value'                
     
    181178               ENDIF 
    182179            ENDIF 
    183          ENDIF 
     180         ENDDO ! End of WHILE LOOP 
    184181     ENDDO 
    185182 
  • branches/2017/dev_METO_2017/NEMOGCM/NEMO/TOP_SRC/PISCES/sms_pisces.F90

    r7646 r8868  
    9797   REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:,:) ::   xfracal    !: ?? 
    9898   REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:,:) ::   nitrfac    !: ?? 
     99   REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:,:) ::   nitrfac2   !: ?? 
    99100   REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:,:) ::   orem       !: ?? 
    100101   REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:,:) ::   xdiss      !: ?? 
     
    159160         ! 
    160161         !*  SMS for the organic matter 
    161          ALLOCATE( xfracal (jpi,jpj,jpk), nitrfac(jpi,jpj,jpk) ,    & 
    162             &      orem    (jpi,jpj,jpk),                           & 
    163             &      prodcal(jpi,jpj,jpk), xdiss   (jpi,jpj,jpk),    & 
     162         ALLOCATE( xfracal (jpi,jpj,jpk), orem(jpi,jpj,jpk)    ,    & 
     163            &      nitrfac(jpi,jpj,jpk), nitrfac2(jpi,jpj,jpk) ,    & 
     164            &      prodcal(jpi,jpj,jpk) , xdiss   (jpi,jpj,jpk),    & 
    164165            &      prodpoc(jpi,jpj,jpk) , conspoc(jpi,jpj,jpk) ,    & 
    165166            &      prodgoc(jpi,jpj,jpk) , consgoc(jpi,jpj,jpk) ,  STAT=ierr(4) ) 
  • branches/2017/dev_METO_2017/NEMOGCM/NEMO/TOP_SRC/TRP/trcnxt.F90

    r7881 r8868  
    106106      ENDIF 
    107107      !                                ! Leap-Frog + Asselin filter time stepping 
    108       IF( neuler == 0 .AND. kt == nittrc000 ) THEN    ! Euler time-stepping at first time-step (only swap) 
     108      IF( (neuler == 0 .AND. kt == nittrc000) .OR. ln_top_euler ) THEN    ! Euler time-stepping (only swap) 
    109109         DO jn = 1, jptra 
    110110            DO jk = 1, jpkm1 
    111111               trn(:,:,jk,jn) = tra(:,:,jk,jn) 
     112               trb(:,:,jk,jn) = trn(:,:,jk,jn)   
    112113            END DO 
    113114         END DO 
  • branches/2017/dev_METO_2017/NEMOGCM/NEMO/TOP_SRC/trc.F90

    r7881 r8868  
    1010   USE par_oce 
    1111   USE par_trc 
    12    USE bdy_oce, only: ln_bdy, nb_bdy, OBC_DATA 
     12   USE bdy_oce, only: jp_bdy, ln_bdy, nb_bdy, OBC_DATA 
    1313    
    1414   IMPLICIT NONE 
     
    169169# endif 
    170170   ! 
    171    CHARACTER(len=20), PUBLIC, ALLOCATABLE,  SAVE,  DIMENSION(:)   ::  cn_trc_dflt          ! Default OBC condition for all tracers 
    172    CHARACTER(len=20), PUBLIC, ALLOCATABLE,  SAVE,  DIMENSION(:)   ::  cn_trc               ! Choice of boundary condition for tracers 
    173    INTEGER,           PUBLIC, ALLOCATABLE,  SAVE,  DIMENSION(:)   ::  nn_trcdmp_bdy        !: =T Tracer damping 
     171   CHARACTER(len=20), PUBLIC, DIMENSION(jp_bdy) :: cn_trc_dflt   ! Default OBC condition for all tracers 
     172   CHARACTER(len=20), PUBLIC, DIMENSION(jp_bdy) :: cn_trc        ! Choice of boundary condition for tracers 
     173   INTEGER,           PUBLIC, DIMENSION(jp_bdy) :: nn_trcdmp_bdy !: =T Tracer damping 
     174!$AGRIF_DO_NOT_TREAT 
    174175   ! External data structure of BDY for TOP. Available elements: cn_obc, ll_trc, trcnow, dmp 
    175176   TYPE(OBC_DATA),    PUBLIC, ALLOCATABLE, DIMENSION(:,:), TARGET ::  trcdta_bdy           !: bdy external data (local process) 
    176    ! 
    177  
     177!$AGRIF_END_DO_NOT_TREAT 
    178178   !!---------------------------------------------------------------------- 
    179179   !! NEMO/TOP 3.3.1 , NEMO Consortium (2010) 
     
    206206      ! 
    207207      IF ( ln_bdy ) THEN 
    208          ALLOCATE( cn_trc_dflt(nb_bdy)   , cn_trc(nb_bdy)     , nn_trcdmp_bdy(nb_bdy) ,       & 
    209          &      trcdta_bdy(jptra,nb_bdy)                                              ,       & 
    210          &      STAT = ierr(2)  ) 
     208         ALLOCATE( trcdta_bdy(jptra, jp_bdy), STAT = ierr(2) ) 
    211209      ENDIF 
    212210      ! 
  • branches/2017/dev_METO_2017/NEMOGCM/NEMO/TOP_SRC/trcdta.F90

    r7646 r8868  
    227227                     ik = mbkt(ji,jj)  
    228228                     IF( ik > 1 ) THEN 
    229                         zl = ( gdept_1d(ik) - gdept_n(ji,jj,ik) ) / ( gdept_1d(ik) - gdept_1d(ik-1) ) 
     229                        zl = ( gdept_1d(ik) - gdept_0(ji,jj,ik) ) / ( gdept_1d(ik) - gdept_1d(ik-1) ) 
    230230                        ptrcdta(ji,jj,ik) = (1.-zl) * ptrcdta(ji,jj,ik) + zl * ptrcdta(ji,jj,ik-1) 
    231231                     ENDIF 
    232232                     ik = mikt(ji,jj) 
    233233                     IF( ik > 1 ) THEN 
    234                         zl = ( gdept_n(ji,jj,ik) - gdept_1d(ik) ) / ( gdept_1d(ik+1) - gdept_1d(ik) ) 
     234                        zl = ( gdept_0(ji,jj,ik) - gdept_1d(ik) ) / ( gdept_1d(ik+1) - gdept_1d(ik) ) 
    235235                        ptrcdta(ji,jj,ik) = (1.-zl) * ptrcdta(ji,jj,ik) + zl * ptrcdta(ji,jj,ik+1) 
    236236                     ENDIF 
  • branches/2017/dev_METO_2017/NEMOGCM/NEMO/TOP_SRC/trcini.F90

    r7753 r8868  
    6969      ! 
    7070      CALL trc_ini_sms   ! SMS 
    71       CALL trc_ini_inv   ! Inventories 
    7271      CALL trc_ini_trp   ! passive tracers transport 
    7372      CALL trc_ice_ini   ! Tracers in sea ice 
     
    7877      IF( nn_dttrc /= 1 ) & 
    7978      CALL trc_sub_ini    ! Initialize variables for substepping passive tracers 
     79      ! 
     80      CALL trc_ini_inv   ! Inventories 
    8081      ! 
    8182      IF( nn_timing == 1 )   CALL timing_stop('trc_init') 
Note: See TracChangeset for help on using the changeset viewer.