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/SBC/sbcwave.F90 – NEMO

Ignore:
Timestamp:
2020-02-12T15:39:06+01:00 (3 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/SBC/sbcwave.F90

    r11536 r12377  
    7272 
    7373   !! * Substitutions 
    74 #  include "vectopt_loop_substitute.h90" 
     74#  include "do_loop_substitute.h90" 
    7575   !!---------------------------------------------------------------------- 
    7676   !! NEMO/OCE 4.0 , NEMO Consortium (2018) 
    … …  
    8080CONTAINS 
    8181 
    82    SUBROUTINE sbc_stokes( ) 
     82   SUBROUTINE sbc_stokes( Kmm ) 
    8383      !!--------------------------------------------------------------------- 
    8484      !!                     ***  ROUTINE sbc_stokes  *** 
    … …  
    9292      !! ** action   
    9393      !!--------------------------------------------------------------------- 
     94      INTEGER, INTENT(in) :: Kmm ! ocean time level index 
    9495      INTEGER  ::   jj, ji, jk   ! dummy loop argument 
    9596      INTEGER  ::   ik           ! local integer  
    … …  
    111112      IF( ll_st_bv_li ) THEN   ! (Eq. (19) in Breivik et al. (2014) ) 
    112113         zfac = 2.0_wp * rpi / 16.0_wp 
    113          DO jj = 1, jpj 
    114             DO ji = 1, jpi 
    115                ! Stokes drift velocity estimated from Hs and Tmean 
    116                ztransp = zfac * hsw(ji,jj)*hsw(ji,jj) / MAX( wmp(ji,jj), 0.0000001_wp ) 
    117                ! Stokes surface speed 
    118                tsd2d(ji,jj) = SQRT( ut0sd(ji,jj)*ut0sd(ji,jj) + vt0sd(ji,jj)*vt0sd(ji,jj)) 
    119                ! Wavenumber scale 
    120                zk_t(ji,jj) = ABS( tsd2d(ji,jj) ) / MAX( ABS( 5.97_wp*ztransp ), 0.0000001_wp ) 
    121             END DO 
    122          END DO 
    123          DO jj = 1, jpjm1              ! exp. wave number & Stokes drift velocity at u- & v-points 
    124             DO ji = 1, jpim1 
    125                zk_u(ji,jj) = 0.5_wp * ( zk_t(ji,jj) + zk_t(ji+1,jj) ) 
    126                zk_v(ji,jj) = 0.5_wp * ( zk_t(ji,jj) + zk_t(ji,jj+1) ) 
    127                ! 
    128                zu0_sd(ji,jj) = 0.5_wp * ( ut0sd(ji,jj) + ut0sd(ji+1,jj) ) 
    129                zv0_sd(ji,jj) = 0.5_wp * ( vt0sd(ji,jj) + vt0sd(ji,jj+1) ) 
    130             END DO 
    131          END DO 
     114         DO_2D_11_11 
     115            ! Stokes drift velocity estimated from Hs and Tmean 
     116            ztransp = zfac * hsw(ji,jj)*hsw(ji,jj) / MAX( wmp(ji,jj), 0.0000001_wp ) 
     117            ! Stokes surface speed 
     118            tsd2d(ji,jj) = SQRT( ut0sd(ji,jj)*ut0sd(ji,jj) + vt0sd(ji,jj)*vt0sd(ji,jj)) 
     119            ! Wavenumber scale 
     120            zk_t(ji,jj) = ABS( tsd2d(ji,jj) ) / MAX( ABS( 5.97_wp*ztransp ), 0.0000001_wp ) 
     121         END_2D 
     122         DO_2D_10_10 
     123            zk_u(ji,jj) = 0.5_wp * ( zk_t(ji,jj) + zk_t(ji+1,jj) ) 
     124            zk_v(ji,jj) = 0.5_wp * ( zk_t(ji,jj) + zk_t(ji,jj+1) ) 
     125            ! 
     126            zu0_sd(ji,jj) = 0.5_wp * ( ut0sd(ji,jj) + ut0sd(ji+1,jj) ) 
     127            zv0_sd(ji,jj) = 0.5_wp * ( vt0sd(ji,jj) + vt0sd(ji,jj+1) ) 
     128         END_2D 
    132129      ELSE IF( ll_st_peakfr ) THEN    ! peak wave number calculated from the peak frequency received by the wave model 
    133          DO jj = 1, jpj 
    134             DO ji = 1, jpi 
    135                zk_t(ji,jj) = ( 2.0_wp * rpi * wfreq(ji,jj) ) * ( 2.0_wp * rpi * wfreq(ji,jj) ) / grav 
    136             END DO 
    137          END DO 
    138          DO jj = 1, jpjm1 
    139             DO ji = 1, jpim1 
    140                zk_u(ji,jj) = 0.5_wp * ( zk_t(ji,jj) + zk_t(ji+1,jj) ) 
    141                zk_v(ji,jj) = 0.5_wp * ( zk_t(ji,jj) + zk_t(ji,jj+1) ) 
    142                ! 
    143                zu0_sd(ji,jj) = 0.5_wp * ( ut0sd(ji,jj) + ut0sd(ji+1,jj) ) 
    144                zv0_sd(ji,jj) = 0.5_wp * ( vt0sd(ji,jj) + vt0sd(ji,jj+1) ) 
    145             END DO 
    146          END DO 
     130         DO_2D_11_11 
     131            zk_t(ji,jj) = ( 2.0_wp * rpi * wfreq(ji,jj) ) * ( 2.0_wp * rpi * wfreq(ji,jj) ) / grav 
     132         END_2D 
     133         DO_2D_10_10 
     134            zk_u(ji,jj) = 0.5_wp * ( zk_t(ji,jj) + zk_t(ji+1,jj) ) 
     135            zk_v(ji,jj) = 0.5_wp * ( zk_t(ji,jj) + zk_t(ji,jj+1) ) 
     136            ! 
     137            zu0_sd(ji,jj) = 0.5_wp * ( ut0sd(ji,jj) + ut0sd(ji+1,jj) ) 
     138            zv0_sd(ji,jj) = 0.5_wp * ( vt0sd(ji,jj) + vt0sd(ji,jj+1) ) 
     139         END_2D 
    147140      ENDIF 
    148141      ! 
    149142      !                       !==  horizontal Stokes Drift 3D velocity  ==! 
    150143      IF( ll_st_bv2014 ) THEN 
    151          DO jk = 1, jpkm1 
    152             DO jj = 2, jpjm1 
    153                DO ji = 2, jpim1 
    154                   zdep_u = 0.5_wp * ( gdept_n(ji,jj,jk) + gdept_n(ji+1,jj,jk) ) 
    155                   zdep_v = 0.5_wp * ( gdept_n(ji,jj,jk) + gdept_n(ji,jj+1,jk) ) 
    156                   !                           
    157                   zkh_u = zk_u(ji,jj) * zdep_u     ! k * depth 
    158                   zkh_v = zk_v(ji,jj) * zdep_v 
    159                   !                                ! Depth attenuation 
    160                   zda_u = EXP( -2.0_wp*zkh_u ) / ( 1.0_wp + 8.0_wp*zkh_u ) 
    161                   zda_v = EXP( -2.0_wp*zkh_v ) / ( 1.0_wp + 8.0_wp*zkh_v ) 
    162                   ! 
    163                   usd(ji,jj,jk) = zda_u * zu0_sd(ji,jj) * umask(ji,jj,jk) 
    164                   vsd(ji,jj,jk) = zda_v * zv0_sd(ji,jj) * vmask(ji,jj,jk) 
    165                END DO 
    166             END DO 
    167          END DO 
     144         DO_3D_00_00( 1, jpkm1 ) 
     145            zdep_u = 0.5_wp * ( gdept(ji,jj,jk,Kmm) + gdept(ji+1,jj,jk,Kmm) ) 
     146            zdep_v = 0.5_wp * ( gdept(ji,jj,jk,Kmm) + gdept(ji,jj+1,jk,Kmm) ) 
     147            !                           
     148            zkh_u = zk_u(ji,jj) * zdep_u     ! k * depth 
     149            zkh_v = zk_v(ji,jj) * zdep_v 
     150            !                                ! Depth attenuation 
     151            zda_u = EXP( -2.0_wp*zkh_u ) / ( 1.0_wp + 8.0_wp*zkh_u ) 
     152            zda_v = EXP( -2.0_wp*zkh_v ) / ( 1.0_wp + 8.0_wp*zkh_v ) 
     153            ! 
     154            usd(ji,jj,jk) = zda_u * zu0_sd(ji,jj) * umask(ji,jj,jk) 
     155            vsd(ji,jj,jk) = zda_v * zv0_sd(ji,jj) * vmask(ji,jj,jk) 
     156         END_3D 
    168157      ELSE IF( ll_st_li2017 .OR. ll_st_peakfr ) THEN 
    169158         ALLOCATE( zstokes_psi_u_top(jpi,jpj), zstokes_psi_v_top(jpi,jpj) ) 
    170          DO jj = 1, jpjm1              ! exp. wave number & Stokes drift velocity at u- & v-points 
    171             DO ji = 1, jpim1 
    172                zstokes_psi_u_top(ji,jj) = 0._wp 
    173                zstokes_psi_v_top(ji,jj) = 0._wp 
    174             END DO 
    175          END DO 
     159         DO_2D_10_10 
     160            zstokes_psi_u_top(ji,jj) = 0._wp 
     161            zstokes_psi_v_top(ji,jj) = 0._wp 
     162         END_2D 
    176163         zsqrtpi = SQRT(rpi) 
    177164         z_two_thirds = 2.0_wp / 3.0_wp 
    178          DO jk = 1, jpkm1 
    179             DO jj = 2, jpjm1 
    180                DO ji = 2, jpim1 
    181                   zbot_u = ( gdepw_n(ji,jj,jk+1) + gdepw_n(ji+1,jj,jk+1) )  ! 2 * bottom depth 
    182                   zbot_v = ( gdepw_n(ji,jj,jk+1) + gdepw_n(ji,jj+1,jk+1) )  ! 2 * bottom depth 
    183                   zkb_u  = zk_u(ji,jj) * zbot_u                             ! 2 * k * bottom depth 
    184                   zkb_v  = zk_v(ji,jj) * zbot_v                             ! 2 * k * bottom depth 
    185                   ! 
    186                   zke3_u = MAX(1.e-8_wp, 2.0_wp * zk_u(ji,jj) * e3u_n(ji,jj,jk))     ! 2k * thickness 
    187                   zke3_v = MAX(1.e-8_wp, 2.0_wp * zk_v(ji,jj) * e3v_n(ji,jj,jk))     ! 2k * thickness 
    188  
    189                   ! Depth attenuation .... do u component first.. 
    190                   zdepth      = zkb_u 
    191                   zsqrt_depth = SQRT(zdepth) 
    192                   zexp_depth  = EXP(-zdepth) 
    193                   zstokes_psi_u_bot = 1.0_wp - zexp_depth  & 
    194                        &              - z_two_thirds * ( zsqrtpi*zsqrt_depth*zdepth*ERFC(zsqrt_depth) & 
    195                        &              + 1.0_wp - (1.0_wp + zdepth)*zexp_depth ) 
    196                   zda_u                    = ( zstokes_psi_u_bot - zstokes_psi_u_top(ji,jj) ) / zke3_u 
    197                   zstokes_psi_u_top(ji,jj) =   zstokes_psi_u_bot 
    198  
    199                   !         ... and then v component 
    200                   zdepth      =zkb_v 
    201                   zsqrt_depth = SQRT(zdepth) 
    202                   zexp_depth  = EXP(-zdepth) 
    203                   zstokes_psi_v_bot = 1.0_wp - zexp_depth  & 
    204                        &              - z_two_thirds * ( zsqrtpi*zsqrt_depth*zdepth*ERFC(zsqrt_depth) & 
    205                        &              + 1.0_wp - (1.0_wp + zdepth)*zexp_depth ) 
    206                   zda_v                    = ( zstokes_psi_v_bot - zstokes_psi_v_top(ji,jj) ) / zke3_v 
    207                   zstokes_psi_v_top(ji,jj) =   zstokes_psi_v_bot 
    208                   ! 
    209                   usd(ji,jj,jk) = zda_u * zu0_sd(ji,jj) * umask(ji,jj,jk) 
    210                   vsd(ji,jj,jk) = zda_v * zv0_sd(ji,jj) * vmask(ji,jj,jk) 
    211                END DO 
    212             END DO 
    213          END DO 
     165         DO_3D_00_00( 1, jpkm1 ) 
     166            zbot_u = ( gdepw(ji,jj,jk+1,Kmm) + gdepw(ji+1,jj,jk+1,Kmm) )  ! 2 * bottom depth 
     167            zbot_v = ( gdepw(ji,jj,jk+1,Kmm) + gdepw(ji,jj+1,jk+1,Kmm) )  ! 2 * bottom depth 
     168            zkb_u  = zk_u(ji,jj) * zbot_u                             ! 2 * k * bottom depth 
     169            zkb_v  = zk_v(ji,jj) * zbot_v                             ! 2 * k * bottom depth 
     170            ! 
     171            zke3_u = MAX(1.e-8_wp, 2.0_wp * zk_u(ji,jj) * e3u(ji,jj,jk,Kmm))     ! 2k * thickness 
     172            zke3_v = MAX(1.e-8_wp, 2.0_wp * zk_v(ji,jj) * e3v(ji,jj,jk,Kmm))     ! 2k * thickness 
     173 
     174            ! Depth attenuation .... do u component first.. 
     175            zdepth      = zkb_u 
     176            zsqrt_depth = SQRT(zdepth) 
     177            zexp_depth  = EXP(-zdepth) 
     178            zstokes_psi_u_bot = 1.0_wp - zexp_depth  & 
     179                 &              - z_two_thirds * ( zsqrtpi*zsqrt_depth*zdepth*ERFC(zsqrt_depth) & 
     180                 &              + 1.0_wp - (1.0_wp + zdepth)*zexp_depth ) 
     181            zda_u                    = ( zstokes_psi_u_bot - zstokes_psi_u_top(ji,jj) ) / zke3_u 
     182            zstokes_psi_u_top(ji,jj) =   zstokes_psi_u_bot 
     183 
     184            !         ... and then v component 
     185            zdepth      =zkb_v 
     186            zsqrt_depth = SQRT(zdepth) 
     187            zexp_depth  = EXP(-zdepth) 
     188            zstokes_psi_v_bot = 1.0_wp - zexp_depth  & 
     189                 &              - z_two_thirds * ( zsqrtpi*zsqrt_depth*zdepth*ERFC(zsqrt_depth) & 
     190                 &              + 1.0_wp - (1.0_wp + zdepth)*zexp_depth ) 
     191            zda_v                    = ( zstokes_psi_v_bot - zstokes_psi_v_top(ji,jj) ) / zke3_v 
     192            zstokes_psi_v_top(ji,jj) =   zstokes_psi_v_bot 
     193            ! 
     194            usd(ji,jj,jk) = zda_u * zu0_sd(ji,jj) * umask(ji,jj,jk) 
     195            vsd(ji,jj,jk) = zda_v * zv0_sd(ji,jj) * vmask(ji,jj,jk) 
     196         END_3D 
    214197         DEALLOCATE( zstokes_psi_u_top, zstokes_psi_v_top ) 
    215198      ENDIF 
    … …  
    220203      !                       !==  vertical Stokes Drift 3D velocity  ==! 
    221204      ! 
    222       DO jk = 1, jpkm1               ! Horizontal e3*divergence 
    223          DO jj = 2, jpj 
    224             DO ji = fs_2, jpi 
    225                ze3divh(ji,jj,jk) = (  e2u(ji  ,jj) * e3u_n(ji  ,jj,jk) * usd(ji  ,jj,jk)    & 
    226                   &                 - e2u(ji-1,jj) * e3u_n(ji-1,jj,jk) * usd(ji-1,jj,jk)    & 
    227                   &                 + e1v(ji,jj  ) * e3v_n(ji,jj  ,jk) * vsd(ji,jj  ,jk)    & 
    228                   &                 - e1v(ji,jj-1) * e3v_n(ji,jj-1,jk) * vsd(ji,jj-1,jk)  ) * r1_e1e2t(ji,jj) 
    229             END DO 
    230          END DO 
    231       END DO 
     205      DO_3D_01_01( 1, jpkm1 ) 
     206         ze3divh(ji,jj,jk) = (  e2u(ji  ,jj) * e3u(ji  ,jj,jk,Kmm) * usd(ji  ,jj,jk)    & 
     207            &                 - e2u(ji-1,jj) * e3u(ji-1,jj,jk,Kmm) * usd(ji-1,jj,jk)    & 
     208            &                 + e1v(ji,jj  ) * e3v(ji,jj  ,jk,Kmm) * vsd(ji,jj  ,jk)    & 
     209            &                 - e1v(ji,jj-1) * e3v(ji,jj-1,jk,Kmm) * vsd(ji,jj-1,jk)  ) * r1_e1e2t(ji,jj) 
     210      END_3D 
    232211      ! 
    233212#if defined key_agrif 
    … …  
    291270      ! 
    292271      IF( ln_tauw ) THEN 
    293          DO jj = 1, jpjm1 
    294             DO ji = 1, jpim1 
    295                ! Stress components at u- & v-points 
    296                utau(ji,jj) = 0.5_wp * ( tauw_x(ji,jj) + tauw_x(ji+1,jj) ) 
    297                vtau(ji,jj) = 0.5_wp * ( tauw_y(ji,jj) + tauw_y(ji,jj+1) ) 
    298                ! 
    299                ! Stress module at t points 
    300                taum(ji,jj) = SQRT( tauw_x(ji,jj)*tauw_x(ji,jj) + tauw_y(ji,jj)*tauw_y(ji,jj) ) 
    301             END DO 
    302          END DO 
     272         DO_2D_10_10 
     273            ! Stress components at u- & v-points 
     274            utau(ji,jj) = 0.5_wp * ( tauw_x(ji,jj) + tauw_x(ji+1,jj) ) 
     275            vtau(ji,jj) = 0.5_wp * ( tauw_y(ji,jj) + tauw_y(ji,jj+1) ) 
     276            ! 
     277            ! Stress module at t points 
     278            taum(ji,jj) = SQRT( tauw_x(ji,jj)*tauw_x(ji,jj) + tauw_y(ji,jj)*tauw_y(ji,jj) ) 
     279         END_2D 
    303280         CALL lbc_lnk_multi( 'sbcwave', utau(:,:), 'U', -1. , vtau(:,:), 'V', -1. , taum(:,:) , 'T', -1. ) 
    304281      ENDIF 
    … …  
    307284 
    308285 
    309    SUBROUTINE sbc_wave( kt ) 
     286   SUBROUTINE sbc_wave( kt, Kmm ) 
    310287      !!--------------------------------------------------------------------- 
    311288      !!                     ***  ROUTINE sbc_wave  *** 
    … …  
    322299      !!--------------------------------------------------------------------- 
    323300      INTEGER, INTENT(in   ) ::   kt   ! ocean time step 
     301      INTEGER, INTENT(in   ) ::   Kmm  ! ocean time index 
    324302      !!--------------------------------------------------------------------- 
    325303      ! 
    … …  
    361339         ! 
    362340         IF( ( ll_st_bv_li   .AND. jp_hsw>0 .AND. jp_wmp>0 .AND. jp_usd>0 .AND. jp_vsd>0 ) .OR. & 
    363            & ( ll_st_peakfr  .AND. jp_wfr>0 .AND. jp_usd>0 .AND. jp_vsd>0                ) ) CALL sbc_stokes() 
     341           & ( ll_st_peakfr  .AND. jp_wfr>0 .AND. jp_usd>0 .AND. jp_vsd>0                ) ) CALL sbc_stokes( Kmm ) 
    364342         ! 
    365343      ENDIF 
    … …  
    395373      !!--------------------------------------------------------------------- 
    396374      ! 
    397       REWIND( numnam_ref )              ! Namelist namsbc_wave in reference namelist : File for drag coeff. from wave model 
    398375      READ  ( numnam_ref, namsbc_wave, IOSTAT = ios, ERR = 901) 
    399376901   IF( ios /= 0 )   CALL ctl_nam ( ios , 'namsbc_wave in reference namelist' ) 
    400377          
    401       REWIND( numnam_cfg )              ! Namelist namsbc_wave in configuration namelist : File for drag coeff. from wave model 
    402378      READ  ( numnam_cfg, namsbc_wave, IOSTAT = ios, ERR = 902 ) 
    403379902   IF( ios >  0 )   CALL ctl_nam ( ios , 'namsbc_wave in configuration namelist' ) 
Note: See TracChangeset for help on using the changeset viewer.