- Timestamp:
- 2015-01-20T15:26:13+01:00 (9 years ago)
- File:
-
- 1 edited
Legend:
- Unmodified
- Added
- Removed
-
branches/2014/dev_r4621_NOC4_BDY_VERT_INTERP/NEMOGCM/NEMO/OPA_SRC/DYN/dynnxt.F90
r4370 r5038 18 18 !! 3.3 ! 2011-03 (P. Oddo) Bug fix for time-splitting+(BDY-OBC) and not VVL 19 19 !! 3.5 ! 2013-07 (J. Chanut) Compliant with time splitting changes 20 !! 3.7 ! 2014-04 (G. Madec) add the diagnostic of the time filter trends 20 21 !!------------------------------------------------------------------------- 21 22 … … 34 35 USE bdydyn ! ocean open boundary conditions 35 36 USE bdyvol ! ocean open boundary condition (bdy_vol routines) 37 USE trd_oce ! trends: ocean variables 38 USE trddyn ! trend manager: dynamics 39 USE trdken ! trend manager: kinetic energy 40 ! 36 41 USE in_out_manager ! I/O manager 42 USE iom ! I/O manager library 37 43 USE lbclnk ! lateral boundary condition (or mpp link) 38 44 USE lib_mpp ! MPP library 39 45 USE wrk_nemo ! Memory Allocation 40 46 USE prtctl ! Print control 41 47 USE timing ! Timing 42 48 #if defined key_agrif 43 49 USE agrif_opa_interp 44 50 #endif 45 USE timing ! Timing46 51 47 52 IMPLICIT NONE … … 79 84 !! at the local domain boundaries through lbc_lnk call, 80 85 !! at the one-way open boundaries (lk_bdy=T), 81 !! at the AGRIF zoom 86 !! at the AGRIF zoom boundaries (lk_agrif=T) 82 87 !! 83 88 !! * Apply the time filter applied and swap of the dynamics … … 99 104 REAL(wp) :: z2dt ! temporary scalar 100 105 #endif 101 REAL(wp) :: zue3a, zue3n, zue3b, zuf, zec ! local scalars102 REAL(wp) :: zve3a, zve3n, zve3b, zvf 103 REAL(wp), POINTER, DIMENSION(:,:) :: zu a, zva104 REAL(wp), POINTER, DIMENSION(:,:,:) :: ze3u_f, ze3v_f 106 REAL(wp) :: zue3a, zue3n, zue3b, zuf, zec ! local scalars 107 REAL(wp) :: zve3a, zve3n, zve3b, zvf, z1_2dt ! - - 108 REAL(wp), POINTER, DIMENSION(:,:) :: zue, zve 109 REAL(wp), POINTER, DIMENSION(:,:,:) :: ze3u_f, ze3v_f, zua, zva 105 110 !!---------------------------------------------------------------------- 106 111 ! 107 IF( nn_timing == 1 ) CALL timing_start('dyn_nxt')108 ! 109 CALL wrk_alloc( jpi,jpj,jpk, ze3u_f, ze3v_f)110 IF ( lk_dynspg_ts ) CALL wrk_alloc( jpi,jpj, zua, zva)112 IF( nn_timing == 1 ) CALL timing_start('dyn_nxt') 113 ! 114 CALL wrk_alloc( jpi,jpj,jpk, ze3u_f, ze3v_f, zua, zva ) 115 IF( lk_dynspg_ts ) CALL wrk_alloc( jpi,jpj, zue, zve ) 111 116 ! 112 117 IF( kt == nit000 ) THEN … … 152 157 153 158 # if defined key_dynspg_ts 159 !!gm IF ( lk_dynspg_ts ) THEN .... 154 160 ! Ensure below that barotropic velocities match time splitting estimate 155 161 ! Compute actual transport and replace it with ts estimate at "after" time step 156 zu a(:,:) = 0._wp157 zv a(:,:) = 0._wp158 DO jk = 1, jpkm1159 zu a(:,:) = zua(:,:) + fse3u_a(:,:,jk) * ua(:,:,jk) * umask(:,:,jk)160 zv a(:,:) = zva(:,:) + fse3v_a(:,:,jk) * va(:,:,jk) * vmask(:,:,jk)162 zue(:,:) = fse3u_a(:,:,1) * ua(:,:,1) * umask(:,:,1) 163 zve(:,:) = fse3v_a(:,:,1) * va(:,:,1) * vmask(:,:,1) 164 DO jk = 2, jpkm1 165 zue(:,:) = zue(:,:) + fse3u_a(:,:,jk) * ua(:,:,jk) * umask(:,:,jk) 166 zve(:,:) = zve(:,:) + fse3v_a(:,:,jk) * va(:,:,jk) * vmask(:,:,jk) 161 167 END DO 162 168 DO jk = 1, jpkm1 163 ua(:,:,jk) = ( ua(:,:,jk) - zu a(:,:) * hur_a(:,:) + ua_b(:,:) ) * umask(:,:,jk)164 va(:,:,jk) = ( va(:,:,jk) - zv a(:,:) * hvr_a(:,:) + va_b(:,:) ) * vmask(:,:,jk)169 ua(:,:,jk) = ( ua(:,:,jk) - zue(:,:) * hur_a(:,:) + ua_b(:,:) ) * umask(:,:,jk) 170 va(:,:,jk) = ( va(:,:,jk) - zve(:,:) * hvr_a(:,:) + va_b(:,:) ) * vmask(:,:,jk) 165 171 END DO 166 172 … … 175 181 END DO 176 182 ENDIF 183 !!gm ENDIF 177 184 # endif 178 185 … … 195 202 # endif 196 203 #endif 204 205 IF( l_trddyn ) THEN ! prepare the atf trend computation + some diagnostics 206 z1_2dt = 1._wp / (2. * rdt) ! Euler or leap-frog time step 207 IF( neuler == 0 .AND. kt == nit000 ) z1_2dt = 1._wp / rdt 208 ! 209 ! ! Kinetic energy and Conversion 210 IF( ln_KE_trd ) CALL trd_dyn( ua, va, jpdyn_ken, kt ) 211 ! 212 IF( ln_dyn_trd ) THEN ! 3D output: total momentum trends 213 zua(:,:,:) = ( ua(:,:,:) - ub(:,:,:) ) * z1_2dt 214 zva(:,:,:) = ( va(:,:,:) - vb(:,:,:) ) * z1_2dt 215 CALL iom_put( "utrd_tot", zua ) ! total momentum trends, except the asselin time filter 216 CALL iom_put( "vtrd_tot", zva ) 217 ENDIF 218 ! 219 zua(:,:,:) = un(:,:,:) ! save the now velocity before the asselin filter 220 zva(:,:,:) = vn(:,:,:) ! (caution: there will be a shift by 1 timestep in the 221 ! ! computation of the asselin filter trends) 222 ENDIF 197 223 198 224 ! Time filter and swap of dynamics arrays … … 217 243 DO jj = 1, jpj 218 244 DO ji = 1, jpi 219 zuf = un(ji,jj,jk) + atfp * ( ub(ji,jj,jk) - 2. e0_wp * un(ji,jj,jk) + ua(ji,jj,jk) )220 zvf = vn(ji,jj,jk) + atfp * ( vb(ji,jj,jk) - 2. e0_wp * vn(ji,jj,jk) + va(ji,jj,jk) )245 zuf = un(ji,jj,jk) + atfp * ( ub(ji,jj,jk) - 2._wp * un(ji,jj,jk) + ua(ji,jj,jk) ) 246 zvf = vn(ji,jj,jk) + atfp * ( vb(ji,jj,jk) - 2._wp * vn(ji,jj,jk) + va(ji,jj,jk) ) 221 247 ! 222 248 ub(ji,jj,jk) = zuf ! ub <-- filtered velocity … … 301 327 ! Revert "before" velocities to time split estimate 302 328 ! Doing it here also means that asselin filter contribution is removed 303 zu a(:,:) = 0._wp304 zv a(:,:) = 0._wp305 DO jk = 1, jpkm1306 zu a(:,:) = zua(:,:) + fse3u_b(:,:,jk) * ub(:,:,jk) * umask(:,:,jk)307 zv a(:,:) = zva(:,:) + fse3v_b(:,:,jk) * vb(:,:,jk) * vmask(:,:,jk)329 zue(:,:) = fse3u_b(:,:,1) * ub(:,:,1) * umask(:,:,1) 330 zve(:,:) = fse3v_b(:,:,1) * vb(:,:,1) * vmask(:,:,1) 331 DO jk = 2, jpkm1 332 zue(:,:) = zue(:,:) + fse3u_b(:,:,jk) * ub(:,:,jk) * umask(:,:,jk) 333 zve(:,:) = zve(:,:) + fse3v_b(:,:,jk) * vb(:,:,jk) * vmask(:,:,jk) 308 334 END DO 309 335 DO jk = 1, jpkm1 310 ub(:,:,jk) = ub(:,:,jk) - (zu a(:,:) * hur(:,:) - un_b(:,:)) * umask(:,:,jk)311 vb(:,:,jk) = vb(:,:,jk) - (zv a(:,:) * hvr(:,:) - vn_b(:,:)) * vmask(:,:,jk)336 ub(:,:,jk) = ub(:,:,jk) - (zue(:,:) * hur(:,:) - un_b(:,:)) * umask(:,:,jk) 337 vb(:,:,jk) = vb(:,:,jk) - (zve(:,:) * hvr(:,:) - vn_b(:,:)) * vmask(:,:,jk) 312 338 END DO 313 339 ENDIF … … 327 353 hv_b(:,:) = hv_b(:,:) + fse3v_b(:,:,jk) * vmask(:,:,jk) 328 354 END DO 329 hur_b(:,:) = umask (:,:,1) / ( hu_b(:,:) + 1. - umask(:,:,1) )330 hvr_b(:,:) = vmask (:,:,1) / ( hv_b(:,:) + 1. - vmask(:,:,1) )355 hur_b(:,:) = umask_i(:,:) / ( hu_b(:,:) + 1._wp - umask_i(:,:) ) 356 hvr_b(:,:) = vmask_i(:,:) / ( hv_b(:,:) + 1._wp - vmask_i(:,:) ) 331 357 ENDIF 332 358 ! … … 335 361 ! 336 362 DO jk = 1, jpkm1 337 #if defined key_vectopt_loop338 DO jj = 1, 1 !Vector opt. => forced unrolling339 DO ji = 1, jpij340 #else341 363 DO jj = 1, jpj 342 364 DO ji = 1, jpi 343 #endif344 365 un_b(ji,jj) = un_b(ji,jj) + fse3u_a(ji,jj,jk) * un(ji,jj,jk) * umask(ji,jj,jk) 345 366 vn_b(ji,jj) = vn_b(ji,jj) + fse3v_a(ji,jj,jk) * vn(ji,jj,jk) * vmask(ji,jj,jk) … … 358 379 ! 359 380 ! 381 382 IF( l_trddyn ) THEN ! 3D output: asselin filter trends on momentum 383 zua(:,:,:) = ( ub(:,:,:) - zua(:,:,:) ) * z1_2dt 384 zva(:,:,:) = ( vb(:,:,:) - zva(:,:,:) ) * z1_2dt 385 CALL trd_dyn( zua, zva, jpdyn_atf, kt ) 386 ENDIF 387 ! 360 388 IF(ln_ctl) CALL prt_ctl( tab3d_1=un, clinfo1=' nxt - Un: ', mask1=umask, & 361 389 & tab3d_2=vn, clinfo2=' Vn: ' , mask2=vmask ) 362 390 ! 363 CALL wrk_dealloc( jpi,jpj,jpk, ze3u_f, ze3v_f)364 IF ( lk_dynspg_ts ) CALL wrk_dealloc( jpi,jpj, zua, zva)391 CALL wrk_dealloc( jpi,jpj,jpk, ze3u_f, ze3v_f, zua, zva ) 392 IF( lk_dynspg_ts ) CALL wrk_dealloc( jpi,jpj, zue, zve ) 365 393 ! 366 394 IF( nn_timing == 1 ) CALL timing_stop('dyn_nxt') … … 370 398 !!========================================================================= 371 399 END MODULE dynnxt 372
Note: See TracChangeset
for help on using the changeset viewer.