Changeset 13337 for NEMO/branches/2020/dev_r13312_AGRIF-03-04_jchanut_vinterp_tstep/src/NST/agrif_oce_update.F90
- Timestamp:
- 2020-07-24T16:01:24+02:00 (4 years ago)
- File:
-
- 1 edited
Legend:
- Unmodified
- Added
- Removed
-
NEMO/branches/2020/dev_r13312_AGRIF-03-04_jchanut_vinterp_tstep/src/NST/agrif_oce_update.F90
r13286 r13337 50 50 IF (lwp.AND.lk_agrif_debug) Write(*,*) 'Update tracers from grid Number',Agrif_Fixed() 51 51 52 #if defined key_vertical 53 ! Effect of this has to be carrefully checked 54 ! depending on what the nesting tools ensure for 55 ! volume conservation: 56 Agrif_UseSpecialValueInUpdate = .FALSE. 57 #else 58 Agrif_UseSpecialValueInUpdate = .TRUE. 59 #endif 52 Agrif_UseSpecialValueInUpdate = .NOT.l_vremap 60 53 Agrif_SpecialValueFineGrid = 0._wp 54 l_vremap = ln_vremap 61 55 ! 62 56 # if ! defined DECAL_FEEDBACK … … 71 65 ! 72 66 Agrif_UseSpecialValueInUpdate = .FALSE. 67 l_vremap = .FALSE. 73 68 ! 74 69 ! … … 86 81 Agrif_UseSpecialValueInUpdate = .FALSE. 87 82 Agrif_SpecialValueFineGrid = 0._wp 88 89 use_sign_north = .TRUE.90 sign_north = -1._wp83 l_vremap = ln_vremap 84 use_sign_north = .TRUE. 85 sign_north = -1._wp 91 86 92 87 ! … … 133 128 ! 134 129 use_sign_north = .FALSE. 130 ln_vremap = .FALSE. 135 131 ! 136 132 END SUBROUTINE Agrif_Update_Dyn … … 291 287 END SUBROUTINE dom_vvl_update_UVF 292 288 293 #if defined key_vertical294 289 295 290 SUBROUTINE updateTS( tabres, i1, i2, j1, j2, k1, k2, n1, n2, before ) … … 311 306 ! 312 307 IF (before) THEN 313 !jc_alt 314 ! AGRIF_SpecialValue = -999._wp 315 DO jn = n1,n2-1 308 IF ( l_vremap ) THEN 309 DO jn = n1,n2-1 310 DO jk=k1,k2 311 DO jj=j1,j2 312 DO ji=i1,i2 313 tabres(ji,jj,jk,jn) = ts(ji,jj,jk,jn,Kmm_a) * e3t(ji,jj,jk,Kmm_a) 314 END DO 315 END DO 316 END DO 317 END DO 316 318 DO jk=k1,k2 317 319 DO jj=j1,j2 318 320 DO ji=i1,i2 319 !jc_alt 320 ! tabres(ji,jj,jk,jn) = (ts(ji,jj,jk,jn,Kmm_a) * e3t(ji,jj,jk,Kmm_a) ) & 321 ! & * tmask(ji,jj,jk) + (tmask(ji,jj,jk)-1._wp) * 999._wp 322 tabres(ji,jj,jk,jn) = ts(ji,jj,jk,jn,Kmm_a) * e3t(ji,jj,jk,Kmm_a) 321 tabres(ji,jj,jk,n2) = tmask(ji,jj,jk) * e3t(ji,jj,jk,Kmm_a) 323 322 END DO 324 323 END DO 325 324 END DO 326 END DO 327 DO jk=k1,k2 325 ELSE 326 DO jn = 1,jpts 327 DO jk=k1,k2 328 DO jj=j1,j2 329 DO ji=i1,i2 330 tabres(ji,jj,jk,jn) = ts(ji,jj,jk,jn,Kmm_a) * e3t(ji,jj,jk,Kmm_a) / e3t_0(ji,jj,jk) 331 END DO 332 END DO 333 END DO 334 END DO 335 336 ENDIF 337 ELSE 338 IF ( l_vremap ) THEN 339 tabres_child(:,:,:,:) = 0._wp 340 AGRIF_SpecialValue = 0._wp 328 341 DO jj=j1,j2 329 342 DO ji=i1,i2 330 !jc_alt 331 ! tabres(ji,jj,jk,n2) = tmask(ji,jj,jk) * e3t(ji,jj,jk,Kmm_a) & 332 ! & + (tmask(ji,jj,jk) - 1._wp) * 999._wp 333 tabres(ji,jj,jk,n2) = tmask(ji,jj,jk) * e3t(ji,jj,jk,Kmm_a) 334 END DO 335 END DO 336 END DO 337 ELSE 338 tabres_child(:,:,:,:) = 0._wp 339 AGRIF_SpecialValue = 0._wp 340 DO jj=j1,j2 341 DO ji=i1,i2 342 N_in = 0 343 DO jk=k1,k2 !k2 = jpk of child grid 344 ! jc_alt 345 ! IF (tabres(ji,jj,jk,n2) < -900._wp ) EXIT 346 IF (tabres(ji,jj,jk,n2) == 0._wp ) EXIT 347 N_in = N_in + 1 348 tabin(jk,:) = tabres(ji,jj,jk,n1:n2-1)/tabres(ji,jj,jk,n2) 349 h_in(N_in) = tabres(ji,jj,jk,n2) 343 N_in = 0 344 DO jk=k1,k2 !k2 = jpk of child grid 345 IF (tabres(ji,jj,jk,n2) == 0._wp ) EXIT 346 N_in = N_in + 1 347 tabin(jk,:) = tabres(ji,jj,jk,n1:n2-1)/tabres(ji,jj,jk,n2) 348 h_in(N_in) = tabres(ji,jj,jk,n2) 349 ENDDO 350 N_out = 0 351 DO jk=1,jpk ! jpk of parent grid 352 IF (tmask(ji,jj,jk) == 0 ) EXIT ! TODO: Will not work with ISF 353 N_out = N_out + 1 354 h_out(N_out) = e3t(ji,jj,jk,Kmm_a) 355 ENDDO 356 IF (N_in*N_out > 0) THEN !Remove this? 357 CALL reconstructandremap(tabin(1:N_in,1:jpts),h_in(1:N_in),tabres_child(ji,jj,1:N_out,1:jpts),h_out(1:N_out),N_in,N_out,jpts) 358 ENDIF 350 359 ENDDO 351 N_out = 0352 DO jk=1,jpk ! jpk of parent grid353 IF (tmask(ji,jj,jk) == 0 ) EXIT ! TODO: Will not work with ISF354 N_out = N_out + 1355 h_out(N_out) = e3t(ji,jj,jk,Kmm_a)356 ENDDO357 IF (N_in*N_out > 0) THEN !Remove this?358 CALL reconstructandremap(tabin(1:N_in,1:jpts),h_in(1:N_in),tabres_child(ji,jj,1:N_out,1:jpts),h_out(1:N_out),N_in,N_out,jpts)359 ENDIF360 360 ENDDO 361 ENDDO 362 363 IF (.NOT.(lk_agrif_fstep.AND.(l_1st_euler))) THEN 364 ! Add asselin part 361 362 IF (.NOT.(lk_agrif_fstep.AND.(l_1st_euler))) THEN 363 ! Add asselin part 364 DO jn = 1,jpts 365 DO jk = 1, jpkm1 366 DO jj = j1, j2 367 DO ji = i1, i2 368 IF( tabres_child(ji,jj,jk,jn) /= 0._wp ) THEN 369 ztb = ts(ji,jj,jk,jn,Kbb_a) * e3t(ji,jj,jk,Kbb_a) ! fse3t_b prior update should be used 370 ztnu = tabres_child(ji,jj,jk,jn) * e3t(ji,jj,jk,Kmm_a) 371 ztno = ts(ji,jj,jk,jn,Kmm_a) * e3t(ji,jj,jk,Krhs_a) 372 ts(ji,jj,jk,jn,Kbb_a) = ( ztb + rn_atfp * ( ztnu - ztno) ) & 373 & * tmask(ji,jj,jk) / e3t(ji,jj,jk,Kbb_a) 374 ENDIF 375 END DO 376 END DO 377 END DO 378 END DO 379 ENDIF 365 380 DO jn = 1,jpts 366 381 DO jk = 1, jpkm1 367 382 DO jj = j1, j2 368 383 DO ji = i1, i2 369 IF( tabres_child(ji,jj,jk,jn) /= 0._wp ) THEN 370 ztb = ts(ji,jj,jk,jn,Kbb_a) * e3t(ji,jj,jk,Kbb_a) ! fse3t_b prior update should be used 371 ztnu = tabres_child(ji,jj,jk,jn) * e3t(ji,jj,jk,Kmm_a) 372 ztno = ts(ji,jj,jk,jn,Kmm_a) * e3t(ji,jj,jk,Krhs_a) 373 ts(ji,jj,jk,jn,Kbb_a) = ( ztb + rn_atfp * ( ztnu - ztno) ) & 374 & * tmask(ji,jj,jk) / e3t(ji,jj,jk,Kbb_a) 375 ENDIF 384 IF( tabres_child(ji,jj,jk,jn) /= 0._wp ) THEN 385 ts(ji,jj,jk,jn,Kmm_a) = tabres_child(ji,jj,jk,jn) 386 END IF 376 387 END DO 377 388 END DO 378 389 END DO 379 390 END DO 380 ENDIF 381 DO jn = 1,jpts 382 DO jk = 1, jpkm1 383 DO jj = j1, j2 384 DO ji = i1, i2 385 IF( tabres_child(ji,jj,jk,jn) /= 0._wp ) THEN 386 ts(ji,jj,jk,jn,Kmm_a) = tabres_child(ji,jj,jk,jn) 387 END IF 391 ELSE 392 DO jn = 1,jpts 393 tabres(i1:i2,j1:j2,k1:k2,jn) = tabres(i1:i2,j1:j2,k1:k2,jn) * e3t_0(i1:i2,j1:j2,k1:k2) & 394 & * tmask(i1:i2,j1:j2,k1:k2) 395 ENDDO 396 397 IF (.NOT.(lk_agrif_fstep.AND.(l_1st_euler))) THEN 398 ! Add asselin part 399 DO jn = 1,jpts 400 DO jk = k1, k2 401 DO jj = j1, j2 402 DO ji = i1, i2 403 IF( tabres(ji,jj,jk,jn) /= 0._wp ) THEN 404 ztb = ts(ji,jj,jk,jn,Kbb_a) * e3t(ji,jj,jk,Kbb_a) ! fse3t_b prior update should be used 405 ztnu = tabres(ji,jj,jk,jn) 406 ztno = ts(ji,jj,jk,jn,Kmm_a) * e3t(ji,jj,jk,Krhs_a) 407 ts(ji,jj,jk,jn,Kbb_a) = ( ztb + rn_atfp * ( ztnu - ztno) ) & 408 & * tmask(ji,jj,jk) / e3t(ji,jj,jk,Kbb_a) 409 ENDIF 410 END DO 411 END DO 388 412 END DO 389 413 END DO 390 END DO 391 END DO 392 ! 414 ENDIF 415 DO jn = 1,jpts 416 DO jk=k1,k2 417 DO jj=j1,j2 418 DO ji=i1,i2 419 IF( tabres(ji,jj,jk,jn) /= 0._wp ) THEN 420 ts(ji,jj,jk,jn,Kmm_a) = tabres(ji,jj,jk,jn) / e3t(ji,jj,jk,Kmm_a) 421 END IF 422 END DO 423 END DO 424 END DO 425 END DO 426 ! 427 ENDIF 393 428 IF ((l_1st_euler).AND.(Agrif_Nb_Step()==0) ) THEN 394 429 ts(i1:i2,j1:j2,1:jpkm1,1:jpts,Kbb_a) = ts(i1:i2,j1:j2,1:jpkm1,1:jpts,Kmm_a) … … 398 433 END SUBROUTINE updateTS 399 434 400 # else401 402 SUBROUTINE updateTS( tabres, i1, i2, j1, j2, k1, k2, n1, n2, before )403 !!---------------------------------------------404 !! *** ROUTINE updateT ***405 !!---------------------------------------------406 INTEGER, INTENT(in) :: i1,i2,j1,j2,k1,k2,n1,n2407 REAL(wp),DIMENSION(i1:i2,j1:j2,k1:k2,n1:n2), INTENT(inout) :: tabres408 LOGICAL, INTENT(in) :: before409 !!410 INTEGER :: ji,jj,jk,jn411 REAL(wp) :: ztb, ztnu, ztno412 !!---------------------------------------------413 !414 IF (before) THEN415 DO jn = 1,jpts416 DO jk=k1,k2417 DO jj=j1,j2418 DO ji=i1,i2419 !> jc tmp420 tabres(ji,jj,jk,jn) = ts(ji,jj,jk,jn,Kmm_a) * e3t(ji,jj,jk,Kmm_a) / e3t_0(ji,jj,jk)421 ! tabres(ji,jj,jk,jn) = ts(ji,jj,jk,jn,Kmm_a) * e3t(ji,jj,jk,Kmm_a)422 !< jc tmp423 END DO424 END DO425 END DO426 END DO427 ELSE428 !> jc tmp429 DO jn = 1,jpts430 tabres(i1:i2,j1:j2,k1:k2,jn) = tabres(i1:i2,j1:j2,k1:k2,jn) * e3t_0(i1:i2,j1:j2,k1:k2) &431 & * tmask(i1:i2,j1:j2,k1:k2)432 ENDDO433 !< jc tmp434 IF (.NOT.(lk_agrif_fstep.AND.(l_1st_euler))) THEN435 ! Add asselin part436 DO jn = 1,jpts437 DO jk = k1, k2438 DO jj = j1, j2439 DO ji = i1, i2440 IF( tabres(ji,jj,jk,jn) /= 0._wp ) THEN441 ztb = ts(ji,jj,jk,jn,Kbb_a) * e3t(ji,jj,jk,Kbb_a) ! fse3t_b prior update should be used442 ztnu = tabres(ji,jj,jk,jn)443 ztno = ts(ji,jj,jk,jn,Kmm_a) * e3t(ji,jj,jk,Krhs_a)444 ts(ji,jj,jk,jn,Kbb_a) = ( ztb + rn_atfp * ( ztnu - ztno) ) &445 & * tmask(ji,jj,jk) / e3t(ji,jj,jk,Kbb_a)446 ENDIF447 END DO448 END DO449 END DO450 END DO451 ENDIF452 DO jn = 1,jpts453 DO jk=k1,k2454 DO jj=j1,j2455 DO ji=i1,i2456 IF( tabres(ji,jj,jk,jn) /= 0._wp ) THEN457 ts(ji,jj,jk,jn,Kmm_a) = tabres(ji,jj,jk,jn) / e3t(ji,jj,jk,Kmm_a)458 END IF459 END DO460 END DO461 END DO462 END DO463 !464 IF ((l_1st_euler).AND.(Agrif_Nb_Step()==0) ) THEN465 ts(i1:i2,j1:j2,k1:k2,1:jpts,Kbb_a) = ts(i1:i2,j1:j2,k1:k2,1:jpts,Kmm_a)466 ENDIF467 !468 ENDIF469 !470 END SUBROUTINE updateTS471 472 # endif473 474 # if defined key_vertical475 435 476 436 SUBROUTINE updateu( tabres, i1, i2, j1, j2, k1, k2, n1, n2, before ) … … 496 456 IF( before ) THEN 497 457 zrhoy = Agrif_Rhoy() 498 !jc_alt499 ! AGRIF_SpecialValue = -999._wp500 458 DO jk=k1,k2 459 tabres(i1:i2,j1:j2,jk,1) = zrhoy * e2u(i1:i2,j1:j2) * e3u(i1:i2,j1:j2,jk,Kmm_a) & 460 & * umask(i1:i2,j1:j2,jk) * uu(i1:i2,j1:j2,jk,Kmm_a) 461 END DO 462 463 IF ( l_vremap ) THEN 464 DO jk=k1,k2 465 tabres(i1:i2,j1:j2,jk,2) = zrhoy * umask(i1:i2,j1:j2,jk) * e2u(i1:i2,j1:j2) * e3u(i1:i2,j1:j2,jk,Kmm_a) 466 END DO 467 ENDIF 468 469 ELSE 470 471 tabres_child(:,:,:) = 0._wp 472 AGRIF_SpecialValue = 0._wp 473 474 IF ( l_vremap ) THEN 475 501 476 DO jj=j1,j2 502 477 DO ji=i1,i2 503 !jc_alt 504 ! tabres(ji,jj,jk,1) = zrhoy * e2u(ji,jj) * e3u(ji,jj,jk,Kmm_a) * umask(ji,jj,jk) * uu(ji,jj,jk,Kmm_a) & 505 ! & + (umask(ji,jj,jk)-1._wp)*999._wp 506 tabres(ji,jj,jk,1) = zrhoy * e2u(ji,jj) * e3u(ji,jj,jk,Kmm_a) * umask(ji,jj,jk) * uu(ji,jj,jk,Kmm_a) 507 !jc_alt 508 ! tabres(ji,jj,jk,2) = zrhoy * umask(ji,jj,jk) * e2u(ji,jj) * e3u(ji,jj,jk,Kmm_a) & 509 ! & + (umask(ji,jj,jk)-1._wp)*999._wp 510 tabres(ji,jj,jk,2) = zrhoy * umask(ji,jj,jk) * e2u(ji,jj) * e3u(ji,jj,jk,Kmm_a) 511 END DO 512 END DO 513 END DO 514 ELSE 515 tabres_child(:,:,:) = 0. 516 AGRIF_SpecialValue = 0._wp 517 DO jj=j1,j2 518 DO ji=i1,i2 519 N_in = 0 520 h_in(:) = 0._wp 521 tabin(:) = 0._wp 522 DO jk=k1,k2 !k2=jpk of child grid 523 !jc_alt 524 ! IF( tabres(ji,jj,jk,2) < -900._wp) EXIT 525 IF( tabres(ji,jj,jk,2) == 0.) EXIT 526 N_in = N_in + 1 527 tabin(jk) = tabres(ji,jj,jk,1)/tabres(ji,jj,jk,2) 528 h_in(N_in) = tabres(ji,jj,jk,2)/e2u(ji,jj) 478 N_in = 0 479 h_in(:) = 0._wp 480 tabin(:) = 0._wp 481 DO jk=k1,k2 !k2=jpk of child grid 482 IF( tabres(ji,jj,jk,2) == 0.) EXIT 483 N_in = N_in + 1 484 tabin(jk) = tabres(ji,jj,jk,1)/tabres(ji,jj,jk,2) 485 h_in(N_in) = tabres(ji,jj,jk,2) * r1_e2u(ji,jj) 486 ENDDO 487 N_out = 0 488 DO jk=1,jpk 489 IF (umask(ji,jj,jk) == 0) EXIT 490 N_out = N_out + 1 491 h_out(N_out) = e3u(ji,jj,jk,Kmm_a) 492 ENDDO 493 IF (N_in * N_out > 0) THEN 494 h_diff = sum(h_out(1:N_out))-sum(h_in(1:N_in)) 495 excess = 0._wp 496 IF (h_diff < -1.e-4) THEN 497 DO jk=N_in,1,-1 498 thick = MIN(-1*h_diff, h_in(jk)) 499 excess = excess + tabin(jk)*thick*e2u(ji,jj) 500 tabin(jk) = tabin(jk)*(1. - thick/h_in(jk)) 501 h_diff = h_diff + thick 502 IF ( h_diff == 0) THEN 503 N_in = jk 504 h_in(jk) = h_in(jk) - thick 505 EXIT 506 ENDIF 507 ENDDO 508 ENDIF 509 CALL reconstructandremap(tabin(1:N_in),h_in(1:N_in),tabres_child(ji,jj,1:N_out),h_out(1:N_out),N_in,N_out,1) 510 tabres_child(ji,jj,N_out) = tabres_child(ji,jj,N_out) + excess/(e2u(ji,jj)*h_out(N_out)) 511 ENDIF 529 512 ENDDO 530 N_out = 0531 DO jk=1,jpk532 IF (umask(ji,jj,jk) == 0) EXIT533 N_out = N_out + 1534 h_out(N_out) = e3u(ji,jj,jk,Kmm_a)535 ENDDO536 IF (N_in * N_out > 0) THEN537 h_diff = sum(h_out(1:N_out))-sum(h_in(1:N_in))538 excess = 0._wp539 IF (h_diff < -1.e-4) THEN540 !Even if bathy at T points match it's possible for the U points to be deeper in the child grid.541 !In this case we need to move transport from the child grid cells below bed of parent grid into the bottom cell.542 DO jk=N_in,1,-1543 thick = MIN(-1*h_diff, h_in(jk))544 excess = excess + tabin(jk)*thick*e2u(ji,jj)545 tabin(jk) = tabin(jk)*(1. - thick/h_in(jk))546 h_diff = h_diff + thick547 IF ( h_diff == 0) THEN548 N_in = jk549 h_in(jk) = h_in(jk) - thick550 EXIT551 ENDIF552 ENDDO553 ENDIF554 CALL reconstructandremap(tabin(1:N_in),h_in(1:N_in),tabres_child(ji,jj,1:N_out),h_out(1:N_out),N_in,N_out,1)555 tabres_child(ji,jj,N_out) = tabres_child(ji,jj,N_out) + excess/(e2u(ji,jj)*h_out(N_out))556 ENDIF557 513 ENDDO 558 ENDDO 514 515 ELSE 516 DO jk=1,jpk 517 DO jj=j1,j2 518 DO ji=i1,i2 519 tabres_child(ji,jj,jk) = tabres(ji,jj,jk,1) * r1_e2u(ji,jj) / e3u(ji,jj,jk,Kmm_a) 520 END DO 521 END DO 522 END DO 523 ENDIF 559 524 ! 560 525 DO jk=1,jpk … … 582 547 END SUBROUTINE updateu 583 548 584 #else585 586 SUBROUTINE updateu( tabres, i1, i2, j1, j2, k1, k2, n1, n2, before )587 !!---------------------------------------------588 !! *** ROUTINE updateu ***589 !!---------------------------------------------590 INTEGER , INTENT(in ) :: i1, i2, j1, j2, k1, k2, n1, n2591 REAL(wp), DIMENSION(i1:i2,j1:j2,k1:k2,n1:n2), INTENT(inout) :: tabres592 LOGICAL , INTENT(in ) :: before593 !594 INTEGER :: ji, jj, jk595 REAL(wp) :: zrhoy, zub, zunu, zuno596 !!---------------------------------------------597 !598 IF( before ) THEN599 zrhoy = Agrif_Rhoy()600 DO jk = k1, k2601 tabres(i1:i2,j1:j2,jk,1) = zrhoy * e2u(i1:i2,j1:j2) * e3u(i1:i2,j1:j2,jk,Kmm_a) * uu(i1:i2,j1:j2,jk,Kmm_a)602 END DO603 ELSE604 DO jk=k1,k2605 DO jj=j1,j2606 DO ji=i1,i2607 tabres(ji,jj,jk,1) = tabres(ji,jj,jk,1) * r1_e2u(ji,jj)608 !609 IF (.NOT.(lk_agrif_fstep.AND.(l_1st_euler))) THEN ! Add asselin part610 zub = uu(ji,jj,jk,Kbb_a) * e3u(ji,jj,jk,Kbb_a) ! fse3t_b prior update should be used611 zuno = uu(ji,jj,jk,Kmm_a) * e3u(ji,jj,jk,Krhs_a)612 zunu = tabres(ji,jj,jk,1)613 uu(ji,jj,jk,Kbb_a) = ( zub + rn_atfp * ( zunu - zuno) ) &614 & * umask(ji,jj,jk) / e3u(ji,jj,jk,Kbb_a)615 ENDIF616 !617 uu(ji,jj,jk,Kmm_a) = tabres(ji,jj,jk,1) * umask(ji,jj,jk) / e3u(ji,jj,jk,Kmm_a)618 END DO619 END DO620 END DO621 !622 IF ((l_1st_euler).AND.(Agrif_Nb_Step()==0) ) THEN623 uu(i1:i2,j1:j2,k1:k2,Kbb_a) = uu(i1:i2,j1:j2,k1:k2,Kmm_a)624 ENDIF625 !626 ENDIF627 !628 END SUBROUTINE updateu629 630 # endif631 632 549 SUBROUTINE correct_u_bdy( tabres, i1, i2, j1, j2, k1, k2, n1, n2, before, nb, ndir ) 633 550 !!--------------------------------------------- … … 674 591 END SUBROUTINE correct_u_bdy 675 592 676 # if defined key_vertical677 593 678 594 SUBROUTINE updatev( tabres, i1, i2, j1, j2, k1, k2, n1, n2, before ) … … 698 614 IF( before ) THEN 699 615 zrhox = Agrif_Rhox() 700 !jc_alt701 ! AGRIF_SpecialValue = -999._wp702 616 DO jk=k1,k2 617 tabres(i1:i2,j1:j2,jk,1) = zrhox * e1v(i1:i2,j1:j2) * e3v(i1:i2,j1:j2,jk,Kmm_a) & 618 & * vmask(i1:i2,j1:j2,jk) * vv(i1:i2,j1:j2,jk,Kmm_a) 619 END DO 620 621 IF ( l_vremap ) THEN 622 DO jk=k1,k2 623 tabres(i1:i2,j1:j2,jk,2) = zrhox * e1v(i1:i2,j1:j2) * e3v(i1:i2,j1:j2,jk,Kmm_a) * vmask(i1:i2,j1:j2,jk) 624 END DO 625 ENDIF 626 627 ELSE 628 629 tabres_child(:,:,:) = 0._wp 630 AGRIF_SpecialValue = 0._wp 631 632 IF ( l_vremap ) THEN 633 703 634 DO jj=j1,j2 704 635 DO ji=i1,i2 705 !jc_alt 706 ! tabres(ji,jj,jk,1) = zrhox * e1v(ji,jj) * e3v(ji,jj,jk,Kmm_a) * vmask(ji,jj,jk) * vv(ji,jj,jk,Kmm_a) & 707 ! & + (vmask(ji,jj,jk)-1._wp) * 999._wp 708 tabres(ji,jj,jk,1) = zrhox * e1v(ji,jj) * e3v(ji,jj,jk,Kmm_a) * vmask(ji,jj,jk) * vv(ji,jj,jk,Kmm_a) 709 !jc_alt 710 ! tabres(ji,jj,jk,2) = zrhox * e1v(ji,jj) * e3v(ji,jj,jk,Kmm_a) * vmask(ji,jj,jk) & 711 ! & + (vmask(ji,jj,jk)-1._wp) * 999._wp 712 tabres(ji,jj,jk,2) = zrhox * e1v(ji,jj) * e3v(ji,jj,jk,Kmm_a) * vmask(ji,jj,jk) 713 END DO 714 END DO 715 END DO 716 ELSE 717 tabres_child(:,:,:) = 0. 718 AGRIF_SpecialValue = 0._wp 719 DO jj=j1,j2 720 DO ji=i1,i2 721 N_in = 0 722 DO jk=k1,k2 723 !jc_alt 724 ! IF (tabres(ji,jj,jk,2) < -900._wp) EXIT 725 IF (tabres(ji,jj,jk,2) == 0) EXIT 726 N_in = N_in + 1 727 tabin(jk) = tabres(ji,jj,jk,1)/tabres(ji,jj,jk,2) 728 h_in(N_in) = tabres(ji,jj,jk,2)/e1v(ji,jj) 729 ENDDO 730 N_out = 0 731 DO jk=1,jpk 732 IF (vmask(ji,jj,jk) == 0) EXIT 733 N_out = N_out + 1 734 h_out(N_out) = e3v(ji,jj,jk,Kmm_a) 735 ENDDO 736 IF (N_in * N_out > 0) THEN 737 h_diff = sum(h_out(1:N_out))-sum(h_in(1:N_in)) 738 excess = 0._wp 739 IF (h_diff < -1.e-4) then 636 N_in = 0 637 DO jk=k1,k2 638 IF (tabres(ji,jj,jk,2) == 0) EXIT 639 N_in = N_in + 1 640 tabin(jk) = tabres(ji,jj,jk,1)/tabres(ji,jj,jk,2) 641 h_in(N_in) = tabres(ji,jj,jk,2) * r1_e1v(ji,jj) 642 ENDDO 643 N_out = 0 644 DO jk=1,jpk 645 IF (vmask(ji,jj,jk) == 0) EXIT 646 N_out = N_out + 1 647 h_out(N_out) = e3v(ji,jj,jk,Kmm_a) 648 ENDDO 649 IF (N_in * N_out > 0) THEN 650 h_diff = sum(h_out(1:N_out))-sum(h_in(1:N_in)) 651 excess = 0._wp 652 IF (h_diff < -1.e-4) then 740 653 !Even if bathy at T points match it's possible for the V points to be deeper in the child grid. 741 654 !In this case we need to move transport from the child grid cells below bed of parent grid into the bottom cell. 742 DO jk=N_in,1,-1 743 thick = MIN(-1*h_diff, h_in(jk)) 744 excess = excess + tabin(jk)*thick*e2u(ji,jj) 745 tabin(jk) = tabin(jk)*(1. - thick/h_in(jk)) 746 h_diff = h_diff + thick 747 IF ( h_diff == 0) THEN 748 N_in = jk 749 h_in(jk) = h_in(jk) - thick 750 EXIT 751 ENDIF 752 ENDDO 655 DO jk=N_in,1,-1 656 thick = MIN(-1*h_diff, h_in(jk)) 657 excess = excess + tabin(jk)*thick*e2u(ji,jj) 658 tabin(jk) = tabin(jk)*(1. - thick/h_in(jk)) 659 h_diff = h_diff + thick 660 IF ( h_diff == 0) THEN 661 N_in = jk 662 h_in(jk) = h_in(jk) - thick 663 EXIT 664 ENDIF 665 ENDDO 666 ENDIF 667 CALL reconstructandremap(tabin(1:N_in),h_in(1:N_in),tabres_child(ji,jj,1:N_out),h_out(1:N_out),N_in,N_out,1) 668 tabres_child(ji,jj,N_out) = tabres_child(ji,jj,N_out) + excess/(e1v(ji,jj)*h_out(N_out)) 753 669 ENDIF 754 CALL reconstructandremap(tabin(1:N_in),h_in(1:N_in),tabres_child(ji,jj,1:N_out),h_out(1:N_out),N_in,N_out,1) 755 tabres_child(ji,jj,N_out) = tabres_child(ji,jj,N_out) + excess/(e1v(ji,jj)*h_out(N_out)) 756 ENDIF 670 ENDDO 757 671 ENDDO 758 ENDDO 672 673 ELSE 674 DO jk=1,jpk 675 DO jj=j1,j2 676 DO ji=i1,i2 677 tabres_child(ji,jj,jk) = tabres(ji,jj,jk,1) * r1_e1v(ji,jj) / e3v(ji,jj,jk,Kmm_a) 678 END DO 679 END DO 680 END DO 681 ENDIF 759 682 ! 760 683 DO jk=1,jpkm1 … … 782 705 END SUBROUTINE updatev 783 706 784 # else785 786 SUBROUTINE updatev( tabres, i1, i2, j1, j2, k1, k2, n1, n2, before)787 !!---------------------------------------------788 !! *** ROUTINE updatev ***789 !!---------------------------------------------790 INTEGER , INTENT(in ) :: i1, i2, j1, j2, k1, k2, n1, n2791 REAL(wp), DIMENSION(i1:i2,j1:j2,k1:k2,n1:n2), INTENT(inout) :: tabres792 LOGICAL , INTENT(in ) :: before793 !794 INTEGER :: ji, jj, jk795 REAL(wp) :: zrhox, zvb, zvnu, zvno796 !!---------------------------------------------797 !798 IF (before) THEN799 zrhox = Agrif_Rhox()800 DO jk=k1,k2801 DO jj=j1,j2802 DO ji=i1,i2803 tabres(ji,jj,jk,1) = zrhox * e1v(ji,jj) * e3v(ji,jj,jk,Kmm_a) * vv(ji,jj,jk,Kmm_a)804 END DO805 END DO806 END DO807 ELSE808 DO jk=k1,k2809 DO jj=j1,j2810 DO ji=i1,i2811 tabres(ji,jj,jk,1) = tabres(ji,jj,jk,1) * r1_e1v(ji,jj)812 !813 IF (.NOT.(lk_agrif_fstep.AND.(l_1st_euler))) THEN ! Add asselin part814 zvb = vv(ji,jj,jk,Kbb_a) * e3v(ji,jj,jk,Kbb_a) ! fse3t_b prior update should be used815 zvno = vv(ji,jj,jk,Kmm_a) * e3v(ji,jj,jk,Krhs_a)816 zvnu = tabres(ji,jj,jk,1)817 vv(ji,jj,jk,Kbb_a) = ( zvb + rn_atfp * ( zvnu - zvno) ) &818 & * vmask(ji,jj,jk) / e3v(ji,jj,jk,Kbb_a)819 ENDIF820 !821 vv(ji,jj,jk,Kmm_a) = tabres(ji,jj,jk,1) * vmask(ji,jj,jk) / e3v(ji,jj,jk,Kmm_a)822 END DO823 END DO824 END DO825 !826 IF ((l_1st_euler).AND.(Agrif_Nb_Step()==0) ) THEN827 vv(i1:i2,j1:j2,k1:k2,Kbb_a) = vv(i1:i2,j1:j2,k1:k2,Kmm_a)828 ENDIF829 !830 ENDIF831 !832 END SUBROUTINE updatev833 834 # endif835 707 836 708 SUBROUTINE correct_v_bdy( tabres, i1, i2, j1, j2, k1, k2, n1, n2, before, nb, ndir )
Note: See TracChangeset
for help on using the changeset viewer.
