Ignore:
Timestamp:
2011-09-28T10:18:52+02:00 (9 years ago)
Author:
hliu
Message:

1) added the semi-implicit bottom friction code (3 in DYN, 1 in ZDF), 2) added par_AMM7/12.h90, 3) added two arch files for NOCL mobius and ubuntu linux

File:
1 edited

Legend:

Unmodified
Added
Removed
  • branches/2011/dev_r2802_NOCL_bfrimp/NEMOGCM/NEMO/OPA_SRC/DYN/dynspg_ts.F90

    r2724 r2872  
    119119      INTEGER  ::   ji, jj, jk, jn   ! dummy loop indices 
    120120      INTEGER  ::   icycle           ! local scalar 
     121      INTEGER  ::   ikbu, ikbv       ! local scalar 
    121122      REAL(wp) ::   zraur, zcoef, z2dt_e, z2dt_b     ! local scalars 
    122123      REAL(wp) ::   z1_8, zx1, zy1                   !   -      - 
     
    124125      REAL(wp) ::   zu_spg, zu_cor, zu_sld, zu_asp   !   -      - 
    125126      REAL(wp) ::   zv_spg, zv_cor, zv_sld, zv_asp   !   -      - 
     127      REAL(wp) ::   ua_btm, va_btm                   !   -      - 
    126128      !!---------------------------------------------------------------------- 
    127129 
     
    147149         hvr_e (:,:) = hvr  (:,:) 
    148150         IF( ln_dynvor_een ) THEN 
    149             ftne(1,:) = 0.e0   ;   ftnw(1,:) = 0.e0   ;   ftse(1,:) = 0.e0   ;   ftsw(1,:) = 0.e0 
     151            ftne(1,:) = 0._wp   ;   ftnw(1,:) = 0._wp   ;   ftse(1,:) = 0._wp   ;   ftsw(1,:) = 0._wp 
    150152            DO jj = 2, jpj 
    151153               DO ji = fs_2, jpi   ! vector opt. 
    152                   ftne(ji,jj) = ( ff(ji-1,jj  ) + ff(ji  ,jj  ) + ff(ji  ,jj-1) ) / 3. 
    153                   ftnw(ji,jj) = ( ff(ji-1,jj-1) + ff(ji-1,jj  ) + ff(ji  ,jj  ) ) / 3. 
    154                   ftse(ji,jj) = ( ff(ji  ,jj  ) + ff(ji  ,jj-1) + ff(ji-1,jj-1) ) / 3. 
    155                   ftsw(ji,jj) = ( ff(ji  ,jj-1) + ff(ji-1,jj-1) + ff(ji-1,jj  ) ) / 3. 
     154                  ftne(ji,jj) = ( ff(ji-1,jj  ) + ff(ji  ,jj  ) + ff(ji  ,jj-1) ) / 3._wp 
     155                  ftnw(ji,jj) = ( ff(ji-1,jj-1) + ff(ji-1,jj  ) + ff(ji  ,jj  ) ) / 3._wp 
     156                  ftse(ji,jj) = ( ff(ji  ,jj  ) + ff(ji  ,jj-1) + ff(ji-1,jj-1) ) / 3._wp 
     157                  ftsw(ji,jj) = ( ff(ji  ,jj-1) + ff(ji-1,jj-1) + ff(ji-1,jj  ) ) / 3._wp 
    156158               END DO 
    157159            END DO 
     
    161163 
    162164      !                                   !* Local constant initialization 
    163       z2dt_b = 2.0 * rdt                                    ! baroclinic time step 
    164       z1_8 = 0.5 * 0.25                                     ! coefficient for vorticity estimates 
    165       z1_4 = 0.5 * 0.5 
    166       zraur  = 1. / rau0                                    ! 1 / volumic mass 
    167       ! 
    168       zhdiv(:,:) = 0.e0                                     ! barotropic divergence 
    169       zu_sld = 0.e0   ;   zu_asp = 0.e0                     ! tides trends (lk_tide=F) 
    170       zv_sld = 0.e0   ;   zv_asp = 0.e0 
     165      z2dt_b = 2.0_wp * rdt                                 ! baroclinic time step 
     166      z1_8   = 0.5_wp * 0.25_wp                             ! coefficient for vorticity estimates 
     167      z1_4   = 0.5_wp * 0.5_wp 
     168      zraur  = 1._wp / rau0                                 ! 1 / volumic mass 
     169      ! 
     170      zhdiv(:,:) = 0._wp                                    ! barotropic divergence 
     171      zu_sld = 0._wp   ;   zu_asp = 0._wp                   ! tides trends (lk_tide=F) 
     172      zv_sld = 0._wp   ;   zv_asp = 0._wp 
    171173 
    172174      ! ----------------------------------------------------------------------------- 
     
    176178      !                                   !* e3*d/dt(Ua), e3*Ub, e3*Vn (Vertically integrated) 
    177179      !                                   ! -------------------------- 
    178       zua(:,:) = 0.e0   ;   zun(:,:) = 0.e0   ;   ub_b(:,:) = 0.e0 
    179       zva(:,:) = 0.e0   ;   zvn(:,:) = 0.e0   ;   vb_b(:,:) = 0.e0 
     180      zua(:,:) = 0._wp   ;   zun(:,:) = 0._wp   ;   ub_b(:,:) = 0._wp 
     181      zva(:,:) = 0._wp   ;   zvn(:,:) = 0._wp   ;   vb_b(:,:) = 0._wp 
    180182      ! 
    181183      DO jk = 1, jpkm1 
     
    205207      END DO 
    206208 
    207       !                                   !* baroclinic momentum trend (remove the vertical mean trend) 
    208       DO jk = 1, jpkm1                    ! -------------------------- 
    209          DO jj = 2, jpjm1 
    210             DO ji = fs_2, fs_jpim1   ! vector opt. 
    211                ua(ji,jj,jk) = ua(ji,jj,jk) - zua(ji,jj) * hur(ji,jj) 
    212                va(ji,jj,jk) = va(ji,jj,jk) - zva(ji,jj) * hvr(ji,jj) 
    213             END DO 
    214          END DO 
    215       END DO 
    216209 
    217210      !                                   !* barotropic Coriolis trends * H (vorticity scheme dependent) 
     
    260253         DO jj = 2, jpjm1  
    261254            DO ji = fs_2, fs_jpim1   ! vector opt. 
    262                zcu(ji,jj) = zcu(ji,jj) - grav * (  ( rhd(ji+1,jj  ,1) + 1 ) * sshn(ji+1,jj  )    & 
    263                   &                              - ( rhd(ji  ,jj  ,1) + 1 ) * sshn(ji  ,jj  )  ) * hu(ji,jj) / e1u(ji,jj) 
    264                zcv(ji,jj) = zcv(ji,jj) - grav * (  ( rhd(ji  ,jj+1,1) + 1 ) * sshn(ji  ,jj+1)    & 
    265                   &                              - ( rhd(ji  ,jj  ,1) + 1 ) * sshn(ji  ,jj  )  ) * hv(ji,jj) / e2v(ji,jj) 
     255      !     remove the rhd term according to J. Chanut's suggestion             
     256               zcu(ji,jj) = zcu(ji,jj) - grav * (  sshn(ji+1,jj  ) - sshn(ji  ,jj  )  ) * hu(ji,jj) / e1u(ji,jj) 
     257               zcv(ji,jj) = zcv(ji,jj) - grav * (  sshn(ji  ,jj+1) - sshn(ji  ,jj  )  ) * hv(ji,jj) / e2v(ji,jj) 
    266258            END DO 
    267259         END DO 
    268260      ENDIF 
     261 
     262     IF(ln_bfrimp.AND.lk_vvl) THEN 
     263!    Remove the bottom stress trend from 3-D sea surface level gradient and Coriolis force            
     264!    H. Liu, Sept. 2011 
     265      DO jj = 2, jpjm1          
     266         DO ji = fs_2, fs_jpim1 
     267            ikbu = mbku(ji,jj) 
     268            ikbv = mbkv(ji,jj) 
     269            ua_btm = zcu(ji,jj) * z2dt_b * hur(ji,jj) * umask (ji,jj,ikbu) 
     270            va_btm = zcv(ji,jj) * z2dt_b * hvr(ji,jj) * vmask (ji,jj,ikbv) 
     271 
     272            ua(ji,jj,ikbu) = ua(ji,jj,ikbu) - bfrua(ji,jj) * ua_btm / fse3u(ji,jj,ikbu) 
     273            va(ji,jj,ikbv) = va(ji,jj,ikbv) - bfrva(ji,jj) * va_btm / fse3v(ji,jj,ikbv) 
     274             
     275            zua(ji,jj) = zua(ji,jj) - bfrua(ji,jj) * ua_btm 
     276            zva(ji,jj) = zva(ji,jj) - bfrva(ji,jj) * va_btm 
     277         END DO 
     278      END DO 
     279     END IF 
     280      !                                   !* baroclinic momentum trend (remove the vertical mean trend) 
     281      DO jk = 1, jpkm1                    ! -------------------------- 
     282         DO jj = 2, jpjm1 
     283            DO ji = fs_2, fs_jpim1   ! vector opt. 
     284               ua(ji,jj,jk) = ua(ji,jj,jk) - zua(ji,jj) * hur(ji,jj) 
     285               va(ji,jj,jk) = va(ji,jj,jk) - zva(ji,jj) * hvr(ji,jj) 
     286            END DO 
     287         END DO 
     288      END DO 
    269289 
    270290      DO jj = 2, jpjm1                             ! Remove coriolis term (and possibly spg) from barotropic trend 
     
    279299      !                                             ! from the barotropic transport trend 
    280300      zcoef = -1. / z2dt_b 
     301      IF(.NOT.ln_bfrimp) THEN 
    281302# if defined key_vectopt_loop 
    282303      DO jj = 1, 1 
     
    311332         END DO 
    312333      ENDIF 
     334      END IF 
    313335 
    314336      !                                   !* d/dt(Ua), Ub, Vn (Vertical mean velocity) 
     
    410432               DO ji = fs_2, fs_jpim1   ! vector opt. 
    411433                  ! surface pressure gradient 
    412                   IF( lk_vvl) THEN 
    413                      zu_spg = -grav * (  ( rhd(ji+1,jj ,1) + 1 ) * sshn_e(ji+1,jj  )    & 
    414                         &              - ( rhd(ji  ,jj ,1) + 1 ) * sshn_e(ji  ,jj  )  ) / e1u(ji,jj) 
    415                      zv_spg = -grav * (  ( rhd(ji ,jj+1,1) + 1 ) * sshn_e(ji  ,jj+1)    & 
    416                         &              - ( rhd(ji ,jj  ,1) + 1 ) * sshn_e(ji  ,jj  )  ) / e2v(ji,jj) 
    417                   ELSE 
    418                      zu_spg = -grav * ( sshn_e(ji+1,jj) - sshn_e(ji,jj) ) / e1u(ji,jj) 
    419                      zv_spg = -grav * ( sshn_e(ji,jj+1) - sshn_e(ji,jj) ) / e2v(ji,jj) 
    420                   ENDIF 
     434                  zu_spg = -grav * ( sshn_e(ji+1,jj) - sshn_e(ji,jj) ) / e1u(ji,jj) 
     435                  zv_spg = -grav * ( sshn_e(ji,jj+1) - sshn_e(ji,jj) ) / e2v(ji,jj) 
    421436                  ! energy conserving formulation for planetary vorticity term 
    422437                  zy1 = ( zwy(ji  ,jj-1) + zwy(ji+1,jj-1) ) / e1u(ji,jj) 
     
    427442                  zv_cor =-z1_4 * ( ff(ji-1,jj  ) * zx1 + ff(ji,jj) * zx2 ) * hvr_e(ji,jj) 
    428443                  ! after velocities with implicit bottom friction 
     444 
     445                  IF(ln_bfrimp) THEN 
     446                  !   A new method to implement the implicit bottom friction.  
     447                  !   H. Liu 
     448                  !   May 2010 
     449                  ua_e(ji,jj) =  zub_e(ji,jj) + z2dt_e * ( zu_cor + zu_spg + zu_sld + zu_asp )  * umask(ji,jj,1)   & 
     450                     &         / ( 1._wp      - z2dt_e * bfrua(ji,jj) * hur_e(ji,jj) ) 
     451                  ua_e(ji,jj) = ( ua_e(ji,jj) + z2dt_e *   zua(ji,jj)  ) * umask(ji,jj,1)    
     452                  va_e(ji,jj) =  zvb_e(ji,jj) + z2dt_e * ( zv_cor + zv_spg + zv_sld + zv_asp )  * vmask(ji,jj,1)   & 
     453                     &         / ( 1._wp      - z2dt_e * bfrva(ji,jj) * hvr_e(ji,jj) ) 
     454                  va_e(ji,jj) = ( va_e(ji,jj) + z2dt_e *   zva(ji,jj)  ) * vmask(ji,jj,1)    
     455                  ELSE 
    429456                  ua_e(ji,jj) = ( zub_e(ji,jj) + z2dt_e * ( zu_cor + zu_spg + zu_sld + zu_asp + zua(ji,jj) ) ) * umask(ji,jj,1)   & 
    430457                     &         / ( 1.e0         - z2dt_e * bfrua(ji,jj) * hur_e(ji,jj) ) 
    431458                  va_e(ji,jj) = ( zvb_e(ji,jj) + z2dt_e * ( zv_cor + zv_spg + zv_sld + zv_asp + zva(ji,jj) ) ) * vmask(ji,jj,1)   & 
    432459                     &         / ( 1.e0         - z2dt_e * bfrva(ji,jj) * hvr_e(ji,jj) ) 
     460                  ENDIF 
    433461               END DO 
    434462            END DO 
     
    438466               DO ji = fs_2, fs_jpim1   ! vector opt. 
    439467                   ! surface pressure gradient 
    440                   IF( lk_vvl) THEN 
    441                      zu_spg = -grav * (  ( rhd(ji+1,jj ,1) + 1 ) * sshn_e(ji+1,jj  )    & 
    442                         &              - ( rhd(ji  ,jj ,1) + 1 ) * sshn_e(ji  ,jj  )  ) / e1u(ji,jj) 
    443                      zv_spg = -grav * (  ( rhd(ji ,jj+1,1) + 1 ) * sshn_e(ji  ,jj+1)    & 
    444                         &              - ( rhd(ji ,jj  ,1) + 1 ) * sshn_e(ji  ,jj  )  ) / e2v(ji,jj) 
    445                   ELSE 
    446                      zu_spg = -grav * ( sshn_e(ji+1,jj) - sshn_e(ji,jj) ) / e1u(ji,jj) 
    447                      zv_spg = -grav * ( sshn_e(ji,jj+1) - sshn_e(ji,jj) ) / e2v(ji,jj) 
    448                   ENDIF 
     468                  zu_spg = -grav * ( sshn_e(ji+1,jj) - sshn_e(ji,jj) ) / e1u(ji,jj) 
     469                  zv_spg = -grav * ( sshn_e(ji,jj+1) - sshn_e(ji,jj) ) / e2v(ji,jj) 
    449470                  ! enstrophy conserving formulation for planetary vorticity term 
    450471                  zy1 =   z1_8 * ( zwy(ji  ,jj-1) + zwy(ji+1,jj-1) + zwy(ji,jj) + zwy(ji+1,jj  ) ) / e1u(ji,jj) 
     
    453474                  zv_cor  = zx1 * ( ff(ji-1,jj  ) + ff(ji,jj) ) * hvr_e(ji,jj) 
    454475                  ! after velocities with implicit bottom friction 
     476                  IF(ln_bfrimp) THEN 
     477                  !   A new method to implement the implicit bottom friction.  
     478                  !   H. Liu 
     479                  !   May 2010 
     480                  ua_e(ji,jj) =  zub_e(ji,jj) + z2dt_e * ( zu_cor + zu_spg + zu_sld + zu_asp )  * umask(ji,jj,1)   & 
     481                     &         / ( 1._wp      - z2dt_e * bfrua(ji,jj) * hur_e(ji,jj) ) 
     482                  ua_e(ji,jj) = ( ua_e(ji,jj) + z2dt_e *   zua(ji,jj)  ) * umask(ji,jj,1)    
     483                  va_e(ji,jj) =  zvb_e(ji,jj) + z2dt_e * ( zv_cor + zv_spg + zv_sld + zv_asp )  * vmask(ji,jj,1)   & 
     484                     &         / ( 1._wp      - z2dt_e * bfrva(ji,jj) * hvr_e(ji,jj) ) 
     485                  va_e(ji,jj) = ( va_e(ji,jj) + z2dt_e *   zva(ji,jj)  ) * vmask(ji,jj,1)    
     486                  ELSE 
    455487                  ua_e(ji,jj) = ( zub_e(ji,jj) + z2dt_e * ( zu_cor + zu_spg + zu_sld + zu_asp + zua(ji,jj) ) ) * umask(ji,jj,1)   & 
    456488                     &         / ( 1.e0         - z2dt_e * bfrua(ji,jj) * hur_e(ji,jj) ) 
    457489                  va_e(ji,jj) = ( zvb_e(ji,jj) + z2dt_e * ( zv_cor + zv_spg + zv_sld + zv_asp + zva(ji,jj) ) ) * vmask(ji,jj,1)   & 
    458490                     &         / ( 1.e0         - z2dt_e * bfrva(ji,jj) * hvr_e(ji,jj) ) 
     491                  ENDIF 
    459492               END DO 
    460493            END DO 
     
    464497               DO ji = fs_2, fs_jpim1   ! vector opt. 
    465498                  ! surface pressure gradient 
    466                   IF( lk_vvl) THEN 
    467                      zu_spg = -grav * (  ( rhd(ji+1,jj  ,1) + 1 ) * sshn_e(ji+1,jj  )    & 
    468                         &              - ( rhd(ji  ,jj  ,1) + 1 ) * sshn_e(ji  ,jj  )  ) / e1u(ji,jj) 
    469                      zv_spg = -grav * (  ( rhd(ji  ,jj+1,1) + 1 ) * sshn_e(ji  ,jj+1)    & 
    470                         &              - ( rhd(ji  ,jj  ,1) + 1 ) * sshn_e(ji  ,jj  )  ) / e2v(ji,jj) 
    471                   ELSE 
    472                      zu_spg = -grav * ( sshn_e(ji+1,jj) - sshn_e(ji,jj) ) / e1u(ji,jj) 
    473                      zv_spg = -grav * ( sshn_e(ji,jj+1) - sshn_e(ji,jj) ) / e2v(ji,jj) 
    474                   ENDIF 
     499                  zu_spg = -grav * ( sshn_e(ji+1,jj) - sshn_e(ji,jj) ) / e1u(ji,jj) 
     500                  zv_spg = -grav * ( sshn_e(ji,jj+1) - sshn_e(ji,jj) ) / e2v(ji,jj) 
    475501                  ! energy/enstrophy conserving formulation for planetary vorticity term 
    476502                  zu_cor = + z1_4 / e1u(ji,jj) * (  ftne(ji,jj  ) * zwy(ji  ,jj  ) + ftnw(ji+1,jj) * zwy(ji+1,jj  )   & 
     
    479505                     &                           + ftnw(ji,jj  ) * zwx(ji-1,jj  ) + ftne(ji,jj  ) * zwx(ji  ,jj  ) ) * hvr_e(ji,jj) 
    480506                  ! after velocities with implicit bottom friction 
     507                  IF(ln_bfrimp) THEN 
     508                  !   A new method to implement the implicit bottom friction.  
     509                  !   H. Liu 
     510                  !   May 2010 
     511                  ua_e(ji,jj) =  zub_e(ji,jj) + z2dt_e * ( zu_cor + zu_spg + zu_sld + zu_asp )  * umask(ji,jj,1)   & 
     512                     &         / ( 1._wp      - z2dt_e * bfrua(ji,jj) * hur_e(ji,jj) ) 
     513                  ua_e(ji,jj) = ( ua_e(ji,jj) + z2dt_e *   zua(ji,jj)  ) * umask(ji,jj,1)    
     514                  va_e(ji,jj) =  zvb_e(ji,jj) + z2dt_e * ( zv_cor + zv_spg + zv_sld + zv_asp )  * vmask(ji,jj,1)   & 
     515                     &         / ( 1._wp      - z2dt_e * bfrva(ji,jj) * hvr_e(ji,jj) ) 
     516                  va_e(ji,jj) = ( va_e(ji,jj) + z2dt_e *   zva(ji,jj)  ) * vmask(ji,jj,1)    
     517                  ELSE 
    481518                  ua_e(ji,jj) = ( zub_e(ji,jj) + z2dt_e * ( zu_cor + zu_spg + zu_sld + zu_asp + zua(ji,jj) ) ) * umask(ji,jj,1)   & 
    482519                     &         / ( 1.e0         - z2dt_e * bfrua(ji,jj) * hur_e(ji,jj) ) 
    483520                  va_e(ji,jj) = ( zvb_e(ji,jj) + z2dt_e * ( zv_cor + zv_spg + zv_sld + zv_asp + zva(ji,jj) ) ) * vmask(ji,jj,1)   & 
    484521                     &         / ( 1.e0         - z2dt_e * bfrva(ji,jj) * hvr_e(ji,jj) ) 
     522                  ENDIF 
    485523               END DO 
    486524            END DO 
Note: See TracChangeset for help on using the changeset viewer.