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

Changeset 15627


Ignore:
Timestamp:
2022-01-04T19:30:17+01:00 (2 years ago)
Author:
techene
Message:

#2605 qco r3. ratios management optimised (RK3) and some cleanning (MLF)

Location:
NEMO/branches/2021/dev_r14318_RK3_stage1/src/OCE
Files:
2 edited

Legend:

Unmodified
Added
Removed
  • NEMO/branches/2021/dev_r14318_RK3_stage1/src/OCE/DOM/domqco.F90

    r15574 r15627  
    3939   PUBLIC  dom_qco_init       ! called by domain.F90 
    4040   PUBLIC  dom_qco_zgr        ! called by isfcpl.F90 
    41    PUBLIC  dom_qco_r3c        ! called by steplf.F90 
     41   PUBLIC  dom_qco_r3c        ! called by stpmlf.F90 
     42   PUBLIC  dom_qco_r3c_RK3    ! called by stprk3_stg.F90 
    4243 
    4344   !                                                      !!* Namelist nam_vvl 
     
    118119#if defined key_RK3 
    119120      CALL dom_qco_r3c( ssh(:,:,Kbb), r3t(:,:,Kbb), r3u(:,:,Kbb), r3v(:,:,Kbb), r3f(:,:) ) 
    120       CALL dom_qco_r3c( ssh(:,:,Kmm), r3t(:,:,Kmm), r3u(:,:,Kmm), r3v(:,:,Kmm)           ) 
     121      CALL dom_qco_r3c( ssh(:,:,Kmm), r3t(:,:,Kmm), r3u(:,:,Kmm), r3v(:,:,Kmm)           )  !! needed for agrif ??? 
    121122#else 
    122123      CALL dom_qco_r3c( ssh(:,:,Kbb), r3t(:,:,Kbb), r3u(:,:,Kbb), r3v(:,:,Kbb)           ) 
     
    156157      !                                      !==  ratio at u-,v-point  ==! 
    157158      ! 
    158 !!st      IF( ln_dynadv_vec ) THEN                     !- Vector Form   (thickness weighted averaging) 
    159 #if ! defined key_qcoTest_FluxForm 
    160       !                                ! no 'key_qcoTest_FluxForm' : surface weighted ssh average 
    161159      DO_2D_OVR( nn_hls, nn_hls-1, nn_hls, nn_hls-1 ) 
    162160         pr3u(ji,jj) = 0.5_wp * (  e1e2t(ji  ,jj) * pssh(ji  ,jj)  & 
     
    164162         pr3v(ji,jj) = 0.5_wp * (  e1e2t(ji,jj  ) * pssh(ji,jj  )  & 
    165163            &                    + e1e2t(ji,jj+1) * pssh(ji,jj+1)  ) * r1_hv_0(ji,jj) * r1_e1e2v(ji,jj) 
     164      END_2D       
     165      ! 
     166      IF( .NOT.PRESENT( pr3f ) ) THEN              !- lbc on ratio at u-, v-points only 
     167         IF (nn_hls==1) CALL lbc_lnk( 'dom_qco_r3c', pr3u, 'U', 1._wp, pr3v, 'V', 1._wp ) 
     168         ! 
     169         ! 
     170      ELSE                                   !==  ratio at f-point  ==! 
     171         ! 
     172         DO_2D_OVR( nn_hls, nn_hls-1, nn_hls, nn_hls-1 ) 
     173            ! round brackets added to fix the order of floating point operations 
     174            ! needed to ensure halo 1 - halo 2 compatibility 
     175            pr3f(ji,jj) = 0.25_wp * (   (  e1e2t(ji  ,jj  ) * pssh(ji  ,jj  )      & 
     176               &                         + e1e2t(ji+1,jj  ) * pssh(ji+1,jj  )  )   & ! bracket for halo 1 - halo 2 compatibility 
     177               &                     +  (  e1e2t(ji  ,jj+1) * pssh(ji  ,jj+1)      & 
     178               &                         + e1e2t(ji+1,jj+1) * pssh(ji+1,jj+1)  )   & ! bracket for halo 1 - halo 2 compatibility 
     179               &                    ) * r1_hf_0(ji,jj) * r1_e1e2f(ji,jj) 
     180         END_2D 
     181         !                                                 ! lbc on ratio at u-,v-,f-points 
     182         IF (nn_hls==1) CALL lbc_lnk( 'dom_qco_r3c', pr3u, 'U', 1._wp, pr3v, 'V', 1._wp, pr3f, 'F', 1._wp ) 
     183         ! 
     184      ENDIF 
     185      ! 
     186   END SUBROUTINE dom_qco_r3c 
     187 
     188    
     189   SUBROUTINE dom_qco_r3c_RK3( pssh, pr3t, pr3u, pr3v, pr3f ) 
     190      !!--------------------------------------------------------------------- 
     191      !!                   ***  ROUTINE r3c  *** 
     192      !! 
     193      !! ** Purpose :   compute the filtered ratio ssh/h_0 at t-,u-,v-,f-points 
     194      !! 
     195      !! ** Method  : - compute the ssh at u- and v-points (f-point optional) 
     196      !!                   Vector Form : surface weighted averaging 
     197      !!                   Flux   Form : simple           averaging 
     198      !!              - compute the ratio ssh/h_0 at t-,u-,v-pts, (f-pt optional) 
     199      !!---------------------------------------------------------------------- 
     200      REAL(wp), DIMENSION(:,:)          , INTENT(in   )  ::   pssh               ! sea surface height   [m] 
     201      REAL(wp), DIMENSION(:,:)          , INTENT(  out)  ::   pr3t, pr3u, pr3v   ! ssh/h0 ratio at t-, u-, v-,points  [-] 
     202      REAL(wp), DIMENSION(:,:), OPTIONAL, INTENT(  out)  ::   pr3f               ! ssh/h0 ratio at f-point   [-] 
     203      ! 
     204      INTEGER ::   ji, jj   ! dummy loop indices 
     205      !!---------------------------------------------------------------------- 
     206      ! 
     207      ! 
     208      DO_2D_OVR( nn_hls, nn_hls, nn_hls, nn_hls ) 
     209         pr3t(ji,jj) = pssh(ji,jj) * r1_ht_0(ji,jj)   !==  ratio at t-point  ==! 
     210      END_2D 
     211      ! 
     212      ! 
     213      !                                      !==  ratio at u-,v-point  ==! 
     214      ! 
     215!!st      IF( ln_dynadv_vec ) THEN                     !- Vector Form   (thickness weighted averaging) 
     216#if ! defined key_qcoTest_FluxForm 
     217      !                                ! no 'key_qcoTest_FluxForm' : surface weighted ssh average 
     218      DO_2D( 0, 0, 0, 0 ) 
     219         pr3u(ji,jj) = 0.5_wp * (  e1e2t(ji  ,jj) * pssh(ji  ,jj)  & 
     220            &                    + e1e2t(ji+1,jj) * pssh(ji+1,jj)  ) * r1_hu_0(ji,jj) * r1_e1e2u(ji,jj) 
     221         pr3v(ji,jj) = 0.5_wp * (  e1e2t(ji,jj  ) * pssh(ji,jj  )  & 
     222            &                    + e1e2t(ji,jj+1) * pssh(ji,jj+1)  ) * r1_hv_0(ji,jj) * r1_e1e2v(ji,jj) 
    166223      END_2D 
    167224!!st      ELSE                                         !- Flux Form   (simple averaging) 
    168225#else 
    169       DO_2D_OVR( nn_hls, nn_hls-1, nn_hls, nn_hls-1 ) 
     226      DO_2D( 0, 0, 0, 0 ) 
    170227         pr3u(ji,jj) = 0.5_wp * (  pssh(ji,jj) + pssh(ji+1,jj  )  ) * r1_hu_0(ji,jj) 
    171228         pr3v(ji,jj) = 0.5_wp * (  pssh(ji,jj) + pssh(ji  ,jj+1)  ) * r1_hv_0(ji,jj) 
     
    174231#endif          
    175232      ! 
    176       IF( .NOT.PRESENT( pr3f ) ) THEN              !- lbc on ratio at u-, v-points only 
    177          IF (nn_hls==1) CALL lbc_lnk( 'dom_qco_r3c', pr3u, 'U', 1._wp, pr3v, 'V', 1._wp ) 
    178          ! 
    179          ! 
    180       ELSE                                   !==  ratio at f-point  ==! 
     233      IF( PRESENT( pr3f ) ) THEN             !==  ratio at f-point  ==! 
    181234         ! 
    182235!!st         IF( ln_dynadv_vec )   THEN                !- Vector Form   (thickness weighted averaging) 
     
    184237         !                                ! no 'key_qcoTest_FluxForm' : surface weighted ssh average 
    185238 
    186       DO_2D_OVR( nn_hls, nn_hls-1, nn_hls, nn_hls-1 ) 
     239      DO_2D( 0, 0, 0, 0 ) 
    187240         ! round brackets added to fix the order of floating point operations 
    188241         ! needed to ensure halo 1 - halo 2 compatibility 
    189          pr3f(ji,jj) = 0.25_wp * ( ( e1e2t(ji  ,jj  ) * pssh(ji  ,jj  )   & 
    190             &                      + e1e2t(ji+1,jj  ) * pssh(ji+1,jj  )   & 
    191             &                      )                                      & ! bracket for halo 1 - halo 2 compatibility 
    192             &                     + ( e1e2t(ji  ,jj+1) * pssh(ji  ,jj+1)  & 
    193             &                       + e1e2t(ji+1,jj+1) * pssh(ji+1,jj+1)  & 
    194             &                       )                                     & ! bracket for halo 1 - halo 2 compatibility 
     242         pr3f(ji,jj) = 0.25_wp * (   (  e1e2t(ji  ,jj  ) * pssh(ji  ,jj  )      & 
     243            &                         + e1e2t(ji+1,jj  ) * pssh(ji+1,jj  )  )   & ! bracket for halo 1 - halo 2 compatibility 
     244            &                     +  (  e1e2t(ji  ,jj+1) * pssh(ji  ,jj+1)      & 
     245            &                         + e1e2t(ji+1,jj+1) * pssh(ji+1,jj+1)  )   & ! bracket for halo 1 - halo 2 compatibility 
    195246            &                    ) * r1_hf_0(ji,jj) * r1_e1e2f(ji,jj) 
    196247      END_2D 
    197248!!st         ELSE                                      !- Flux Form   (simple averaging) 
    198249#else 
    199       DO_2D_OVR( nn_hls, nn_hls-1, nn_hls, nn_hls-1 ) 
     250      DO_2D( 0, 0, 0, 0 ) 
    200251         ! round brackets added to fix the order of floating point operations 
    201252         ! needed to ensure halo 1 - halo 2 compatibility 
    202          pr3f(ji,jj) = 0.25_wp * ( ( pssh(ji,jj  ) + pssh(ji+1,jj  ) ) & 
    203             &                    + ( pssh(ji,jj+1) + pssh(ji+1,jj+1)   & 
    204             &                       )                                  & ! bracket for halo 1 - halo 2 compatibility 
     253         pr3f(ji,jj) = 0.25_wp * (   (  pssh(ji,jj  ) + pssh(ji+1,jj  )  )   & 
     254            &                     +  (  pssh(ji,jj+1) + pssh(ji+1,jj+1)  )   & ! bracket for halo 1 - halo 2 compatibility 
    205255            &                    ) * r1_hf_0(ji,jj) 
    206256      END_2D 
    207257!!st         ENDIF 
    208258#endif 
    209          !                                                 ! lbc on ratio at u-,v-,f-points 
    210          IF (nn_hls==1) CALL lbc_lnk( 'dom_qco_r3c', pr3u, 'U', 1._wp, pr3v, 'V', 1._wp, pr3f, 'F', 1._wp ) 
    211259         ! 
    212260      ENDIF 
    213261      ! 
    214    END SUBROUTINE dom_qco_r3c 
    215  
    216  
     262   END SUBROUTINE dom_qco_r3c_RK3 
     263 
     264     
    217265   SUBROUTINE qco_ctl 
    218266      !!--------------------------------------------------------------------- 
  • NEMO/branches/2021/dev_r14318_RK3_stage1/src/OCE/stprk3_stg.F90

    r15621 r15627  
    4949   REAL(wp) ::   r2_3 = 2._wp / 3._wp 
    5050 
    51    REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:) ::   ssha         ! sea-surface height  at N+1 
    52    REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:) ::   ua_b, va_b   ! barotropic velocity at N+1 
     51   REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:) ::   ssha                    ! sea-surface height  at N+1 
     52   REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:) ::   r3ta, r3ua, r3va, r3fa, r3fb   ! ssh/h_0 ratio at t,u,v-column 
     53   REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:) ::   ua_b, va_b              ! barotropic velocity at N+1 
    5354 
    5455   !! * Substitutions 
     
    120121         vv_b(:,:,Kaa) = r2_3 * vv_b(:,:,Kbb) + r1_3 * va_b(:,:) 
    121122         ! 
     123         ! 
     124         !                     !==  ssh/h0 ratio at Kaa  ==!  
     125         ! 
     126         IF( .NOT.lk_linssh ) THEN     ! "after" ssh/h_0 ratio at t,u,v-column 
     127            ! 
     128            ALLOCATE( r3ta(jpi,jpj) , r3ua(jpi,jpj) , r3va(jpi,jpj) , r3fa(jpi,jpj) , r3fb(jpi,jpj) ) 
     129            ! 
     130            r3fb(:,:) = r3f(:,:)       !!st dirty fix check with gm 
     131            CALL dom_qco_r3c_RK3( ssha, r3ta, r3ua, r3va, r3fa ) 
     132            ! 
     133            CALL lbc_lnk( 'stprk3_stg', r3ua, 'U', 1._wp, r3va, 'V', 1._wp, r3fa, 'F', 1._wp ) 
     134            ! 
     135            r3t(:,:,Kaa) = r2_3 * r3t(:,:,Kbb) + r1_3 * r3ta(:,:) 
     136            r3u(:,:,Kaa) = r2_3 * r3u(:,:,Kbb) + r1_3 * r3ua(:,:) 
     137            r3v(:,:,Kaa) = r2_3 * r3v(:,:,Kbb) + r1_3 * r3va(:,:) 
     138            r3f(:,:) = r2_3 * r3fb(:,:) + r1_3 * r3fa(:,:) 
     139         ENDIF 
     140         ! 
    122141         !                 !---------------! 
    123142      CASE ( 2 )           !==  Stage 2  ==!   Kbb = N   ;   Kmm = N+1/3   ;   Kaa = N+1/2 
     
    132151         vv_b(:,:,Kaa) = r1_2 * ( vv_b(:,:,Kbb) + va_b(:,:) ) 
    133152         ! 
     153         IF( .NOT.lk_linssh ) THEN 
     154            r3t(:,:,Kaa) = r1_2 * ( r3t(:,:,Kbb) + r3ta(:,:) ) 
     155            r3u(:,:,Kaa) = r1_2 * ( r3u(:,:,Kbb) + r3ua(:,:) ) 
     156            r3v(:,:,Kaa) = r1_2 * ( r3v(:,:,Kbb) + r3va(:,:) ) 
     157            r3f(:,:) = r1_2 * ( r3fb(:,:) + r3fa(:,:) ) 
     158         ENDIF 
     159         ! 
    134160         !                 !---------------! 
    135161      CASE ( 3 )           !==  Stage 3  ==!   Kbb = N   ;   Kmm = N+1/2   ;   Kaa = N+1 
     
    145171         DEALLOCATE( ssha , ua_b , va_b ) 
    146172         ! 
     173         IF( .NOT.lk_linssh ) THEN 
     174            r3t(:,:,Kaa) = r3ta(:,:) 
     175            r3u(:,:,Kaa) = r3ua(:,:) 
     176            r3v(:,:,Kaa) = r3va(:,:) 
     177            r3f(:,:    ) = r3fa(:,:) 
     178            DEALLOCATE( r3ta, r3ua, r3va, r3fa, r3fb ) 
     179            ! 
     180         ENDIF 
     181         ! 
    147182      END SELECT 
    148       ! 
    149       !                     !==  ssh/h0 ratio at Kaa  ==!  
    150       ! 
    151       IF( .NOT.lk_linssh )   CALL dom_qco_r3c( ssh(:,:,Kaa), r3t(:,:,Kaa), r3u(:,:,Kaa), r3v(:,:,Kaa), r3f(:,:) )   ! "after" ssh/h_0 ratio at t,u,v-column 
    152 !!      SELECT CASE( kstg ) 
    153          ! 
    154 !!      CASE ( 3 ) 
    155          !!st required at each stage for div hor loops 
    156          CALL lbc_lnk( 'stprk3_stg', r3u(:,:,Kaa), 'U', 1._wp, r3v(:,:,Kaa), 'V', 1._wp, r3f(:,:), 'F', 1._wp ) 
    157          ! 
    158 !!      END SELECT 
    159183      ! 
    160184      ! 
Note: See TracChangeset for help on using the changeset viewer.