Changeset 14967


Ignore:
Timestamp:
2021-06-10T09:56:11+02:00 (6 months ago)
Author:
cbricaud
Message:

introduce limit conditions under ice caps in GLS; see #2593 and #2604

Location:
NEMO/trunk/src/OCE/ZDF
Files:
2 edited

Legend:

Unmodified
Added
Removed
  • NEMO/trunk/src/OCE/ZDF/zdfgls.F90

    r14966 r14967  
    354354         END_2D 
    355355         ! 
     356         IF( ln_isfcav) THEN     ! top boundary   (ocean cavity) 
     357            DO_2D_OVR( nn_hls-1, nn_hls-1, nn_hls-1, nn_hls-1 ) 
     358               IF( mikt(ji,jj) > 1 )THEN 
     359                  itop   = mikt(ji,jj)       ! k   top w-point 
     360                  itopp1 = mikt(ji,jj) + 1   ! k+1 1st w-point below the top one 
     361                  !                                                ! mask at the 
     362                  !                                                ocean surface 
     363                  !                                                points 
     364                  z_en = MAX( rc02r * ustar2_top(ji,jj), rn_emin ) * ( 1._wp - tmask(ji,jj,1) ) 
     365                  ! 
     366                  ! Dirichlet condition applied at: 
     367                  !     top level (itop)         &      Just below it (itopp1) 
     368                  zd_lw(ji,jj,itop) = 0._wp   ;   zd_lw(ji,jj,itopp1) = 0._wp 
     369                  zd_up(ji,jj,itop) = 0._wp   ;   zd_up(ji,jj,itopp1) = 0._wp 
     370                  zdiag(ji,jj,itop) = 1._wp   ;   zdiag(ji,jj,itopp1) = 1._wp 
     371                  en   (ji,jj,itop) = z_en    ;   en   (ji,jj,itopp1) = z_en 
     372               ENDIF 
     373            END_2D 
     374         ENDIF 
    356375         ! 
    357376      CASE ( 1 )             ! Neumann boundary condition (set d(e)/dz) 
     
    377396         END_2D 
    378397         ! 
     398         IF( ln_isfcav) THEN     ! top boundary   (ocean cavity) 
     399            DO_2D_OVR( nn_hls-1, nn_hls-1, nn_hls-1, nn_hls-1 ) 
     400               IF( mikt(ji,jj) > 1 )THEN 
     401                  itop   = mikt(ji,jj)       ! k   top w-point 
     402                  itopp1 = mikt(ji,jj) + 1   ! k+1 1st w-point below the top one 
     403                  !                                                ! mask at the 
     404                  !                                                ocean surface 
     405                  !                                                points 
     406                  z_en = MAX( rc02r * ustar2_top(ji,jj), rn_emin ) * ( 1._wp - tmask(ji,jj,1) ) 
     407                  ! 
     408                  ! Bottom level Dirichlet condition: 
     409                  !     Bottom level (ibot)      &      Just above it (ibotm1) 
     410                  !         Dirichlet            !         Neumann 
     411                  zd_lw(ji,jj,itop) = 0._wp   !   ! Remove zd_up from zdiag 
     412                  zdiag(ji,jj,itop) = 1._wp   ;   zdiag(ji,jj,itopp1) = zdiag(ji,jj,itopp1) + zd_up(ji,jj,itopp1) 
     413                  zd_up(ji,jj,itop) = 0._wp   ;   zd_up(ji,jj,itopp1) = 0._wp 
     414                  en   (ji,jj,itop) = z_en 
     415               ENDIF 
     416            END_2D 
     417         ENDIF 
    379418         ! 
    380419      END SELECT 
     
    632671         END_2D 
    633672         ! 
     673         IF( ln_isfcav) THEN     ! top boundary   (ocean cavity) 
     674            DO_2D( nn_hls-1, nn_hls-1, nn_hls-1, nn_hls-1 ) 
     675               IF ( mikt(ji,jj) > 1 ) THEN 
     676                  itop   = mikt(ji,jj)       ! k   top w-point 
     677                  itopp1 = mikt(ji,jj) + 1   ! k+1 1st w-point below the top one 
     678                  ! 
     679                  zdep(ji,jj) = vkarmn * r_z0_top 
     680                  psi (ji,jj,itop) = rc0**rpp * en(ji,jj,itop)**rmm *zdep(ji,jj)**rnn 
     681                  zd_lw(ji,jj,itop) = 0._wp 
     682                  zd_up(ji,jj,itop) = 0._wp 
     683                  zdiag(ji,jj,itop) = 1._wp 
     684                  ! 
     685                  ! Just above last level, Dirichlet condition again (GOTM like) 
     686                  zdep(ji,jj) = vkarmn * ( r_z0_top + e3t(ji,jj,itopp1,Kmm) ) 
     687                  psi (ji,jj,itopp1) = rc0**rpp * en(ji,jj,itop  )**rmm *zdep(ji,jj)**rnn 
     688                  zd_lw(ji,jj,itopp1) = 0._wp 
     689                  zd_up(ji,jj,itopp1) = 0._wp 
     690                  zdiag(ji,jj,itopp1) = 1._wp 
     691               END IF 
     692            END_2D 
     693         END IF 
     694         ! 
    634695      CASE ( 1 )             ! Neumman boundary condition 
    635696         ! 
     
    656717            psi(ji,jj,ibotm1) = psi(ji,jj,ibotm1) + zflxb / e3w(ji,jj,ibotm1,Kmm) 
    657718         END_2D 
     719         ! 
     720         IF( ln_isfcav) THEN     ! top boundary   (ocean cavity) 
     721            DO_2D( nn_hls-1, nn_hls-1, nn_hls-1, nn_hls-1 ) 
     722               IF ( mikt(ji,jj) > 1 ) THEN 
     723                  itop   = mikt(ji,jj)       ! k   top w-point 
     724                  itopp1 = mikt(ji,jj) + 1   ! k+1 1st w-point below the top one 
     725                  ! 
     726                  ! Bottom level Dirichlet condition: 
     727                  zdep(ji,jj) = vkarmn * r_z0_top 
     728                  psi (ji,jj,itop) = rc0**rpp * en(ji,jj,itop)**rmm *zdep(ji,jj)**rnn 
     729                  ! 
     730                  zd_lw(ji,jj,itop) = 0._wp 
     731                  zd_up(ji,jj,itop) = 0._wp 
     732                  zdiag(ji,jj,itop) = 1._wp 
     733                  ! 
     734                  ! Just below cavity level: Neumann condition with flux 
     735                  ! injection 
     736                  zdiag(ji,jj,itopp1) = zdiag(ji,jj,itopp1) + zd_up(ji,jj,itopp1) ! Remove zd_up from zdiag 
     737                  zd_up(ji,jj,itopp1) = 0._wp 
     738                  ! 
     739                  ! Set psi vertical flux below cavity: 
     740                  zdep(ji,jj) = r_z0_top + 0.5_wp*e3t(ji,jj,itopp1,Kmm) 
     741                  zflxb = rsbc_psi2 * ( p_avm(ji,jj,itop) + p_avm(ji,jj,itopp1))   & 
     742                     &  * (0.5_wp*(en(ji,jj,itop)+en(ji,jj,itopp1)))**rmm * zdep(ji,jj)**(rnn-1._wp) 
     743                  psi(ji,jj,itopp1) = psi(ji,jj,itopp1) + zflxb / e3w(ji,jj,itopp1,Kmm) 
     744               END IF 
     745            END_2D 
     746         END IF 
     747 
    658748         ! 
    659749      END SELECT 
  • NEMO/trunk/src/OCE/ZDF/zdfphy.F90

    r14921 r14967  
    213213      IF( ioptio /= 1 )    CALL ctl_stop( 'zdf_phy_init: one and only one vertical diffusion option has to be defined ' ) 
    214214      IF( ln_isfcav ) THEN 
    215       IF( ln_zdfric .OR. ln_zdfgls )    CALL ctl_stop( 'zdf_phy_init: zdfric and zdfgls never tested with ice shelves cavities ' ) 
     215      IF( ln_zdfric )      CALL ctl_stop( 'zdf_phy_init: zdfric never tested with ice shelves cavities ' ) 
    216216      ENDIF 
    217217      !                                ! shear production term flag 
Note: See TracChangeset for help on using the changeset viewer.