Changeset 4246
- Timestamp:
- 2013-11-19T12:46:44+01:00 (10 years ago)
- Location:
- branches/2013/dev_r3858_NOC_ZTC/NEMOGCM/NEMO/OPA_SRC
- Files:
-
- 4 edited
Legend:
- Unmodified
- Added
- Removed
-
branches/2013/dev_r3858_NOC_ZTC/NEMOGCM/NEMO/OPA_SRC/DYN/dynspg_ts.F90
r4241 r4246 901 901 ENDIF 902 902 ! 903 DO jk = 1, jpkm1 904 ! Correct velocities: 905 un(:,:,jk) = ( un(:,:,jk) + un_adv(:,:) - un_b(:,:) )*umask(:,:,jk) 906 vn(:,:,jk) = ( vn(:,:,jk) + vn_adv(:,:) - vn_b(:,:) )*vmask(:,:,jk) 907 ! 908 END DO 909 ! 903 910 ! !* write time-spliting arrays in the restart 904 911 IF(lrst_oce .AND.ln_bt_fw) CALL ts_rst( kt, 'WRITE' ) -
branches/2013/dev_r3858_NOC_ZTC/NEMOGCM/NEMO/OPA_SRC/DYN/sshwzv.F90
r4178 r4246 50 50 51 51 PUBLIC ssh_nxt ! called by step.F90 52 PUBLIC wzv_1 ! called by step.F90 53 PUBLIC wzv_2 ! called by step.F90 52 PUBLIC wzv ! called by step.F90 54 53 PUBLIC ssh_swp ! called by step.F90 55 54 … … 152 151 153 152 154 SUBROUTINE wzv _1( kt )155 !!---------------------------------------------------------------------- 156 !! *** ROUTINE wzv _1***153 SUBROUTINE wzv( kt ) 154 !!---------------------------------------------------------------------- 155 !! *** ROUTINE wzv *** 157 156 !! 158 157 !! ** Purpose : compute the now vertical velocity … … 169 168 ! 170 169 INTEGER, INTENT(in) :: kt ! time step 171 REAL(wp), POINTER, DIMENSION(:,:,:) :: zhdiv 170 REAL(wp), POINTER, DIMENSION(:,: ) :: z2d 171 REAL(wp), POINTER, DIMENSION(:,:,:) :: z3d, zhdiv 172 172 ! 173 173 INTEGER :: ji, jj, jk ! dummy loop indices … … 175 175 !!---------------------------------------------------------------------- 176 176 177 IF( nn_timing == 1 ) CALL timing_start('wzv_1') 178 ! 179 CALL wrk_alloc( jpi, jpj, jpk, zhdiv ) 177 IF( nn_timing == 1 ) CALL timing_start('wzv') 180 178 ! 181 179 IF( kt == nit000 ) THEN 182 180 ! 183 181 IF(lwp) WRITE(numout,*) 184 IF(lwp) WRITE(numout,*) 'wzv _1: now vertical velocity '182 IF(lwp) WRITE(numout,*) 'wzv : now vertical velocity ' 185 183 IF(lwp) WRITE(numout,*) '~~~~~ ' 186 184 ! … … 195 193 ! 196 194 IF( ln_vvl_ztilde .OR. ln_vvl_layer ) THEN ! z_tilde and layer cases 195 CALL wrk_alloc( jpi, jpj, jpk, zhdiv ) 196 ! 197 197 DO jk = 1, jpkm1 198 198 ! horizontal divergence of thickness diffusion transport ( velocity multiplied by e3t) … … 213 213 END DO 214 214 ! IF( ln_vvl_layer ) wn(:,:,:) = 0.e0 215 CALL wrk_dealloc( jpi, jpj, jpk, zhdiv ) 215 216 ELSE ! z_star and linear free surface cases 216 217 DO jk = jpkm1, 1, -1 ! integrate from the bottom the hor. divergence … … 224 225 wn(:,:,jk) = wn(:,:,jk) * bdytmask(:,:) 225 226 #endif 226 227 ! 228 CALL wrk_dealloc( jpi, jpj, jpk, zhdiv ) 229 ! 230 IF( nn_timing == 1 ) CALL timing_stop('wzv_1') 231 232 233 END SUBROUTINE wzv_1 234 235 SUBROUTINE wzv_2( kt ) 236 !!---------------------------------------------------------------------- 237 !! *** ROUTINE wzv_2 *** 238 !! 239 !! ** Purpose : compute the now vertical velocity 240 !! 241 !! ** Method : - Using the incompressibility hypothesis, the vertical 242 !! velocity is computed by integrating the horizontal divergence 243 !! from the bottom to the surface minus the scale factor evolution. 244 !! The boundary conditions are w=0 at the bottom (no flux) and. 245 !! 246 !! ** action : wn : now vertical velocity 247 !! 248 !! Reference : Leclair, M., and G. Madec, 2009, Ocean Modelling. 249 !!---------------------------------------------------------------------- 250 ! 251 INTEGER, INTENT(in) :: kt ! time step 252 REAL(wp), POINTER, DIMENSION(:,: ) :: z2d 253 REAL(wp), POINTER, DIMENSION(:,:,:) :: z3d, zhdiv 254 ! 255 INTEGER :: ji, jj, jk ! dummy loop indices 256 REAL(wp) :: z1_2dt ! local scalars 257 !!---------------------------------------------------------------------- 258 259 IF( nn_timing == 1 ) CALL timing_start('wzv_2') 260 ! 261 CALL wrk_alloc( jpi, jpj, jpk, z3d, zhdiv ) 262 ! 263 IF( kt == nit000 ) THEN 264 ! 265 IF(lwp) WRITE(numout,*) 266 IF(lwp) WRITE(numout,*) 'wzv_2 : now vertical velocity ' 267 IF(lwp) WRITE(numout,*) '~~~~~ ' 268 ! 269 ENDIF 270 ! !------------------------------! 271 ! ! Now Vertical Velocity ! 272 ! !------------------------------! 273 z1_2dt = 1. / ( 2. * rdt ) ! set time step size (Euler/Leapfrog) 274 IF( neuler == 0 .AND. kt == nit000 ) z1_2dt = 1. / rdt 275 276 IF( lk_dynspg_ts ) THEN 277 ! At this stage: 278 ! 1) vertical scale factors have been updated. 279 ! 2) Time averaged barotropic velocities around step kt+1 (ua_b, va_b) as well as 280 ! "advective" barotropic velocities (un_adv, vn_adv) at step kt in agreement 281 ! with continuity equation are available. 282 ! 3) 3d velocities (ua, va) hold ua_b, va_b issued from time splitting as barotropic components. 283 ! Below => Update "now" velocities, divergence, then vertical velocity 284 ! NB: hdivn is NOT updated such that hdivb is not updated also. This means that horizontal 285 ! momentum diffusion is still performed on Instantaneous barotropic arrays. I experencied 286 ! some issues with UBS with the current method. See "divup" comments in development branch to 287 ! update the divergence fully if necessary (2013/dev_r3867_MERCATOR1_DYN). 288 ! 289 DO jk = 1, jpkm1 290 ! Correct velocities: 291 un(:,:,jk) = ( un(:,:,jk) + un_adv(:,:) - un_b(:,:) )*umask(:,:,jk) 292 vn(:,:,jk) = ( vn(:,:,jk) + vn_adv(:,:) - vn_b(:,:) )*vmask(:,:,jk) 293 ! 294 ! Compute divergence: 295 DO jj = 2, jpjm1 296 DO ji = fs_2, fs_jpim1 ! vector opt. 297 z3d(ji,jj,jk) = & 298 ( e2u(ji,jj)*fse3u(ji,jj,jk) * un(ji,jj,jk) - e2u(ji-1,jj)*fse3u(ji-1,jj,jk) * un(ji-1,jj,jk) & 299 + e1v(ji,jj)*fse3v(ji,jj,jk) * vn(ji,jj,jk) - e1v(ji,jj-1)*fse3v(ji,jj-1,jk) * vn(ji,jj-1,jk) ) & 300 / ( e1t(ji,jj) * e2t(ji,jj) * fse3t(ji,jj,jk) ) 301 END DO 302 END DO 303 END DO 304 ! 305 IF( ln_rnf ) CALL sbc_rnf_div( z3d ) ! runoffs (update divergence) 306 ELSE 307 z3d(:,:,:) = hdivn(:,:,:) 308 ENDIF 309 ! 310 IF( ln_vvl_ztilde .OR. ln_vvl_layer ) THEN ! z_tilde and layer cases 311 DO jk = 1, jpkm1 312 ! horizontal divergence of thickness diffusion transport ( velocity multiplied by e3t) 313 ! - ML - note: computation allready done in dom_vvl_sf_nxt. Could be optimized (not critical and clearer this way) 314 DO jj = 2, jpjm1 315 DO ji = fs_2, fs_jpim1 ! vector opt. 316 zhdiv(ji,jj,jk) = r1_e12t(ji,jj) * ( un_td(ji,jj,jk) - un_td(ji-1,jj,jk) + vn_td(ji,jj,jk) - vn_td(ji,jj-1,jk) ) 317 END DO 318 END DO 319 END DO 320 CALL lbc_lnk(zhdiv, 'T', 1.) ! - ML - Perhaps not necessary: not used for horizontal "connexions" 321 ! ! Is it problematic to have a wrong vertical velocity in boundary cells? 322 ! ! Same question holds for hdivn. Perhaps just for security 323 DO jk = jpkm1, 1, -1 ! integrate from the bottom the hor. divergence 324 ! computation of w 325 wn(:,:,jk) = wn(:,:,jk+1) - ( fse3t_n(:,:,jk) * z3d(:,:,jk) + zhdiv(:,:,jk) & 326 & + z1_2dt * ( fse3t_a(:,:,jk) - fse3t_b(:,:,jk) ) ) * tmask(:,:,jk) 327 END DO 328 ! IF( ln_vvl_layer ) wn(:,:,:) = 0.e0 329 ELSE ! z_star and linear free surface cases 330 DO jk = jpkm1, 1, -1 ! integrate from the bottom the hor. divergence 331 ! computation of w 332 wn(:,:,jk) = wn(:,:,jk+1) - ( fse3t_n(:,:,jk) * z3d(:,:,jk) & 333 & + z1_2dt * ( fse3t_a(:,:,jk) - fse3t_b(:,:,jk) ) ) * tmask(:,:,jk) 334 END DO 335 ENDIF 336 337 #if defined key_bdy 338 wn(:,:,jk) = wn(:,:,jk) * bdytmask(:,:) 339 #endif 340 227 ! 341 228 ! !------------------------------! 342 229 ! ! outputs ! … … 345 232 IF( lk_diaar5 ) THEN ! vertical mass transport & its square value 346 233 CALL wrk_alloc( jpi, jpj, z2d ) 234 CALL wrk_alloc( jpi, jpj, jpk, z3d ) 347 235 ! Caution: in the VVL case, it only correponds to the baroclinic mass transport. 348 236 z2d(:,:) = rau0 * e12t(:,:) … … 353 241 CALL iom_put( "w_masstr2", z3d(:,:,:) * z3d(:,:,:) ) 354 242 CALL wrk_dealloc( jpi, jpj, z2d ) 355 ENDIF 356 ! 357 CALL wrk_dealloc( jpi, jpj, jpk, z3d, zhdiv ) 358 ! 359 IF( nn_timing == 1 ) CALL timing_stop('wzv_2') 360 361 362 END SUBROUTINE wzv_2 243 CALL wrk_dealloc( jpi, jpj, jpk, z3d ) 244 ENDIF 245 ! 246 IF( nn_timing == 1 ) CALL timing_stop('wzv') 247 248 249 END SUBROUTINE wzv 363 250 364 251 SUBROUTINE ssh_swp( kt ) -
branches/2013/dev_r3858_NOC_ZTC/NEMOGCM/NEMO/OPA_SRC/step.F90
r4178 r4246 98 98 ! Ocean dynamics : hdiv, rot, ssh, e3, wn 99 99 !<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<< 100 CALL ssh_nxt ( kstp ) ! after ssh 100 CALL ssh_nxt ( kstp ) ! after ssh (includes call to div_cur) 101 101 IF( lk_dynspg_ts ) THEN 102 CALL wzv _1( kstp ) ! now cross-level velocity102 CALL wzv ( kstp ) ! now cross-level velocity 103 103 ! In case the time splitting case, update almost all momentum trends here: 104 104 ! Note that the computation of vertical velocity above, hence "after" sea level … … 124 124 CALL dyn_hpg( kstp ) ! horizontal gradient of Hydrostatic pressure 125 125 CALL dyn_spg( kstp, indic ) ! surface pressure gradient 126 CALL div_cur( kstp ) ! Horizontal divergence & Relative vorticity (2nd call in time-split case) 126 127 127 128 ua_sv(:,:,:) = ua(:,:,:) ! save next velocities (not trends !) … … 129 130 ENDIF 130 131 IF( lk_vvl ) CALL dom_vvl_sf_nxt( kstp ) ! after vertical scale factors 131 CALL wzv _2( kstp ) ! now cross-level velocity (original)132 CALL wzv ( kstp ) ! now cross-level velocity (original) 132 133 133 134 !>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>> -
branches/2013/dev_r3858_NOC_ZTC/NEMOGCM/NEMO/OPA_SRC/step_oce.F90
r3865 r4246 11 11 USE ldftra_oce ! ocean tracer - trends 12 12 USE ldfdyn_oce ! ocean dynamics - trends 13 USE divcur ! hor. divergence and curl (div & cur routines) 13 14 USE in_out_manager ! I/O manager 14 15 USE iom !
Note: See TracChangeset
for help on using the changeset viewer.