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 15173 – NEMO

Changeset 15173


Ignore:
Timestamp:
2021-08-04T16:43:12+02:00 (3 years ago)
Author:
deazer
Message:

Bug fixes for WAD and Zenv

Location:
NEMO/branches/NERC/NEMO_4.0.4_CO9_package_tides/src/OCE
Files:
4 edited

Legend:

Unmodified
Added
Removed
  • NEMO/branches/NERC/NEMO_4.0.4_CO9_package_tides/src/OCE/DOM/dom_oce.F90

    r14075 r15173  
    152152   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) ::           r1_hv_b , r1_hv_n , r1_hv_a   !: inverse of v-depth [1/m] 
    153153 
     154!CEOD Scale of water column down to shallowest of neighbourinbg points over total 
     155!water depth 
     156   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) :: scaled_e3t_0_ik   , scaled_e3t_0_jk 
     157   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) :: scaled_e3t_0_ip1k , scaled_e3t_0_jp1k 
     158!CEOD 
    154159   INTEGER, PUBLIC ::   nla10              !: deepest    W level Above  ~10m (nlb10 - 1) 
    155160   INTEGER, PUBLIC ::   nlb10              !: shallowest W level Bellow ~10m (nla10 + 1)  
     
    294299         ! 
    295300      ALLOCATE( wmask(jpi,jpj,jpk) , wumask(jpi,jpj,jpk), wvmask(jpi,jpj,jpk) , STAT=ierr(12) ) 
     301!CEOD  
     302      ALLOCATE(  scaled_e3t_0_ik  (jpi,jpj) , scaled_e3t_0_jk  (jpi,jpj) ,      &  
     303         &       scaled_e3t_0_ip1k(jpi,jpj) , scaled_e3t_0_jp1k(jpi,jpj) ,  STAT=ierr(13) ) 
     304         ! 
    296305      ! 
    297306      dom_oce_alloc = MAXVAL(ierr) 
  • NEMO/branches/NERC/NEMO_4.0.4_CO9_package_tides/src/OCE/DOM/domain.F90

    r14078 r15173  
    7979      INTEGER , DIMENSION(jpi,jpj) ::   ik_top , ik_bot       ! top and bottom ocean level 
    8080      REAL(wp), DIMENSION(jpi,jpj) ::   z1_hu_0, z1_hv_0 
     81!CEOD 
     82      REAL(wp), DIMENSION(jpi,jpj) ::  k_bot_i_min, k_bot_j_min 
     83      REAL(wp), DIMENSION(jpi,jpj) :: sum_e3t_0_min_ik 
     84      REAL(wp), DIMENSION(jpi,jpj) :: sum_e3t_0_min_jk 
     85      REAL(wp), DIMENSION(jpi,jpj) :: sum_e3t_0_min_ip1k 
     86      REAL(wp), DIMENSION(jpi,jpj) :: sum_e3t_0_min_jp1k 
    8187      !!---------------------------------------------------------------------- 
    8288      ! 
     
    155161         hv_0(:,:) = hv_0(:,:) + e3v_0(:,:,jk) * vmask(:,:,jk) 
    156162      END DO 
     163!CEOD Here we need to work out which is the shallowest point in terms of k 
     164! levels at neighbouring points in i and j directions. 
     165! Then we need to work out the proportion of the water column down to ht_0 that 
     166! that level is for both i points. 
     167! e.g. point i could say go down to level 10 but point i+1 only level 5 
     168! then work out depth of bottom of level 5 at point i over depth at point i to 
     169! level 10 
     170! This is needed later to work out where a u point depth should be when using 
     171! enveloping coordinates, as it will be 1/2 depth at t point at level 5 for i,i+1 
     172! points Because level 5 at i can change in time, thus u bed can change in time! 
     173! What about under ice shelves? Need rto return to that later 
     174 
     175! 1) get the shallowest k level at neighbouring i and j points store in ! k_bot_i_min 
     176 
     177         DO jj = 1, jpjm1 
     178            DO ji = 1, jpim1                ! SPG with the application of W/D gravity filters 
     179                 k_bot_i_min(ji,jj) = MIN( ik_bot(ji,jj), ik_bot(ji+1,jj  )) 
     180                 k_bot_j_min(ji,jj) = MIN( ik_bot(ji,jj), ik_bot(ji,  jj+1)) 
     181            ENDDO 
     182        ENDDO 
     183! 2) Work out the depth down to the shallowest k level both at point i and i+1 
     184! (likewise for j)  
     185 
     186        sum_e3t_0_min_ik(:,:) = 0 
     187        sum_e3t_0_min_jk(:,:) = 0 
     188        sum_e3t_0_min_ip1k(:,:) = 0 
     189        sum_e3t_0_min_jp1k(:,:) = 0 
     190 
     191        !Get sum of e3t_0s down to local min 
     192        DO jj = 1, jpjm1 
     193           DO ji = 1, jpim1                
     194             DO jk = 1, k_bot_i_min(ji,jj) 
     195  
     196                sum_e3t_0_min_ik  (ji,jj) =  sum_e3t_0_min_ik  (ji,jj) + e3t_0(ji  ,jj,jk)*tmask(ji,jj,jk) 
     197                sum_e3t_0_min_ip1k(ji,jj) =  sum_e3t_0_min_ip1k(ji,jj) + e3t_0(ji+1,jj,jk)*tmask(ji+1,jj,jk) 
     198             ENDDO 
     199            
     200             DO jk = 1, k_bot_j_min(ji,jj) 
     201                sum_e3t_0_min_jk  (ji,jj) =  sum_e3t_0_min_jk  (ji,jj) + e3t_0(ji,jj  ,jk)*tmask(ji,jj,jk) 
     202                sum_e3t_0_min_jp1k(ji,jj) =  sum_e3t_0_min_jp1k(ji,jj) + e3t_0(ji,jj+1,jk)*tmask(ji,jj+1,jk) 
     203             ENDDO 
     204! 3)  Then work out what fraction of the at rest water column that is, we later 
     205! multiply the now water depth by this scale to work out the bottom of the kth 
     206! level at time now in dynspg_ts  
     207             scaled_e3t_0_ik  (ji,jj) = sum_e3t_0_min_ik  (ji,jj)/( ht_0(ji  ,jj  ) + 1._wp - ssmask(ji,jj)) 
     208             scaled_e3t_0_jk  (ji,jj) = sum_e3t_0_min_jk  (ji,jj)/( ht_0(ji  ,jj  ) + 1._wp - ssmask(ji,jj)) 
     209             scaled_e3t_0_ip1k(ji,jj) = sum_e3t_0_min_ip1k(ji,jj)/( ht_0(ji+1,jj  ) + 1._wp - ssmask(ji+1,jj)) 
     210             scaled_e3t_0_jp1k(ji,jj) = sum_e3t_0_min_jp1k(ji,jj)/( ht_0(ji  ,jj+1) + 1._wp - ssmask(ji,jj+1)) 
     211          ENDDO 
     212       ENDDO 
    157213      ! 
    158214      !           !==  time varying part of coordinate system  ==! 
  • NEMO/branches/NERC/NEMO_4.0.4_CO9_package_tides/src/OCE/DOM/dtatsd.F90

    r14075 r15173  
    179179      ptsd(:,:,:,jp_sal) = sf_tsd(jp_sal)%fnow(:,:,:)  
    180180      ! 
    181       IF( ln_sco ) THEN                   !==   s- or mixed s-zps-coordinate   ==! 
    182          ! 
    183          IF( kt == nit000 .AND. lwp )THEN 
    184             WRITE(numout,*) 
    185             WRITE(numout,*) 'dta_tsd: interpolates T & S data onto the s- or mixed s-z-coordinate mesh' 
    186          ENDIF 
    187          ! 
    188          DO jj = 1, jpj                         ! vertical interpolation of T & S 
    189             DO ji = 1, jpi 
    190                DO jk = 1, jpk                        ! determines the intepolated T-S profiles at each (i,j) points 
    191                   zl = gdept_0(ji,jj,jk) 
    192                   IF(     zl < gdept_1d(1  ) ) THEN          ! above the first level of data 
    193                      ztp(jk) =  ptsd(ji,jj,1    ,jp_tem) 
    194                      zsp(jk) =  ptsd(ji,jj,1    ,jp_sal) 
    195                   ELSEIF( zl > gdept_1d(jpk) ) THEN          ! below the last level of data 
    196                      ztp(jk) =  ptsd(ji,jj,jpkm1,jp_tem) 
    197                      zsp(jk) =  ptsd(ji,jj,jpkm1,jp_sal) 
    198                   ELSE                                      ! inbetween : vertical interpolation between jkk & jkk+1 
    199                      DO jkk = 1, jpkm1                                  ! when  gdept(jkk) < zl < gdept(jkk+1) 
    200                         IF( (zl-gdept_1d(jkk)) * (zl-gdept_1d(jkk+1)) <= 0._wp ) THEN 
    201                            zi = ( zl - gdept_1d(jkk) ) / (gdept_1d(jkk+1)-gdept_1d(jkk)) 
    202                            ztp(jk) = ptsd(ji,jj,jkk,jp_tem) + ( ptsd(ji,jj,jkk+1,jp_tem) - ptsd(ji,jj,jkk,jp_tem) ) * zi  
    203                            zsp(jk) = ptsd(ji,jj,jkk,jp_sal) + ( ptsd(ji,jj,jkk+1,jp_sal) - ptsd(ji,jj,jkk,jp_sal) ) * zi 
    204                         ENDIF 
    205                      END DO 
    206                   ENDIF 
    207                END DO 
    208                DO jk = 1, jpkm1 
    209                   ptsd(ji,jj,jk,jp_tem) = ztp(jk) * tmask(ji,jj,jk)     ! mask required for mixed zps-s-coord 
    210                   ptsd(ji,jj,jk,jp_sal) = zsp(jk) * tmask(ji,jj,jk) 
    211                END DO 
    212                ptsd(ji,jj,jpk,jp_tem) = 0._wp 
    213                ptsd(ji,jj,jpk,jp_sal) = 0._wp 
    214             END DO 
    215          END DO 
    216          !  
    217       ELSE                                !==   z- or zps- coordinate   ==! 
     181!CEOD, We think this is incorrect we dont want to do this. 
    218182         !                              
    219183         ptsd(:,:,:,jp_tem) = ptsd(:,:,:,jp_tem) * tmask(:,:,:)    ! Mask 
     
    239203         ENDIF 
    240204         ! 
    241       ENDIF 
     205!CEOD 
    242206      ! 
    243207      IF( .NOT.ln_tsd_dmp ) THEN                   !==   deallocate T & S structure   ==!  
  • NEMO/branches/NERC/NEMO_4.0.4_CO9_package_tides/src/OCE/DYN/dynspg_ts.F90

    r14075 r15173  
    150150      REAL(wp) ::   r1_2dt_b, z1_hu, z1_hv          ! local scalars 
    151151      REAL(wp) ::   za0, za1, za2, za3              !   -      - 
     152      REAL(wp), DIMENSION(jpi,jpj) :: zdep_u, zdep_v 
    152153      REAL(wp) ::   zztmp, zldg               !   -      - 
    153154      REAL(wp) ::   zhu_bck, zhv_bck, zhdiv         !   -      - 
     
    459460            DO jj = 1, jpj 
    460461               DO ji = 1, jpim1   ! not jpi-column 
    461                   zhup2_e(ji,jj) = hu_0(ji,jj) + r1_2 * r1_e1e2u(ji,jj)                        & 
    462                        &                              * (  e1e2t(ji  ,jj) * zsshp2_e(ji  ,jj)  & 
    463                        &                                 + e1e2t(ji+1,jj) * zsshp2_e(ji+1,jj)  ) * ssumask(ji,jj) 
     462                  zhup2_e(ji,jj) = r1_2 * r1_e1e2u(ji,jj) * ( e1e2t(ji  ,jj) * zhtp2_e(ji  ,jj)*scaled_e3t_0_ik  (ji,jj) & 
     463                                                 &         +  e1e2t(ji+1,jj) * zhtp2_e(ji+1,jj)*scaled_e3t_0_ip1k(ji,jj) ) *ssumask(ji,jj) 
    464464               END DO 
    465465            END DO 
    466466            DO jj = 1, jpjm1        ! not jpj-row 
    467467               DO ji = 1, jpi 
    468                   zhvp2_e(ji,jj) = hv_0(ji,jj) + r1_2 * r1_e1e2v(ji,jj)                        & 
    469                        &                              * (  e1e2t(ji,jj  ) * zsshp2_e(ji,jj  )  & 
    470                        &                                 + e1e2t(ji,jj+1) * zsshp2_e(ji,jj+1)  ) * ssvmask(ji,jj) 
     468                  zhvp2_e(ji,jj) = r1_2 * r1_e1e2v(ji,jj) * ( e1e2t(ji,jj  ) * zhtp2_e(ji,  jj)*scaled_e3t_0_jk  (ji,jj) & 
     469                                                     &     +  e1e2t(ji,jj+1) * zhtp2_e(ji,jj+1)*scaled_e3t_0_jp1k(ji,jj) ) *ssvmask(ji,jj) 
    471470               END DO 
    472471            END DO 
     
    647646                  !                    ! hu_e, hv_e hold depth at jn,  zhup2_e, zhvp2_e hold extrapolated depth at jn+1/2 
    648647                  !                    ! backward interpolated depth used in spg terms at jn+1/2 
    649                   zhu_bck = hu_0(ji,jj) + r1_2*r1_e1e2u(ji,jj) * (  e1e2t(ji  ,jj) * zsshp2_e(ji  ,jj)    & 
    650                        &                                          + e1e2t(ji+1,jj) * zsshp2_e(ji+1,jj)  ) * ssumask(ji,jj) 
    651                   zhv_bck = hv_0(ji,jj) + r1_2*r1_e1e2v(ji,jj) * (  e1e2t(ji,jj  ) * zsshp2_e(ji,jj  )    & 
    652                        &                                          + e1e2t(ji,jj+1) * zsshp2_e(ji,jj+1)  ) * ssvmask(ji,jj) 
    653                   !                    ! inverse depth at jn+1 
    654                   z1_hu = ssumask(ji,jj) / ( hu_0(ji,jj) + zsshu_a(ji,jj) + 1._wp - ssumask(ji,jj) ) 
    655                   z1_hv = ssvmask(ji,jj) / ( hv_0(ji,jj) + zsshv_a(ji,jj) + 1._wp - ssvmask(ji,jj) ) 
    656                   ! 
     648                  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) & 
     649                           &                         + e1e2t(ji+1,jj)*(ht_0(ji+1,jj) + zsshp2_e(ji+1,jj))*scaled_e3t_0_ip1k(ji,jj) & 
     650                           &        ) * ssumask(ji,jj)  
     651 
     652                  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) & 
     653                           &                         + e1e2t(ji,jj+1)*(ht_0(ji,jj+1) + zsshp2_e(ji,jj+1))*scaled_e3t_0_jp1k(ji,jj) & 
     654                           &        ) * ssvmask(ji,jj)  
     655 
     656                  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) & 
     657                                            &      +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) ) 
     658 
     659                  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) &  
     660                                            &      + 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) ) 
     661 
     662 
    657663                  ua_e(ji,jj) = (               hu_e  (ji,jj) *   un_e (ji,jj)      &  
    658664                       &            + rdtbt * (  zhu_bck        * zu_spg (ji,jj)  &   ! 
     
    678684        
    679685         IF( .NOT.ln_linssh ) THEN   !* Update ocean depth (variable volume case only) 
    680             hu_e (2:jpim1,2:jpjm1) = hu_0(2:jpim1,2:jpjm1) + zsshu_a(2:jpim1,2:jpjm1) 
     686         DO jj = 1, jpjm1 
     687            DO ji = 1, jpim1      ! NO Vector Opt. 
     688               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) & 
     689                                                 &     +  e1e2t(ji+1,jj) * (ssha_e(ji+1,jj)+ht_0(ji+1,jj))*scaled_e3t_0_ip1k(ji,jj) )*ssumask(ji,jj) 
     690               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  ) & 
     691                                                 &     +  e1e2t(ji,jj+1) * (ssha_e(ji,jj+1)+ht_0(ji,jj+1))*scaled_e3t_0_jp1k(ji,jj) )*ssvmask(ji,jj) 
     692            ENDDO 
     693        ENDDO 
     694         CALL lbc_lnk_multi( 'dynspg_ts', zdep_u, 'U', -1._wp ) 
     695         CALL lbc_lnk_multi( 'dynspg_ts', zdep_v, 'V', -1._wp ) 
     696 
     697            hu_e (2:jpim1,2:jpjm1) = zdep_u(2:jpim1,2:jpjm1)  
    681698            hur_e(2:jpim1,2:jpjm1) = ssumask(2:jpim1,2:jpjm1) / ( hu_e(2:jpim1,2:jpjm1) + 1._wp - ssumask(2:jpim1,2:jpjm1) ) 
    682             hv_e (2:jpim1,2:jpjm1) = hv_0(2:jpim1,2:jpjm1) + zsshv_a(2:jpim1,2:jpjm1) 
     699            hv_e (2:jpim1,2:jpjm1) = zdep_v(2:jpim1,2:jpjm1) 
    683700            hvr_e(2:jpim1,2:jpjm1) = ssvmask(2:jpim1,2:jpjm1) / ( hv_e(2:jpim1,2:jpjm1) + 1._wp - ssvmask(2:jpim1,2:jpjm1) ) 
    684701            CALL lbc_lnk_multi( 'dynspg_ts', ua_e , 'U', -1._wp, va_e , 'V', -1._wp  & 
     
    778795                  &              * ( e1e2t(ji,jj  ) * ssha(ji,jj  )      & 
    779796                  &              +   e1e2t(ji,jj+1) * ssha(ji,jj+1) ) 
    780             END DO 
    781          END DO 
     797               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) & 
     798                                                 &     +  e1e2t(ji+1,jj) * (ssha(ji+1,jj)+ht_0(ji+1,jj))*scaled_e3t_0_ip1k(ji,jj) ) 
     799               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  ) & 
     800                                                 &     +  e1e2t(ji,jj+1) * (ssha(ji,jj+1)+ht_0(ji,jj+1))*scaled_e3t_0_jp1k(ji,jj) ) 
     801            END DO 
     802         END DO 
     803         CALL lbc_lnk_multi( 'dynspg_ts', zdep_u, 'U', -1._wp ) 
     804         CALL lbc_lnk_multi( 'dynspg_ts', zdep_v, 'V', -1._wp ) 
     805 
    782806         CALL lbc_lnk_multi( 'dynspg_ts', zsshu_a, 'U', 1._wp, zsshv_a, 'V', 1._wp ) ! Boundary conditions 
    783807         ! 
     
    787811         END DO 
    788812         ! Save barotropic velocities not transport: 
    789          ua_b(:,:) =  ua_b(:,:) / ( hu_0(:,:) + zsshu_a(:,:) + 1._wp - ssumask(:,:) ) 
    790          va_b(:,:) =  va_b(:,:) / ( hv_0(:,:) + zsshv_a(:,:) + 1._wp - ssvmask(:,:) ) 
     813         ua_b(:,:) =  ua_b(:,:) / ( zdep_u(:,:) + 1._wp - ssumask(:,:) ) 
     814         va_b(:,:) =  va_b(:,:) / ( zdep_v(:,:) + 1._wp - ssvmask(:,:) ) 
    791815      ENDIF 
    792816 
Note: See TracChangeset for help on using the changeset viewer.