Changeset 7753 for trunk/NEMOGCM/NEMO/OPA_SRC/DYN/dynnxt.F90
- Timestamp:
- 2017-03-03T12:46:59+01:00 (7 years ago)
- File:
-
- 1 edited
-
trunk/NEMOGCM/NEMO/OPA_SRC/DYN/dynnxt.F90 (modified) (12 diffs)
Legend:
- Unmodified
- Added
- Removed
-
trunk/NEMOGCM/NEMO/OPA_SRC/DYN/dynnxt.F90
r7698 r7753 115 115 ! Ensure below that barotropic velocities match time splitting estimate 116 116 ! Compute actual transport and replace it with ts estimate at "after" time step 117 !$OMP PARALLEL 118 !$OMP DO schedule(static) private(jj, ji) 119 DO jj = 1, jpj 120 DO ji = 1, jpi 121 zue(ji,jj) = e3u_a(ji,jj,1) * ua(ji,jj,1) * umask(ji,jj,1) 122 zve(ji,jj) = e3v_a(ji,jj,1) * va(ji,jj,1) * vmask(ji,jj,1) 123 END DO 117 zue(:,:) = e3u_a(:,:,1) * ua(:,:,1) * umask(:,:,1) 118 zve(:,:) = e3v_a(:,:,1) * va(:,:,1) * vmask(:,:,1) 119 DO jk = 2, jpkm1 120 zue(:,:) = zue(:,:) + e3u_a(:,:,jk) * ua(:,:,jk) * umask(:,:,jk) 121 zve(:,:) = zve(:,:) + e3v_a(:,:,jk) * va(:,:,jk) * vmask(:,:,jk) 124 122 END DO 125 DO jk = 2, jpkm1 126 !$OMP DO schedule(static) private(jj,ji) 127 DO jj = 1, jpj 128 DO ji = 1, jpi 129 zue(ji,jj) = zue(ji,jj) + e3u_a(ji,jj,jk) * ua(ji,jj,jk) * umask(ji,jj,jk) 130 zve(ji,jj) = zve(ji,jj) + e3v_a(ji,jj,jk) * va(ji,jj,jk) * vmask(ji,jj,jk) 131 END DO 132 END DO 123 DO jk = 1, jpkm1 124 ua(:,:,jk) = ( ua(:,:,jk) - zue(:,:) * r1_hu_a(:,:) + ua_b(:,:) ) * umask(:,:,jk) 125 va(:,:,jk) = ( va(:,:,jk) - zve(:,:) * r1_hv_a(:,:) + va_b(:,:) ) * vmask(:,:,jk) 133 126 END DO 134 !$OMP DO schedule(static) private(jk,jj,ji)135 DO jk = 1, jpkm1136 DO jj = 1, jpj137 DO ji = 1, jpi138 ua(ji,jj,jk) = ( ua(ji,jj,jk) - zue(ji,jj) * r1_hu_a(ji,jj) + ua_b(ji,jj) ) * umask(ji,jj,jk)139 va(ji,jj,jk) = ( va(ji,jj,jk) - zve(ji,jj) * r1_hv_a(ji,jj) + va_b(ji,jj) ) * vmask(ji,jj,jk)140 END DO141 END DO142 END DO143 !$OMP END PARALLEL144 127 ! 145 128 IF( .NOT.ln_bt_fw ) THEN … … 148 131 ! In the forward case, this is done below after asselin filtering 149 132 ! so that asselin contribution is removed at the same time 150 !$OMP PARALLEL DO schedule(static) private(jk,jj,ji)151 133 DO jk = 1, jpkm1 152 DO jj = 1, jpj 153 DO ji = 1, jpi 154 un(ji,jj,jk) = ( un(ji,jj,jk) - un_adv(ji,jj) + un_b(ji,jj) )*umask(ji,jj,jk) 155 vn(ji,jj,jk) = ( vn(ji,jj,jk) - vn_adv(ji,jj) + vn_b(ji,jj) )*vmask(ji,jj,jk) 156 END DO 157 END DO 158 END DO 159 134 un(:,:,jk) = ( un(:,:,jk) - un_adv(:,:) + un_b(:,:) )*umask(:,:,jk) 135 vn(:,:,jk) = ( vn(:,:,jk) - vn_adv(:,:) + vn_b(:,:) )*vmask(:,:,jk) 136 END DO 160 137 ENDIF 161 138 ENDIF … … 184 161 ! 185 162 IF( ln_dyn_trd ) THEN ! 3D output: total momentum trends 186 !$OMP PARALLEL DO schedule(static) private(jk, jj, ji) 187 DO jk = 1, jpk 188 DO jj = 1, jpj 189 DO ji = 1, jpi 190 zua(ji,jj,jk) = ( ua(ji,jj,jk) - ub(ji,jj,jk) ) * z1_2dt 191 zva(ji,jj,jk) = ( va(ji,jj,jk) - vb(ji,jj,jk) ) * z1_2dt 192 END DO 193 END DO 194 END DO 163 zua(:,:,:) = ( ua(:,:,:) - ub(:,:,:) ) * z1_2dt 164 zva(:,:,:) = ( va(:,:,:) - vb(:,:,:) ) * z1_2dt 195 165 CALL iom_put( "utrd_tot", zua ) ! total momentum trends, except the asselin time filter 196 166 CALL iom_put( "vtrd_tot", zva ) 197 167 ENDIF 198 168 ! 199 !$OMP PARALLEL DO schedule(static) private(jk, jj, ji) 200 DO jk = 1, jpk 201 DO jj = 1, jpj 202 DO ji = 1, jpi 203 zua(ji,jj,jk) = un(ji,jj,jk) ! save the now velocity before the asselin filter 204 zva(ji,jj,jk) = vn(ji,jj,jk) ! (caution: there will be a shift by 1 timestep in the 205 ! ! computation of the asselin filter trends) 206 END DO 207 END DO 208 END DO 169 zua(:,:,:) = un(:,:,:) ! save the now velocity before the asselin filter 170 zva(:,:,:) = vn(:,:,:) ! (caution: there will be a shift by 1 timestep in the 171 ! ! computation of the asselin filter trends) 209 172 ENDIF 210 173 … … 212 175 ! ------------------------------------------ 213 176 IF( neuler == 0 .AND. kt == nit000 ) THEN !* Euler at first time-step: only swap 214 !$OMP PARALLEL215 !$OMP DO schedule(static) private(jk,jj,ji)216 177 DO jk = 1, jpkm1 217 DO jj = 1, jpj 218 DO ji = 1, jpi 219 un(ji,jj,jk) = ua(ji,jj,jk) ! un <-- ua 220 vn(ji,jj,jk) = va(ji,jj,jk) 221 END DO 178 un(:,:,jk) = ua(:,:,jk) ! un <-- ua 179 vn(:,:,jk) = va(:,:,jk) 180 END DO 181 IF(.NOT.ln_linssh ) THEN 182 DO jk = 1, jpkm1 183 e3t_b(:,:,jk) = e3t_n(:,:,jk) 184 e3u_b(:,:,jk) = e3u_n(:,:,jk) 185 e3v_b(:,:,jk) = e3v_n(:,:,jk) 222 186 END DO 223 END DO 224 !$OMP END DO NOWAIT 225 IF(.NOT.ln_linssh ) THEN 226 !$OMP DO schedule(static) private(jk,jj,ji) 227 DO jk = 1, jpkm1 228 DO jj = 1, jpj 229 DO ji = 1, jpi 230 e3t_b(ji,jj,jk) = e3t_n(ji,jj,jk) 231 e3u_b(ji,jj,jk) = e3u_n(ji,jj,jk) 232 e3v_b(ji,jj,jk) = e3v_n(ji,jj,jk) 233 END DO 234 END DO 235 END DO 236 ENDIF 237 !$OMP END PARALLEL 187 ENDIF 238 188 ELSE !* Leap-Frog : Asselin filter and swap 239 189 ! ! =============! 240 190 IF( ln_linssh ) THEN ! Fixed volume ! 241 191 ! ! =============! 242 !$OMP PARALLEL DO schedule(static) private(jk, jj, ji, zuf, zvf)243 192 DO jk = 1, jpkm1 244 193 DO jj = 1, jpj … … 261 210 ! ---------------------------------------------------- 262 211 IF( ln_dynspg_ts .AND. ln_bt_fw ) THEN ! No asselin filtering on thicknesses if forward time splitting 263 !$OMP PARALLEL DO schedule(static) private(jj,ji) 264 DO jj = 1, jpj 265 DO ji = 1, jpi 266 e3t_b(ji,jj,1:jpkm1) = e3t_n(ji,jj,1:jpkm1) 267 END DO 268 END DO 212 e3t_b(:,:,1:jpkm1) = e3t_n(:,:,1:jpkm1) 269 213 ELSE 270 !$OMP PARALLEL DO schedule(static) private(jk,jj,ji)271 214 DO jk = 1, jpkm1 272 DO jj = 1, jpj 273 DO ji = 1, jpi 274 e3t_b(ji,jj,jk) = e3t_n(ji,jj,jk) + atfp * ( e3t_b(ji,jj,jk) - 2._wp * e3t_n(ji,jj,jk) + e3t_a(ji,jj,jk) ) 275 END DO 276 END DO 215 e3t_b(:,:,jk) = e3t_n(:,:,jk) + atfp * ( e3t_b(:,:,jk) - 2._wp * e3t_n(:,:,jk) + e3t_a(:,:,jk) ) 277 216 END DO 278 217 ! Add volume filter correction: compatibility with tracer advection scheme … … 280 219 zcoef = atfp * rdt * r1_rau0 281 220 IF ( .NOT. ln_isf ) THEN ! if no ice shelf melting 282 !$OMP PARALLEL DO schedule(static) private(jj,ji) 283 DO jj = 1, jpj 284 DO ji = 1, jpi 285 e3t_b(ji,jj,1) = e3t_b(ji,jj,1) - zcoef * ( emp_b(ji,jj) - emp(ji,jj) & 286 & - rnf_b(ji,jj) + rnf(ji,jj) ) * tmask(ji,jj,1) 287 END DO 288 END DO 221 e3t_b(:,:,1) = e3t_b(:,:,1) - zcoef * ( emp_b(:,:) - emp(:,:) & 222 & - rnf_b(:,:) + rnf(:,:) ) * tmask(:,:,1) 289 223 ELSE ! if ice shelf melting 290 !$OMP PARALLEL DO schedule(static) private(jj,ji,ikt)291 224 DO jj = 1, jpj 292 225 DO ji = 1, jpi … … 304 237 CALL dom_vvl_interpol( e3t_b(:,:,:), e3u_b(:,:,:), 'U' ) 305 238 CALL dom_vvl_interpol( e3t_b(:,:,:), e3v_b(:,:,:), 'V' ) 306 !$OMP PARALLEL DO schedule(static) private(jk, jj, ji, zuf, zvf)307 239 DO jk = 1, jpkm1 308 240 DO jj = 1, jpj … … 325 257 CALL dom_vvl_interpol( e3t_b(:,:,:), ze3u_f, 'U' ) 326 258 CALL dom_vvl_interpol( e3t_b(:,:,:), ze3v_f, 'V' ) 327 !$OMP PARALLEL328 !$OMP DO schedule(static) private(jk, jj, ji, zue3a, zve3a, zue3n, zve3n, zue3b, zve3b, zuf, zvf)329 259 DO jk = 1, jpkm1 330 260 DO jj = 1, jpj … … 347 277 END DO 348 278 END DO 349 !$OMP DO schedule(static) private(jj, ji) 350 DO jj = 1, jpj 351 DO ji = 1, jpi 352 e3u_b(ji,jj,1:jpkm1) = ze3u_f(ji,jj,1:jpkm1) ! e3u_b <-- filtered scale factor 353 e3v_b(ji,jj,1:jpkm1) = ze3v_f(ji,jj,1:jpkm1) 354 END DO 355 END DO 356 !$OMP END PARALLEL 279 e3u_b(:,:,1:jpkm1) = ze3u_f(:,:,1:jpkm1) ! e3u_b <-- filtered scale factor 280 e3v_b(:,:,1:jpkm1) = ze3v_f(:,:,1:jpkm1) 357 281 ! 358 282 CALL wrk_dealloc( jpi,jpj,jpk, ze3u_f, ze3v_f ) … … 364 288 ! Revert "before" velocities to time split estimate 365 289 ! Doing it here also means that asselin filter contribution is removed 366 !$OMP PARALLEL 367 !$OMP DO schedule(static) private(jj, ji) 368 DO jj = 1, jpj 369 DO ji = 1, jpi 370 zue(ji,jj) = e3u_b(ji,jj,1) * ub(ji,jj,1) * umask(ji,jj,1) 371 zve(ji,jj) = e3v_b(ji,jj,1) * vb(ji,jj,1) * vmask(ji,jj,1) 372 END DO 290 zue(:,:) = e3u_b(:,:,1) * ub(:,:,1) * umask(:,:,1) 291 zve(:,:) = e3v_b(:,:,1) * vb(:,:,1) * vmask(:,:,1) 292 DO jk = 2, jpkm1 293 zue(:,:) = zue(:,:) + e3u_b(:,:,jk) * ub(:,:,jk) * umask(:,:,jk) 294 zve(:,:) = zve(:,:) + e3v_b(:,:,jk) * vb(:,:,jk) * vmask(:,:,jk) 373 295 END DO 374 DO jk = 2, jpkm1 375 !$OMP DO schedule(static) private(jj, ji) 376 DO jj = 1, jpj 377 DO ji = 1, jpi 378 zue(ji,jj) = zue(ji,jj) + e3u_b(ji,jj,jk) * ub(ji,jj,jk) * umask(ji,jj,jk) 379 zve(ji,jj) = zve(ji,jj) + e3v_b(ji,jj,jk) * vb(ji,jj,jk) * vmask(ji,jj,jk) 380 END DO 381 END DO 296 DO jk = 1, jpkm1 297 ub(:,:,jk) = ub(:,:,jk) - (zue(:,:) * r1_hu_n(:,:) - un_b(:,:)) * umask(:,:,jk) 298 vb(:,:,jk) = vb(:,:,jk) - (zve(:,:) * r1_hv_n(:,:) - vn_b(:,:)) * vmask(:,:,jk) 382 299 END DO 383 !$OMP DO schedule(static) private(jk,jj,ji)384 DO jk = 1, jpkm1385 DO jj = 1, jpj386 DO ji = 1, jpi387 ub(ji,jj,jk) = ub(ji,jj,jk) - (zue(ji,jj) * r1_hu_n(ji,jj) - un_b(ji,jj)) * umask(ji,jj,jk)388 vb(ji,jj,jk) = vb(ji,jj,jk) - (zve(ji,jj) * r1_hv_n(ji,jj) - vn_b(ji,jj)) * vmask(ji,jj,jk)389 END DO390 END DO391 END DO392 !$OMP END PARALLEL393 300 ENDIF 394 301 ! … … 401 308 ! 402 309 IF(.NOT.ln_linssh ) THEN 403 !$OMP PARALLEL 404 !$OMP DO schedule(static) private(jj, ji) 405 DO jj = 1, jpj 406 DO ji = 1, jpi 407 hu_b(ji,jj) = e3u_b(ji,jj,1) * umask(ji,jj,1) 408 hv_b(ji,jj) = e3v_b(ji,jj,1) * vmask(ji,jj,1) 409 END DO 310 hu_b(:,:) = e3u_b(:,:,1) * umask(:,:,1) 311 hv_b(:,:) = e3v_b(:,:,1) * vmask(:,:,1) 312 DO jk = 2, jpkm1 313 hu_b(:,:) = hu_b(:,:) + e3u_b(:,:,jk) * umask(:,:,jk) 314 hv_b(:,:) = hv_b(:,:) + e3v_b(:,:,jk) * vmask(:,:,jk) 410 315 END DO 411 DO jk = 2, jpkm1 412 !$OMP DO schedule(static) private(jj, ji) 413 DO jj = 1, jpj 414 DO ji = 1, jpi 415 hu_b(ji,jj) = hu_b(ji,jj) + e3u_b(ji,jj,jk) * umask(ji,jj,jk) 416 hv_b(ji,jj) = hv_b(ji,jj) + e3v_b(ji,jj,jk) * vmask(ji,jj,jk) 417 END DO 418 END DO 419 END DO 420 !$OMP DO schedule(static) private(jj, ji) 421 DO jj = 1, jpj 422 DO ji = 1, jpi 423 r1_hu_b(ji,jj) = ssumask(ji,jj) / ( hu_b(ji,jj) + 1._wp - ssumask(ji,jj) ) 424 r1_hv_b(ji,jj) = ssvmask(ji,jj) / ( hv_b(ji,jj) + 1._wp - ssvmask(ji,jj) ) 425 END DO 426 END DO 427 !$OMP END PARALLEL 428 ENDIF 429 ! 430 !$OMP PARALLEL 431 !$OMP DO schedule(static) private(jj, ji) 432 DO jj = 1, jpj 433 DO ji = 1, jpi 434 un_b(ji,jj) = e3u_a(ji,jj,1) * un(ji,jj,1) * umask(ji,jj,1) 435 ub_b(ji,jj) = e3u_b(ji,jj,1) * ub(ji,jj,1) * umask(ji,jj,1) 436 vn_b(ji,jj) = e3v_a(ji,jj,1) * vn(ji,jj,1) * vmask(ji,jj,1) 437 vb_b(ji,jj) = e3v_b(ji,jj,1) * vb(ji,jj,1) * vmask(ji,jj,1) 438 END DO 316 r1_hu_b(:,:) = ssumask(:,:) / ( hu_b(:,:) + 1._wp - ssumask(:,:) ) 317 r1_hv_b(:,:) = ssvmask(:,:) / ( hv_b(:,:) + 1._wp - ssvmask(:,:) ) 318 ENDIF 319 ! 320 un_b(:,:) = e3u_a(:,:,1) * un(:,:,1) * umask(:,:,1) 321 ub_b(:,:) = e3u_b(:,:,1) * ub(:,:,1) * umask(:,:,1) 322 vn_b(:,:) = e3v_a(:,:,1) * vn(:,:,1) * vmask(:,:,1) 323 vb_b(:,:) = e3v_b(:,:,1) * vb(:,:,1) * vmask(:,:,1) 324 DO jk = 2, jpkm1 325 un_b(:,:) = un_b(:,:) + e3u_a(:,:,jk) * un(:,:,jk) * umask(:,:,jk) 326 ub_b(:,:) = ub_b(:,:) + e3u_b(:,:,jk) * ub(:,:,jk) * umask(:,:,jk) 327 vn_b(:,:) = vn_b(:,:) + e3v_a(:,:,jk) * vn(:,:,jk) * vmask(:,:,jk) 328 vb_b(:,:) = vb_b(:,:) + e3v_b(:,:,jk) * vb(:,:,jk) * vmask(:,:,jk) 439 329 END DO 440 DO jk = 2, jpkm1 441 !$OMP DO schedule(static) private(jj, ji) 442 DO jj = 1, jpj 443 DO ji = 1, jpi 444 un_b(ji,jj) = un_b(ji,jj) + e3u_a(ji,jj,jk) * un(ji,jj,jk) * umask(ji,jj,jk) 445 ub_b(ji,jj) = ub_b(ji,jj) + e3u_b(ji,jj,jk) * ub(ji,jj,jk) * umask(ji,jj,jk) 446 vn_b(ji,jj) = vn_b(ji,jj) + e3v_a(ji,jj,jk) * vn(ji,jj,jk) * vmask(ji,jj,jk) 447 vb_b(ji,jj) = vb_b(ji,jj) + e3v_b(ji,jj,jk) * vb(ji,jj,jk) * vmask(ji,jj,jk) 448 END DO 449 END DO 450 END DO 451 !$OMP DO schedule(static) private(jj, ji) 452 DO jj = 1, jpj 453 DO ji = 1, jpi 454 un_b(ji,jj) = un_b(ji,jj) * r1_hu_a(ji,jj) 455 vn_b(ji,jj) = vn_b(ji,jj) * r1_hv_a(ji,jj) 456 ub_b(ji,jj) = ub_b(ji,jj) * r1_hu_b(ji,jj) 457 vb_b(ji,jj) = vb_b(ji,jj) * r1_hv_b(ji,jj) 458 END DO 459 END DO 460 !$OMP END PARALLEL 330 un_b(:,:) = un_b(:,:) * r1_hu_a(:,:) 331 vn_b(:,:) = vn_b(:,:) * r1_hv_a(:,:) 332 ub_b(:,:) = ub_b(:,:) * r1_hu_b(:,:) 333 vb_b(:,:) = vb_b(:,:) * r1_hv_b(:,:) 461 334 ! 462 335 IF( .NOT.ln_dynspg_ts ) THEN ! output the barotropic currents … … 465 338 ENDIF 466 339 IF( l_trddyn ) THEN ! 3D output: asselin filter trends on momentum 467 !$OMP PARALLEL DO schedule(static) private(jk, jj, ji) 468 DO jk = 1, jpkm1 469 DO jj = 1, jpj 470 DO ji = 1, jpi 471 zua(ji,jj,jk) = ( ub(ji,jj,jk) - zua(ji,jj,jk) ) * z1_2dt 472 zva(ji,jj,jk) = ( vb(ji,jj,jk) - zva(ji,jj,jk) ) * z1_2dt 473 END DO 474 END DO 475 END DO 340 zua(:,:,:) = ( ub(:,:,:) - zua(:,:,:) ) * z1_2dt 341 zva(:,:,:) = ( vb(:,:,:) - zva(:,:,:) ) * z1_2dt 476 342 CALL trd_dyn( zua, zva, jpdyn_atf, kt ) 477 343 ENDIF
Note: See TracChangeset
for help on using the changeset viewer.
