Ignore:
Timestamp:
2020-01-16T07:56:47+01:00 (17 months ago)
Author:
deazer
Message:

Adjusted to account for flux form advection
should now work with both types of advection and enveloping bathymetry
moved one off caluculations of depth scales to domain.F90

File:
1 edited

Legend:

Unmodified
Added
Removed
  • NEMO/branches/UKMO/NEMO_4.0.1_MIRROR_WAD_ZENV/src/OCE/DYN/dynspg_ts.F90

    r12095 r12326  
    152152      REAL(wp) ::   r1_2dt_b, z1_hu, z1_hv          ! local scalars 
    153153      REAL(wp) ::   za0, za1, za2, za3              !   -      - 
    154 !CEOD 
    155       REAL(wp), DIMENSION(jpi,jpj) :: e3u_0_tot, e3v_0_tot 
    156       REAL(wp), DIMENSION(jpi,jpj) :: sum_e3t_0_min_ik 
    157       REAL(wp), DIMENSION(jpi,jpj) :: sum_e3t_0_min_jk 
    158       REAL(wp), DIMENSION(jpi,jpj) :: scaled_e3u_0_tot, scaled_e3v_0_tot 
    159       REAL(wp), DIMENSION(jpi,jpj) :: k_bot 
    160       REAL(wp), DIMENSION(jpi,jpj) :: k_bot_i_min 
    161       REAL(wp), DIMENSION(jpi,jpj) :: k_bot_j_min 
    162154      REAL(wp), DIMENSION(jpi,jpj) :: zdep_u, zdep_v 
    163 !CEOD 
    164155      REAL(wp) ::   zmdi, zztmp, zldg               !   -      - 
    165156      REAL(wp) ::   zhu_bck, zhv_bck, zhdiv         !   -      - 
     
    191182!     zwdramp = 1._wp / (rn_wdmin2 - rn_wdmin1) ! more general ramp 
    192183      !                                         ! inverse of baroclinic time step  
    193 !CEOD 
    194          e3u_0_tot(:,:)=0.0 
    195          e3v_0_tot(:,:)=0.0 
    196 !CEOD 
    197          DO jk=1,jpk 
    198             e3u_0_tot(:,:) = e3u_0_tot(:,:) + e3u_0(:,:,jk) 
    199             e3v_0_tot(:,:) = e3v_0_tot(:,:) + e3v_0(:,:,jk) 
    200          ENDDO  
    201          k_bot(:,:) = 1. 
    202          DO jk=1,jpkm1 
    203           DO jj = 1, jpj 
    204             DO ji = 1, jpi             ! SPG with the application of W/D gravity filters 
    205               IF ( tmask(ji,jj,jk) .ne. 0 ) THEN ! If its a wet point we keep going down 
    206                  k_bot(ji,jj) = jk 
    207               ENDIF 
    208             ENDDO 
    209            ENDDO 
    210          ENDDO  
    211 !work out k_bot min 
    212          DO jj = 1, jpjm1 
    213             DO ji = 1, jpim1                ! SPG with the application of W/D gravity filters 
    214                  k_bot_i_min(ji,jj) = MIN( k_bot(ji,jj), k_bot(ji+1,jj  )) 
    215                  k_bot_j_min(ji,jj) = MIN( k_bot(ji,jj), k_bot(ji,  jj+1)) 
    216             ENDDO 
    217         ENDDO 
    218         sum_e3t_0_min_ik(:,:) = 0 
    219         sum_e3t_0_min_jk(:,:) = 0 
    220  
    221         !Get sum of e3t_0s sown to local min 
    222         DO jj = 1, jpjm1 
    223            DO ji = 1, jpim1                
    224              DO jk = 1, k_bot_i_min(ji,jj) 
    225                 sum_e3t_0_min_ik(ji,jj) =  sum_e3t_0_min_ik(ji,jj) + e3t_0(ji,jj,jk) 
    226              ENDDO 
    227             
    228              DO jk = 1, k_bot_j_min(ji,jj) 
    229                 sum_e3t_0_min_jk(ji,jj) =  sum_e3t_0_min_jk(ji,jj) + e3t_0(ji,jj,jk) 
    230              ENDDO 
    231              !CEOD Now scale that against ht_0 will be one in the case it matches of course(Potential for ht_0 = zero) 
    232              scaled_e3u_0_tot(ji,jj) = sum_e3t_0_min_ik(ji,jj)/( ht_0(ji,jj) + 1._wp - ssmask(ji,jj)) 
    233              scaled_e3v_0_tot(ji,jj) = sum_e3t_0_min_jk(ji,jj)/( ht_0(ji,jj) + 1._wp - ssmask(ji,jj)) 
    234 !CEOD             WRITE(numout,*)'scaled values', ji,jj ,scaled_e3u_0_tot(ji,jj) , sum_e3t_0_min_ik(ji,jj), ht_0(ji,jj) + 1._wp - ssmask(ji,jj) 
    235           ENDDO 
    236        ENDDO 
    237 !CEO to join things up 
    238          CALL lbc_lnk_multi( 'dynspg_ts', scaled_e3u_0_tot, 'U', -1._wp ) 
    239          CALL lbc_lnk_multi( 'dynspg_ts', scaled_e3v_0_tot, 'V', -1._wp ) 
    240                
    241  
    242184      IF( kt == nit000 .AND. neuler == 0 ) THEN   ;   r1_2dt_b = 1._wp / (         rdt ) 
    243185      ELSE                                        ;   r1_2dt_b = 1._wp / ( 2._wp * rdt ) 
     
    522464            DO jj = 1, jpj 
    523465               DO ji = 1, jpim1   ! not jpi-column 
    524 !CEOD                  zhup2_e(ji,jj) = hu_0(ji,jj) + r1_2 * r1_e1e2u(ji,jj)                        & 
    525 !CEOD                       &                              * (  e1e2t(ji  ,jj) * zsshp2_e(ji  ,jj)  & 
    526 !CEOD                       &                                 + e1e2t(ji+1,jj) * zsshp2_e(ji+1,jj)  ) * ssumask(ji,jj) 
    527                   zhup2_e(ji,jj) = r1_2 * r1_e1e2u(ji,jj) * ( e1e2t(ji  ,jj) * zhtp2_e(ji  ,jj)*scaled_e3u_0_tot(ji  ,jj) & 
    528                                                  &         +  e1e2t(ji+1,jj) * zhtp2_e(ji+1,jj)*scaled_e3u_0_tot(ji+1,jj) ) *ssumask(ji,jj) 
     466                  zhup2_e(ji,jj) = r1_2 * r1_e1e2u(ji,jj) * ( e1e2t(ji  ,jj) * zhtp2_e(ji  ,jj)*scaled_e3t_0_ik  (ji,jj) & 
     467                                                 &         +  e1e2t(ji+1,jj) * zhtp2_e(ji+1,jj)*scaled_e3t_0_ip1k(ji,jj) ) *ssumask(ji,jj) 
    529468               END DO 
    530469            END DO 
    531470            DO jj = 1, jpjm1        ! not jpj-row 
    532471               DO ji = 1, jpi 
    533 !CEOD                  zhvp2_e(ji,jj) = hv_0(ji,jj) + r1_2 * r1_e1e2v(ji,jj)                        & 
    534 !CEOD                       &                              * (  e1e2t(ji,jj  ) * zsshp2_e(ji,jj  )  & 
    535 !CEOD                       &                                 + e1e2t(ji,jj+1) * zsshp2_e(ji,jj+1)  ) * ssvmask(ji,jj) 
    536                   zhvp2_e(ji,jj) = r1_2 * r1_e1e2v(ji,jj) * ( e1e2t(ji,jj  ) * zhtp2_e(ji,  jj)*scaled_e3v_0_tot(ji,jj  ) & 
    537                                                      &     +  e1e2t(ji,jj+1) * zhtp2_e(ji,jj+1)*scaled_e3v_0_tot(ji,jj+1) ) *ssvmask(ji,jj) 
     472                  zhvp2_e(ji,jj) = r1_2 * r1_e1e2v(ji,jj) * ( e1e2t(ji,jj  ) * zhtp2_e(ji,  jj)*scaled_e3t_0_jk  (ji,jj) & 
     473                                                     &     +  e1e2t(ji,jj+1) * zhtp2_e(ji,jj+1)*scaled_e3t_0_jp1k(ji,jj) ) *ssvmask(ji,jj) 
    538474               END DO 
    539475            END DO 
     
    714650                  !                    ! hu_e, hv_e hold depth at jn,  zhup2_e, zhvp2_e hold extrapolated depth at jn+1/2 
    715651                  !                    ! backward interpolated depth used in spg terms at jn+1/2 
    716                   !CEODzhu_bck = hu_0(ji,jj) + r1_2*r1_e1e2u(ji,jj) * (  e1e2t(ji  ,jj) * zsshp2_e(ji  ,jj)    & 
    717                   zhu_bck = (hu_0(ji,jj)/e3u_0_tot(ji,jj))*(e3u_0_tot(ji,jj)+ r1_2*r1_e1e2u(ji,jj) * (  e1e2t(ji  ,jj) * zsshp2_e(ji  ,jj)    & 
    718                        &                                          + e1e2t(ji+1,jj) * zsshp2_e(ji+1,jj)  ) * ssumask(ji,jj) ) 
    719                   !CEODzhv_bck = hv_0(ji,jj) + r1_2*r1_e1e2v(ji,jj) * (  e1e2t(ji,jj  ) * zsshp2_e(ji,jj  )    & 
    720                   zhv_bck = (hv_0(ji,jj)/e3v_0_tot(ji,jj))*(e3v_0_tot(ji,jj) + r1_2*r1_e1e2v(ji,jj) * (  e1e2t(ji,jj  ) * zsshp2_e(ji,jj  )    & 
    721                        &                                          + e1e2t(ji,jj+1) * zsshp2_e(ji,jj+1)  ) * ssvmask(ji,jj) ) 
    722                   !                    ! inverse depth at jn+1 
    723                   !CEODz1_hu = ssumask(ji,jj) / ( hu_0(ji,jj) + zsshu_a(ji,jj) + 1._wp - ssumask(ji,jj) ) 
    724                   !CEODz1_hv = ssvmask(ji,jj) / ( hv_0(ji,jj) + zsshv_a(ji,jj) + 1._wp - ssvmask(ji,jj) ) 
    725                   z1_hu = ssumask(ji,jj) / ( (hu_0(ji,jj)/e3u_0_tot(ji,jj))*(e3u_0_tot(ji,jj) + zsshu_a(ji,jj) + 1._wp - ssumask(ji,jj)) )  
    726                   z1_hv = ssvmask(ji,jj) / ( (hv_0(ji,jj)/e3v_0_tot(ji,jj))*(e3v_0_tot(ji,jj) + zsshv_a(ji,jj) + 1._wp - ssvmask(ji,jj)) ) 
    727                   ! 
     652                  zhu_bck = r1_2 * r1_e1e2u(ji,jj) * ( e1e2t(ji  ,jj)*(ht_0(ji  ,jj) + zsshp2_e(ji  ,jj))*scaled_e3t_0_ik  (ji,jj) & 
     653                           &                         + e1e2t(ji+1,jj)*(ht_0(ji+1,jj) + zsshp2_e(ji+1,jj))*scaled_e3t_0_ip1k(ji,jj) & 
     654                           &        ) * ssumask(ji,jj)  
     655 
     656                  zhv_bck = r1_2 * r1_e1e2v(ji,jj) * ( e1e2t(ji,jj  )*(ht_0(ji,jj  ) + zsshp2_e(ji,jj  ))*scaled_e3t_0_jk  (ji,jj) & 
     657                           &                         + e1e2t(ji,jj+1)*(ht_0(ji,jj+1) + zsshp2_e(ji,jj+1))*scaled_e3t_0_jp1k(ji,jj) & 
     658                           &        ) * ssvmask(ji,jj)  
     659 
     660                  z1_hu = ssumask(ji,jj) / ( r1_2 * r1_e1e2u(ji,jj) *(e1e2t(ji  ,jj)*(ht_0(ji,jj)   + ssha_e(ji,jj)  )*scaled_e3t_0_ik(ji  ,jj) & 
     661                                            &      +e1e2t(ji+1,jj)*(ht_0(ji+1,jj) + ssha_e(ji+1,jj))*scaled_e3t_0_ip1k(ji,jj)) + 1._wp - ssumask(ji,jj) ) 
     662 
     663                  z1_hv = ssvmask(ji,jj) / ( r1_2 * r1_e1e2v(ji,jj)* ( e1e2t(ji,jj  )*(ht_0(ji,jj)   + ssha_e(ji,jj)  )*scaled_e3t_0_jk(ji  ,jj) &  
     664                                            &      + e1e2t(ji,jj+1)* (ht_0(ji,jj+1) + ssha_e(ji,jj+1))*scaled_e3t_0_jp1k(ji,jj)) + 1._wp - ssvmask(ji,jj) ) 
     665 
     666 
    728667                  ua_e(ji,jj) = (               hu_e  (ji,jj) *   un_e (ji,jj)      &  
    729668                       &            + rdtbt * (  zhu_bck        * zu_spg (ji,jj)  &   ! 
     
    749688        
    750689         IF( .NOT.ln_linssh ) THEN   !* Update ocean depth (variable volume case only) 
    751 !CEOD TBD 
    752690         DO jj = 1, jpjm1 
    753691            DO ji = 1, jpim1      ! NO Vector Opt. 
    754                zdep_u(ji,jj) = r1_2 * r1_e1e2u(ji,jj) * ( e1e2t(ji,jj  ) * (ssha_e(ji,  jj)+ht_0(ji,  jj))*scaled_e3u_0_tot(ji  ,jj) & 
    755                                                  &     +  e1e2t(ji+1,jj) * (ssha_e(ji+1,jj)+ht_0(ji+1,jj))*scaled_e3u_0_tot(ji+1,jj) )*ssumask(ji,jj) 
    756                zdep_v(ji,jj) = r1_2 * r1_e1e2v(ji,jj) * ( e1e2t(ji,jj  ) * (ssha_e(ji,  jj)+ht_0(ji,  jj))*scaled_e3v_0_tot(ji,jj  ) & 
    757                                                  &     +  e1e2t(ji,jj+1) * (ssha_e(ji,jj+1)+ht_0(ji,jj+1))*scaled_e3v_0_tot(ji,jj+1) )*ssvmask(ji,jj) 
    758  
    759 !CEOD             WRITE(numout,*)'zdep_u', ji,jj ,zdep_u(ji,jj) ,scaled_e3u_0_tot(ji,jj),ssha_e(ji,jj),jn 
     692               zdep_u(ji,jj) = r1_2 * r1_e1e2u(ji,jj) * ( e1e2t(ji,jj  ) * (ssha_e(ji,  jj)+ht_0(ji,  jj))*scaled_e3t_0_ik  (ji,jj) & 
     693                                                 &     +  e1e2t(ji+1,jj) * (ssha_e(ji+1,jj)+ht_0(ji+1,jj))*scaled_e3t_0_ip1k(ji,jj) )*ssumask(ji,jj) 
     694               zdep_v(ji,jj) = r1_2 * r1_e1e2v(ji,jj) * ( e1e2t(ji,jj  ) * (ssha_e(ji,  jj)+ht_0(ji,  jj))*scaled_e3t_0_jk  (ji,jj  ) & 
     695                                                 &     +  e1e2t(ji,jj+1) * (ssha_e(ji,jj+1)+ht_0(ji,jj+1))*scaled_e3t_0_jp1k(ji,jj) )*ssvmask(ji,jj) 
    760696            ENDDO 
    761697        ENDDO 
     
    763699         CALL lbc_lnk_multi( 'dynspg_ts', zdep_v, 'V', -1._wp ) 
    764700 
    765             !CEODhu_e (2:jpim1,2:jpjm1) = hu_0(2:jpim1,2:jpjm1) + zsshu_a(2:jpim1,2:jpjm1) 
    766701            hu_e (2:jpim1,2:jpjm1) = zdep_u(2:jpim1,2:jpjm1)  
    767702            hur_e(2:jpim1,2:jpjm1) = ssumask(2:jpim1,2:jpjm1) / ( hu_e(2:jpim1,2:jpjm1) + 1._wp - ssumask(2:jpim1,2:jpjm1) ) 
    768             !CEODhv_e (2:jpim1,2:jpjm1) = hv_0(2:jpim1,2:jpjm1) + zsshv_a(2:jpim1,2:jpjm1) 
    769703            hv_e (2:jpim1,2:jpjm1) = zdep_v(2:jpim1,2:jpjm1) 
    770704            hvr_e(2:jpim1,2:jpjm1) = ssvmask(2:jpim1,2:jpjm1) / ( hv_e(2:jpim1,2:jpjm1) + 1._wp - ssvmask(2:jpim1,2:jpjm1) ) 
     
    865799                  &              * ( e1e2t(ji,jj  ) * ssha(ji,jj  )      & 
    866800                  &              +   e1e2t(ji,jj+1) * ssha(ji,jj+1) ) 
    867                zdep_u(ji,jj) = r1_2 * r1_e1e2u(ji,jj) * ( e1e2t(ji,jj  ) * (ssha(ji,  jj)+ht_0(ji,  jj))*scaled_e3u_0_tot(ji  ,jj) & 
    868                                                  &     +  e1e2t(ji+1,jj) * (ssha(ji+1,jj)+ht_0(ji+1,jj))*scaled_e3u_0_tot(ji+1,jj) ) 
    869                zdep_v(ji,jj) = r1_2 * r1_e1e2v(ji,jj) * ( e1e2t(ji,jj  ) * (ssha(ji,  jj)+ht_0(ji,  jj))*scaled_e3v_0_tot(ji,jj  ) & 
    870                                                  &     +  e1e2t(ji,jj+1) * (ssha(ji,jj+1)+ht_0(ji,jj+1))*scaled_e3v_0_tot(ji,jj+1) ) 
     801               zdep_u(ji,jj) = r1_2 * r1_e1e2u(ji,jj) * ( e1e2t(ji,jj  ) * (ssha(ji,  jj)+ht_0(ji,  jj))*scaled_e3t_0_ik(ji  ,jj) & 
     802                                                 &     +  e1e2t(ji+1,jj) * (ssha(ji+1,jj)+ht_0(ji+1,jj))*scaled_e3t_0_ip1k(ji,jj) ) 
     803               zdep_v(ji,jj) = r1_2 * r1_e1e2v(ji,jj) * ( e1e2t(ji,jj  ) * (ssha(ji,  jj)+ht_0(ji,  jj))*scaled_e3t_0_jk(ji,jj  ) & 
     804                                                 &     +  e1e2t(ji,jj+1) * (ssha(ji,jj+1)+ht_0(ji,jj+1))*scaled_e3t_0_jp1k(ji,jj) ) 
    871805            END DO 
    872806         END DO 
     
    881815         END DO 
    882816         ! Save barotropic velocities not transport: 
    883 !CEOD         ua_b(:,:) =  ua_b(:,:) / ( hu_0(:,:) + zsshu_a(:,:) + 1._wp - ssumask(:,:) ) 
    884 !CEOD         va_b(:,:) =  va_b(:,:) / ( hv_0(:,:) + zsshv_a(:,:) + 1._wp - ssvmask(:,:) ) 
    885817         ua_b(:,:) =  ua_b(:,:) / ( zdep_u(:,:) + 1._wp - ssumask(:,:) ) 
    886818         va_b(:,:) =  va_b(:,:) / ( zdep_v(:,:) + 1._wp - ssvmask(:,:) ) 
    887       ENDIF 
    888  
     819      ENDIF  
    889820 
    890821      ! Correct velocities so that the barotropic velocity equals (un_adv, vn_adv) (in all cases)   
Note: See TracChangeset for help on using the changeset viewer.