Changeset 14419 for NEMO/branches/2021
- Timestamp:
- 2021-02-09T13:22:16+01:00 (3 years ago)
- Location:
- NEMO/branches/2021/dev_r14318_RK3_stage1/src/OCE/DYN
- Files:
-
- 3 edited
-
dynadv.F90 (modified) (3 diffs)
-
dynadv_cen2.F90 (modified) (4 diffs)
-
dynadv_ubs.F90 (modified) (6 diffs)
Legend:
- Unmodified
- Added
- Removed
-
NEMO/branches/2021/dev_r14318_RK3_stage1/src/OCE/DYN/dynadv.F90
r14381 r14419 51 51 CONTAINS 52 52 53 SUBROUTINE dyn_adv( kt, Kbb, Kmm, puu, pvv, Krhs, pau, pav, paw )53 SUBROUTINE dyn_adv( kt, Kbb, Kmm, puu, pvv, Krhs, pau, pav, paw, no_zad ) 54 54 !!--------------------------------------------------------------------- 55 55 !! *** ROUTINE dyn_adv *** … … 64 64 !! (see dynvor module). 65 65 !!---------------------------------------------------------------------- 66 INTEGER , INTENT(in ) :: kt ! ocean time-step index67 INTEGER , INTENT(in ) :: Kbb, Kmm, Krhs ! ocean time level indices68 REAL(wp), DIMENSION(:,:,:) , OPTIONAL, TARGET, INTENT(in ) :: pau, pav, paw! advective velocity69 REAL(wp), DIMENSION(jpi,jpj,jpk,jpt) , TARGET, INTENT(inout) :: puu, pvv! ocean velocities and RHS of momentum Eq.66 INTEGER , INTENT(in ) :: kt , Kbb, Kmm, Krhs ! ocean time step and level indices 67 INTEGER , OPTIONAL , INTENT(in ) :: no_zad ! no vertical advection compotation 68 REAL(wp), DIMENSION(:,:,:), OPTIONAL, TARGET, INTENT(in ) :: pau, pav, paw ! advective velocity 69 REAL(wp), DIMENSION(jpi,jpj,jpk,jpt), TARGET, INTENT(inout) :: puu, pvv ! ocean velocities and RHS of momentum Eq. 70 70 !!---------------------------------------------------------------------- 71 71 ! … … 73 73 ! 74 74 SELECT CASE( n_dynadv ) !== compute advection trend and add it to general trend ==! 75 CASE( np_VEC_c2 ) 76 CALL dyn_keg ( kt, nn_dynkeg , Kmm, puu, pvv, Krhs ) ! vector form :horizontal gradient of kinetic energy77 CALL dyn_zad ( kt , Kmm, puu, pvv, Krhs ) ! vector form :vertical advection78 CASE( np_FLX_c2 ) 79 CALL dyn_adv_cen2( kt , Kmm, puu, pvv, Krhs, pau, pav, paw ) ! 2nd order centered scheme75 CASE( np_VEC_c2 ) != vector form =! 76 CALL dyn_keg ( kt, nn_dynkeg , Kmm, puu, pvv, Krhs ) ! horizontal gradient of kinetic energy 77 CALL dyn_zad ( kt , Kmm, puu, pvv, Krhs ) ! vertical advection 78 CASE( np_FLX_c2 ) != flux form =! 79 CALL dyn_adv_cen2( kt , Kmm, puu, pvv, Krhs, pau, pav, paw, no_zad ) ! 2nd order centered scheme 80 80 CASE( np_FLX_ubs ) 81 CALL dyn_adv_ubs ( kt , Kbb, Kmm, puu, pvv, Krhs, pau, pav, paw ) ! 3rd order UBS scheme (UP3)81 CALL dyn_adv_ubs ( kt , Kbb, Kmm, puu, pvv, Krhs, pau, pav, paw, no_zad ) ! 3rd order UBS scheme (UP3) 82 82 END SELECT 83 83 ! -
NEMO/branches/2021/dev_r14318_RK3_stage1/src/OCE/DYN/dynadv_cen2.F90
r14381 r14419 36 36 CONTAINS 37 37 38 SUBROUTINE dyn_adv_cen2( kt, Kmm, puu, pvv, Krhs, pau, pav, paw )38 SUBROUTINE dyn_adv_cen2( kt, Kmm, puu, pvv, Krhs, pau, pav, paw, no_zad ) 39 39 !!---------------------------------------------------------------------- 40 40 !! *** ROUTINE dyn_adv_cen2 *** 41 41 !! 42 !! ** Purpose : Compute the nowmomentum advection trend in flux form42 !! ** Purpose : Compute the momentum advection trend in flux form 43 43 !! and the general trend of the momentum equation. 44 44 !! 45 !! ** Method : Trend evaluated using now fields (centered in time) 45 !! ** Method : Trend evaluated with a 2nd order centered scheme 46 !! using fields at Kmm time-level. 47 !! In RK3 time stepping case, the optional arguments (pau,pav,paw) 48 !! are present. They are used as advective velocity while 49 !! the advected velocity remains (puu,pvv). 46 50 !! 47 !! ** Action : (puu (:,:,:,Krhs),pvv(:,:,:,Krhs)) updated with the now vorticity termtrend51 !! ** Action : (puu,pvv)(:,:,:,Krhs) updated with the advective trend 48 52 !!---------------------------------------------------------------------- 49 INTEGER , INTENT(in ) :: kt ! ocean time-step index50 INTEGER , INTENT(in ) :: Kmm, Krhs ! ocean time level indices51 REAL(wp), DIMENSION(jpi,jpj,jpk,jpt) , TARGET, INTENT(inout) :: puu, pvv! ocean velocities and RHS of momentum equation52 REAL(wp), DIMENSION(:,:,:) , OPTIONAL, TARGET, INTENT(in ) :: pau, pav, paw! advective velocity53 INTEGER , INTENT(in ) :: kt , Kmm, Krhs ! ocean time-step and level indices 54 INTEGER , OPTIONAL , INTENT(in ) :: no_zad ! no vertical advection computation 55 REAL(wp), DIMENSION(jpi,jpj,jpk,jpt), TARGET, INTENT(inout) :: puu, pvv ! ocean velocities and RHS of momentum equation 56 REAL(wp), DIMENSION(:,:,:), OPTIONAL, TARGET, INTENT(in ) :: pau, pav, paw ! advective velocity 53 57 ! 54 58 INTEGER :: ji, jj, jk ! dummy loop indices 59 REAL(wp) :: zzu, zzv ! local scalars 55 60 REAL(wp), DIMENSION(jpi,jpj,jpk) :: zfu_t, zfu_f, zfu_uw, zfu 56 61 REAL(wp), DIMENSION(jpi,jpj,jpk) :: zfv_t, zfv_f, zfv_vw, zfv, zfw 57 REAL(wp), DIMENSION(:,:,:), POINTER :: zpt u, zptv, zptw62 REAL(wp), DIMENSION(:,:,:), POINTER :: zpt_u, zpt_v, zpt_w 58 63 !!---------------------------------------------------------------------- 59 64 ! … … 70 75 ! 71 76 IF( PRESENT( pau ) ) THEN ! RK3: advective velocity (pau,pav,paw) /= advected velocity (puu,pvv,ww) 72 zpt u => pau(:,:,:)73 zpt v => pav(:,:,:)74 zpt w => paw(:,:,:)77 zpt_u => pau(:,:,:) 78 zpt_v => pav(:,:,:) 79 zpt_w => paw(:,:,:) 75 80 ELSE ! MLF: advective velocity = (puu,pvv,ww) 76 zpt u => puu(:,:,:,Kmm)77 zpt v => pvv(:,:,:,Kmm)78 zpt w => ww (:,:,: )81 zpt_u => puu(:,:,:,Kmm) 82 zpt_v => pvv(:,:,:,Kmm) 83 zpt_w => ww (:,:,: ) 79 84 ENDIF 80 85 ! … … 82 87 ! 83 88 DO jk = 1, jpkm1 ! horizontal transport 84 zfu(:,:,jk) = 0.25_wp * e2u(:,:) * e3u(:,:,jk,Kmm) * zpt u(:,:,jk)85 zfv(:,:,jk) = 0.25_wp * e1v(:,:) * e3v(:,:,jk,Kmm) * zpt v(:,:,jk)89 zfu(:,:,jk) = 0.25_wp * e2u(:,:) * e3u(:,:,jk,Kmm) * zpt_u(:,:,jk) 90 zfv(:,:,jk) = 0.25_wp * e1v(:,:) * e3v(:,:,jk,Kmm) * zpt_v(:,:,jk) 86 91 DO_2D( 1, 0, 1, 0 ) ! horizontal momentum fluxes (at T- and F-point) 87 92 zfu_t(ji+1,jj ,jk) = ( zfu(ji,jj,jk) + zfu(ji+1,jj,jk) ) * ( puu(ji,jj,jk,Kmm) + puu(ji+1,jj ,jk,Kmm) ) … … 108 113 ENDIF 109 114 ! 110 ! !== Vertical advection ==! 111 ! 112 DO_2D( 0, 0, 0, 0 ) ! surface/bottom advective fluxes set to zero 113 zfu_uw(ji,jj,jpk) = 0._wp ; zfv_vw(ji,jj,jpk) = 0._wp 114 zfu_uw(ji,jj, 1 ) = 0._wp ; zfv_vw(ji,jj, 1 ) = 0._wp 115 END_2D 116 IF( ln_linssh ) THEN ! linear free surface: advection through the surface 117 DO_2D( 0, 0, 0, 0 ) 118 zfu_uw(ji,jj,1) = 0.5_wp * ( e1e2t(ji,jj) * zptw(ji,jj,1) + e1e2t(ji+1,jj) * zptw(ji+1,jj,1) ) * puu(ji,jj,1,Kmm) 119 zfv_vw(ji,jj,1) = 0.5_wp * ( e1e2t(ji,jj) * zptw(ji,jj,1) + e1e2t(ji,jj+1) * zptw(ji,jj+1,1) ) * pvv(ji,jj,1,Kmm) 115 IF( PRESENT( no_zad ) ) THEN !== No vertical advection ==! (except if linear free surface) 116 ! == 117 IF( ln_linssh ) THEN ! linear free surface: advection through the surface z=0 118 DO_2D( 0, 0, 0, 0 ) 119 zzu = 0.5_wp * ( e1e2t(ji,jj) * zpt_w(ji,jj,1) + e1e2t(ji+1,jj) * zpt_w(ji+1,jj,1) ) * puu(ji,jj,1,Kmm) 120 zzv = 0.5_wp * ( e1e2t(ji,jj) * zpt_w(ji,jj,1) + e1e2t(ji,jj+1) * zpt_w(ji,jj+1,1) ) * pvv(ji,jj,1,Kmm) 121 puu(ji,jj,1,Krhs) = puu(ji,jj,1,Krhs) - zzu * r1_e1e2u(ji,jj) & 122 & / e3u(ji,jj,1,Kmm) 123 pvv(ji,jj,1,Krhs) = pvv(ji,jj,1,Krhs) - zzv * r1_e1e2v(ji,jj) & 124 & / e3v(ji,jj,1,Kmm) 125 END_2D 126 ENDIF 127 ! 128 ELSE !== Vertical advection ==! 129 ! 130 DO_2D( 0, 0, 0, 0 ) ! surface/bottom advective fluxes set to zero 131 zfu_uw(ji,jj,jpk) = 0._wp ; zfv_vw(ji,jj,jpk) = 0._wp 132 zfu_uw(ji,jj, 1 ) = 0._wp ; zfv_vw(ji,jj, 1 ) = 0._wp 120 133 END_2D 134 IF( ln_linssh ) THEN ! linear free surface: advection through the surface z=0 135 DO_2D( 0, 0, 0, 0 ) 136 zfu_uw(ji,jj,1) = 0.5_wp * ( e1e2t(ji,jj) * zpt_w(ji,jj,1) + e1e2t(ji+1,jj) * zpt_w(ji+1,jj,1) ) * puu(ji,jj,1,Kmm) 137 zfv_vw(ji,jj,1) = 0.5_wp * ( e1e2t(ji,jj) * zpt_w(ji,jj,1) + e1e2t(ji,jj+1) * zpt_w(ji,jj+1,1) ) * pvv(ji,jj,1,Kmm) 138 END_2D 139 ENDIF 140 DO jk = 2, jpkm1 ! interior advective fluxes 141 DO_2D( 0, 1, 0, 1 ) ! 1/4 * Vertical transport 142 zfw(ji,jj,jk) = 0.25_wp * e1e2t(ji,jj) * zpt_w(ji,jj,jk) 143 END_2D 144 DO_2D( 0, 0, 0, 0 ) 145 zfu_uw(ji,jj,jk) = ( zfw(ji,jj,jk) + zfw(ji+1,jj ,jk) ) * ( puu(ji,jj,jk,Kmm) + puu(ji,jj,jk-1,Kmm) ) 146 zfv_vw(ji,jj,jk) = ( zfw(ji,jj,jk) + zfw(ji ,jj+1,jk) ) * ( pvv(ji,jj,jk,Kmm) + pvv(ji,jj,jk-1,Kmm) ) 147 END_2D 148 END DO 149 DO_3D( 0, 0, 0, 0, 1, jpkm1 ) ! divergence of vertical momentum flux divergence 150 puu(ji,jj,jk,Krhs) = puu(ji,jj,jk,Krhs) - ( zfu_uw(ji,jj,jk) - zfu_uw(ji,jj,jk+1) ) * r1_e1e2u(ji,jj) & 151 & / e3u(ji,jj,jk,Kmm) 152 pvv(ji,jj,jk,Krhs) = pvv(ji,jj,jk,Krhs) - ( zfv_vw(ji,jj,jk) - zfv_vw(ji,jj,jk+1) ) * r1_e1e2v(ji,jj) & 153 & / e3v(ji,jj,jk,Kmm) 154 END_3D 155 ! 156 IF( l_trddyn ) THEN ! trends: send trend to trddyn for diagnostic 157 zfu_t(:,:,:) = puu(:,:,:,Krhs) - zfu_t(:,:,:) 158 zfv_t(:,:,:) = pvv(:,:,:,Krhs) - zfv_t(:,:,:) 159 CALL trd_dyn( zfu_t, zfv_t, jpdyn_zad, kt, Kmm ) 160 ENDIF 161 ! ! Control print 162 IF(sn_cfctl%l_prtctl) CALL prt_ctl( tab3d_1=puu(:,:,:,Krhs), clinfo1=' cen2 adv - Ua: ', mask1=umask, & 163 & tab3d_2=pvv(:,:,:,Krhs), clinfo2= ' Va: ', mask2=vmask, clinfo3='dyn' ) 164 ! 121 165 ENDIF 122 DO jk = 2, jpkm1 ! interior advective fluxes123 DO_2D( 0, 1, 0, 1 ) ! 1/4 * Vertical transport124 zfw(ji,jj,jk) = 0.25_wp * e1e2t(ji,jj) * zptw(ji,jj,jk)125 END_2D126 DO_2D( 0, 0, 0, 0 )127 zfu_uw(ji,jj,jk) = ( zfw(ji,jj,jk) + zfw(ji+1,jj ,jk) ) * ( puu(ji,jj,jk,Kmm) + puu(ji,jj,jk-1,Kmm) )128 zfv_vw(ji,jj,jk) = ( zfw(ji,jj,jk) + zfw(ji ,jj+1,jk) ) * ( pvv(ji,jj,jk,Kmm) + pvv(ji,jj,jk-1,Kmm) )129 END_2D130 END DO131 DO_3D( 0, 0, 0, 0, 1, jpkm1 ) ! divergence of vertical momentum flux divergence132 puu(ji,jj,jk,Krhs) = puu(ji,jj,jk,Krhs) - ( zfu_uw(ji,jj,jk) - zfu_uw(ji,jj,jk+1) ) * r1_e1e2u(ji,jj) &133 & / e3u(ji,jj,jk,Kmm)134 pvv(ji,jj,jk,Krhs) = pvv(ji,jj,jk,Krhs) - ( zfv_vw(ji,jj,jk) - zfv_vw(ji,jj,jk+1) ) * r1_e1e2v(ji,jj) &135 & / e3v(ji,jj,jk,Kmm)136 END_3D137 !138 IF( l_trddyn ) THEN ! trends: send trend to trddyn for diagnostic139 zfu_t(:,:,:) = puu(:,:,:,Krhs) - zfu_t(:,:,:)140 zfv_t(:,:,:) = pvv(:,:,:,Krhs) - zfv_t(:,:,:)141 CALL trd_dyn( zfu_t, zfv_t, jpdyn_zad, kt, Kmm )142 ENDIF143 ! ! Control print144 IF(sn_cfctl%l_prtctl) CALL prt_ctl( tab3d_1=puu(:,:,:,Krhs), clinfo1=' cen2 adv - Ua: ', mask1=umask, &145 & tab3d_2=pvv(:,:,:,Krhs), clinfo2= ' Va: ', mask2=vmask, clinfo3='dyn' )146 166 ! 147 167 END SUBROUTINE dyn_adv_cen2 -
NEMO/branches/2021/dev_r14318_RK3_stage1/src/OCE/DYN/dynadv_ubs.F90
r14381 r14419 42 42 CONTAINS 43 43 44 SUBROUTINE dyn_adv_ubs( kt, Kbb, Kmm, puu, pvv, Krhs, pau, pav, paw )44 SUBROUTINE dyn_adv_ubs( kt, Kbb, Kmm, puu, pvv, Krhs, pau, pav, paw, no_zad ) 45 45 !!---------------------------------------------------------------------- 46 46 !! *** ROUTINE dyn_adv_ubs *** … … 65 65 !! gamma1=1/3 and gamma2=1/32. 66 66 !! 67 !! ** Action : - (puu(:,:,:,Krhs),pvv(:,:,:,Krhs)) updated with the 3D advective momentum trends 67 !! In RK3 time stepping case, the optional arguments 68 !! (pau,pav,paw) are present. They are used as advective velocity 69 !! while the advected velocity remains (puu,pvv). 70 !! 71 !! ** Action : (puu,pvv)(:,:,:,Krhs) updated with the advective trend 68 72 !! 69 73 !! Reference : Shchepetkin & McWilliams, 2005, Ocean Modelling. 70 74 !!---------------------------------------------------------------------- 71 INTEGER , INTENT(in ) :: kt ! ocean time-step index72 INTEGER , INTENT(in ) :: Kbb, Kmm, Krhs ! ocean time level indices73 REAL(wp), DIMENSION(jpi,jpj,jpk,jpt) , TARGET, INTENT(inout) :: puu, pvv! ocean velocities and RHS of momentum equation74 REAL(wp), DIMENSION(:,:,:) , OPTIONAL, TARGET, INTENT(in ) :: pau, pav, paw! advective velocity75 INTEGER , INTENT(in ) :: kt , Kbb, Kmm, Krhs ! ocean time-step and level indices 76 INTEGER , OPTIONAL , INTENT(in ) :: no_zad ! no vertical advection compotation 77 REAL(wp), DIMENSION(jpi,jpj,jpk,jpt), TARGET, INTENT(inout) :: puu, pvv ! ocean velocities and RHS of momentum equation 78 REAL(wp), DIMENSION(:,:,:), OPTIONAL, TARGET, INTENT(in ) :: pau, pav, paw ! advective velocity 75 79 ! 76 80 INTEGER :: ji, jj, jk ! dummy loop indices 77 REAL(wp) :: zui, zvj, zfuj, zfvi, zl_u, zl_v ! local scalars81 REAL(wp) :: zui, zvj, zfuj, zfvi, zl_u, zl_v, zzu, zzv ! local scalars 78 82 REAL(wp), DIMENSION(jpi,jpj,jpk) :: zfu_t, zfu_f, zfu_uw, zfu 79 83 REAL(wp), DIMENSION(jpi,jpj,jpk) :: zfv_t, zfv_f, zfv_vw, zfv, zfw 80 84 REAL(wp), DIMENSION(jpi,jpj,jpk,2) :: zlu_uu, zlu_uv 81 85 REAL(wp), DIMENSION(jpi,jpj,jpk,2) :: zlv_vv, zlv_vu 82 REAL(wp), DIMENSION(:,:,:), POINTER :: zpt u, zptv, zptw86 REAL(wp), DIMENSION(:,:,:), POINTER :: zpt_u, zpt_v, zpt_w 83 87 !!---------------------------------------------------------------------- 84 88 ! … … 105 109 ! 106 110 IF( PRESENT( pau ) ) THEN ! RK3: advective velocity (pau,pav,paw) /= advected velocity (puu,pvv,ww) 107 zpt u => pau(:,:,:)108 zpt v => pav(:,:,:)109 zpt w => paw(:,:,:)111 zpt_u => pau(:,:,:) 112 zpt_v => pav(:,:,:) 113 zpt_w => paw(:,:,:) 110 114 ELSE ! MLF: advective velocity = (puu,pvv,ww) 111 zpt u => puu(:,:,:,Kmm)112 zpt v => pvv(:,:,:,Kmm)113 zpt w => ww (:,:,: )115 zpt_u => puu(:,:,:,Kmm) 116 zpt_v => pvv(:,:,:,Kmm) 117 zpt_w => ww (:,:,: ) 114 118 ENDIF 115 119 ! … … 118 122 ! ! =========================== ! 119 123 ! ! horizontal volume fluxes 120 zfu(:,:,jk) = e2u(:,:) * e3u(:,:,jk,Kmm) * zpt u(:,:,jk)121 zfv(:,:,jk) = e1v(:,:) * e3v(:,:,jk,Kmm) * zpt v(:,:,jk)124 zfu(:,:,jk) = e2u(:,:) * e3u(:,:,jk,Kmm) * zpt_u(:,:,jk) 125 zfv(:,:,jk) = e1v(:,:) * e3v(:,:,jk,Kmm) * zpt_v(:,:,jk) 122 126 ! 123 127 DO_2D( 0, 0, 0, 0 ) ! laplacian … … 146 150 DO jk = 1, jpkm1 ! ====================== ! 147 151 ! ! horizontal volume fluxes 148 zfu(:,:,jk) = 0.25_wp * e2u(:,:) * e3u(:,:,jk,Kmm) * zpt u(:,:,jk)149 zfv(:,:,jk) = 0.25_wp * e1v(:,:) * e3v(:,:,jk,Kmm) * zpt v(:,:,jk)152 zfu(:,:,jk) = 0.25_wp * e2u(:,:) * e3u(:,:,jk,Kmm) * zpt_u(:,:,jk) 153 zfv(:,:,jk) = 0.25_wp * e1v(:,:) * e3v(:,:,jk,Kmm) * zpt_v(:,:,jk) 150 154 ! 151 155 DO_2D( 1, 0, 1, 0 ) ! horizontal momentum fluxes at T- and F-point … … 200 204 ! ! Vertical advection ! 201 205 ! ! ==================== ! 202 DO_2D( 0, 0, 0, 0 ) ! surface/bottom advective fluxes set to zero 203 zfu_uw(ji,jj,jpk) = 0._wp 204 zfv_vw(ji,jj,jpk) = 0._wp 205 zfu_uw(ji,jj, 1 ) = 0._wp 206 zfv_vw(ji,jj, 1 ) = 0._wp 207 END_2D 208 IF( ln_linssh ) THEN ! constant volume : advection through the surface 209 DO_2D( 0, 0, 0, 0 ) 210 zfu_uw(ji,jj,1) = 0.5_wp * ( e1e2t(ji,jj) * zptw(ji,jj,1) + e1e2t(ji+1,jj) * zptw(ji+1,jj,1) ) * puu(ji,jj,1,Kmm) 211 zfv_vw(ji,jj,1) = 0.5_wp * ( e1e2t(ji,jj) * zptw(ji,jj,1) + e1e2t(ji,jj+1) * zptw(ji,jj+1,1) ) * pvv(ji,jj,1,Kmm) 212 END_2D 213 ENDIF 214 DO jk = 2, jpkm1 ! interior fluxes 215 DO_2D( 0, 1, 0, 1 ) 216 zfw(ji,jj,jk) = 0.25_wp * e1e2t(ji,jj) * zptw(ji,jj,jk) 217 END_2D 218 DO_2D( 0, 0, 0, 0 ) 219 zfu_uw(ji,jj,jk) = ( zfw(ji,jj,jk)+ zfw(ji+1,jj,jk) ) * ( puu(ji,jj,jk,Kmm) + puu(ji,jj,jk-1,Kmm) ) 220 zfv_vw(ji,jj,jk) = ( zfw(ji,jj,jk)+ zfw(ji,jj+1,jk) ) * ( pvv(ji,jj,jk,Kmm) + pvv(ji,jj,jk-1,Kmm) ) 221 END_2D 222 END DO 223 DO_3D( 0, 0, 0, 0, 1, jpkm1 ) ! divergence of vertical momentum flux divergence 224 puu(ji,jj,jk,Krhs) = puu(ji,jj,jk,Krhs) - ( zfu_uw(ji,jj,jk) - zfu_uw(ji,jj,jk+1) ) * r1_e1e2u(ji,jj) & 225 & / e3u(ji,jj,jk,Kmm) 226 pvv(ji,jj,jk,Krhs) = pvv(ji,jj,jk,Krhs) - ( zfv_vw(ji,jj,jk) - zfv_vw(ji,jj,jk+1) ) * r1_e1e2v(ji,jj) & 227 & / e3v(ji,jj,jk,Kmm) 228 END_3D 229 ! 230 IF( l_trddyn ) THEN ! save the vertical advection trend for diagnostic 231 zfu_t(:,:,:) = puu(:,:,:,Krhs) - zfu_t(:,:,:) 232 zfv_t(:,:,:) = pvv(:,:,:,Krhs) - zfv_t(:,:,:) 233 CALL trd_dyn( zfu_t, zfv_t, jpdyn_zad, kt, Kmm ) 234 ENDIF 235 ! ! Control print 236 IF(sn_cfctl%l_prtctl) CALL prt_ctl( tab3d_1=puu(:,:,:,Krhs), clinfo1=' ubs2 adv - Ua: ', mask1=umask, & 237 & tab3d_2=pvv(:,:,:,Krhs), clinfo2= ' Va: ', mask2=vmask, clinfo3='dyn' ) 206 ! 207 ! ! ======================== ! 208 IF( PRESENT( no_zad ) ) THEN ! No vertical advection ! (except if linear free surface) 209 ! ! ======================== ! ------ 210 ! 211 IF( ln_linssh ) THEN ! linear free surface: advection through the surface z=0 212 DO_2D( 0, 0, 0, 0 ) 213 zzu = 0.5_wp * ( e1e2t(ji,jj) * zpt_w(ji,jj,1) + e1e2t(ji+1,jj) * zpt_w(ji+1,jj,1) ) * puu(ji,jj,1,Kmm) 214 zzv = 0.5_wp * ( e1e2t(ji,jj) * zpt_w(ji,jj,1) + e1e2t(ji,jj+1) * zpt_w(ji,jj+1,1) ) * pvv(ji,jj,1,Kmm) 215 puu(ji,jj,1,Krhs) = puu(ji,jj,1,Krhs) - zzu * r1_e1e2u(ji,jj) & 216 & / e3u(ji,jj,1,Kmm) 217 pvv(ji,jj,1,Krhs) = pvv(ji,jj,1,Krhs) - zzv * r1_e1e2v(ji,jj) & 218 & / e3v(ji,jj,1,Kmm) 219 END_2D 220 ENDIF 221 ! ! =================== ! 222 ELSE ! Vertical advection ! 223 ! ! =================== ! 224 DO_2D( 0, 0, 0, 0 ) ! surface/bottom advective fluxes set to zero 225 zfu_uw(ji,jj,jpk) = 0._wp 226 zfv_vw(ji,jj,jpk) = 0._wp 227 zfu_uw(ji,jj, 1 ) = 0._wp 228 zfv_vw(ji,jj, 1 ) = 0._wp 229 END_2D 230 IF( ln_linssh ) THEN ! constant volume : advection through the surface 231 DO_2D( 0, 0, 0, 0 ) 232 zfu_uw(ji,jj,1) = 0.5_wp * ( e1e2t(ji,jj) * zpt_w(ji,jj,1) + e1e2t(ji+1,jj) * zpt_w(ji+1,jj,1) ) * puu(ji,jj,1,Kmm) 233 zfv_vw(ji,jj,1) = 0.5_wp * ( e1e2t(ji,jj) * zpt_w(ji,jj,1) + e1e2t(ji,jj+1) * zpt_w(ji,jj+1,1) ) * pvv(ji,jj,1,Kmm) 234 END_2D 235 ENDIF 236 DO jk = 2, jpkm1 ! interior fluxes 237 DO_2D( 0, 1, 0, 1 ) 238 zfw(ji,jj,jk) = 0.25_wp * e1e2t(ji,jj) * zpt_w(ji,jj,jk) 239 END_2D 240 DO_2D( 0, 0, 0, 0 ) 241 zfu_uw(ji,jj,jk) = ( zfw(ji,jj,jk)+ zfw(ji+1,jj,jk) ) * ( puu(ji,jj,jk,Kmm) + puu(ji,jj,jk-1,Kmm) ) 242 zfv_vw(ji,jj,jk) = ( zfw(ji,jj,jk)+ zfw(ji,jj+1,jk) ) * ( pvv(ji,jj,jk,Kmm) + pvv(ji,jj,jk-1,Kmm) ) 243 END_2D 244 END DO 245 DO_3D( 0, 0, 0, 0, 1, jpkm1 ) ! divergence of vertical momentum flux divergence 246 puu(ji,jj,jk,Krhs) = puu(ji,jj,jk,Krhs) - ( zfu_uw(ji,jj,jk) - zfu_uw(ji,jj,jk+1) ) * r1_e1e2u(ji,jj) & 247 & / e3u(ji,jj,jk,Kmm) 248 pvv(ji,jj,jk,Krhs) = pvv(ji,jj,jk,Krhs) - ( zfv_vw(ji,jj,jk) - zfv_vw(ji,jj,jk+1) ) * r1_e1e2v(ji,jj) & 249 & / e3v(ji,jj,jk,Kmm) 250 END_3D 251 ! 252 IF( l_trddyn ) THEN ! save the vertical advection trend for diagnostic 253 zfu_t(:,:,:) = puu(:,:,:,Krhs) - zfu_t(:,:,:) 254 zfv_t(:,:,:) = pvv(:,:,:,Krhs) - zfv_t(:,:,:) 255 CALL trd_dyn( zfu_t, zfv_t, jpdyn_zad, kt, Kmm ) 256 ENDIF 257 ! ! Control print 258 IF(sn_cfctl%l_prtctl) CALL prt_ctl( tab3d_1=puu(:,:,:,Krhs), clinfo1=' ubs2 adv - Ua: ', mask1=umask, & 259 & tab3d_2=pvv(:,:,:,Krhs), clinfo2= ' Va: ', mask2=vmask, clinfo3='dyn' ) 260 ! 261 ENDIF 238 262 ! 239 263 END SUBROUTINE dyn_adv_ubs
Note: See TracChangeset
for help on using the changeset viewer.
