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 12377 for NEMO/trunk/src/OCE/ASM/asminc.F90 – NEMO

Ignore:
Timestamp:
2020-02-12T15:39:06+01:00 (4 years ago)
Author:
acc
Message:

The big one. Merging all 2019 developments from the option 1 branch back onto the trunk.

This changeset reproduces 2019/dev_r11943_MERGE_2019 on the trunk using a 2-URL merge
onto a working copy of the trunk. I.e.:

svn merge --ignore-ancestry \

​svn+ssh://acc@forge.ipsl.jussieu.fr/ipsl/forge/projets/nemo/svn/NEMO/trunk \
​svn+ssh://acc@forge.ipsl.jussieu.fr/ipsl/forge/projets/nemo/svn/NEMO/branches/2019/dev_r11943_MERGE_2019 ./

The --ignore-ancestry flag avoids problems that may otherwise arise from the fact that
the merge history been trunk and branch may have been applied in a different order but
care has been taken before this step to ensure that all applicable fixes and updates
are present in the merge branch.

The trunk state just before this step has been branched to releases/release-4.0-HEAD
and that branch has been immediately tagged as releases/release-4.0.2. Any fixes
or additions in response to tickets on 4.0, 4.0.1 or 4.0.2 should be done on
releases/release-4.0-HEAD. From now on future 'point' releases (e.g. 4.0.2) will
remain unchanged with periodic releases as needs demand. Note release-4.0-HEAD is a
transitional naming convention. Future full releases, say 4.2, will have a release-4.2
branch which fulfills this role and the first point release (e.g. 4.2.0) will be made
immediately following the release branch creation.

2020 developments can be started from any trunk revision later than this one.

Location:
NEMO/trunk
Files:
2 edited

Legend:

Unmodified
Added
Removed

  • NEMO/trunk

    • Property svn:externals

      •  

        old new  
        33^/utils/build/mk@HEAD         mk 
        44^/utils/tools@HEAD            tools 
        5 ^/vendors/AGRIF/dev@HEAD      ext/AGRIF 
         5^/vendors/AGRIF/dev_r11615_ENHANCE-04_namelists_as_internalfiles_agrif@HEAD      ext/AGRIF 
        66^/vendors/FCM@HEAD            ext/FCM 
        77^/vendors/IOIPSL@HEAD         ext/IOIPSL 

  • NEMO/trunk/src/OCE/ASM/asminc.F90

    r11536 r12377  
    9494 
    9595   !! * Substitutions 
    96 #  include "vectopt_loop_substitute.h90" 
     96#  include "do_loop_substitute.h90" 
    9797   !!---------------------------------------------------------------------- 
    9898   !! NEMO/OCE 4.0 , NEMO Consortium (2018) 
    … …  
    102102CONTAINS 
    103103 
    104    SUBROUTINE asm_inc_init 
     104   SUBROUTINE asm_inc_init( Kbb, Kmm, Krhs ) 
    105105      !!---------------------------------------------------------------------- 
    106106      !!                    ***  ROUTINE asm_inc_init  *** 
    … …  
    112112      !! ** Action  :  
    113113      !!---------------------------------------------------------------------- 
     114      INTEGER, INTENT(in) ::  Kbb, Kmm, Krhs  ! time level indices 
     115      ! 
    114116      INTEGER :: ji, jj, jk, jt  ! dummy loop indices 
    115117      INTEGER :: imid, inum      ! local integers 
    … …  
    145147      ln_temnofreeze = .FALSE. 
    146148 
    147       REWIND( numnam_ref )              ! Namelist nam_asminc in reference namelist : Assimilation increment 
    148149      READ  ( numnam_ref, nam_asminc, IOSTAT = ios, ERR = 901) 
    149150901   IF( ios /= 0 )   CALL ctl_nam ( ios , 'nam_asminc in reference namelist' ) 
    150       REWIND( numnam_cfg )              ! Namelist nam_asminc in configuration namelist : Assimilation increment 
    151151      READ  ( numnam_cfg, nam_asminc, IOSTAT = ios, ERR = 902 ) 
    152152902   IF( ios >  0 )   CALL ctl_nam ( ios , 'nam_asminc in configuration namelist' ) 
    … …  
    413413            DO jk = 1, jpkm1           ! zhdiv = e1e1 * div 
    414414               zhdiv(:,:) = 0._wp 
    415                DO jj = 2, jpjm1 
    416                   DO ji = fs_2, fs_jpim1   ! vector opt. 
    417                      zhdiv(ji,jj) = (  e2u(ji  ,jj) * e3u_n(ji  ,jj,jk) * u_bkginc(ji  ,jj,jk)    & 
    418                         &            - e2u(ji-1,jj) * e3u_n(ji-1,jj,jk) * u_bkginc(ji-1,jj,jk)    & 
    419                         &            + e1v(ji,jj  ) * e3v_n(ji,jj  ,jk) * v_bkginc(ji,jj  ,jk)    & 
    420                         &            - e1v(ji,jj-1) * e3v_n(ji,jj-1,jk) * v_bkginc(ji,jj-1,jk)  ) / e3t_n(ji,jj,jk) 
    421                   END DO 
    422                END DO 
     415               DO_2D_00_00 
     416                  zhdiv(ji,jj) = (  e2u(ji  ,jj) * e3u(ji  ,jj,jk,Kmm) * u_bkginc(ji  ,jj,jk)    & 
     417                     &            - e2u(ji-1,jj) * e3u(ji-1,jj,jk,Kmm) * u_bkginc(ji-1,jj,jk)    & 
     418                     &            + e1v(ji,jj  ) * e3v(ji,jj  ,jk,Kmm) * v_bkginc(ji,jj  ,jk)    & 
     419                     &            - e1v(ji,jj-1) * e3v(ji,jj-1,jk,Kmm) * v_bkginc(ji,jj-1,jk)  ) / e3t(ji,jj,jk,Kmm) 
     420               END_2D 
    423421               CALL lbc_lnk( 'asminc', zhdiv, 'T', 1. )   ! lateral boundary cond. (no sign change) 
    424422               ! 
    425                DO jj = 2, jpjm1 
    426                   DO ji = fs_2, fs_jpim1   ! vector opt. 
    427                      u_bkginc(ji,jj,jk) = u_bkginc(ji,jj,jk)                         & 
    428                         &               + 0.2_wp * ( zhdiv(ji+1,jj) - zhdiv(ji  ,jj) ) * r1_e1u(ji,jj) * umask(ji,jj,jk) 
    429                      v_bkginc(ji,jj,jk) = v_bkginc(ji,jj,jk)                         & 
    430                         &               + 0.2_wp * ( zhdiv(ji,jj+1) - zhdiv(ji,jj  ) ) * r1_e2v(ji,jj) * vmask(ji,jj,jk)  
    431                   END DO 
    432                END DO 
     423               DO_2D_00_00 
     424                  u_bkginc(ji,jj,jk) = u_bkginc(ji,jj,jk)                         & 
     425                     &               + 0.2_wp * ( zhdiv(ji+1,jj) - zhdiv(ji  ,jj) ) * r1_e1u(ji,jj) * umask(ji,jj,jk) 
     426                  v_bkginc(ji,jj,jk) = v_bkginc(ji,jj,jk)                         & 
     427                     &               + 0.2_wp * ( zhdiv(ji,jj+1) - zhdiv(ji,jj  ) ) * r1_e2v(ji,jj) * vmask(ji,jj,jk)  
     428               END_2D 
    433429            END DO 
    434430            ! 
    … …  
    494490      ! 
    495491      IF( lk_asminc ) THEN                            !==  data assimilation  ==! 
    496          IF( ln_bkgwri )   CALL asm_bkg_wri( nit000 - 1 )      ! Output background fields 
     492         IF( ln_bkgwri )   CALL asm_bkg_wri( nit000 - 1, Kmm )      ! Output background fields 
    497493         IF( ln_asmdin ) THEN                                  ! Direct initialization 
    498             IF( ln_trainc )   CALL tra_asm_inc( nit000 - 1 )      ! Tracers 
    499             IF( ln_dyninc )   CALL dyn_asm_inc( nit000 - 1 )      ! Dynamics 
    500             IF( ln_sshinc )   CALL ssh_asm_inc( nit000 - 1 )      ! SSH 
     494            IF( ln_trainc )   CALL tra_asm_inc( nit000 - 1, Kbb, Kmm, ts    , Krhs )      ! Tracers 
     495            IF( ln_dyninc )   CALL dyn_asm_inc( nit000 - 1, Kbb, Kmm, uu, vv, Krhs )      ! Dynamics 
     496            IF( ln_sshinc )   CALL ssh_asm_inc( nit000 - 1, Kbb, Kmm )                    ! SSH 
    501497         ENDIF 
    502498      ENDIF 
    … …  
    505501    
    506502    
    507    SUBROUTINE tra_asm_inc( kt ) 
     503   SUBROUTINE tra_asm_inc( kt, Kbb, Kmm, pts, Krhs ) 
    508504      !!---------------------------------------------------------------------- 
    509505      !!                    ***  ROUTINE tra_asm_inc  *** 
    … …  
    515511      !! ** Action  :  
    516512      !!---------------------------------------------------------------------- 
    517       INTEGER, INTENT(IN) ::   kt   ! Current time step 
     513      INTEGER                                  , INTENT(in   ) :: kt             ! Current time step 
     514      INTEGER                                  , INTENT(in   ) :: Kbb, Kmm, Krhs ! Time level indices 
     515      REAL(wp), DIMENSION(jpi,jpj,jpk,jpts,jpt), INTENT(inout) :: pts            ! active tracers and RHS of tracer equation 
    518516      ! 
    519517      INTEGER  :: ji, jj, jk 
    … …  
    526524      ! used to prevent the applied increments taking the temperature below the local freezing point  
    527525      DO jk = 1, jpkm1 
    528         CALL eos_fzp( tsn(:,:,jk,jp_sal), fzptnz(:,:,jk), gdept_n(:,:,jk) ) 
     526        CALL eos_fzp( pts(:,:,jk,jp_sal,Kmm), fzptnz(:,:,jk), gdept(:,:,jk,Kmm) ) 
    529527      END DO 
    530528         ! 
    … …  
    549547                  ! Do not apply negative increments if the temperature will fall below freezing 
    550548                  WHERE(t_bkginc(:,:,jk) > 0.0_wp .OR. & 
    551                      &   tsn(:,:,jk,jp_tem) + tsa(:,:,jk,jp_tem) + t_bkginc(:,:,jk) * wgtiau(it) > fzptnz(:,:,jk) )  
    552                      tsa(:,:,jk,jp_tem) = tsa(:,:,jk,jp_tem) + t_bkginc(:,:,jk) * zincwgt   
     549                     &   pts(:,:,jk,jp_tem,Kmm) + pts(:,:,jk,jp_tem,Krhs) + t_bkginc(:,:,jk) * wgtiau(it) > fzptnz(:,:,jk) )  
     550                     pts(:,:,jk,jp_tem,Krhs) = pts(:,:,jk,jp_tem,Krhs) + t_bkginc(:,:,jk) * zincwgt   
    553551                  END WHERE 
    554552               ELSE 
    555                   tsa(:,:,jk,jp_tem) = tsa(:,:,jk,jp_tem) + t_bkginc(:,:,jk) * zincwgt   
     553                  pts(:,:,jk,jp_tem,Krhs) = pts(:,:,jk,jp_tem,Krhs) + t_bkginc(:,:,jk) * zincwgt   
    556554               ENDIF 
    557555               IF (ln_salfix) THEN 
    … …  
    559557                  ! minimum value salfixmin 
    560558                  WHERE(s_bkginc(:,:,jk) > 0.0_wp .OR. & 
    561                      &   tsn(:,:,jk,jp_sal) + tsa(:,:,jk,jp_sal) + s_bkginc(:,:,jk) * wgtiau(it) > salfixmin )  
    562                      tsa(:,:,jk,jp_sal) = tsa(:,:,jk,jp_sal) + s_bkginc(:,:,jk) * zincwgt 
     559                     &   pts(:,:,jk,jp_sal,Kmm) + pts(:,:,jk,jp_sal,Krhs) + s_bkginc(:,:,jk) * wgtiau(it) > salfixmin )  
     560                     pts(:,:,jk,jp_sal,Krhs) = pts(:,:,jk,jp_sal,Krhs) + s_bkginc(:,:,jk) * zincwgt 
    563561                  END WHERE 
    564562               ELSE 
    565                   tsa(:,:,jk,jp_sal) = tsa(:,:,jk,jp_sal) + s_bkginc(:,:,jk) * zincwgt 
     563                  pts(:,:,jk,jp_sal,Krhs) = pts(:,:,jk,jp_sal,Krhs) + s_bkginc(:,:,jk) * zincwgt 
    566564               ENDIF 
    567565            END DO 
    … …  
    584582            IF (ln_temnofreeze) THEN 
    585583               ! Do not apply negative increments if the temperature will fall below freezing 
    586                WHERE( t_bkginc(:,:,:) > 0.0_wp .OR. tsn(:,:,:,jp_tem) + t_bkginc(:,:,:) > fzptnz(:,:,:) )  
    587                   tsn(:,:,:,jp_tem) = t_bkg(:,:,:) + t_bkginc(:,:,:)    
     584               WHERE( t_bkginc(:,:,:) > 0.0_wp .OR. pts(:,:,:,jp_tem,Kmm) + t_bkginc(:,:,:) > fzptnz(:,:,:) )  
     585                  pts(:,:,:,jp_tem,Kmm) = t_bkg(:,:,:) + t_bkginc(:,:,:)    
    588586               END WHERE 
    589587            ELSE 
    590                tsn(:,:,:,jp_tem) = t_bkg(:,:,:) + t_bkginc(:,:,:)    
     588               pts(:,:,:,jp_tem,Kmm) = t_bkg(:,:,:) + t_bkginc(:,:,:)    
    591589            ENDIF 
    592590            IF (ln_salfix) THEN 
    593591               ! Do not apply negative increments if the salinity will fall below a specified 
    594592               ! minimum value salfixmin 
    595                WHERE( s_bkginc(:,:,:) > 0.0_wp .OR. tsn(:,:,:,jp_sal) + s_bkginc(:,:,:) > salfixmin )  
    596                   tsn(:,:,:,jp_sal) = s_bkg(:,:,:) + s_bkginc(:,:,:)    
     593               WHERE( s_bkginc(:,:,:) > 0.0_wp .OR. pts(:,:,:,jp_sal,Kmm) + s_bkginc(:,:,:) > salfixmin )  
     594                  pts(:,:,:,jp_sal,Kmm) = s_bkg(:,:,:) + s_bkginc(:,:,:)    
    597595               END WHERE 
    598596            ELSE 
    599                tsn(:,:,:,jp_sal) = s_bkg(:,:,:) + s_bkginc(:,:,:)    
     597               pts(:,:,:,jp_sal,Kmm) = s_bkg(:,:,:) + s_bkginc(:,:,:)    
    600598            ENDIF 
    601599 
    602             tsb(:,:,:,:) = tsn(:,:,:,:)                 ! Update before fields 
    603  
    604             CALL eos( tsb, rhd, rhop, gdept_0(:,:,:) )  ! Before potential and in situ densities 
     600            pts(:,:,:,:,Kbb) = pts(:,:,:,:,Kmm)                 ! Update before fields 
     601 
     602            CALL eos( pts(:,:,:,:,Kbb), rhd, rhop, gdept_0(:,:,:) )  ! Before potential and in situ densities 
    605603!!gm  fabien 
    606 !            CALL eos( tsb, rhd, rhop )                ! Before potential and in situ densities 
     604!            CALL eos( pts(:,:,:,:,Kbb), rhd, rhop )                ! Before potential and in situ densities 
    607605!!gm 
    608606 
    609             IF( ln_zps .AND. .NOT. lk_c1d .AND. .NOT. ln_isfcav)      & 
    610                &  CALL zps_hde    ( kt, jpts, tsb, gtsu, gtsv,        &  ! Partial steps: before horizontal gradient 
    611                &                              rhd, gru , grv          )  ! of t, s, rd at the last ocean level 
    612             IF( ln_zps .AND. .NOT. lk_c1d .AND.       ln_isfcav)      & 
    613                &  CALL zps_hde_isf( nit000, jpts, tsb, gtsu, gtsv, gtui, gtvi,    &    ! Partial steps for top cell (ISF) 
    614                &                                  rhd, gru , grv , grui, grvi     )    ! of t, s, rd at the last ocean level 
     607            IF( ln_zps .AND. .NOT. lk_c1d .AND. .NOT. ln_isfcav)           & 
     608               &  CALL zps_hde    ( kt, Kmm, jpts, pts(:,:,:,:,Kbb), gtsu, gtsv,        &  ! Partial steps: before horizontal gradient 
     609               &                              rhd, gru , grv               )  ! of t, s, rd at the last ocean level 
     610            IF( ln_zps .AND. .NOT. lk_c1d .AND.       ln_isfcav)                       & 
     611               &  CALL zps_hde_isf( nit000, Kmm, jpts, pts(:,:,:,:,Kbb), gtsu, gtsv, gtui, gtvi,    &  ! Partial steps for top cell (ISF) 
     612               &                                  rhd, gru , grv , grui, grvi          )  ! of t, s, rd at the last ocean level 
    615613 
    616614            DEALLOCATE( t_bkginc ) 
    … …  
    627625 
    628626 
    629    SUBROUTINE dyn_asm_inc( kt ) 
     627   SUBROUTINE dyn_asm_inc( kt, Kbb, Kmm, puu, pvv, Krhs ) 
    630628      !!---------------------------------------------------------------------- 
    631629      !!                    ***  ROUTINE dyn_asm_inc  *** 
    … …  
    637635      !! ** Action  :  
    638636      !!---------------------------------------------------------------------- 
    639       INTEGER, INTENT(IN) :: kt   ! Current time step 
     637      INTEGER                             , INTENT( in )  ::  kt             ! ocean time-step index 
     638      INTEGER                             , INTENT( in )  ::  Kbb, Kmm, Krhs ! ocean time level indices 
     639      REAL(wp), DIMENSION(jpi,jpj,jpk,jpt), INTENT(inout) ::  puu, pvv       ! ocean velocities and RHS of momentum equation 
    640640      ! 
    641641      INTEGER :: jk 
    … …  
    661661            ! Update the dynamic tendencies 
    662662            DO jk = 1, jpkm1 
    663                ua(:,:,jk) = ua(:,:,jk) + u_bkginc(:,:,jk) * zincwgt 
    664                va(:,:,jk) = va(:,:,jk) + v_bkginc(:,:,jk) * zincwgt 
     663               puu(:,:,jk,Krhs) = puu(:,:,jk,Krhs) + u_bkginc(:,:,jk) * zincwgt 
     664               pvv(:,:,jk,Krhs) = pvv(:,:,jk,Krhs) + v_bkginc(:,:,jk) * zincwgt 
    665665            END DO 
    666666            ! 
    … …  
    680680            ! 
    681681            ! Initialize the now fields with the background + increment 
    682             un(:,:,:) = u_bkg(:,:,:) + u_bkginc(:,:,:) 
    683             vn(:,:,:) = v_bkg(:,:,:) + v_bkginc(:,:,:)   
    684             ! 
    685             ub(:,:,:) = un(:,:,:)         ! Update before fields 
    686             vb(:,:,:) = vn(:,:,:) 
     682            puu(:,:,:,Kmm) = u_bkg(:,:,:) + u_bkginc(:,:,:) 
     683            pvv(:,:,:,Kmm) = v_bkg(:,:,:) + v_bkginc(:,:,:)   
     684            ! 
     685            puu(:,:,:,Kbb) = puu(:,:,:,Kmm)         ! Update before fields 
     686            pvv(:,:,:,Kbb) = pvv(:,:,:,Kmm) 
    687687            ! 
    688688            DEALLOCATE( u_bkg    ) 
    … …  
    697697 
    698698 
    699    SUBROUTINE ssh_asm_inc( kt ) 
     699   SUBROUTINE ssh_asm_inc( kt, Kbb, Kmm ) 
    700700      !!---------------------------------------------------------------------- 
    701701      !!                    ***  ROUTINE ssh_asm_inc  *** 
    … …  
    707707      !! ** Action  :  
    708708      !!---------------------------------------------------------------------- 
    709       INTEGER, INTENT(IN) :: kt   ! Current time step 
     709      INTEGER, INTENT(IN) :: kt         ! Current time step 
     710      INTEGER, INTENT(IN) :: Kbb, Kmm   ! Current time step 
    710711      ! 
    711712      INTEGER :: it 
    … …  
    754755            neuler = 0                                   ! Force Euler forward step 
    755756            ! 
    756             sshn(:,:) = ssh_bkg(:,:) + ssh_bkginc(:,:)   ! Initialize the now fields the background + increment 
    757             ! 
    758             sshb(:,:) = sshn(:,:)                        ! Update before fields 
    759             e3t_b(:,:,:) = e3t_n(:,:,:) 
    760 !!gm why not e3u_b, e3v_b, gdept_b ???? 
     757            ssh(:,:,Kmm) = ssh_bkg(:,:) + ssh_bkginc(:,:)   ! Initialize the now fields the background + increment 
     758            ! 
     759            ssh(:,:,Kbb) = ssh(:,:,Kmm)                        ! Update before fields 
     760            e3t(:,:,:,Kbb) = e3t(:,:,:,Kmm) 
     761!!gm why not e3u(:,:,:,Kbb), e3v(:,:,:,Kbb), gdept(:,:,:,Kbb) ???? 
    761762            ! 
    762763            DEALLOCATE( ssh_bkg    ) 
    … …  
    770771 
    771772 
    772    SUBROUTINE ssh_asm_div( kt, phdivn ) 
     773   SUBROUTINE ssh_asm_div( kt, Kbb, Kmm, phdivn ) 
    773774      !!---------------------------------------------------------------------- 
    774775      !!                  ***  ROUTINE ssh_asm_div  *** 
    … …  
    784785      !!---------------------------------------------------------------------- 
    785786      INTEGER, INTENT(IN) :: kt                               ! ocean time-step index 
     787      INTEGER, INTENT(IN) :: Kbb, Kmm                         ! time level indices 
    786788      REAL(wp), DIMENSION(:,:,:), INTENT(inout) ::   phdivn   ! horizontal divergence 
    787789      !! 
    … …  
    791793      !  
    792794#if defined key_asminc 
    793       CALL ssh_asm_inc( kt ) !==   (calculate increments) 
     795      CALL ssh_asm_inc( kt, Kbb, Kmm ) !==   (calculate increments) 
    794796      ! 
    795797      IF( ln_linssh ) THEN  
    796          phdivn(:,:,1) = phdivn(:,:,1) - ssh_iau(:,:) / e3t_n(:,:,1) * tmask(:,:,1) 
     798         phdivn(:,:,1) = phdivn(:,:,1) - ssh_iau(:,:) / e3t(:,:,1,Kmm) * tmask(:,:,1) 
    797799      ELSE  
    798800         ALLOCATE( ztim(jpi,jpj) ) 
    799          ztim(:,:) = ssh_iau(:,:) / ( ht_n(:,:) + 1.0 - ssmask(:,:) ) 
     801         ztim(:,:) = ssh_iau(:,:) / ( ht(:,:) + 1.0 - ssmask(:,:) ) 
    800802         DO jk = 1, jpkm1                                  
    801803            phdivn(:,:,jk) = phdivn(:,:,jk) - ztim(:,:) * tmask(:,:,jk)  
Note: See TracChangeset for help on using the changeset viewer.