- Timestamp:
- 2013-11-20T17:28:04+01:00 (11 years ago)
- File:
-
- 1 edited
Legend:
- Unmodified
- Added
- Removed
-
branches/2013/dev_MERGE_2013/NEMOGCM/NEMO/NST_SRC/agrif_opa_interp.F90
r3294 r4292 27 27 USE agrif_opa_sponge 28 28 USE lib_mpp 29 USE wrk_nemo 29 USE wrk_nemo 30 USE dynspg_oce 30 31 31 32 IMPLICIT NONE 32 33 PRIVATE 34 35 ! Barotropic arrays used to store open boundary data during 36 ! time-splitting loop: 37 REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:) :: ubdy_w, vbdy_w, hbdy_w 38 REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:) :: ubdy_e, vbdy_e, hbdy_e 39 REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:) :: ubdy_n, vbdy_n, hbdy_n 40 REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:) :: ubdy_s, vbdy_s, hbdy_s 33 41 34 PUBLIC Agrif_tra, Agrif_dyn, Agrif_ssh, interpu, interpv 42 PUBLIC Agrif_tra, Agrif_dyn, Agrif_ssh, Agrif_dyn_ts, Agrif_ssh_ts 43 PUBLIC interpu, interpv, interpunb, interpvnb, interpsshn 35 44 36 45 # include "domzgr_substitute.h90" … … 169 178 REAL(wp) :: timeref 170 179 REAL(wp) :: z2dt, znugdt 171 REAL(wp) :: zrhox, rhoy180 REAL(wp) :: zrhox, zrhoy 172 181 REAL(wp), POINTER, DIMENSION(:,:,:) :: zua, zva 173 182 REAL(wp), POINTER, DIMENSION(:,:) :: spgv1, spgu1, zua2d, zva2d … … 180 189 181 190 zrhox = Agrif_Rhox() 182 rhoy = Agrif_Rhoy()191 zrhoy = Agrif_Rhoy() 183 192 184 193 timeref = 1. … … 201 210 zva2d = 0. 202 211 212 #if defined key_dynspg_flt 203 213 Agrif_SpecialValue=0. 204 214 Agrif_UseSpecialValue = ln_spc_dyn 205 215 CALL Agrif_Bc_variable(zua2d,e1u_id,calledweight=1.,procname=interpu2d) 206 216 CALL Agrif_Bc_variable(zva2d,e2v_id,calledweight=1.,procname=interpv2d) 217 #endif 207 218 Agrif_UseSpecialValue = .FALSE. 208 219 … … 210 221 IF((nbondi == -1).OR.(nbondi == 2)) THEN 211 222 223 #if defined key_dynspg_flt 212 224 DO jj=1,jpj 213 laplacu(2,jj) = timeref * (zua2d(2,jj)/(rhoy*e2u(2,jj)))*umask(2,jj,1) 214 END DO 215 216 DO jk=1,jpkm1 217 DO jj=1,jpj 218 ua(1:2,jj,jk) = (zua(1:2,jj,jk)/(rhoy*e2u(1:2,jj))) 225 laplacu(2,jj) = timeref * (zua2d(2,jj)/(zrhoy*e2u(2,jj)))*umask(2,jj,1) 226 END DO 227 #endif 228 229 DO jk=1,jpkm1 230 DO jj=1,jpj 231 ua(1:2,jj,jk) = (zua(1:2,jj,jk)/(zrhoy*e2u(1:2,jj))) 219 232 ua(1:2,jj,jk) = ua(1:2,jj,jk) / fse3u(1:2,jj,jk) 220 233 END DO 221 234 END DO 222 235 236 #if defined key_dynspg_flt 223 237 DO jk=1,jpkm1 224 238 DO jj=1,jpj … … 240 254 ENDIF 241 255 END DO 256 #else 257 spgu(2,:) = ua_b(2,:) 258 #endif 242 259 243 260 DO jk=1,jpkm1 … … 278 295 279 296 IF((nbondi == 1).OR.(nbondi == 2)) THEN 280 297 #if defined key_dynspg_flt 281 298 DO jj=1,jpj 282 laplacu(nlci-2,jj) = timeref * (zua2d(nlci-2,jj)/(rhoy*e2u(nlci-2,jj))) 283 END DO 284 285 DO jk=1,jpkm1 286 DO jj=1,jpj 287 ua(nlci-2:nlci-1,jj,jk) = (zua(nlci-2:nlci-1,jj,jk)/(rhoy*e2u(nlci-2:nlci-1,jj))) 299 laplacu(nlci-2,jj) = timeref * (zua2d(nlci-2,jj)/(zrhoy*e2u(nlci-2,jj))) 300 END DO 301 #endif 302 303 DO jk=1,jpkm1 304 DO jj=1,jpj 305 ua(nlci-2:nlci-1,jj,jk) = (zua(nlci-2:nlci-1,jj,jk)/(zrhoy*e2u(nlci-2:nlci-1,jj))) 288 306 289 307 ua(nlci-2:nlci-1,jj,jk) = ua(nlci-2:nlci-1,jj,jk) / fse3u(nlci-2:nlci-1,jj,jk) … … 292 310 END DO 293 311 312 #if defined key_dynspg_flt 294 313 DO jk=1,jpkm1 295 314 DO jj=1,jpj … … 312 331 ENDIF 313 332 END DO 333 #else 334 spgu(nlci-2,:) = ua_b(nlci-2,:) 335 #endif 314 336 315 337 DO jk=1,jpkm1 … … 353 375 IF((nbondj == -1).OR.(nbondj == 2)) THEN 354 376 377 #if defined key_dynspg_flt 355 378 DO ji=1,jpi 356 379 laplacv(ji,2) = timeref * (zva2d(ji,2)/(zrhox*e1v(ji,2))) 357 380 END DO 381 #endif 358 382 359 383 DO jk=1,jpkm1 … … 364 388 END DO 365 389 390 #if defined key_dynspg_flt 366 391 DO jk=1,jpkm1 367 392 DO ji=1,jpi … … 383 408 ENDIF 384 409 END DO 410 #else 411 spgv(:,2)=va_b(:,2) 412 #endif 385 413 386 414 DO jk=1,jpkm1 … … 413 441 DO jk=1,jpkm1 414 442 DO ji=1,jpi 415 ua(ji,2,jk) = (zua(ji,2,jk)/( rhoy*e2u(ji,2)))*umask(ji,2,jk)443 ua(ji,2,jk) = (zua(ji,2,jk)/(zrhoy*e2u(ji,2)))*umask(ji,2,jk) 416 444 ua(ji,2,jk) = ua(ji,2,jk) / fse3u(ji,2,jk) 417 445 END DO … … 422 450 IF((nbondj == 1).OR.(nbondj == 2)) THEN 423 451 452 #if defined key_dynspg_flt 424 453 DO ji=1,jpi 425 454 laplacv(ji,nlcj-2) = timeref * (zva2d(ji,nlcj-2)/(zrhox*e1v(ji,nlcj-2))) 426 455 END DO 456 #endif 427 457 428 458 DO jk=1,jpkm1 … … 433 463 END DO 434 464 465 #if defined key_dynspg_flt 435 466 DO jk=1,jpkm1 436 467 DO ji=1,jpi … … 438 469 END DO 439 470 END DO 440 441 471 442 472 spgv(:,nlcj-2)=0. … … 453 483 ENDIF 454 484 END DO 485 #else 486 spgv(:,nlcj-2)=va_b(:,nlcj-2) 487 #endif 455 488 456 489 DO jk=1,jpkm1 … … 483 516 DO jk=1,jpkm1 484 517 DO ji=1,jpi 485 ua(ji,nlcj-1,jk) = (zua(ji,nlcj-1,jk)/( rhoy*e2u(ji,nlcj-1)))*umask(ji,nlcj-1,jk)518 ua(ji,nlcj-1,jk) = (zua(ji,nlcj-1,jk)/(zrhoy*e2u(ji,nlcj-1)))*umask(ji,nlcj-1,jk) 486 519 ua(ji,nlcj-1,jk) = ua(ji,nlcj-1,jk) / fse3u(ji,nlcj-1,jk) 487 520 END DO … … 495 528 END SUBROUTINE Agrif_dyn 496 529 530 SUBROUTINE Agrif_dyn_ts( kt, jn ) 531 !!---------------------------------------------------------------------- 532 !! *** ROUTINE Agrif_dyn_ts *** 533 !!---------------------------------------------------------------------- 534 !! 535 INTEGER, INTENT(in) :: kt, jn 536 !! 537 INTEGER :: ji, jj 538 REAL(wp) :: zrhox, zrhoy 539 REAL(wp), POINTER, DIMENSION(:,:) :: spgv1, spgu1 540 REAL(wp), POINTER, DIMENSION(:,:) :: zunb, zvnb, zsshn 541 !!---------------------------------------------------------------------- 542 543 IF( Agrif_Root() ) RETURN 544 545 IF ((kt==nit000).AND.(jn==1)) THEN 546 ALLOCATE( ubdy_w(jpj), vbdy_w(jpj), hbdy_w(jpj)) 547 ALLOCATE( ubdy_e(jpj), vbdy_e(jpj), hbdy_e(jpj)) 548 ALLOCATE( ubdy_n(jpi), vbdy_n(jpi), hbdy_n(jpi)) 549 ALLOCATE( ubdy_s(jpi), vbdy_s(jpi), hbdy_s(jpi)) 550 ENDIF 551 552 IF (jn==1) THEN 553 ! Fill boundary arrays at each baroclinic step 554 ! with Parent grid barotropic fluxes and sea level 555 ! 556 CALL wrk_alloc( jpi, jpj, zunb, zvnb, zsshn ) 557 558 zrhox = Agrif_Rhox() 559 zrhoy = Agrif_Rhoy() 560 561 !alt Agrif_SpecialValue = 0.e0 562 !alt Agrif_UseSpecialValue = .TRUE. 563 !alt CALL Agrif_Bc_variable(zsshn, sshn_id, procname=interpsshn ) 564 !alt Agrif_UseSpecialValue = .FALSE. 565 566 Agrif_SpecialValue=0. 567 Agrif_UseSpecialValue = ln_spc_dyn 568 zunb(:,:) = 0._wp ; zvnb(:,:) = 0._wp 569 CALL Agrif_Bc_variable(zunb,unb_id,procname=interpunb) 570 CALL Agrif_Bc_variable(zvnb,vnb_id,procname=interpvnb) 571 Agrif_UseSpecialValue = .FALSE. 572 573 IF((nbondi == -1).OR.(nbondi == 2)) THEN 574 DO jj=1,jpj 575 ubdy_w(jj) = (zunb(2,jj)/(zrhoy*e2u(2,jj))) 576 vbdy_w(jj) = (zvnb(2,jj)/(zrhox*e1v(2,jj))) 577 hbdy_w(jj) = zsshn(2,jj) 578 END DO 579 ENDIF 580 581 IF((nbondi == 1).OR.(nbondi == 2)) THEN 582 DO jj=1,jpj 583 ubdy_e(jj) = zunb(nlci-2,jj)/(zrhoy*e2u(nlci-2,jj)) 584 vbdy_e(jj) = zvnb(nlci-1,jj)/(zrhox*e1v(nlci-1,jj)) 585 hbdy_e(jj) = zsshn(nlci-1,jj) 586 END DO 587 ENDIF 588 589 IF((nbondj == -1).OR.(nbondj == 2)) THEN 590 DO ji=1,jpi 591 ubdy_s(ji) = zunb(ji,2)/(zrhoy*e2u(ji,2)) 592 vbdy_s(ji) = zvnb(ji,2)/(zrhox*e1v(ji,2)) 593 hbdy_s(ji) = zsshn(ji,2) 594 END DO 595 ENDIF 596 597 IF((nbondj == 1).OR.(nbondj == 2)) THEN 598 DO ji=1,jpi 599 ubdy_n(ji) = zunb(ji,nlcj-1)/(zrhoy*e2u(ji,nlcj-1)) 600 vbdy_n(ji) = zvnb(ji,nlcj-2)/(zrhox*e1v(ji,nlcj-2)) 601 hbdy_n(ji) = zsshn(ji,nlcj-1) 602 END DO 603 ENDIF 604 605 CALL wrk_dealloc( jpi, jpj, zunb, zvnb, zsshn ) 606 ENDIF ! jn==1 607 608 ! Then update velocities at each barotropic time step 609 IF((nbondi == -1).OR.(nbondi == 2)) THEN 610 DO jj=1,jpj 611 va_e(2,jj) = vbdy_w(jj) * hvr_e(2,jj) 612 ! Specified fluxes: 613 ua_e(2,jj) = ubdy_w(jj) * hur_e(2,jj) 614 ! Characteristics method: 615 !alt ua_e(2,jj) = 0.5_wp * ( ubdy_w(jj) * hur_e(2,jj) + ua_e(3,jj) & 616 !alt & - sqrt(grav * hur_e(2,jj)) * (sshn_e(3,jj) - hbdy_w(jj)) ) 617 END DO 618 ENDIF 619 620 IF((nbondi == 1).OR.(nbondi == 2)) THEN 621 DO jj=1,jpj 622 va_e(nlci-1,jj) = vbdy_e(jj) * hvr_e(nlci-1,jj) 623 ! Specified fluxes: 624 ua_e(nlci-2,jj) = ubdy_e(jj) * hur_e(nlci-2,jj) 625 ! Characteristics method: 626 !alt ua_e(nlci-2,jj) = 0.5_wp * ( ubdy_e(jj) * hur_e(nlci-2,jj) + ua_e(nlci-3,jj) & 627 !alt & + sqrt(grav * hur_e(nlci-2,jj)) * (sshn_e(nlci-2,jj) - hbdy_e(jj)) ) 628 END DO 629 ENDIF 630 631 IF((nbondj == -1).OR.(nbondj == 2)) THEN 632 DO ji=1,jpi 633 ua_e(ji,2) = ubdy_s(ji) * hur_e(ji,2) 634 ! Specified fluxes: 635 va_e(ji,2) = vbdy_s(ji) * hvr_e(ji,2) 636 ! Characteristics method: 637 !alt va_e(ji,2) = 0.5_wp * ( vbdy_s(ji) * hvr_e(ji,2) + va_e(ji,3) & 638 !alt & - sqrt(grav * hvr_e(ji,2)) * (sshn_e(ji,3) - hbdy_s(ji)) ) 639 END DO 640 ENDIF 641 642 IF((nbondj == 1).OR.(nbondj == 2)) THEN 643 DO ji=1,jpi 644 ua_e(ji,nlcj-1) = ubdy_n(ji) * hur_e(ji,nlcj-1) 645 ! Specified fluxes: 646 va_e(ji,nlcj-2) = vbdy_n(ji) * hvr_e(ji,nlcj-2) 647 ! Characteristics method: 648 !alt va_e(ji,nlcj-2) = 0.5_wp * ( vbdy_n(ji) * hvr_e(ji,nlcj-2) + va_e(ji,nlcj-3) & 649 !alt & + sqrt(grav * hvr_e(ji,nlcj-2)) * (sshn_e(ji,nlcj-2) - hbdy_n(ji)) ) 650 END DO 651 ENDIF 652 ! 653 END SUBROUTINE Agrif_dyn_ts 497 654 498 655 SUBROUTINE Agrif_ssh( kt ) … … 518 675 519 676 IF((nbondj == -1).OR.(nbondj == 2)) THEN 520 ssha(:,2)=ssh n(:,3)521 sshn(:,2)=ssh b(:,3)677 ssha(:,2)=ssha(:,3) 678 sshn(:,2)=sshn(:,3) 522 679 ENDIF 523 680 524 681 IF((nbondj == 1).OR.(nbondj == 2)) THEN 525 682 ssha(:,nlcj-1)=ssha(:,nlcj-2) 526 ssh a(:,nlcj-1)=sshn(:,nlcj-2)683 sshn(:,nlcj-1)=sshn(:,nlcj-2) 527 684 ENDIF 528 685 529 686 END SUBROUTINE Agrif_ssh 530 687 688 SUBROUTINE Agrif_ssh_ts( kt ) 689 !!---------------------------------------------------------------------- 690 !! *** ROUTINE Agrif_ssh_ts *** 691 !!---------------------------------------------------------------------- 692 INTEGER, INTENT(in) :: kt 693 !! 694 !!---------------------------------------------------------------------- 695 696 IF((nbondi == -1).OR.(nbondi == 2)) THEN 697 ssha_e(2,:) = ssha_e(3,:) 698 ENDIF 699 700 IF((nbondi == 1).OR.(nbondi == 2)) THEN 701 ssha_e(nlci-1,:) = ssha_e(nlci-2,:) 702 ENDIF 703 704 IF((nbondj == -1).OR.(nbondj == 2)) THEN 705 ssha_e(:,2) = ssha_e(:,3) 706 ENDIF 707 708 IF((nbondj == 1).OR.(nbondj == 2)) THEN 709 ssha_e(:,nlcj-1) = ssha_e(:,nlcj-2) 710 ENDIF 711 712 END SUBROUTINE Agrif_ssh_ts 713 714 SUBROUTINE interpsshn(tabres,i1,i2,j1,j2) 715 !!---------------------------------------------------------------------- 716 !! *** ROUTINE interpsshn *** 717 !!---------------------------------------------------------------------- 718 INTEGER, INTENT(in) :: i1,i2,j1,j2 719 REAL(wp), DIMENSION(i1:i2,j1:j2), INTENT(inout) :: tabres 720 !! 721 INTEGER :: ji,jj 722 !!---------------------------------------------------------------------- 723 724 tabres(i1:i2,j1:j2) = sshn(i1:i2,j1:j2) 725 726 END SUBROUTINE interpsshn 531 727 532 728 SUBROUTINE interpu(tabres,i1,i2,j1,j2,k1,k2) … … 611 807 612 808 END SUBROUTINE interpv2d 809 810 SUBROUTINE interpunb(tabres,i1,i2,j1,j2) 811 !!---------------------------------------------------------------------- 812 !! *** ROUTINE interpunb *** 813 !!---------------------------------------------------------------------- 814 INTEGER, INTENT(in) :: i1,i2,j1,j2 815 REAL(wp), DIMENSION(i1:i2,j1:j2), INTENT(inout) :: tabres 816 !! 817 INTEGER :: ji,jj,jk 818 !!---------------------------------------------------------------------- 819 820 tabres(:,:) = 0.e0 821 DO jk=1,jpkm1 822 DO jj=j1,j2 823 DO ji=i1,i2 824 tabres(ji,jj) = tabres(ji,jj) + e2u(ji,jj) * un(ji,jj,jk) & 825 * umask(ji,jj,jk) * fse3u(ji,jj,jk) 826 END DO 827 END DO 828 END DO 829 830 END SUBROUTINE interpunb 831 832 SUBROUTINE interpvnb(tabres,i1,i2,j1,j2) 833 !!---------------------------------------------------------------------- 834 !! *** ROUTINE interpvnb *** 835 !!---------------------------------------------------------------------- 836 INTEGER, INTENT(in) :: i1,i2,j1,j2 837 REAL(wp), DIMENSION(i1:i2,j1:j2), INTENT(inout) :: tabres 838 !! 839 INTEGER :: ji,jj,jk 840 !!---------------------------------------------------------------------- 841 842 tabres(:,:) = 0.e0 843 DO jk=1,jpkm1 844 DO jj=j1,j2 845 DO ji=i1,i2 846 tabres(ji,jj) = tabres(ji,jj) + e1v(ji,jj) * vn(ji,jj,jk) & 847 * vmask(ji,jj,jk) * fse3v(ji,jj,jk) 848 END DO 849 END DO 850 END DO 851 852 END SUBROUTINE interpvnb 613 853 614 854 #else
Note: See TracChangeset
for help on using the changeset viewer.