Changeset 6225 for branches/2014/dev_r4704_NOC5_MPP_BDY_UPDATE/NEMOGCM/NEMO/OPA_SRC/DYN/dynadv_ubs.F90
- Timestamp:
- 2016-01-08T10:35:19+01:00 (8 years ago)
- File:
-
- 1 edited
Legend:
- Unmodified
- Added
- Removed
-
branches/2014/dev_r4704_NOC5_MPP_BDY_UPDATE/NEMOGCM/NEMO/OPA_SRC/DYN/dynadv_ubs.F90
r4153 r6225 16 16 USE oce ! ocean dynamics and tracers 17 17 USE dom_oce ! ocean space and time domain 18 USE trdmod ! ocean dynamics trends 19 USE trdmod_oce ! ocean variables trends 18 USE trd_oce ! trends: ocean variables 19 USE trddyn ! trend manager: dynamics 20 ! 20 21 USE in_out_manager ! I/O manager 21 22 USE prtctl ! Print control 22 23 USE lbclnk ! ocean lateral boundary conditions (or mpp link) 23 24 USE lib_mpp ! MPP library 24 USE wrk_nemo 25 USE timing 25 USE wrk_nemo ! Memory Allocation 26 USE timing ! Timing 26 27 27 28 IMPLICIT NONE … … 34 35 35 36 !! * Substitutions 36 # include "domzgr_substitute.h90"37 37 # include "vectopt_loop_substitute.h90" 38 38 !!---------------------------------------------------------------------- … … 70 70 !! Reference : Shchepetkin & McWilliams, 2005, Ocean Modelling. 71 71 !!---------------------------------------------------------------------- 72 INTEGER, INTENT(in) :: kt ! ocean time-step index 73 ! 74 INTEGER :: ji, jj, jk ! dummy loop indices 75 REAL(wp) :: zbu, zbv ! temporary scalars 76 REAL(wp) :: zui, zvj, zfuj, zfvi, zl_u, zl_v ! temporary scalars 72 INTEGER, INTENT(in) :: kt ! ocean time-step index 73 ! 74 INTEGER :: ji, jj, jk ! dummy loop indices 75 REAL(wp) :: zui, zvj, zfuj, zfvi, zl_u, zl_v ! local scalars 77 76 REAL(wp), POINTER, DIMENSION(:,:,: ) :: zfu, zfv 78 77 REAL(wp), POINTER, DIMENSION(:,:,: ) :: zfu_t, zfv_t, zfu_f, zfv_f, zfu_uw, zfv_vw, zfw … … 82 81 IF( nn_timing == 1 ) CALL timing_start('dyn_adv_ubs') 83 82 ! 84 CALL wrk_alloc( jpi, jpj, jpk,zfu_t , zfv_t , zfu_f , zfv_f, zfu_uw, zfv_vw, zfu, zfv, zfw )85 CALL wrk_alloc( jpi, jpj, jpk, jpts,zlu_uu, zlv_vv, zlu_uv, zlv_vu )83 CALL wrk_alloc( jpi,jpj,jpk, zfu_t , zfv_t , zfu_f , zfv_f, zfu_uw, zfv_vw, zfu, zfv, zfw ) 84 CALL wrk_alloc( jpi,jpj,jpk,jpts, zlu_uu, zlv_vv, zlu_uv, zlv_vu ) 86 85 ! 87 86 IF( kt == nit000 ) THEN … … 100 99 zlu_uv(:,:,:,:) = 0._wp 101 100 zlv_vu(:,:,:,:) = 0._wp 102 103 IF( l_trddyn ) THEN ! Save ua and vatrends101 ! 102 IF( l_trddyn ) THEN ! trends: store the input trends 104 103 zfu_uw(:,:,:) = ua(:,:,:) 105 104 zfv_vw(:,:,:) = va(:,:,:) 106 105 ENDIF 107 108 106 ! ! =========================== ! 109 107 DO jk = 1, jpkm1 ! Laplacian of the velocity ! 110 108 ! ! =========================== ! 111 109 ! ! horizontal volume fluxes 112 zfu(:,:,jk) = e2u(:,:) * fse3u(:,:,jk) * un(:,:,jk)113 zfv(:,:,jk) = e1v(:,:) * fse3v(:,:,jk) * vn(:,:,jk)110 zfu(:,:,jk) = e2u(:,:) * e3u_n(:,:,jk) * un(:,:,jk) 111 zfv(:,:,jk) = e1v(:,:) * e3v_n(:,:,jk) * vn(:,:,jk) 114 112 ! 115 113 DO jj = 2, jpjm1 ! laplacian 116 114 DO ji = fs_2, fs_jpim1 ! vector opt. 117 zlu_uu(ji,jj,jk,1) = ( ub (ji+1,jj,jk)-2.*ub (ji,jj,jk)+ub (ji-1,jj,jk) ) * umask(ji,jj,jk) 118 zlv_vv(ji,jj,jk,1) = ( vb (ji,jj+1,jk)-2.*vb (ji,jj,jk)+vb (ji,jj-1,jk) ) * vmask(ji,jj,jk) 119 zlu_uv(ji,jj,jk,1) = ( ub (ji,jj+1,jk)-2.*ub (ji,jj,jk)+ub (ji,jj-1,jk) ) * umask(ji,jj,jk) 120 zlv_vu(ji,jj,jk,1) = ( vb (ji+1,jj,jk)-2.*vb (ji,jj,jk)+vb (ji-1,jj,jk) ) * vmask(ji,jj,jk) 121 ! 122 zlu_uu(ji,jj,jk,2) = ( zfu(ji+1,jj,jk)-2.*zfu(ji,jj,jk)+zfu(ji-1,jj,jk) ) * umask(ji,jj,jk) 123 zlv_vv(ji,jj,jk,2) = ( zfv(ji,jj+1,jk)-2.*zfv(ji,jj,jk)+zfv(ji,jj-1,jk) ) * vmask(ji,jj,jk) 124 zlu_uv(ji,jj,jk,2) = ( zfu(ji,jj+1,jk)-2.*zfu(ji,jj,jk)+zfu(ji,jj-1,jk) ) * umask(ji,jj,jk) 125 zlv_vu(ji,jj,jk,2) = ( zfv(ji+1,jj,jk)-2.*zfv(ji,jj,jk)+zfv(ji-1,jj,jk) ) * vmask(ji,jj,jk) 126 END DO 127 END DO 128 END DO 129 !!gm BUG !!! just below this should be +1 in all the communications 130 ! CALL lbc_lnk( zlu_uu(:,:,:,1), 'U', -1.) ; CALL lbc_lnk( zlu_uv(:,:,:,1), 'U', -1.) 131 ! CALL lbc_lnk( zlu_uu(:,:,:,2), 'U', -1.) ; CALL lbc_lnk( zlu_uv(:,:,:,2), 'U', -1.) 132 ! CALL lbc_lnk( zlv_vv(:,:,:,1), 'V', -1.) ; CALL lbc_lnk( zlv_vu(:,:,:,1), 'V', -1.) 133 ! CALL lbc_lnk( zlv_vv(:,:,:,2), 'V', -1.) ; CALL lbc_lnk( zlv_vu(:,:,:,2), 'V', -1.) 134 ! 135 !!gm corrected: 115 zlu_uu(ji,jj,jk,1) = ( ub (ji+1,jj ,jk) - 2.*ub (ji,jj,jk) + ub (ji-1,jj ,jk) ) * umask(ji,jj,jk) 116 zlv_vv(ji,jj,jk,1) = ( vb (ji ,jj+1,jk) - 2.*vb (ji,jj,jk) + vb (ji ,jj-1,jk) ) * vmask(ji,jj,jk) 117 zlu_uv(ji,jj,jk,1) = ( ub (ji ,jj+1,jk) - ub (ji ,jj ,jk) ) * fmask(ji ,jj ,jk) & 118 & - ( ub (ji ,jj ,jk) - ub (ji ,jj-1,jk) ) * fmask(ji ,jj-1,jk) 119 zlv_vu(ji,jj,jk,1) = ( vb (ji+1,jj ,jk) - vb (ji ,jj ,jk) ) * fmask(ji ,jj ,jk) & 120 & - ( vb (ji ,jj ,jk) - vb (ji-1,jj ,jk) ) * fmask(ji-1,jj ,jk) 121 ! 122 zlu_uu(ji,jj,jk,2) = ( zfu(ji+1,jj ,jk) - 2.*zfu(ji,jj,jk) + zfu(ji-1,jj ,jk) ) * umask(ji,jj,jk) 123 zlv_vv(ji,jj,jk,2) = ( zfv(ji ,jj+1,jk) - 2.*zfv(ji,jj,jk) + zfv(ji ,jj-1,jk) ) * vmask(ji,jj,jk) 124 zlu_uv(ji,jj,jk,2) = ( zfu(ji ,jj+1,jk) - zfu(ji ,jj ,jk) ) * fmask(ji ,jj ,jk) & 125 & - ( zfu(ji ,jj ,jk) - zfu(ji ,jj-1,jk) ) * fmask(ji ,jj-1,jk) 126 zlv_vu(ji,jj,jk,2) = ( zfv(ji+1,jj ,jk) - zfv(ji ,jj ,jk) ) * fmask(ji ,jj ,jk) & 127 & - ( zfv(ji ,jj ,jk) - zfv(ji-1,jj ,jk) ) * fmask(ji-1,jj ,jk) 128 END DO 129 END DO 130 END DO 136 131 CALL lbc_lnk( zlu_uu(:,:,:,1), 'U', 1. ) ; CALL lbc_lnk( zlu_uv(:,:,:,1), 'U', 1. ) 137 132 CALL lbc_lnk( zlu_uu(:,:,:,2), 'U', 1. ) ; CALL lbc_lnk( zlu_uv(:,:,:,2), 'U', 1. ) 138 133 CALL lbc_lnk( zlv_vv(:,:,:,1), 'V', 1. ) ; CALL lbc_lnk( zlv_vu(:,:,:,1), 'V', 1. ) 139 134 CALL lbc_lnk( zlv_vv(:,:,:,2), 'V', 1. ) ; CALL lbc_lnk( zlv_vu(:,:,:,2), 'V', 1. ) 140 !!gm end 141 135 ! 142 136 ! ! ====================== ! 143 137 ! ! Horizontal advection ! 144 138 DO jk = 1, jpkm1 ! ====================== ! 145 139 ! ! horizontal volume fluxes 146 zfu(:,:,jk) = 0.25 * e2u(:,:) * fse3u(:,:,jk) * un(:,:,jk)147 zfv(:,:,jk) = 0.25 * e1v(:,:) * fse3v(:,:,jk) * vn(:,:,jk)140 zfu(:,:,jk) = 0.25_wp * e2u(:,:) * e3u_n(:,:,jk) * un(:,:,jk) 141 zfv(:,:,jk) = 0.25_wp * e1v(:,:) * e3v_n(:,:,jk) * vn(:,:,jk) 148 142 ! 149 143 DO jj = 1, jpjm1 ! horizontal momentum fluxes at T- and F-point … … 152 146 zvj = ( vn(ji,jj,jk) + vn(ji ,jj+1,jk) ) 153 147 ! 154 IF (zui > 0) THEN ; zl_u = zlu_uu(ji ,jj,jk,1)155 ELSE ; zl_u = zlu_uu(ji+1,jj,jk,1)156 ENDIF 157 IF (zvj > 0) THEN ; zl_v = zlv_vv(ji,jj ,jk,1)158 ELSE ; zl_v = zlv_vv(ji,jj+1,jk,1)148 IF( zui > 0 ) THEN ; zl_u = zlu_uu(ji ,jj,jk,1) 149 ELSE ; zl_u = zlu_uu(ji+1,jj,jk,1) 150 ENDIF 151 IF( zvj > 0 ) THEN ; zl_v = zlv_vv(ji,jj ,jk,1) 152 ELSE ; zl_v = zlv_vv(ji,jj+1,jk,1) 159 153 ENDIF 160 154 ! … … 168 162 zfuj = ( zfu(ji,jj,jk) + zfu(ji ,jj+1,jk) ) 169 163 zfvi = ( zfv(ji,jj,jk) + zfv(ji+1,jj ,jk) ) 170 IF (zfuj > 0) THEN ; zl_v = zlv_vu( ji ,jj ,jk,1)171 ELSE ; zl_v = zlv_vu( ji+1,jj,jk,1)172 ENDIF 173 IF (zfvi > 0) THEN ; zl_u = zlu_uv( ji,jj ,jk,1)174 ELSE ; zl_u = zlu_uv( ji,jj+1,jk,1)164 IF( zfuj > 0 ) THEN ; zl_v = zlv_vu( ji ,jj ,jk,1) 165 ELSE ; zl_v = zlv_vu( ji+1,jj,jk,1) 166 ENDIF 167 IF( zfvi > 0 ) THEN ; zl_u = zlu_uv( ji,jj ,jk,1) 168 ELSE ; zl_u = zlu_uv( ji,jj+1,jk,1) 175 169 ENDIF 176 170 ! … … 183 177 DO jj = 2, jpjm1 ! divergence of horizontal momentum fluxes 184 178 DO ji = fs_2, fs_jpim1 ! vector opt. 185 zbu = e1u(ji,jj) * e2u(ji,jj) * fse3u(ji,jj,jk) 186 zbv = e1v(ji,jj) * e2v(ji,jj) * fse3v(ji,jj,jk) 187 ! 188 ua(ji,jj,jk) = ua(ji,jj,jk) - ( zfu_t(ji+1,jj ,jk) - zfu_t(ji ,jj ,jk) & 189 & + zfv_f(ji ,jj ,jk) - zfv_f(ji ,jj-1,jk) ) / zbu 190 va(ji,jj,jk) = va(ji,jj,jk) - ( zfu_f(ji ,jj ,jk) - zfu_f(ji-1,jj ,jk) & 191 & + zfv_t(ji ,jj+1,jk) - zfv_t(ji ,jj ,jk) ) / zbv 192 END DO 193 END DO 194 END DO 195 IF( l_trddyn ) THEN ! save the horizontal advection trend for diagnostic 179 ua(ji,jj,jk) = ua(ji,jj,jk) - ( zfu_t(ji+1,jj,jk) - zfu_t(ji,jj ,jk) & 180 & + zfv_f(ji ,jj,jk) - zfv_f(ji,jj-1,jk) ) * r1_e1e2u(ji,jj) / e3u_n(ji,jj,jk) 181 va(ji,jj,jk) = va(ji,jj,jk) - ( zfu_f(ji,jj ,jk) - zfu_f(ji-1,jj,jk) & 182 & + zfv_t(ji,jj+1,jk) - zfv_t(ji ,jj,jk) ) * r1_e1e2v(ji,jj) / e3v_n(ji,jj,jk) 183 END DO 184 END DO 185 END DO 186 IF( l_trddyn ) THEN ! trends: send trends to trddyn for diagnostic 196 187 zfu_uw(:,:,:) = ua(:,:,:) - zfu_uw(:,:,:) 197 188 zfv_vw(:,:,:) = va(:,:,:) - zfv_vw(:,:,:) 198 CALL trd_ mod( zfu_uw, zfv_vw, jpdyn_trd_had, 'DYN', kt )189 CALL trd_dyn( zfu_uw, zfv_vw, jpdyn_keg, kt ) 199 190 zfu_t(:,:,:) = ua(:,:,:) 200 191 zfv_t(:,:,:) = va(:,:,:) 201 192 ENDIF 202 203 193 ! ! ==================== ! 204 194 ! ! Vertical advection ! 205 DO jk = 1, jpkm1! ==================== !206 ! ! Vertical volume fluxesÊ207 zfw(:,:,jk) = 0.25 * e1t(:,:) * e2t(:,:) * wn(:,:,jk)208 !209 IF( jk == 1 ) THEN ! surface/bottom advective fluxes210 zfu_uw( :,:,jpk) = 0.e0 ! Bottom value : flux set to zero211 zfv_vw( :,:,jpk) = 0.e0212 ! ! Surface value :213 IF( lk_vvl ) THEN ! variable volume : flux set to zero214 zfu_uw(:,:, 1 ) = 0.e0215 zfv_vw(:,:, 1 ) = 0.e0216 ELSE ! constant volume : advection through the surface217 DO jj = 2, jpjm1218 DO ji = fs_2, fs_jpim1219 zfu_uw(ji,jj, 1 ) = 2.e0 * ( zfw(ji,jj,1) + zfw(ji+1,jj ,1) ) * un(ji,jj,1)220 zfv_vw(ji,jj, 1 ) = 2.e0 * ( zfw(ji,jj,1) + zfw(ji ,jj+1,1) ) * vn(ji,jj,1)221 END DO222 END DO223 ENDIF224 ELSE ! interior fluxes225 DO jj = 2, jpjm1226 DO ji = fs_2, fs_jpim1 ! vector opt.227 zfu_uw(ji,jj,jk) = ( zfw(ji,jj,jk)+ zfw(ji+1,jj ,jk) ) * ( un(ji,jj,jk) + un(ji,jj,jk-1) )228 zfv_vw(ji,jj,jk) = ( zfw(ji,jj,jk)+ zfw(ji ,jj+1,jk) ) * ( vn(ji,jj,jk) + vn(ji,jj,jk-1) )229 END DO230 END DO231 ENDIF232 END DO233 DO jk = 1, jpkm1 ! divergence of vertical momentum flux divergence234 DO jj = 2, jpjm1235 DO ji = fs_2, fs_jpim1 ! vector opt.236 ua(ji,jj,jk) = ua(ji,jj,jk) - ( zfu_uw(ji,jj,jk) - zfu_uw(ji,jj,jk+1) ) &237 & / ( e1u(ji,jj) * e2u(ji,jj) * fse3u(ji,jj,jk) )238 va(ji,jj,jk) = va(ji,jj,jk) - ( zfv_vw(ji,jj,jk) - zfv_vw(ji,jj,jk+1) ) &239 & / ( e1v(ji,jj) * e2v(ji,jj) * fse3v(ji,jj,jk))240 END DO 241 END DO 242 END DO 243 ! 244 IF( l_trddyn ) THEN 195 ! ! ==================== ! 196 DO jj = 2, jpjm1 ! surface/bottom advective fluxes set to zero 197 DO ji = fs_2, fs_jpim1 198 zfu_uw(ji,jj,jpk) = 0._wp 199 zfv_vw(ji,jj,jpk) = 0._wp 200 zfu_uw(ji,jj, 1 ) = 0._wp 201 zfv_vw(ji,jj, 1 ) = 0._wp 202 END DO 203 END DO 204 IF( ln_linssh ) THEN ! constant volume : advection through the surface 205 DO jj = 2, jpjm1 206 DO ji = fs_2, fs_jpim1 207 zfu_uw(ji,jj,1) = 0.5_wp * ( e1e2t(ji,jj) * wn(ji,jj,1) + e1e2t(ji+1,jj) * wn(ji+1,jj,1) ) * un(ji,jj,1) 208 zfv_vw(ji,jj,1) = 0.5_wp * ( e1e2t(ji,jj) * wn(ji,jj,1) + e1e2t(ji,jj+1) * wn(ji,jj+1,1) ) * vn(ji,jj,1) 209 END DO 210 END DO 211 ENDIF 212 DO jk = 2, jpkm1 ! interior fluxes 213 DO jj = 2, jpjm1 214 DO ji = fs_2, fs_jpim1 215 zfw(ji,jj,jk) = 0.25_wp * e1e2t(ji,jj) * wn(ji,jj,jk) 216 END DO 217 END DO 218 DO jj = 2, jpjm1 219 DO ji = fs_2, fs_jpim1 ! vector opt. 220 zfu_uw(ji,jj,jk) = ( zfw(ji,jj,jk)+ zfw(ji+1,jj,jk) ) * ( un(ji,jj,jk) + un(ji,jj,jk-1) ) 221 zfv_vw(ji,jj,jk) = ( zfw(ji,jj,jk)+ zfw(ji,jj+1,jk) ) * ( vn(ji,jj,jk) + vn(ji,jj,jk-1) ) 222 END DO 223 END DO 224 END DO 225 DO jk = 1, jpkm1 ! divergence of vertical momentum flux divergence 226 DO jj = 2, jpjm1 227 DO ji = fs_2, fs_jpim1 ! vector opt. 228 ua(ji,jj,jk) = ua(ji,jj,jk) - ( zfu_uw(ji,jj,jk) - zfu_uw(ji,jj,jk+1) ) * r1_e1e2u(ji,jj) / e3u_n(ji,jj,jk) 229 va(ji,jj,jk) = va(ji,jj,jk) - ( zfv_vw(ji,jj,jk) - zfv_vw(ji,jj,jk+1) ) * r1_e1e2v(ji,jj) / e3v_n(ji,jj,jk) 230 END DO 231 END DO 232 END DO 233 ! 234 IF( l_trddyn ) THEN ! save the vertical advection trend for diagnostic 245 235 zfu_t(:,:,:) = ua(:,:,:) - zfu_t(:,:,:) 246 236 zfv_t(:,:,:) = va(:,:,:) - zfv_t(:,:,:) 247 CALL trd_ mod( zfu_t, zfv_t, jpdyn_trd_zad, 'DYN', kt )248 ENDIF 249 ! 237 CALL trd_dyn( zfu_t, zfv_t, jpdyn_zad, kt ) 238 ENDIF 239 ! ! Control print 250 240 IF(ln_ctl) CALL prt_ctl( tab3d_1=ua, clinfo1=' ubs2 adv - Ua: ', mask1=umask, & 251 241 & tab3d_2=va, clinfo2= ' Va: ', mask2=vmask, clinfo3='dyn' ) 252 242 ! 253 CALL wrk_dealloc( jpi, jpj, jpk,zfu_t , zfv_t , zfu_f , zfv_f, zfu_uw, zfv_vw, zfu, zfv, zfw )254 CALL wrk_dealloc( jpi, jpj, jpk, jpts,zlu_uu, zlv_vv, zlu_uv, zlv_vu )243 CALL wrk_dealloc( jpi,jpj,jpk, zfu_t , zfv_t , zfu_f , zfv_f, zfu_uw, zfv_vw, zfu, zfv, zfw ) 244 CALL wrk_dealloc( jpi,jpj,jpk,jpts, zlu_uu, zlv_vv, zlu_uv, zlv_vu ) 255 245 ! 256 246 IF( nn_timing == 1 ) CALL timing_stop('dyn_adv_ubs')
Note: See TracChangeset
for help on using the changeset viewer.
