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 3294 for trunk/NEMOGCM/NEMO/OPA_SRC/DYN/dynnxt.F90 – NEMO

Ignore:
Timestamp:
2012-01-28T17:44:18+01:00 (12 years ago)
Author:
rblod
Message:

Merge of 3.4beta into the trunk

File:
1 edited

Legend:

Unmodified
Added
Removed
  • trunk/NEMOGCM/NEMO/OPA_SRC/DYN/dynnxt.F90

    r2779 r3294  
    3333   USE obcdyn_bt       ! 2D open boundary condition for momentum (obc_dyn_bt routine) 
    3434   USE obcvol          ! ocean open boundary condition (obc_vol routines) 
    35    USE bdy_oce         ! unstructured open boundary conditions 
    36    USE bdydta          ! unstructured open boundary conditions 
    37    USE bdydyn          ! unstructured open boundary conditions 
     35   USE bdy_oce         ! ocean open boundary conditions 
     36   USE bdydta          ! ocean open boundary conditions 
     37   USE bdydyn          ! ocean open boundary conditions 
     38   USE bdyvol          ! ocean open boundary condition (bdy_vol routines) 
    3839   USE in_out_manager  ! I/O manager 
    3940   USE lbclnk          ! lateral boundary condition (or mpp link) 
    4041   USE lib_mpp         ! MPP library 
     42   USE wrk_nemo        ! Memory Allocation 
    4143   USE prtctl          ! Print control 
    4244#if defined key_agrif 
    4345   USE agrif_opa_interp 
    4446#endif 
     47   USE timing          ! Timing 
    4548 
    4649   IMPLICIT NONE 
     
    7780      !!              * Apply lateral boundary conditions on after velocity  
    7881      !!             at the local domain boundaries through lbc_lnk call, 
    79       !!             at the radiative open boundaries (lk_obc=T), 
    80       !!             at the relaxed   open boundaries (lk_bdy=T), and 
     82      !!             at the one-way open boundaries (lk_obc=T), 
    8183      !!             at the AGRIF zoom     boundaries (lk_agrif=T) 
    8284      !! 
     
    9294      !!               un,vn   now horizontal velocity of next time-step 
    9395      !!---------------------------------------------------------------------- 
    94       USE wrk_nemo, ONLY:   wrk_in_use, wrk_not_released 
    95       USE oce     , ONLY:   ze3u_f => ta       , ze3v_f => sa       ! (ta,sa) used as 3D workspace 
    96       USE wrk_nemo, ONLY:   zs_t   => wrk_2d_1 , zs_u_1 => wrk_2d_2 , zs_v_1 => wrk_2d_3 
    97       ! 
    9896      INTEGER, INTENT( in ) ::   kt      ! ocean time-step index 
    9997      ! 
    10098      INTEGER  ::   ji, jj, jk   ! dummy loop indices 
     99      INTEGER  ::   iku, ikv     ! local integers 
    101100#if ! defined key_dynspg_flt 
    102101      REAL(wp) ::   z2dt         ! temporary scalar 
    103102#endif 
    104       REAL(wp) ::   zue3a, zue3n, zue3b, zuf    ! local scalars 
    105       REAL(wp) ::   zve3a, zve3n, zve3b, zvf    !   -      - 
    106       REAL(wp) ::   zec, zv_t_ij, zv_t_ip1j, zv_t_ijp1 
     103      REAL(wp) ::   zue3a, zue3n, zue3b, zuf, zec   ! local scalars 
     104      REAL(wp) ::   zve3a, zve3n, zve3b, zvf        !   -      - 
     105      REAL(wp), POINTER, DIMENSION(:,:,:) ::  ze3u_f, ze3v_f  
    107106      !!---------------------------------------------------------------------- 
    108  
    109       IF( wrk_in_use(2, 1,2,3) ) THEN 
    110          CALL ctl_stop('dyn_nxt: requested workspace arrays unavailable')   ;   RETURN 
    111       ENDIF 
    112  
     107      ! 
     108      IF( nn_timing == 1 )  CALL timing_start('dyn_nxt') 
     109      ! 
     110      CALL wrk_alloc( jpi,jpj,jpk, ze3u_f, ze3v_f ) 
     111      ! 
    113112      IF( kt == nit000 ) THEN 
    114113         IF(lwp) WRITE(numout,*) 
     
    174173      ENDIF 
    175174      ! 
    176 # elif defined key_bdy  
     175# elif defined key_bdy 
    177176      !                                !* BDY open boundaries 
    178       IF( .NOT. lk_dynspg_flt ) THEN 
    179          CALL bdy_dyn_frs( kt ) 
    180 #  if ! defined key_vvl 
    181          ua_e(:,:) = 0.e0 
    182          va_e(:,:) = 0.e0 
    183          ! Set these variables for use in bdy_dyn_fla 
    184          hur_e(:,:) = hur(:,:) 
    185          hvr_e(:,:) = hvr(:,:) 
    186          DO jk = 1, jpkm1   !! Vertically integrated momentum trends 
    187             ua_e(:,:) = ua_e(:,:) + fse3u(:,:,jk) * umask(:,:,jk) * ua(:,:,jk) 
    188             va_e(:,:) = va_e(:,:) + fse3v(:,:,jk) * vmask(:,:,jk) * va(:,:,jk) 
    189          END DO 
    190          ua_e(:,:) = ua_e(:,:) * hur(:,:) 
    191          va_e(:,:) = va_e(:,:) * hvr(:,:) 
    192          DO jk = 1 , jpkm1 
    193             ua(:,:,jk) = ua(:,:,jk) - ua_e(:,:) 
    194             va(:,:,jk) = va(:,:,jk) - va_e(:,:) 
    195          END DO 
    196          CALL bdy_dta_fla( kt+1, 0,2*nn_baro) 
    197          CALL bdy_dyn_fla( sshn_b ) 
    198          CALL lbc_lnk( ua_e, 'U', -1. )   ! Boundary points should be updated 
    199          CALL lbc_lnk( va_e, 'V', -1. )   ! 
    200          DO jk = 1 , jpkm1 
    201             ua(:,:,jk) = ( ua(:,:,jk) + ua_e(:,:) ) * umask(:,:,jk) 
    202             va(:,:,jk) = ( va(:,:,jk) + va_e(:,:) ) * vmask(:,:,jk) 
    203          END DO 
    204 #  endif 
    205       ENDIF 
     177      IF( lk_dynspg_exp ) CALL bdy_dyn( kt ) 
     178      IF( lk_dynspg_ts )  CALL bdy_dyn( kt, dyn3d_only=.true. ) 
     179 
     180!!$   Do we need a call to bdy_vol here?? 
     181      ! 
    206182# endif 
    207183      ! 
     
    238214         ELSE                             ! Variable volume ! 
    239215            !                             ! ================! 
    240             ! Before scale factor at t-points 
    241             ! ------------------------------- 
    242             DO jk = 1, jpkm1 
     216            ! 
     217            DO jk = 1, jpkm1                 ! Before scale factor at t-points 
    243218               fse3t_b(:,:,jk) = fse3t_n(:,:,jk)                                   & 
    244219                  &              + atfp * (  fse3t_b(:,:,jk) + fse3t_a(:,:,jk)     & 
    245                   &                         - 2.e0 * fse3t_n(:,:,jk)            ) 
    246             ENDDO 
    247             ! Add volume filter correction only at the first level of t-point scale factors 
    248             zec = atfp * rdt / rau0 
     220                  &                         - 2._wp * fse3t_n(:,:,jk)            ) 
     221            END DO 
     222            zec = atfp * rdt / rau0          ! Add filter correction only at the 1st level of t-point scale factors 
    249223            fse3t_b(:,:,1) = fse3t_b(:,:,1) - zec * ( emp_b(:,:) - emp(:,:) ) * tmask(:,:,1) 
    250             ! surface at t-points and inverse surface at (u/v)-points used in surface averaging computations 
    251             zs_t  (:,:) =       e1t(:,:) * e2t(:,:) 
    252             zs_u_1(:,:) = 0.5 / ( e1u(:,:) * e2u(:,:) ) 
    253             zs_v_1(:,:) = 0.5 / ( e1v(:,:) * e2v(:,:) ) 
    254224            ! 
    255             IF( ln_dynadv_vec ) THEN 
    256                ! Before scale factor at (u/v)-points 
    257                ! ----------------------------------- 
    258                ! Scale factor anomaly at (u/v)-points: surface averaging of scale factor at t-points 
    259                DO jk = 1, jpkm1 
    260                   DO jj = 1, jpjm1 
    261                      DO ji = 1, jpim1 
    262                         zv_t_ij           = zs_t(ji  ,jj  ) * fse3t_b(ji  ,jj  ,jk) 
    263                         zv_t_ip1j         = zs_t(ji+1,jj  ) * fse3t_b(ji+1,jj  ,jk) 
    264                         zv_t_ijp1         = zs_t(ji  ,jj+1) * fse3t_b(ji  ,jj+1,jk) 
    265                         fse3u_b(ji,jj,jk) = umask(ji,jj,jk) * ( zs_u_1(ji,jj) * ( zv_t_ij + zv_t_ip1j ) - fse3u_0(ji,jj,jk) ) 
    266                         fse3v_b(ji,jj,jk) = vmask(ji,jj,jk) * ( zs_v_1(ji,jj) * ( zv_t_ij + zv_t_ijp1 ) - fse3v_0(ji,jj,jk) ) 
    267                      END DO 
    268                   END DO 
    269                END DO 
    270                ! lateral boundary conditions 
    271                CALL lbc_lnk( fse3u_b(:,:,:), 'U', 1. ) 
    272                CALL lbc_lnk( fse3v_b(:,:,:), 'V', 1. ) 
    273                ! Add initial scale factor to scale factor anomaly 
    274                fse3u_b(:,:,:) = fse3u_b(:,:,:) + fse3u_0(:,:,:) 
    275                fse3v_b(:,:,:) = fse3v_b(:,:,:) + fse3v_0(:,:,:) 
    276                ! Leap-Frog - Asselin filter and swap: applied on velocity 
    277                ! ----------------------------------- 
    278                DO jk = 1, jpkm1 
    279                   DO jj = 1, jpj 
     225            IF( ln_dynadv_vec ) THEN         ! vector invariant form (no thickness weighted calulation) 
     226               ! 
     227               !                                      ! before scale factors at u- & v-pts (computed from fse3t_b) 
     228               CALL dom_vvl_2( kt, fse3u_b(:,:,:), fse3v_b(:,:,:) ) 
     229               ! 
     230               DO jk = 1, jpkm1                       ! Leap-Frog - Asselin filter and swap: applied on velocity 
     231                  DO jj = 1, jpj                      !                                                 -------- 
    280232                     DO ji = 1, jpi 
    281233                        zuf = un(ji,jj,jk) + atfp * ( ub(ji,jj,jk) - 2.e0 * un(ji,jj,jk) + ua(ji,jj,jk) ) 
     
    290242               END DO 
    291243               ! 
    292             ELSE 
    293                ! Temporary filered scale factor at (u/v)-points (will become before scale factor) 
    294                !----------------------------------------------- 
    295                ! Scale factor anomaly at (u/v)-points: surface averaging of scale factor at t-points 
    296                DO jk = 1, jpkm1 
    297                   DO jj = 1, jpjm1 
    298                      DO ji = 1, jpim1 
    299                         zv_t_ij          = zs_t(ji  ,jj  ) * fse3t_b(ji  ,jj  ,jk) 
    300                         zv_t_ip1j        = zs_t(ji+1,jj  ) * fse3t_b(ji+1,jj  ,jk) 
    301                         zv_t_ijp1        = zs_t(ji  ,jj+1) * fse3t_b(ji  ,jj+1,jk) 
    302                         ze3u_f(ji,jj,jk) = umask(ji,jj,jk) * ( zs_u_1(ji,jj) * ( zv_t_ij + zv_t_ip1j ) - fse3u_0(ji,jj,jk) ) 
    303                         ze3v_f(ji,jj,jk) = vmask(ji,jj,jk) * ( zs_v_1(ji,jj) * ( zv_t_ij + zv_t_ijp1 ) - fse3v_0(ji,jj,jk) ) 
    304                      END DO 
    305                   END DO 
    306                END DO 
    307                ! lateral boundary conditions 
    308                CALL lbc_lnk( ze3u_f, 'U', 1. ) 
    309                CALL lbc_lnk( ze3v_f, 'V', 1. ) 
    310                ! Add initial scale factor to scale factor anomaly 
    311                ze3u_f(:,:,:) = ze3u_f(:,:,:) + fse3u_0(:,:,:) 
    312                ze3v_f(:,:,:) = ze3v_f(:,:,:) + fse3v_0(:,:,:) 
    313                ! Leap-Frog - Asselin filter and swap: applied on thickness weighted velocity 
    314                ! -----------------------------------             =========================== 
    315                DO jk = 1, jpkm1 
    316                   DO jj = 1, jpj 
    317                      DO ji = 1, jpim1 
     244            ELSE                             ! flux form (thickness weighted calulation) 
     245               ! 
     246               CALL dom_vvl_2( kt, ze3u_f, ze3v_f )   ! before scale factors at u- & v-pts (computed from fse3t_b) 
     247               ! 
     248               DO jk = 1, jpkm1                       ! Leap-Frog - Asselin filter and swap:  
     249                  DO jj = 1, jpj                      !                   applied on thickness weighted velocity 
     250                     DO ji = 1, jpim1                 !                              --------------------------- 
    318251                        zue3a = ua(ji,jj,jk) * fse3u_a(ji,jj,jk) 
    319252                        zve3a = va(ji,jj,jk) * fse3v_a(ji,jj,jk) 
     
    323256                        zve3b = vb(ji,jj,jk) * fse3v_b(ji,jj,jk) 
    324257                        ! 
    325                         zuf  = ( zue3n + atfp * ( zue3b - 2.e0 * zue3n  + zue3a ) ) / ze3u_f(ji,jj,jk) 
    326                         zvf  = ( zve3n + atfp * ( zve3b - 2.e0 * zve3n  + zve3a ) ) / ze3v_f(ji,jj,jk) 
     258                        zuf = ( zue3n + atfp * ( zue3b - 2._wp * zue3n  + zue3a ) ) / ze3u_f(ji,jj,jk) 
     259                        zvf = ( zve3n + atfp * ( zve3b - 2._wp * zve3n  + zve3a ) ) / ze3v_f(ji,jj,jk) 
    327260                        ! 
    328                         ub(ji,jj,jk) = zuf                      ! ub <-- filtered velocity 
     261                        ub(ji,jj,jk) = zuf                     ! ub <-- filtered velocity 
    329262                        vb(ji,jj,jk) = zvf 
    330                         un(ji,jj,jk) = ua(ji,jj,jk)             ! un <-- ua 
     263                        un(ji,jj,jk) = ua(ji,jj,jk)            ! un <-- ua 
    331264                        vn(ji,jj,jk) = va(ji,jj,jk) 
    332265                     END DO 
    333266                  END DO 
    334267               END DO 
    335                fse3u_b(:,:,:) = ze3u_f(:,:,:)                   ! e3u_b <-- filtered scale factor 
    336                fse3v_b(:,:,:) = ze3v_f(:,:,:) 
    337                CALL lbc_lnk( ub, 'U', -1. )                     ! lateral boundary conditions 
     268               fse3u_b(:,:,1:jpkm1) = ze3u_f(:,:,1:jpkm1)      ! e3u_b <-- filtered scale factor 
     269               fse3v_b(:,:,1:jpkm1) = ze3v_f(:,:,1:jpkm1) 
     270               CALL lbc_lnk( ub, 'U', -1. )                    ! lateral boundary conditions 
    338271               CALL lbc_lnk( vb, 'V', -1. ) 
    339272            ENDIF 
     
    346279         &                       tab3d_2=vn, clinfo2=' Vn: '       , mask2=vmask ) 
    347280      !  
    348       IF( wrk_not_released(2, 1,2,3) )   CALL ctl_stop('dyn_nxt: failed to release workspace arrays') 
     281      CALL wrk_dealloc( jpi,jpj,jpk, ze3u_f, ze3v_f ) 
     282      ! 
     283      IF( nn_timing == 1 )  CALL timing_stop('dyn_nxt') 
    349284      ! 
    350285   END SUBROUTINE dyn_nxt 
Note: See TracChangeset for help on using the changeset viewer.