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 11082 for NEMO/branches/UKMO/NEMO_4.0_GO8_package/src/OCE/DYN – NEMO

Ignore:
Timestamp:
2019-06-06T16:21:52+02:00 (5 years ago)
Author:
davestorkey
Message:

UKMO/NEMO_4.0_GO8_package : update to be relative to 11081 of NEMO_4.0_mirror.

Location:
NEMO/branches/UKMO/NEMO_4.0_GO8_package/src/OCE/DYN
Files:
2 edited

Legend:

Unmodified
Added
Removed
  • NEMO/branches/UKMO/NEMO_4.0_GO8_package/src/OCE/DYN/dynkeg.F90

    r10888 r11082  
    7474      INTEGER, INTENT( in ) ::   kscheme   ! =0/1   type of KEG scheme  
    7575      ! 
    76       INTEGER  ::   ji, jj, jk, jb    ! dummy loop indices 
    77       INTEGER  ::   ii, ifu, ib_bdy   ! local integers 
    78       INTEGER  ::   ij, ifv, igrd     !   -       - 
    79       REAL(wp) ::   zu, zv            ! local scalars 
     76      INTEGER  ::   ji, jj, jk, jb           ! dummy loop indices 
     77      INTEGER  ::   ifu, ifv, igrd, ib_bdy   ! local integers 
     78      REAL(wp) ::   zu, zv                   ! local scalars 
    8079      REAL(wp), DIMENSION(jpi,jpj,jpk)        ::   zhke 
    8180      REAL(wp), ALLOCATABLE, DIMENSION(:,:,:) ::   ztrdu, ztrdv  
     81      REAL(wp)  :: zweightu, zweightv 
    8282      !!---------------------------------------------------------------------- 
    8383      ! 
     
    9797       
    9898      zhke(:,:,jpk) = 0._wp 
    99        
    100       IF (ln_bdy) THEN 
    101          ! Maria Luneva & Fred Wobus: July-2016 
    102          ! compensate for lack of turbulent kinetic energy on liquid bdy points 
    103          DO ib_bdy = 1, nb_bdy 
    104             IF( cn_dyn3d(ib_bdy) /= 'none' ) THEN 
    105                igrd = 2           ! Copying normal velocity into points outside bdy 
    106                DO jb = 1, idx_bdy(ib_bdy)%nblenrim(igrd) 
    107                   DO jk = 1, jpkm1 
    108                      ii   = idx_bdy(ib_bdy)%nbi(jb,igrd) 
    109                      ij   = idx_bdy(ib_bdy)%nbj(jb,igrd) 
    110                      ifu   = NINT( idx_bdy(ib_bdy)%flagu(jb,igrd) ) 
    111                      un(ii-ifu,ij,jk) = un(ii,ij,jk) * umask(ii,ij,jk) 
    112                   END DO 
    113                END DO 
    114                ! 
    115                igrd = 3           ! Copying normal velocity into points outside bdy 
    116                DO jb = 1, idx_bdy(ib_bdy)%nblenrim(igrd) 
    117                   DO jk = 1, jpkm1 
    118                      ii   = idx_bdy(ib_bdy)%nbi(jb,igrd) 
    119                      ij   = idx_bdy(ib_bdy)%nbj(jb,igrd) 
    120                      ifv   = NINT( idx_bdy(ib_bdy)%flagv(jb,igrd) ) 
    121                      vn(ii,ij-ifv,jk) = vn(ii,ij,jk) * vmask(ii,ij,jk) 
    122                   END DO 
    123                END DO 
    124             ENDIF 
    125          ENDDO   
    126       ENDIF  
    12799 
    128100      SELECT CASE ( kscheme )             !== Horizontal kinetic energy at T-point  ==! 
     
    140112            END DO 
    141113         END DO 
     114         ! 
     115         IF (ln_bdy) THEN 
     116            ! Maria Luneva & Fred Wobus: July-2016 
     117            ! compensate for lack of turbulent kinetic energy on liquid bdy points 
     118            DO ib_bdy = 1, nb_bdy 
     119               IF( cn_dyn3d(ib_bdy) /= 'none' ) THEN 
     120                  igrd = 1           ! compensating null velocity on the bdy 
     121                  DO jb = 1, idx_bdy(ib_bdy)%nblenrim(igrd) 
     122                     ji   = idx_bdy(ib_bdy)%nbi(jb,igrd)   ! maximum extent : from 2 to jpi-1 
     123                     jj   = idx_bdy(ib_bdy)%nbj(jb,igrd)   ! maximum extent : from 2 to jpj-1 
     124                     DO jk = 1, jpkm1 
     125                        zhke(ji,jj,jk) = 0._wp 
     126                        zweightu = umask(ji-1,jj  ,jk) + umask(ji,jj,jk) 
     127                        zweightv = vmask(ji  ,jj-1,jk) + vmask(ji,jj,jk) 
     128                        zu = un(ji-1,jj  ,jk) * un(ji-1,jj  ,jk)  +  un(ji  ,jj  ,jk) * un(ji  ,jj  ,jk) 
     129                        zv = vn(ji  ,jj-1,jk) * vn(ji  ,jj-1,jk)  +  vn(ji  ,jj  ,jk) * vn(ji  ,jj  ,jk) 
     130                        IF( zweightu > 0._wp )   zhke(ji,jj,jk) =  zhke(ji,jj,jk) + zu / (2._wp * zweightu)  
     131                        IF( zweightv > 0._wp )   zhke(ji,jj,jk) =  zhke(ji,jj,jk) + zv / (2._wp * zweightv)  
     132                     END DO 
     133                  END DO 
     134               END IF 
     135               CALL lbc_bdy_lnk( 'dynkeg', zhke, 'T', 1., ib_bdy )   ! send 2 and recv jpi, jpj used in the computation of the speed tendencies 
     136            END DO 
     137         END IF 
    142138         ! 
    143139      CASE ( nkeg_HW )                          !--  Hollingsworth scheme  --! 
     
    158154            END DO 
    159155         END DO 
     156         IF (ln_bdy) THEN 
     157            ! Maria Luneva & Fred Wobus: July-2016 
     158            ! compensate for lack of turbulent kinetic energy on liquid bdy points 
     159            DO ib_bdy = 1, nb_bdy 
     160               IF( cn_dyn3d(ib_bdy) /= 'none' ) THEN 
     161                  igrd = 1           ! compensation null velocity on land at the bdy 
     162                  DO jb = 1, idx_bdy(ib_bdy)%nblenrim(igrd) 
     163                     ji   = idx_bdy(ib_bdy)%nbi(jb,igrd)   ! maximum extent : from 2 to jpi-1 
     164                     jj   = idx_bdy(ib_bdy)%nbj(jb,igrd)   ! maximum extent : from 2 to jpj-1 
     165                     DO jk = 1, jpkm1 
     166                        zhke(ji,jj,jk) = 0._wp 
     167                        zweightu = 8._wp * ( umask(ji-1,jj  ,jk) + umask(ji  ,jj  ,jk) ) & 
     168                             &   + 2._wp * ( umask(ji-1,jj-1,jk) + umask(ji-1,jj+1,jk) + umask(ji  ,jj-1,jk) + umask(ji  ,jj+1,jk) )  
     169                        zweightv = 8._wp * ( vmask(ji  ,jj-1,jk) + vmask(ji  ,jj-1,jk) ) & 
     170                             &   + 2._wp * ( vmask(ji-1,jj-1,jk) + vmask(ji+1,jj-1,jk) + vmask(ji-1,jj  ,jk) + vmask(ji+1,jj  ,jk) ) 
     171                        zu = 8._wp * ( un(ji-1,jj  ,jk) * un(ji-1,jj  ,jk)    & 
     172                           &         + un(ji  ,jj  ,jk) * un(ji  ,jj  ,jk) )  & 
     173                           &   +     ( un(ji-1,jj-1,jk) + un(ji-1,jj+1,jk) ) * ( un(ji-1,jj-1,jk) + un(ji-1,jj+1,jk) )   & 
     174                           &   +     ( un(ji  ,jj-1,jk) + un(ji  ,jj+1,jk) ) * ( un(ji  ,jj-1,jk) + un(ji  ,jj+1,jk) ) 
     175                        zv = 8._wp * ( vn(ji  ,jj-1,jk) * vn(ji  ,jj-1,jk)    & 
     176                           &         + vn(ji  ,jj  ,jk) * vn(ji  ,jj  ,jk) )  & 
     177                           &  +      ( vn(ji-1,jj-1,jk) + vn(ji+1,jj-1,jk) ) * ( vn(ji-1,jj-1,jk) + vn(ji+1,jj-1,jk) )   & 
     178                           &  +      ( vn(ji-1,jj  ,jk) + vn(ji+1,jj  ,jk) ) * ( vn(ji-1,jj  ,jk) + vn(ji+1,jj  ,jk) ) 
     179                        IF( zweightu > 0._wp )   zhke(ji,jj,jk) =  zhke(ji,jj,jk) +  zu / ( 2._wp * zweightu ) 
     180                        IF( zweightv > 0._wp )   zhke(ji,jj,jk) =  zhke(ji,jj,jk) +  zv / ( 2._wp * zweightv ) 
     181                     END DO 
     182                  END DO 
     183               END IF 
     184            END DO 
     185         END IF 
    160186         CALL lbc_lnk( 'dynkeg', zhke, 'T', 1. ) 
    161187         ! 
    162       END SELECT 
    163  
    164       IF (ln_bdy) THEN 
    165          ! restore velocity masks at points outside boundary 
    166          un(:,:,:) = un(:,:,:) * umask(:,:,:) 
    167          vn(:,:,:) = vn(:,:,:) * vmask(:,:,:) 
    168       ENDIF       
    169  
     188      END SELECT  
    170189      ! 
    171190      DO jk = 1, jpkm1                    !==  grad( KE ) added to the general momentum trends  ==! 
  • NEMO/branches/UKMO/NEMO_4.0_GO8_package/src/OCE/DYN/sshwzv.F90

    r10888 r11082  
    297297         IF(lwp) WRITE(numout,*) 'wAimp : Courant number-based partitioning of now vertical velocity ' 
    298298         IF(lwp) WRITE(numout,*) '~~~~~ ' 
    299          ! 
    300          Cu_adv(:,:,jpk) = 0._wp              ! bottom value : Cu_adv=0 (set once for all) 
    301       ENDIF 
     299      ENDIF 
     300      ! 
    302301      ! 
    303302      DO jk = 1, jpkm1            ! calculate Courant numbers 
     
    305304            DO ji = 2, fs_jpim1   ! vector opt. 
    306305               z1_e3w = 1._wp / e3w_n(ji,jj,jk) 
    307                Cu_adv(ji,jj,jk) = r2dt * ( ( MAX( wn(ji,jj,jk) , 0._wp ) - MIN( wn(ji,jj,jk+1) , 0._wp ) )    & 
    308                   &                      + ( MAX( e2u(ji  ,jj)*e3uw_n(ji  ,jj,jk)*un(ji  ,jj,jk), 0._wp ) -   & 
    309                   &                          MIN( e2u(ji-1,jj)*e3uw_n(ji-1,jj,jk)*un(ji-1,jj,jk), 0._wp ) )   & 
    310                   &                        * r1_e1e2t(ji,jj)                                                  & 
    311                   &                      + ( MAX( e1v(ji,jj  )*e3vw_n(ji,jj  ,jk)*vn(ji,jj  ,jk), 0._wp ) -   & 
    312                   &                          MIN( e1v(ji,jj-1)*e3vw_n(ji,jj-1,jk)*vn(ji,jj-1,jk), 0._wp ) )   & 
    313                   &                        * r1_e1e2t(ji,jj)                                                  & 
    314                   &                      ) * z1_e3w 
     306               Cu_adv(ji,jj,jk) = 2._wp * rdt * ( ( MAX( wn(ji,jj,jk) , 0._wp ) - MIN( wn(ji,jj,jk+1) , 0._wp ) )    &  ! 2*rdt and not r2dt (for restartability) 
     307                  &                             + ( MAX( e2u(ji  ,jj)*e3uw_n(ji  ,jj,jk)*un(ji  ,jj,jk), 0._wp ) -   & 
     308                  &                                 MIN( e2u(ji-1,jj)*e3uw_n(ji-1,jj,jk)*un(ji-1,jj,jk), 0._wp ) )   & 
     309                  &                               * r1_e1e2t(ji,jj)                                                  & 
     310                  &                             + ( MAX( e1v(ji,jj  )*e3vw_n(ji,jj  ,jk)*vn(ji,jj  ,jk), 0._wp ) -   & 
     311                  &                                 MIN( e1v(ji,jj-1)*e3vw_n(ji,jj-1,jk)*vn(ji,jj-1,jk), 0._wp ) )   & 
     312                  &                               * r1_e1e2t(ji,jj)                                                  & 
     313                  &                             ) * z1_e3w 
    315314            END DO 
    316315         END DO 
    317316      END DO 
     317      CALL lbc_lnk( 'sshwzv', Cu_adv, 'T', 1. ) 
    318318      ! 
    319319      CALL iom_put("Courant",Cu_adv) 
    320320      ! 
    321       wi(:,:,:) = 0._wp                                 ! Includes top and bottom values set to zero 
    322321      IF( MAXVAL( Cu_adv(:,:,:) ) > Cu_min ) THEN       ! Quick check if any breaches anywhere 
    323322         DO jk = 1, jpkm1                               ! or scan Courant criterion and partition 
    324             DO jj = 2, jpjm1                            ! w where necessary 
    325                DO ji = 2, fs_jpim1   ! vector opt. 
     323            DO jj = 1, jpj                              ! w where necessary 
     324               DO ji = 1, jpi 
    326325                  ! 
    327326                  zCu = MAX( Cu_adv(ji,jj,jk) , Cu_adv(ji,jj,jk+1) ) 
    328327                  ! 
    329                   IF( zCu < Cu_min ) THEN               !<-- Fully explicit 
     328                  IF( zCu <= Cu_min ) THEN              !<-- Fully explicit 
    330329                     zcff = 0._wp 
    331330                  ELSEIF( zCu < Cu_cut ) THEN           !<-- Mixed explicit 
     
    346345      ELSE 
    347346         ! Fully explicit everywhere 
    348          Cu_adv = 0.0_wp                                ! Reuse array to output coefficient 
     347         Cu_adv(:,:,:) = 0._wp                          ! Reuse array to output coefficient 
     348         wi    (:,:,:) = 0._wp 
    349349      ENDIF 
    350350      CALL iom_put("wimp",wi)  
Note: See TracChangeset for help on using the changeset viewer.