Changeset 7340


Ignore:
Timestamp:
2016-11-25T16:41:40+01:00 (4 years ago)
Author:
emanuelaclementi
Message:

#1643 Correction after review in development branch 2015/dev_r5936_INGV1_WAVE

Location:
branches/2015/dev_r5936_INGV1_WAVE/NEMOGCM/NEMO/OPA_SRC
Files:
5 edited

Legend:

Unmodified
Added
Removed
  • branches/2015/dev_r5936_INGV1_WAVE/NEMOGCM/NEMO/OPA_SRC/DYN/sshwzv.F90

    r7171 r7340  
    3939   USE wrk_nemo        ! Memory Allocation 
    4040   USE timing          ! Timing 
    41    USE sbcwave,  ONLY:  usd3dt, vsd3dt,wsd3d 
     41   USE sbcwave         ! Stokes velocities 
    4242 
    4343   IMPLICIT NONE 
     
    162162      REAL(wp), POINTER, DIMENSION(:,:  ) ::  z2d 
    163163      REAL(wp), POINTER, DIMENSION(:,:,:) ::  z3d, zhdiv 
    164       REAL(wp) ::  dsshnu, dsshnv 
     164      REAL(wp) ::  zdsshu, zdsshv 
    165165      !!---------------------------------------------------------------------- 
    166166      ! 
     
    209209         END DO 
    210210      ENDIF 
    211 ! 
    212 !     In case ln_wave and ln_sdw, the surface vertical velocity is modified 
    213 !     accounting for Sokes Drift velocity 
    214 ! 
    215       IF ( ln_wave .AND. ln_sdw )  THEN 
    216       ! Compute d(ssh)/dx  and d(ssh)/dy   
    217       ! Compute the surface vertical velocity accounting for the Stokes Drift 
    218       !--------------------------------- 
    219        DO jj = 2 , jpjm1 
    220         DO ji = 2 , jpim1 
    221           dsshnu = ( ssha(ji+1,jj) - ssha(ji,jj) ) / e1u(ji,jj) 
    222           dsshnv = ( ssha(ji,jj+1) - ssha(ji,jj) ) / e2v(ji,jj) 
    223           wn(ji,jj,1) = wn(ji,jj,1) +( usd3dt(ji,jj,1) * dsshnu     & 
    224                 &     + vsd3dt(ji,jj,1) * dsshnv                    & 
    225                 &     - ( wsd3d (ji,jj,1) )) * tmask(ji,jj,1) 
    226         ENDDO 
    227        ENDDO 
     211 
     212      IF( ln_wave .AND. ln_sdw ) THEN 
     213         ! Compute d(ssh)/dx  and d(ssh)/dy   
     214         ! Compute the surface vertical velocity accounting for the Stokes Drift 
     215         DO jj = 1 , jpjm1 
     216            DO ji = 1 , fs_jpim1 
     217              zdsshu = ( sshn(ji+1,jj) - sshn(ji,jj) ) / e1u(ji,jj) 
     218              zdsshv = ( sshn(ji,jj+1) - sshn(ji,jj) ) / e2v(ji,jj) 
     219              wn(ji,jj,1) = wn(ji,jj,1) + ( usd3dt(ji,jj,1) * zdsshu   & 
     220                 &                        + vsd3dt(ji,jj,1) * zdsshv   & 
     221                 &        - wsd3d (ji,jj,1) ) * tmask(ji,jj,1) 
     222            END DO 
     223         END DO 
    228224      ENDIF 
    229225      ! 
  • branches/2015/dev_r5936_INGV1_WAVE/NEMOGCM/NEMO/OPA_SRC/SBC/sbcwave.F90

    r7194 r7340  
    44   !! Wave module  
    55   !!====================================================================== 
    6    !! History :  3.3  !   2011-09  (Adani M)  Original code: Drag Coefficient  
    7    !!         :  3.4  !   2012-10  (Adani M)                 Stokes Drift  
    8    !!            3.6  !   2014-09  (Clementi E, Oddo P)New Stokes Drift Computation 
     6   !! History :  3.3  !   2011-09  (M. Adani)  Original code: Drag Coefficient  
     7   !!         :  3.4  !   2012-10  (M. Adani)  Stokes Drift  
     8   !!            3.6  !   2014-09  (E. Clementi,P. Oddo) New Stokes Drift Computation 
    99   !!---------------------------------------------------------------------- 
    1010 
     
    1313   !!---------------------------------------------------------------------- 
    1414   USE oce            !  
    15    USE sbc_oce       ! Surface boundary condition: ocean fields 
     15   USE sbc_oce        ! Surface boundary condition: ocean fields 
    1616   USE bdy_oce        ! 
    1717   USE domvvl         ! 
    18    ! 
    1918   USE iom            ! I/O manager library 
    2019   USE in_out_manager ! I/O manager 
    2120   USE lib_mpp        ! distribued memory computing library 
    22    USE fldread       ! read input fields 
     21   USE fldread        ! read input fields 
    2322   USE wrk_nemo       ! 
    2423   USE phycst         ! physical constants  
     
    2928   PUBLIC   sbc_wave    ! routine called in sbcmod 
    3029    
    31    INTEGER , PARAMETER ::   jpfld  = 4           ! number of files to read for stokes drift 
    32    INTEGER , PARAMETER ::   jp_usd = 1           ! index of stokes drift  (i-component) (m/s)    at T-point 
    33    INTEGER , PARAMETER ::   jp_vsd = 2           ! index of stokes drift  (j-component) (m/s)    at T-point 
    34    INTEGER , PARAMETER ::   jp_swh = 3           ! index of significant wave hight      (m)      at T-point 
    35    INTEGER , PARAMETER ::   jp_wmp = 4           ! index of mean wave period            (s)      at T-point 
     30   INTEGER, PARAMETER ::   jpfld  = 4           ! number of files to read for stokes drift 
     31   INTEGER, PARAMETER ::   jp_usd = 1           ! index of stokes drift  (i-component) (m/s)    at T-point 
     32   INTEGER, PARAMETER ::   jp_vsd = 2           ! index of stokes drift  (j-component) (m/s)    at T-point 
     33   INTEGER, PARAMETER ::   jp_swh = 3           ! index of significant wave hight      (m)      at T-point 
     34   INTEGER, PARAMETER ::   jp_wmp = 4           ! index of mean wave period            (s)      at T-point 
    3635 
    3736   TYPE(FLD), ALLOCATABLE, DIMENSION(:)  :: sf_cd    ! structure of input fields (file informations, fields read) Drag Coefficient 
     
    4342   REAL(wp), PUBLIC, ALLOCATABLE, DIMENSION(:,:)       :: tauoc_wave 
    4443   REAL(wp), PUBLIC, ALLOCATABLE, DIMENSION(:,:)       :: tsd2d 
    45    REAL(wp), PUBLIC, ALLOCATABLE, DIMENSION(:,:)       :: usd2d, vsd2d 
    46    REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:)     :: usd3d, vsd3d, wsd3d  
    47    REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:)     :: usd3dt, vsd3dt 
    48    REAL(wp), PUBLIC, ALLOCATABLE, DIMENSION(:,:)       :: zusd2dt, zvsd2dt 
     44   REAL(wp), PUBLIC, ALLOCATABLE, DIMENSION(:,:,:)     :: usd3d, vsd3d, wsd3d  
     45   REAL(wp), PUBLIC, ALLOCATABLE, DIMENSION(:,:,:)     :: usd3dt , vsd3dt  
    4946 
    5047   !! * Substitutions 
     
    7471      USE zdf_oce,  ONLY : ln_zdfqiao 
    7572 
    76       INTEGER, INTENT( in  ) ::   kt       ! ocean time step 
    77       ! 
    78       INTEGER                ::   ierror   ! return error code 
    79       INTEGER                ::   ifpr, jj,ji,jk  
    80       INTEGER                ::   ios      ! Local integer output status for namelist read 
    81       ! 
    82       CHARACTER(len=100)     ::  cn_dir                          ! Root directory for location of drag coefficient files 
    83       REAL(wp)                       ::  ztransp, zsp0, zk, zus, zvs 
    84       REAL(wp), DIMENSION(jpi,jpj)   ::  zfac  
    85       REAL(wp), DIMENSION(:,:,:), POINTER :: ze3hdiv   ! 3D workspace 
    86       TYPE(FLD_N), DIMENSION(jpfld) ::   slf_i     ! array of namelist informations on the fields to read 
    87       TYPE(FLD_N)            ::  sn_cdg, sn_usd, sn_vsd,  & 
    88                              &   sn_swh, sn_wmp, sn_wnum, sn_tauoc      ! informations about the fields to be read 
    89       !! 
     73      INTEGER, INTENT( in  ) :: kt       ! ocean time step 
     74      ! 
     75      INTEGER            ::  ierror           ! return error code 
     76      INTEGER            ::  ifpr, jj,ji,jk   ! dummy loop indice  
     77      INTEGER            ::  ios              ! Local integer output status for namelist read 
     78      CHARACTER(len=100) ::  cn_dir           ! Root directory for location of drag coefficient files 
     79      REAL(wp)           ::  ztransp, zfac, zsp0, zk, zus, zvs 
     80      REAL(wp), DIMENSION(jpi,jpj) :: zusd2dt, zvsd2dt 
     81      REAL(wp), DIMENSION(jpi,jpj,jpk) :: ze3hdiv   ! 3D workspace 
     82      TYPE(FLD_N) ::  sn_cdg, sn_usd, sn_vsd            ! informations about the fields to be read 
     83      TYPE(FLD_N) ::  sn_swh, sn_wmp, sn_wnum, sn_tauoc !     "          "    "     "    "  "   " 
     84      TYPE(FLD_N), DIMENSION(jpfld) ::   slf_i          ! array of namelist informations on the fields to read 
    9085      NAMELIST/namsbc_wave/  sn_cdg, cn_dir, sn_usd, sn_vsd, sn_swh, sn_wmp, sn_wnum, sn_tauoc 
    9186      !!--------------------------------------------------------------------- 
     
    10398         IF(lwm) WRITE ( numond, namsbc_wave ) 
    10499         ! 
    105          IF ( ln_cdgw ) THEN 
    106             ALLOCATE( sf_cd(1), STAT=ierror )           !* allocate and fill sf_wave with sn_cdg 
     100         IF( ln_cdgw ) THEN 
     101            ALLOCATE( sf_cd(1), STAT=ierror )           ! allocate and fill sf_wave with sn_cdg 
    107102            IF( ierror > 0 )   CALL ctl_stop( 'STOP', 'sbc_wave: unable to allocate sf_wave structure' ) 
    108103            ! 
    109                                    ALLOCATE( sf_cd(1)%fnow(jpi,jpj,1)  ) 
    110             IF( sn_cdg%ln_tint )   ALLOCATE( sf_cd(1)%fdta(jpi,jpj,1,2) ) 
     104            ALLOCATE( sf_cd(1)%fnow(jpi,jpj,1) ) 
     105            IF( sn_cdg%ln_tint ) ALLOCATE( sf_cd(1)%fdta(jpi,jpj,1,2) ) 
    111106            CALL fld_fill( sf_cd, (/ sn_cdg /), cn_dir, 'sbc_wave', 'Wave module ', 'namsbc_wave' ) 
    112107            ALLOCATE( cdn_wave(jpi,jpj) ) 
    113             cdn_wave(:,:) = 0.0 
    114          ENDIF 
    115  
    116          IF ( ln_tauoc ) THEN 
    117             ALLOCATE( sf_tauoc(1), STAT=ierror )           !* allocate and fill sf_wave with sn_tauoc 
     108         ENDIF 
     109 
     110         IF( ln_tauoc ) THEN 
     111            ALLOCATE( sf_tauoc(1), STAT=ierror )           ! allocate and fill sf_wave with sn_tauoc 
    118112            IF( ierror > 0 )   CALL ctl_stop( 'STOP', 'sbc_wave: unable to allocate sf_wave structure' ) 
    119113            ! 
    120                                    ALLOCATE( sf_tauoc(1)%fnow(jpi,jpj,1)  ) 
    121             IF( sn_tauoc%ln_tint )   ALLOCATE( sf_tauoc(1)%fdta(jpi,jpj,1,2) ) 
     114            ALLOCATE( sf_tauoc(1)%fnow(jpi,jpj,1) ) 
     115            IF( sn_tauoc%ln_tint ) ALLOCATE( sf_tauoc(1)%fdta(jpi,jpj,1,2) ) 
    122116            CALL fld_fill( sf_tauoc, (/ sn_tauoc /), cn_dir, 'sbc_wave', 'Wave module', 'namsbc_wave' ) 
    123117            ALLOCATE( tauoc_wave(jpi,jpj) ) 
    124118            tauoc_wave(:,:) = 0.0 
    125         ENDIF 
    126  
    127          IF ( ln_sdw ) THEN 
     119         ENDIF 
     120 
     121         IF( ln_sdw ) THEN 
    128122            slf_i(jp_usd) = sn_usd ; slf_i(jp_vsd) = sn_vsd; 
    129123            slf_i(jp_swh) = sn_swh ; slf_i(jp_wmp) = sn_wmp; 
    130             ALLOCATE( sf_sd(jpfld), STAT=ierror )           !* allocate and fill sf_sd with stokes drift 
     124            ALLOCATE( sf_sd(jpfld), STAT=ierror )           ! allocate and fill sf_sd with stokes drift 
    131125            IF( ierror > 0 )   CALL ctl_stop( 'STOP', 'sbc_wave: unable to allocate sf_wave structure' ) 
    132126            ! 
    133             DO ifpr= 1, jpfld 
     127            DO ifpr = 1, jpfld 
    134128               ALLOCATE( sf_sd(ifpr)%fnow(jpi,jpj,1) ) 
    135129               IF( slf_i(ifpr)%ln_tint )   ALLOCATE( sf_sd(ifpr)%fdta(jpi,jpj,1,2) ) 
     
    137131 
    138132            CALL fld_fill( sf_sd, slf_i, cn_dir, 'sbc_wave', 'Wave module ', 'namsbc_wave' ) 
    139             ALLOCATE( usd2d(jpi,jpj),vsd2d(jpi,jpj) ) 
    140             ALLOCATE( usd3d(jpi,jpj,jpk),vsd3d(jpi,jpj,jpk),wsd3d(jpi,jpj,jpk) ) 
    141             ALLOCATE( usd3dt(jpi,jpj,jpk),vsd3dt(jpi,jpj,jpk) ) 
     133            ALLOCATE( usd3dt(jpi,jpj,jpk), vsd3dt(jpi,jpj,jpk), wsd3d(jpi,jpj,jpk) ) 
     134            ALLOCATE( usd3d (jpi,jpj,jpk), vsd3d (jpi,jpj,jpk) ) 
    142135            ALLOCATE( swh(jpi,jpj), wmp(jpi,jpj) ) 
    143             ALLOCATE( zusd2dt(jpi,jpj), zvsd2dt(jpi,jpj) ) 
    144             usd3d(:,:,:) = 0._wp   ;   usd2d(:,:) = 0._wp   ;    
    145             vsd3d(:,:,:) = 0._wp   ;   vsd2d(:,:) = 0._wp   ;  
    146             wsd3d(:,:,:) = 0._wp   ; 
    147             usd3dt(:,:,:) = 0._wp  ;   vsd3dt(:,:,:) = 0._wp   ; 
    148             swh  (:,:)   = 0._wp   ;    wmp (:,:) = 0._wp ; 
    149             IF ( ln_zdfqiao ) THEN     !==  Vertical mixing enhancement using Qiao,2010  ==! 
     136            usd3d(:,:,:) = 0._wp         
     137            vsd3d(:,:,:) = 0._wp      
     138            wsd3d(:,:,:) = 0._wp    
     139            IF( ln_zdfqiao ) THEN     !==  Vertical mixing enhancement using Qiao,2010  ==! 
    150140               ALLOCATE( sf_wn(1), STAT=ierror )           !* allocate and fill sf_wave with sn_wnum 
    151                IF( ierror > 0 )   CALL ctl_stop( 'STOP', 'sbc_wave: unable toallocate sf_wave structure' ) 
    152                                       ALLOCATE( sf_wn(1)%fnow(jpi,jpj,1)  ) 
     141               IF( ierror > 0 )   CALL ctl_stop( 'STOP', 'sbc_wave: unable to allocate sf_wave structure' ) 
     142               ALLOCATE( sf_wn(1)%fnow(jpi,jpj,1) ) 
    153143               IF( sn_wnum%ln_tint )  ALLOCATE( sf_wn(1)%fdta(jpi,jpj,1,2) ) 
    154144               CALL fld_fill( sf_wn, (/ sn_wnum /), cn_dir, 'sbc_wave', 'Wave module', 'namsbc_wave' ) 
    155145               ALLOCATE( wnum(jpi,jpj),tsd2d(jpi,jpj) ) 
    156                wnum(:,:) = 0._wp ; tsd2d(:,:) = 0._wp 
    157146            ENDIF 
    158147         ENDIF 
    159148      ENDIF 
    160149      ! 
    161       IF ( ln_cdgw ) THEN              !==  Neutral drag coefficient  ==! 
    162          CALL fld_read( kt, nn_fsbc, sf_cd )      ! read from external forcing 
     150      IF( ln_cdgw ) THEN              !==  Neutral drag coefficient  ==! 
     151         CALL fld_read( kt, nn_fsbc, sf_cd )        ! read from external forcing 
    163152         cdn_wave(:,:) = sf_cd(1)%fnow(:,:,1) 
    164153      ENDIF 
    165154 
    166       IF ( ln_tauoc ) THEN             !==  Wave induced stress  ==! 
    167          CALL fld_read( kt, nn_fsbc, sf_tauoc )      !* read wave norm stress from external forcing 
     155      IF( ln_tauoc ) THEN             !==  Wave induced stress  ==! 
     156         CALL fld_read( kt, nn_fsbc, sf_tauoc )     ! read wave norm stress from external forcing 
    168157         tauoc_wave(:,:) = sf_tauoc(1)%fnow(:,:,1) 
    169158      ENDIF 
    170159 
    171       IF ( ln_sdw )  THEN                         !==  Computation of the 3d Stokes Drift  ==!  
    172          ! 
    173          CALL fld_read( kt, nn_fsbc, sf_sd )      !* read wave parameters from external forcing 
     160      IF( ln_sdw )  THEN              !==  Computation of the 3d Stokes Drift  ==!  
     161         CALL fld_read( kt, nn_fsbc, sf_sd )        ! read wave parameters from external forcing 
    174162         swh(:,:)     = sf_sd(jp_swh)%fnow(:,:,1)   ! significant wave height 
    175163         wmp(:,:)     = sf_sd(jp_wmp)%fnow(:,:,1)   ! wave mean period 
     
    180168         !(DOI: 10.1175/JPO-D-14-0020.1)==!  
    181169         ! 
    182          CALL wrk_alloc( jpi,jpj,jpk, ze3hdiv ) 
    183170         DO jk = 1, jpk 
    184171            DO jj = 1, jpj 
    185172               DO ji = 1, jpi 
    186                ! On T grid 
    187                ! Stokes transport speed estimated from Hs and Tmean 
    188                ztransp = 2.0_wp*rpi*swh(ji,jj)**2.0_wp/(16.0_wp*MAX(wmp(ji,jj),0.0000001_wp)) 
    189                ! Stokes surface speed 
    190                zsp0 = SQRT( zusd2dt(ji,jj)**2 + zvsd2dt(ji,jj)**2) 
    191                ! Wavenumber scale 
    192                zk = ABS(zsp0)/MAX(ABS(5.97_wp*ztransp),0.0000001_wp) 
    193                ! Depth attenuation 
    194                zfac(ji,jj) = EXP(-2.0_wp*zk*fsdept(ji,jj,jk))/(1.0_wp+8.0_wp*zk*fsdept(ji,jj,jk)) 
    195                END DO 
    196             END DO 
    197             ! 
    198             DO jj = 1, jpj 
    199                DO ji = 1, jpi 
    200                  usd3dt(ji,jj,jk) = zfac(ji,jj) * zusd2dt (ji,jj) * tmask(ji,jj,jk) 
    201                  vsd3dt(ji,jj,jk) = zfac(ji,jj) * zvsd2dt (ji,jj) * tmask(ji,jj,jk) 
     173                  ! On T grid 
     174                  ! Stokes transport speed estimated from Hs and Tmean 
     175                  ztransp = 2.0_wp*rpi*swh(ji,jj)**2.0_wp/(16.0_wp*MAX(wmp(ji,jj),0.0000001_wp)) 
     176                  ! Stokes surface speed 
     177                  zsp0 = SQRT( zusd2dt(ji,jj)**2 + zvsd2dt(ji,jj)**2 ) 
     178                  ! Wavenumber scale 
     179                  zk = ABS(zsp0)/MAX(ABS(5.97_wp*ztransp),0.0000001_wp) 
     180                  ! Depth attenuation 
     181                  zfac = EXP(-2.0_wp*zk*fsdept(ji,jj,jk))/(1.0_wp+8.0_wp*zk*fsdept(ji,jj,jk)) 
     182                  ! 
     183                  usd3dt(ji,jj,jk) = zfac * zusd2dt(ji,jj) * tmask(ji,jj,jk) 
     184                  vsd3dt(ji,jj,jk) = zfac * zvsd2dt(ji,jj) * tmask(ji,jj,jk) 
    202185               END DO 
    203186            END DO 
     
    207190            DO jj = 1, jpjm1 
    208191               DO ji = 1, jpim1 
    209                usd3d(ji,jj,jk) =   0.5 *  umask(ji,jj,jk)  *        & 
    210                                &  (usd3dt(ji,jj,jk) + usd3dt(ji+1,jj,jk)) 
    211                vsd3d(ji,jj,jk) =  0.5 *  vmask(ji,jj,jk)  *           & 
    212                                &  (vsd3dt(ji,jj,jk) + vsd3dt(ji,jj+1,jk)) 
     192                  usd3d(ji,jj,jk) = 0.5 *  umask(ji,jj,jk) *   & 
     193                                  &  ( usd3dt(ji,jj,jk) + usd3dt(ji+1,jj,jk) ) 
     194                  vsd3d(ji,jj,jk) = 0.5 *  vmask(ji,jj,jk) *   & 
     195                                  &  ( vsd3dt(ji,jj,jk) + vsd3dt(ji,jj+1,jk) ) 
    213196               END DO 
    214197            END DO 
     
    218201         CALL lbc_lnk( vsd3d(:,:,:), 'V', -1. ) 
    219202         ! 
    220          DO jk = 1, jpkm1                       ! e3t * Horizontal divergence 
    221             DO jj = 2, jpjm1 
    222                DO ji = fs_2, fs_jpim1   ! vector opt. 
    223                   ze3hdiv(ji,jj,jk) = (  e2u(ji  ,jj) * fse3u_n(ji  ,jj,jk) * usd3d(ji  ,jj,jk)     & 
    224                      &                 - e2u(ji-1,jj) * fse3u_n(ji-1,jj,jk) * usd3d(ji-1,jj,jk)     & 
    225                      &                 + e1v(ji,jj  ) * fse3v_n(ji,jj  ,jk) * vsd3d(ji,jj  ,jk)     & 
    226                      &                 - e1v(ji,jj-1) * fse3v_n(ji,jj-1,jk) * vsd3d(ji,jj-1,jk)   ) * r1_e1e2t(ji,jj) 
     203         DO jk = 1, jpkm1               ! Horizontal divergence 
     204            DO jj = 2, jpj 
     205               DO ji = fs_2, jpi    
     206                  ze3hdiv(ji,jj,jk) = (  e2u(ji  ,jj) * usd3d(ji  ,jj,jk)     & 
     207                     &                 - e2u(ji-1,jj) * usd3d(ji-1,jj,jk)     & 
     208                     &                 + e1v(ji,jj  ) * vsd3d(ji,jj  ,jk)     & 
     209                     &                 - e1v(ji,jj-1) * vsd3d(ji,jj-1,jk)   ) * r1_e1e2t(ji,jj) 
    227210               END DO   
    228211            END DO   
    229             IF( .NOT. AGRIF_Root() ) THEN 
    230                IF( nbondi ==  1 .OR. nbondi == 2 )   ze3hdiv(nlci-1,   :  ,jk) = 0._wp      ! east 
    231                IF( nbondi == -1 .OR. nbondi == 2 )   ze3hdiv(  2   ,   :  ,jk) = 0._wp      ! west 
    232                IF( nbondj ==  1 .OR. nbondj == 2 )   ze3hdiv(  :   ,nlcj-1,jk) = 0._wp      ! north 
    233                IF( nbondj == -1 .OR. nbondj == 2 )   ze3hdiv(  :   ,  2   ,jk) = 0._wp      ! south 
    234             ENDIF 
    235          END DO 
     212         END DO 
     213         ! 
     214         IF( .NOT. AGRIF_Root() ) THEN 
     215            IF( nbondi ==  1 .OR. nbondi == 2 )   ze3hdiv(nlci-1,   :  ,:) = 0._wp      ! east 
     216            IF( nbondi == -1 .OR. nbondi == 2 )   ze3hdiv(  2   ,   :  ,:) = 0._wp      ! west 
     217            IF( nbondj ==  1 .OR. nbondj == 2 )   ze3hdiv(  :   ,nlcj-1,:) = 0._wp      ! north 
     218            IF( nbondj == -1 .OR. nbondj == 2 )   ze3hdiv(  :   ,  2   ,:) = 0._wp      ! south 
     219         ENDIF 
     220         ! 
    236221         CALL lbc_lnk( ze3hdiv, 'T', 1. )  
    237222         ! 
    238          DO jk = jpkm1, 1, -1                   !* integrate from the bottom the e3t * hor. divergence 
    239             wsd3d(:,:,jk) = wsd3d(:,:,jk+1) - ze3hdiv(:,:,jk) 
     223         DO jk = jpkm1, 1, -1                   ! integrate from the bottom the e3t * hor. divergence 
     224            wsd3d(:,:,jk) = wsd3d(:,:,jk+1) - fse3t_n(:,:,jk) * ze3hdiv(:,:,jk) 
    240225         END DO 
    241226#if defined key_bdy 
     
    246231         ENDIF 
    247232#endif 
    248          CALL wrk_dealloc( jpi,jpj,jpk, ze3hdiv ) 
    249233 
    250234         IF ( ln_zdfqiao ) THEN 
    251             CALL fld_read( kt, nn_fsbc, sf_wn )      !* read wave parameters from external forcing 
     235            CALL fld_read( kt, nn_fsbc, sf_wn )      ! read wave parameters from external forcing 
    252236            wnum(:,:) = sf_wn(1)%fnow(:,:,1) 
    253237            
     
    256240            DO jj = 1, jpj 
    257241               DO ji = 1, jpi 
    258                    tsd2d(ji,jj) = SQRT( zusd2dt(ji,jj)**2 + zvsd2dt(ji,jj)**2) 
     242                  tsd2d(ji,jj) = SQRT( zusd2dt(ji,jj)**2 + zvsd2dt(ji,jj)**2 )  
    259243               END DO 
    260244            END DO 
  • branches/2015/dev_r5936_INGV1_WAVE/NEMOGCM/NEMO/OPA_SRC/TRA/traadv.F90

    r5983 r7340  
    3737   USE sbcwave        ! wave module 
    3838   USE sbc_oce        ! surface boundary condition: ocean 
    39  
    4039   USE diaptr         ! Poleward heat transport  
    41    USE sbcwave        ! wave module 
    42    USE sbc_oce        ! surface boundary condition: ocean 
    4340 
    4441   IMPLICIT NONE 
     
    111108      ! 
    112109      !                                         !==  effective transport  ==! 
    113       IF (ln_wave .AND. ln_sdw) THEN 
     110      IF( ln_wave .AND. ln_sdw ) THEN 
    114111         DO jk = 1, jpkm1 
    115112            zun(:,:,jk) = e2u(:,:) * fse3u(:,:,jk) *      & 
    116                         &  ( un(:,:,jk) + usd3d(:,:,jk) )                     !eulerian transport + Stokes Drift 
     113                        &  ( un(:,:,jk) + usd3d(:,:,jk) )                       ! eulerian transport + Stokes Drift 
    117114            zvn(:,:,jk) = e1v(:,:) * fse3v(:,:,jk) *      & 
    118115                        &  ( vn(:,:,jk) + vsd3d(:,:,jk) ) 
     
    121118         END DO 
    122119      ELSE 
    123       DO jk = 1, jpkm1 
    124          zun(:,:,jk) = e2u  (:,:) * fse3u(:,:,jk) * un(:,:,jk)                  ! eulerian transport only 
    125          zvn(:,:,jk) = e1v  (:,:) * fse3v(:,:,jk) * vn(:,:,jk) 
    126          zwn(:,:,jk) = e1e2t(:,:)                 * wn(:,:,jk) 
    127       END DO 
     120         DO jk = 1, jpkm1 
     121            zun(:,:,jk) = e2u  (:,:) * fse3u(:,:,jk) * un(:,:,jk)               ! eulerian transport only 
     122            zvn(:,:,jk) = e1v  (:,:) * fse3v(:,:,jk) * vn(:,:,jk) 
     123            zwn(:,:,jk) = e1e2t(:,:)                 * wn(:,:,jk) 
     124         END DO 
    128125      ENDIF 
    129126      ! 
  • branches/2015/dev_r5936_INGV1_WAVE/NEMOGCM/NEMO/OPA_SRC/ZDF/zdfqiao.F90

    r7171 r7340  
    66   !! History :  3.6  !  2014-10  (E. Clementi)  Original code 
    77   !!---------------------------------------------------------------------- 
    8    !!---------------------------------------------------------------------- 
    9    !!   qiao_init       
    108   !!   zdf_qiao        : compute Qiao parameters 
    119   !!---------------------------------------------------------------------- 
    1210 
    13    USE iom             ! I/O manager library 
    1411   USE in_out_manager  ! I/O manager 
    1512   USE lib_mpp         ! distribued memory computing library 
     
    1815   USE sbcwave         ! wave module 
    1916   USE dom_oce 
     17   USE lbclnk          ! ocean lateral boundary conditions (or mpp link)   
    2018    
    21    !!---------------------------------------------------------------------- 
    22    !!   qiao_init       : compute QBv: Qiao terms to be added to vertical eddy 
    23    !!                     diffusivity and viscosity coefficients  
    24    !!---------------------------------------------------------------------- 
    25  
    2619   IMPLICIT NONE 
    2720   PRIVATE 
    2821 
    29    PUBLIC   zdf_qiao    ! routine called in zdf_ric 
     22   PUBLIC zdf_qiao    ! routine called in step 
    3023 
    31    REAL(wp),PUBLIC,ALLOCATABLE,DIMENSION (:,:,:)     :: QBv, QBvu, QBvv 
     24   REAL(wp), PUBLIC, ALLOCATABLE, DIMENSION(:,:,:) :: qbv, qbvu, qbvv 
    3225 
    3326   !! * Substitutions 
    3427#  include "domzgr_substitute.h90" 
     28#  include "vectopt_loop_substitute.h90" 
    3529   !!---------------------------------------------------------------------- 
    3630   !! NEMO/OPA 4.0 , NEMO Consortium (2011)  
     
    4539      !!                     ***  ROUTINE zdf_qiao *** 
    4640      !! 
    47       !! ** Purpose :Compute the Qiao term (QBv) to be added to 
     41      !! ** Purpose :Compute the Qiao term (qbv) to be added to 
    4842      !!             vertical viscosity and diffusivity coeffs.   
    4943      !! 
    50       !! ** Method  :QBv = alpha * A * Us(0) * exp (3 * k * z) 
     44      !! ** Method  :qbv = alpha * A * Us(0) * exp (3 * k * z) 
    5145      !!              
    5246      !! ** action  :Compute the Qiao wave dependent term  
     
    5650      INTEGER, INTENT( in  ) ::  kt   ! ocean time step 
    5751      ! 
    58       INTEGER                ::  jj, ji, jk 
     52      INTEGER :: jj, ji, jk   ! dummy loop indices 
    5953      !!--------------------------------------------------------------------- 
    6054      ! 
    61       ! 
    62       !                                         ! -------------------- ! 
    6355      IF( kt == nit000 ) THEN                   ! First call kt=nit000 ! 
    64          ALLOCATE(QBv(jpi,jpj,jpk))             ! -------------------- ! 
    65          ALLOCATE(QBvu(jpi,jpj,jpk)) 
    66          ALLOCATE(QBvv(jpi,jpj,jpk)) 
     56         IF( .NOT. ( ln_wave .AND. ln_sdw ) )   & 
     57            &   CALL ctl_stop ( 'Ask for wave Qiao enhanced turbulence but ln_wave   & 
     58            &                    and ln_sdw have to be activated') 
     59         IF( zdf_qiao_alloc() /= 0 )   & 
     60            &   CALL ctl_stop( 'STOP', 'zdf_qiao : unable to allocate arrays' ) 
    6761      ENDIF 
    6862 
    69       QBv (:,:,:) = 0.0 
    70       QBvu(:,:,:) = 0.0 
    71       QBvv(:,:,:) = 0.0 
    72  
    7363      ! 
    74       ! Compute the Qiao term Bv (QBv) to be added to 
     64      ! Compute the Qiao term Bv (qbv) to be added to 
    7565      ! vertical viscosity and diffusivity 
    76       ! QBv = alpha * A * Us(0) * exp (3 * k * z) 
     66      ! qbv = alpha * A * Us(0) * exp (3 * k * z) 
    7767      ! alpha here is set to 1 
    7868      !--------------------------------------------------------------------------------- 
    7969      ! 
    80       IF ( ln_wave ) THEN 
    81          DO jk = 1, jpk 
    82             DO jj = 1, jpjm1 
    83                DO ji = 1, jpim1 
    84                   QBv(ji,jj,jk) = 1.0 * 0.353553 * swh(ji,jj) * tsd2d(ji,jj) *       & 
    85                &              exp(3.0 * wnum(ji,jj) *                                &                      
    86                &              (-MIN( fsdept(ji  ,jj  ,jk) , fsdept(ji+1,jj  ,jk),    & 
    87                &                     fsdept(ji  ,jj+1,jk) , fsdept(ji+1,jj+1,jk)))) 
    88                END DO 
     70      DO jk = 1, jpk 
     71         DO jj = 1, jpjm1 
     72            DO ji = 1, fs_jpim1 
     73               qbv(ji,jj,jk) = 1.0 * 0.353553 * swh(ji,jj) * tsd2d(ji,jj) *           & 
     74            &                  EXP(3.0 * wnum(ji,jj) *                                &                      
     75            &                  (-MIN( fsdepw(ji  ,jj  ,jk), fsdepw(ji+1,jj  ,jk),     & 
     76            &                         fsdepw(ji  ,jj+1,jk), fsdepw(ji+1,jj+1,jk))))   & 
     77            &                          * wmask(ji,jj,jk) 
    8978            END DO 
    9079         END DO 
    91  
    92          QBv(jpi,:,:)=QBv(jpim1,:,:) 
    93          QBv(:,jpj,:)=QBv(:,jpjm1,:) 
    94  
    95          ! 
    96          ! Interpolate Qiao parameter QBv into the grid_U and grid_V 
    97          !------------------------------------------------- 
    98          ! 
    99          DO jk = 1, jpk 
    100             DO jj = 1, jpjm1 
    101                DO ji = 1, jpim1 
    102                   QBvu(ji,jj,jk) = 0.5 *  umask(ji,jj,jk)  *               & 
    103                &           ( QBv(ji  ,jj,jk) * tmask(ji  ,jj,jk)           & 
    104                &           + QBv(ji+1,jj,jk) * tmask(ji+1,jj,jk) ) 
    105                   QBvv(ji,jj,jk) = 0.5 *  vmask(ji,jj,jk)  *               & 
    106                &           ( QBv(ji,jj  ,jk) * tmask(ji,jj  ,jk)           & 
    107                &           + QBv(ji,jj+1,jk) * tmask(ji,jj+1,jk) ) 
    108                END DO 
     80      END DO 
     81      ! 
     82      CALL lbc_lnk( qbv, 'W', 1. )   ! Lateral boundary conditions 
     83          
     84      ! 
     85      ! Interpolate Qiao parameter qbv into the grid_U and grid_V 
     86      !---------------------------------------------------------- 
     87      ! 
     88      DO jk = 1, jpk 
     89         DO jj = 1, jpjm1 
     90            DO ji = 1, fs_jpim1 
     91               qbvu(ji,jj,jk) = 0.5 * wumask(ji,jj,jk)  *              &   
     92            &                  ( qbv(ji,jj,jk) + qbv(ji+1,jj  ,jk) ) 
     93               qbvv(ji,jj,jk) = 0.5 * wvmask(ji,jj,jk)  *              & 
     94            &                  ( qbv(ji,jj,jk) + qbv(ji  ,jj+1,jk) ) 
    10995            END DO 
    11096         END DO 
    111          !  
    112          QBvu(jpi,:,:)=QBvu(jpim1,:,:) 
    113          QBvu(:,jpj,:)=QBvu(:,jpjm1,:) 
    114          QBvv(jpi,:,:)=QBvv(jpim1,:,:) 
    115          QBvv(:,jpj,:)=QBvv(:,jpjm1,:) 
     97      END DO 
     98      !  
     99      CALL lbc_lnk( qbvu, 'U', 1. ) ; CALL lbc_lnk( qbvv, 'V', 1. )   ! Lateral boundary conditions 
    116100 
    117         ELSE 
    118            CALL ctl_stop( 'STOP', 'To use Qiao formulation you have to set: ln_wave=.true.') 
    119         ENDIF 
    120         ! 
     101      ! Enhance vertical mixing coeff.          
     102      !------------------------------- 
     103      ! 
     104      DO jk = 1, jpkm1 
     105         DO jj = 1, jpj 
     106            DO ji = 1, jpi 
     107               avmu(ji,jj,jk) = ( avmu(ji,jj,jk) + qbvu(ji,jj,jk) ) * umask(ji,jj,jk) 
     108               avmv(ji,jj,jk) = ( avmv(ji,jj,jk) + qbvv(ji,jj,jk) ) * vmask(ji,jj,jk) 
     109               avt (ji,jj,jk) = ( avt (ji,jj,jk) + qbv (ji,jj,jk) ) * tmask(ji,jj,jk) 
     110            END DO 
     111         END DO 
     112      END DO 
     113      ! 
    121114   END SUBROUTINE zdf_qiao 
     115 
     116   INTEGER FUNCTION zdf_qiao_alloc() 
     117      !!---------------------------------------------------------------------- 
     118      !!                ***  FUNCTION zdf_qiao_alloc  *** 
     119      !!---------------------------------------------------------------------- 
     120      ALLOCATE( qbv(jpi,jpj,jpk), qbvu(jpi,jpj,jpk), qbvv(jpi,jpj,jpk),   & 
     121         &      STAT = zdf_qiao_alloc ) 
     122      ! 
     123      IF( lk_mpp             )  CALL mpp_sum ( zdf_qiao_alloc ) 
     124      IF( zdf_qiao_alloc > 0 )  CALL ctl_warn('zdf_qiao_alloc: allocation of arrays failed.') 
     125      ! 
     126   END FUNCTION zdf_qiao_alloc 
    122127       
    123128   !!====================================================================== 
  • branches/2015/dev_r5936_INGV1_WAVE/NEMOGCM/NEMO/OPA_SRC/step.F90

    r7221 r7340  
    130130                         CALL zdf_bfr( kstp )         ! bottom friction (if quadratic) 
    131131      !                                               ! Vertical eddy viscosity and diffusivity coefficients 
    132       IF( lk_zdfric  )   CALL zdf_ric( kstp )            ! Richardson number dependent Kz 
    133       IF( lk_zdftke  )   CALL zdf_tke( kstp )            ! TKE closure scheme for Kz 
    134       IF( lk_zdfgls  )   CALL zdf_gls( kstp )            ! GLS closure scheme for Kz 
    135       IF( ln_zdfqiao )   THEN 
    136       !Activated Qiao enhanced turbulence but neither ln_wave or ln_sdw are activated 
    137        IF ( .NOT.( ln_wave .AND. ln_sdw ) )   THEN 
    138         CALL ctl_stop ( 'Ask for wave Qiao enhanced turbulence but ln_wave and ln_sdw have to be activated') 
    139        ELSE 
    140         CALL zdf_qiao(kstp )                             ! Qiao vertical mixing  
    141          DO jk = 1, jpkm1 
    142           DO jj = 1, jpj 
    143             DO ji = 1, jpi 
    144                avmu(ji,jj,jk) = (avmu(ji,jj,jk) + QBvu(ji,jj,jk)) * umask(ji,jj,jk) 
    145                avmv(ji,jj,jk) = (avmv(ji,jj,jk) + QBvv(ji,jj,jk)) * vmask(ji,jj,jk) 
    146                avt( ji,jj,jk) = (avt( ji,jj,jk) + QBv(ji,jj,jk))  * tmask(ji,jj,jk) 
    147             END DO 
    148           END DO 
    149          END DO 
    150        ENDIF  
    151       ENDIF 
    152       ! 
    153       IF( lk_zdfcst  ) THEN                              ! Constant Kz (reset avt, avm[uv] to the background value) 
     132      IF( lk_zdfric  )   CALL zdf_ric ( kstp )             ! Richardson number dependent Kz 
     133      IF( lk_zdftke  )   CALL zdf_tke ( kstp )             ! TKE closure scheme for Kz 
     134      IF( lk_zdfgls  )   CALL zdf_gls ( kstp )             ! GLS closure scheme for Kz 
     135      IF( ln_zdfqiao )   CALL zdf_qiao( kstp )             ! Qiao vertical mixing  
     136      ! 
     137      IF( lk_zdfcst  ) THEN                                ! Constant Kz (reset avt, avm[uv] to the background value) 
    154138         avt (:,:,:) = rn_avt0 * wmask (:,:,:) 
    155139         avmu(:,:,:) = rn_avm0 * wumask(:,:,:) 
Note: See TracChangeset for help on using the changeset viewer.