- Timestamp:
- 2010-10-11T18:24:37+02:00 (14 years ago)
- File:
-
- 1 edited
Legend:
- Unmodified
- Added
- Removed
-
branches/DEV_r2191_3partymerge2010/NEMO/OPA_SRC/LDF/ldfslp.F90
r1515 r2205 24 24 USE phycst ! physical constants 25 25 USE zdfmxl ! mixed layer depth 26 USE eosbn2 26 27 USE lbclnk ! ocean lateral boundary conditions (or mpp link) 27 28 USE in_out_manager ! I/O manager … … 32 33 33 34 PUBLIC ldf_slp ! routine called by step.F90 35 PUBLIC ldf_slp_grif ! " 34 36 35 37 LOGICAL , PUBLIC, PARAMETER :: lk_ldfslp = .TRUE. !: slopes flag 36 38 REAL(wp), PUBLIC, DIMENSION(jpi,jpj,jpk) :: uslp, wslpi !: i_slope at U- and W-points 37 39 REAL(wp), PUBLIC, DIMENSION(jpi,jpj,jpk) :: vslp, wslpj !: j-slope at V- and W-points 40 REAL(wp), PUBLIC, DIMENSION(jpi,jpj,jpk) :: wslp2 !: wslp**2 from Griffies quarter cells 41 REAL(wp), PUBLIC, DIMENSION(jpi,jpj,jpk) :: alpha, beta !: alpha,beta at T points 42 REAL(wp), PUBLIC, DIMENSION(jpi,jpj,jpk) :: tfw,sfw,ftu,fsu,ftv,fsv,ftud,fsud,ftvd,fsvd 43 REAL(wp), PUBLIC, DIMENSION(jpi,jpj,jpk) :: psix_eiv 44 REAL(wp), PUBLIC, DIMENSION(jpi,jpj,jpk) :: psiy_eiv 38 45 39 46 REAL(wp), DIMENSION(jpi,jpj,jpk) :: omlmask ! mask of the surface mixed layer at T-pt … … 43 50 !! * Substitutions 44 51 # include "domzgr_substitute.h90" 52 # include "ldftra_substitute.h90" 53 # include "ldfeiv_substitute.h90" 45 54 # include "vectopt_loop_substitute.h90" 46 55 !!---------------------------------------------------------------------- … … 305 314 END DO 306 315 307 308 316 ! III. Specific grid points 309 317 ! =========================== … … 340 348 ! 341 349 END SUBROUTINE ldf_slp 350 351 SUBROUTINE ldf_slp_grif ( kt ) 352 !!---------------------------------------------------------------------- 353 !! *** ROUTINE ldf_slp_grif *** 354 !! 355 !! ** Purpose : Compute the squared slopes of neutral surfaces (slope 356 !! of iso-pycnal surfaces referenced locally) ('key_traldfiso') 357 !! at W-points using the Griffies quarter-cells. Also calculates 358 !! alpha and beta at T-points for use in the Griffies isopycnal 359 !! scheme. 360 !! 361 !! ** Method : The slope in the i-direction is computed at U- and 362 !! W-points (uslp, wslpi) and the slope in the j-direction is 363 !! computed at V- and W-points (vslp, wslpj). 364 !! 365 !! ** Action : - alpha, beta 366 !! wslp2 squared slope of neutral surfaces at w-points. 367 !! 368 !! History : 369 !! 9.0 ! 06-10 (C. Harris) New subroutine 370 !!---------------------------------------------------------------------- 371 !! * Modules used 372 USE oce , zdit => ua, & ! use ua as workspace 373 zdis => va, & ! use va as workspace 374 zdjt => ta, & ! use ta as workspace 375 zdjs => sa ! use sa as workspace 376 !! * Arguments 377 INTEGER, INTENT( in ) :: kt ! ocean time-step index 378 379 !! * Local declarations 380 INTEGER :: ji, jj, jk, ip, jp, kp ! dummy loop indices 381 INTEGER :: iku, ikv ! temporary integer 382 REAL(wp) :: & 383 zt, zs, zh, zt2, zsp5, zp1t1, & ! temporary scalars 384 zdenr, zrhotmp, zdndt, zdddt, & ! " " 385 zdnds, zddds, znum, zden, & ! " " 386 zslope, za_sxe, zslopec, zdsloper,& ! " " 387 zfact, zepsln, zatempw,zatempu,zatempv, & ! " " 388 ze1ur,ze2vr,ze3wr,zdxt,zdxs,zdyt,zdys,zdzt,zdzs,zvolf,& 389 zr_slpmax,zdxrho,zdyrho,zabs_dzrho 390 REAL(wp), DIMENSION(jpi,jpj,jpk,0:1,0:1) :: & 391 zsx,zsy 392 REAL(wp), DIMENSION(jpi,jpj,0:1,0:1) :: & 393 zsx_ml_base,zsy_ml_base 394 REAL(wp), DIMENSION(jpi,jpj,jpk) :: & 395 zdkt,zdks 396 REAL(wp), DIMENSION(jpi,jpj) :: & 397 zr_ml_basew 398 !!---------------------------------------------------------------------- 399 400 ! 0. Local constant initialization 401 ! -------------------------------- 402 zr_slpmax = 1.0_wp/slpmax 403 404 ! zslopec=0.004 405 ! zdsloper=1000.0 406 zepsln=1e-25 407 408 SELECT CASE ( nn_eos ) 409 410 CASE ( 0 ) ! Jackett and McDougall (1994) formulation 411 412 ! ! =============== 413 DO jk = 1, jpk ! Horizontal slab 414 ! ! =============== 415 DO jj = 1, jpjm1 416 DO ji = 1, fs_jpim1 417 zt = tb(ji,jj,jk) ! potential temperature 418 zs = sb(ji,jj,jk) - 35.0 ! salinity anomaly (s-35) 419 zh = fsdept(ji,jj,jk) ! depth in meters 420 421 beta(ji,jj,jk) = ( ( -0.415613e-09 * zt + 0.555579e-07 ) * zt & 422 & - 0.301985e-05 ) * zt & 423 & + 0.785567e-03 & 424 & + ( 0.515032e-08 * zs & 425 & + 0.788212e-08 * zt - 0.356603e-06 ) * zs & 426 & +( ( 0.121551e-17 * zh & 427 & - 0.602281e-15 * zs & 428 & - 0.175379e-14 * zt + 0.176621e-12 ) * zh & 429 & + 0.408195e-10 * zs & 430 & + ( - 0.213127e-11 * zt + 0.192867e-09 ) * zt & 431 & - 0.121555e-07 ) * zh 432 433 alpha(ji,jj,jk) = - beta(ji,jj,jk) * & 434 & (((( - 0.255019e-07 * zt + 0.298357e-05 ) * zt & 435 & - 0.203814e-03 ) * zt & 436 & + 0.170907e-01 ) * zt & 437 & + 0.665157e-01 & 438 & + ( - 0.678662e-05 * zs & 439 & - 0.846960e-04 * zt + 0.378110e-02 ) * zs & 440 & + ( ( - 0.302285e-13 * zh & 441 & - 0.251520e-11 * zs & 442 & + 0.512857e-12 * zt * zt ) * zh & 443 & - 0.164759e-06 * zs & 444 & +( 0.791325e-08 * zt - 0.933746e-06 ) * zt & 445 & + 0.380374e-04 ) * zh) 446 447 ENDDO 448 ENDDO 449 ENDDO 450 451 CASE ( 1 ) 452 453 alpha(:,:,:)=-rn_alpha 454 beta(:,:,:)=0.0 455 456 CASE ( 2 ) 457 458 alpha(:,:,:)=-rn_alpha 459 beta (:,:,:)=rn_beta 460 461 CASE ( 3 ) 462 463 DO jk =1, jpk 464 DO jj = 1, jpjm1 465 DO ji = 1, fs_jpim1 466 zt = tb(ji,jj,jk) 467 zs = sb(ji,jj,jk) 468 zh = fsdept(ji,jj,jk) 469 zt2 = zt**2 470 zsp5 = sqrt(ABS(zs)) 471 zp1t1=zh*zt 472 znum=9.99843699e+02+zt*(7.35212840e+00+zt*(-5.45928211e-02+3.98476704e-04*zt)) & 473 +zs*(2.96938239e+00-7.23268813e-03*zt+2.12382341e-03*zs) & 474 +zh*(1.04004591e-02+1.03970529e-07*zt2+5.18761880e-06*zs+ & 475 zh*(-3.24041825e-08-1.23869360e-11*zt2)) 476 zden=1.00000000e+00+zt*(7.28606739e-03+zt*(-4.60835542e-05+zt*(3.68390573e-07+zt*1.80809186e-10))) & 477 +zs*(2.14691708e-03+zt*(-9.27062484e-06-1.78343643e-10*zt2)+zsp5*(4.76534122e-06+1.63410736e-09*zt2)) & 478 + zh*(5.30848875e-06+zh*zt*(-3.03175128e-16*zt2-1.27934137e-17*zh)) 479 zdenr=1.0/zden 480 zrhotmp=znum*zdenr 481 zdndt=7.35212840e+00+zt*(-1.091856422e-01+1.195430112e-03*zt)-7.23268813e-03*zs & 482 +zp1t1*(2.07941058e-07-2.4773872e-11*zh) 483 zdddt=7.28606739e-03+zt*(-9.21671084e-05+zt*(1.105171719e-06+7.23236744e-10*zt)) & 484 +zs*(-9.27062484e-06+zt*(-5.35030929e-10*zt+3.26821472e-09*zsp5)) & 485 +zh*zh*(-9.09525384e-16*zt2-1.27934137e-17*zh) 486 zdnds=2.96938239e+00-7.23268813e-03*zt+2*2.12382341e-03*zs+5.18761880e-06*zh 487 zddds=2.14691708e-03+zt*(-9.27062484e-06-1.78343643e-10*zt2)+zsp5*(7.14801183e-06+2.45116104e-09*zt2) 488 alpha(ji,jj,jk)=(zdndt-zrhotmp*zdddt)*zdenr 489 beta(ji,jj,jk)=zdenr*(zdnds-zrhotmp*zddds) 490 491 END DO 492 END DO 493 END DO 494 495 CASE DEFAULT 496 497 IF(lwp) WRITE(numout,cform_err) 498 IF(lwp) WRITE(numout,*) ' bad flag value for nn_eos = ', nn_eos 499 nstop = nstop + 1 500 501 END SELECT 502 503 CALL lbc_lnk( alpha, 'T', 1. ) 504 CALL lbc_lnk( beta, 'T', 1. ) 505 506 507 ! --------------------------------------- 508 ! 1. Horizontal tracer gradients at T-level jk 509 ! --------------------------------------- 510 DO jk = 1, jpkm1 511 DO jj = 1, jpjm1 512 DO ji = 1, fs_jpim1 ! vector opt. 513 ! i-gradient of T and S at jj 514 zdit (ji,jj,jk) = ( tb(ji+1,jj,jk)-tb(ji,jj,jk) ) * umask(ji,jj,jk) 515 zdis (ji,jj,jk) = ( sb(ji+1,jj,jk)-sb(ji,jj,jk) ) * umask(ji,jj,jk) 516 ! j-gradient of T and S at jj 517 zdjt (ji,jj,jk) = ( tb(ji,jj+1,jk)-tb(ji,jj,jk) ) * vmask(ji,jj,jk) 518 zdjs (ji,jj,jk) = ( sb(ji,jj+1,jk)-sb(ji,jj,jk) ) * vmask(ji,jj,jk) 519 END DO 520 END DO 521 END DO 522 523 IF( ln_zps ) THEN ! partial steps correction at the last level 524 # if defined key_vectopt_loop && ! defined key_mpp_omp 525 jj = 1 526 DO ji = 1, jpij-jpi ! vector opt. (forced unrolling) 527 # else 528 DO jj = 1, jpjm1 529 DO ji = 1, jpim1 530 # endif 531 ! last ocean level 532 iku = MIN( mbathy(ji,jj), mbathy(ji+1,jj ) ) - 1 533 ikv = MIN( mbathy(ji,jj), mbathy(ji ,jj+1) ) - 1 534 ! i-gradient of T and S 535 zdit (ji,jj,iku) = gtu(ji,jj) 536 zdis (ji,jj,iku) = gsu(ji,jj) 537 ! j-gradient of T and S 538 zdjt (ji,jj,ikv) = gtv(ji,jj) 539 zdjs (ji,jj,ikv) = gsv(ji,jj) 540 # if ! defined key_vectopt_loop || defined key_mpp_omp 541 END DO 542 # endif 543 END DO 544 ENDIF 545 546 ! --------------------------------------- 547 ! 1. Vertical tracer gradient at w-level jk 548 ! --------------------------------------- 549 DO jk = 2, jpk 550 zdkt(:,:,jk) = ( tb(:,:,jk-1) - tb(:,:,jk) ) * tmask(:,:,jk) 551 zdks(:,:,jk) = ( sb(:,:,jk-1) - sb(:,:,jk) ) * tmask(:,:,jk) 552 END DO 553 554 zdkt(:,:,1) = 0.0 555 zdks(:,:,1) = 0.0 556 557 ! --------------------------------------- 558 ! Depth of the w-point below ML base 559 ! --------------------------------------- 560 DO jj = 1, jpj 561 DO ji = 1, jpi 562 jk = nmln(ji,jj) 563 zr_ml_basew(ji,jj)=1.0/fsdepw(ji,jj,jk+1) 564 END DO 565 END DO 566 567 568 wslp2(:,:,:)=0.0 569 tfw(:,:,:) = 0.0 570 sfw(:,:,:) = 0.0 571 ftu(:,:,:) = 0.0 572 fsu(:,:,:) = 0.0 573 ftv(:,:,:) = 0.0 574 fsv(:,:,:) = 0.0 575 ftud(:,:,:) = 0.0 576 fsud(:,:,:) = 0.0 577 ftvd(:,:,:) = 0.0 578 fsvd(:,:,:) = 0.0 579 psix_eiv(:,:,:) = 0.0 580 psiy_eiv(:,:,:) = 0.0 581 582 ! ---------------------------------------------------------------------- 583 ! x-z plane 584 ! ---------------------------------------------------------------------- 585 586 ! calculate limited triad x-slopes zsx in interior (1=<jk=<jpk-1) 587 DO ip=0,1 588 DO kp=0,1 589 590 DO jk = 1, jpkm1 591 DO jj = 1, jpjm1 592 DO ji = 1, fs_jpim1 ! vector opt. 593 594 ze1ur=1.0/e1u(ji,jj) 595 zdxt=zdit(ji,jj,jk)*ze1ur 596 zdxs=zdis(ji,jj,jk)*ze1ur 597 598 ze3wr=1.0/fse3w(ji+ip,jj,jk+kp) 599 zdzt=zdkt(ji+ip,jj,jk+kp)*ze3wr 600 zdzs=zdks(ji+ip,jj,jk+kp)*ze3wr 601 ! Calculate the horizontal and vertical density gradient 602 zdxrho = alpha(ji+ip,jj,jk)*zdxt+beta(ji+ip,jj,jk)*zdxs 603 zabs_dzrho = ABS(alpha(ji+ip,jj,jk)*zdzt+beta(ji+ip,jj,jk)*zdzs)+zepsln 604 ! Limit by slpmax, and mask by psi-point 605 zsx(ji+ip,jj,jk,1-ip,kp) = umask(ji,jj,jk+kp) & 606 & *zdxrho/MAX(zabs_dzrho,zr_slpmax*ABS(zdxrho)) 607 END DO 608 END DO 609 END DO 610 611 END DO 612 END DO 613 614 ! calculate slope of triad immediately below mixed-layer base 615 DO ip=0,1 616 DO kp=0,1 617 DO jj = 1, jpjm1 618 DO ji = 1, fs_jpim1 619 jk = nmln(ji+ip,jj) 620 zsx_ml_base(ji+ip,jj,1-ip,kp)=zsx(ji+ip,jj,jk+1-kp,1-ip,kp) 621 END DO 622 END DO 623 END DO 624 END DO 625 626 ! Below ML use limited zsx as is 627 ! Inside ML replace by linearly reducing zsx_ml_base towards surface 628 DO ip=0,1 629 DO kp=0,1 630 631 DO jk = 1, jpkm1 632 633 DO jj = 1, jpjm1 634 635 DO ji = 1, fs_jpim1 ! vector opt. 636 ! k index of uppermost point(s) of triad is jk+kp-1 637 ! must be .ge. nmln(ji,jj) for zfact=1. 638 ! otherwise zfact=0. 639 zfact = 1 - 1/(1 + (jk+kp-1)/nmln(ji+ip,jj)) 640 zsx(ji+ip,jj,jk,1-ip,kp) = zfact*zsx(ji+ip,jj,jk,1-ip,kp) + & 641 & (1.0_wp-zfact)*(fsdepw(ji+ip,jj,jk+kp)*zr_ml_basew(ji+ip,jj))*zsx_ml_base(ji+ip,jj,1-ip,kp) 642 END DO 643 644 END DO 645 646 END DO 647 END DO 648 END DO 649 650 ! Use zsx to calculate fluxes and save averaged slopes psix_eiv at psi-points 651 DO ip=0,1 652 DO kp=0,1 653 654 DO jk = 1, jpkm1 655 656 DO jj = 1, jpjm1 657 658 DO ji = 1, fs_jpim1 ! vector opt. 659 660 ze1ur=1.0/e1u(ji,jj) 661 zdxt=zdit(ji,jj,jk)*ze1ur 662 zdxs=zdis(ji,jj,jk)*ze1ur 663 664 ze3wr=1.0/fse3w(ji+ip,jj,jk+kp) 665 zdzt=zdkt(ji+ip,jj,jk+kp)*ze3wr 666 zdzs=zdks(ji+ip,jj,jk+kp)*ze3wr 667 zslope=zsx(ji+ip,jj,jk,1-ip,kp) 668 669 zvolf = 0.25_wp*e1u(ji,jj)*e2u(ji,jj)*fse3u(ji,jj,jk) 670 671 ftu(ji,jj,jk)= ftu(ji,jj,jk)+zslope*zdzt*zvolf*ze1ur 672 fsu(ji,jj,jk)= fsu(ji,jj,jk)+zslope*zdzs*zvolf*ze1ur 673 ftud(ji,jj,jk)=ftud(ji,jj,jk)+fsahtu(ji,jj,jk)*zdxt*zvolf*ze1ur 674 fsud(ji,jj,jk)=fsud(ji,jj,jk)+fsahtu(ji,jj,jk)*zdxs*zvolf*ze1ur 675 tfw(ji+ip,jj,jk+kp)=tfw(ji+ip,jj,jk+kp)+(zvolf*ze3wr)*zslope*zdxt 676 sfw(ji+ip,jj,jk+kp)=sfw(ji+ip,jj,jk+kp)+(zvolf*ze3wr)*zslope*zdxs 677 wslp2(ji+ip,jj,jk+kp)=wslp2(ji+ip,jj,jk+kp)+ & 678 & ((zvolf*ze3wr)/(e1t(ji+ip,jj)*e2t(ji+ip,jj)))*(zslope)**2 679 psix_eiv(ji,jj,jk+kp) = psix_eiv(ji,jj,jk+kp) + 0.25_wp*zslope 680 681 END DO 682 END DO 683 684 END DO 685 END DO 686 END DO 687 688 ! ---------------------------------------------------------------------- 689 ! y-z plane 690 ! ---------------------------------------------------------------------- 691 692 ! calculate limited triad y-slopes zsy in interior (1=<jk=<jpk-1) 693 DO jp=0,1 694 DO kp=0,1 695 696 DO jk = 1, jpkm1 697 DO jj = 1, jpjm1 698 DO ji = 1, fs_jpim1 ! vector opt. 699 700 ze2vr=1.0/e2v(ji,jj) 701 zdyt=zdjt(ji,jj,jk)*ze2vr 702 zdys=zdjs(ji,jj,jk)*ze2vr 703 704 ze3wr=1.0/fse3w(ji,jj+jp,jk+kp) 705 zdzt=zdkt(ji,jj+jp,jk+kp)*ze3wr 706 zdzs=zdks(ji,jj+jp,jk+kp)*ze3wr 707 ! Calculate the horizontal and vertical density gradient 708 zdyrho = alpha(ji,jj+jp,jk)*zdyt+beta(ji,jj+jp,jk)*zdys 709 zabs_dzrho = ABS(alpha(ji,jj+jp,jk)*zdzt+beta(ji,jj+jp,jk)*zdzs)+zepsln 710 ! Limit by slpmax, and mask by psi-point 711 zsy(ji,jj+jp,jk,1-jp,kp) = vmask(ji,jj,jk+kp) & 712 & *zdyrho/MAX(zabs_dzrho,zr_slpmax*ABS(zdyrho)) 713 END DO 714 END DO 715 END DO 716 717 END DO 718 END DO 719 720 ! calculate slope of triad immediately below mixed-layer base 721 DO jp=0,1 722 DO kp=0,1 723 DO jj = 1, jpjm1 724 DO ji = 1, fs_jpim1 725 jk = nmln(ji,jj+jp) 726 zsy_ml_base(ji,jj+jp,1-jp,kp)=zsy(ji,jj+jp,jk+1-kp,1-jp,kp) 727 END DO 728 END DO 729 END DO 730 END DO 731 732 ! Below ML use limited zsy as is 733 ! Inside ML replace by linearly reducing zsy_ml_base towards surface 734 DO jp=0,1 735 DO kp=0,1 736 737 DO jk = 1, jpkm1 738 739 DO jj = 1, jpjm1 740 741 DO ji = 1, fs_jpim1 ! vector opt. 742 ! k index of uppermost point(s) of triad is jk+kp-1 743 ! must be .ge. nmln(ji,jj) for zfact=1. 744 ! otherwise zfact=0. 745 zfact = 1 - 1/(1 + (jk+kp-1)/nmln(ji,jj+jp)) 746 zsy(ji,jj+jp,jk,1-jp,kp) = zfact*zsy(ji,jj+jp,jk,1-jp,kp) + & 747 & (1.0_wp-zfact)*(fsdepw(ji,jj+jp,jk+kp)*zr_ml_basew(ji,jj+jp))*zsy_ml_base(ji,jj+jp,1-jp,kp) 748 END DO 749 750 END DO 751 752 END DO 753 END DO 754 END DO 755 756 ! Use zsy to calculate fluxes and save averaged slopes psiy_eiv at psi-points 757 DO jp=0,1 758 DO kp=0,1 759 760 DO jk = 1, jpkm1 761 762 DO jj = 1, jpjm1 763 764 DO ji = 1, fs_jpim1 ! vector opt. 765 766 ze2vr=1.0/e2v(ji,jj) 767 zdyt=zdjt(ji,jj,jk)*ze2vr 768 zdys=zdjs(ji,jj,jk)*ze2vr 769 770 ze3wr=1.0/fse3w(ji,jj+jp,jk+kp) 771 zdzt=zdkt(ji,jj+jp,jk+kp)*ze3wr 772 zdzs=zdks(ji,jj+jp,jk+kp)*ze3wr 773 zslope=zsy(ji,jj+jp,jk,1-jp,kp) 774 775 zvolf = 0.25_wp*e1v(ji,jj)*e2v(ji,jj)*fse3v(ji,jj,jk) 776 777 ftv(ji,jj,jk)= ftv(ji,jj,jk)+zslope*zdzt*zvolf*ze2vr 778 fsv(ji,jj,jk)= fsv(ji,jj,jk)+zslope*zdzs*zvolf*ze2vr 779 ftvd(ji,jj,jk)=ftvd(ji,jj,jk)+fsahtv(ji,jj,jk)*zdyt*zvolf*ze2vr 780 fsvd(ji,jj,jk)=fsvd(ji,jj,jk)+fsahtv(ji,jj,jk)*zdys*zvolf*ze2vr 781 tfw(ji,jj+jp,jk+kp)=tfw(ji,jj+jp,jk+kp)+(zvolf*ze3wr)*zslope*zdyt 782 sfw(ji,jj+jp,jk+kp)=sfw(ji,jj+jp,jk+kp)+(zvolf*ze3wr)*zslope*zdys 783 wslp2(ji,jj+jp,jk+kp)=wslp2(ji,jj+jp,jk+kp)+ & 784 & ((zvolf*ze3wr)/(e1t(ji,jj+jp)*e2t(ji,jj+jp)))*(zslope)**2 785 psiy_eiv(ji,jj,jk+kp) = psiy_eiv(ji,jj,jk+kp) + 0.25_wp*zslope 786 787 END DO 788 END DO 789 790 END DO 791 END DO 792 END DO 793 794 tfw(:,:,1)=0.0 795 sfw(:,:,1)=0.0 796 wslp2(:,:,1)=0.0 797 798 CALL lbc_lnk( wslp2, 'W', 1. ) 799 CALL lbc_lnk( tfw, 'W', 1. ) 800 CALL lbc_lnk( sfw, 'W', 1. ) 801 CALL lbc_lnk( ftu, 'U', -1. ) 802 CALL lbc_lnk( fsu, 'U', -1. ) 803 CALL lbc_lnk( ftv, 'V', -1. ) 804 CALL lbc_lnk( fsv, 'V', -1. ) 805 CALL lbc_lnk( ftud, 'U', -1. ) 806 CALL lbc_lnk( fsud, 'U', -1. ) 807 CALL lbc_lnk( ftvd, 'V', -1. ) 808 CALL lbc_lnk( fsvd, 'V', -1. ) 809 CALL lbc_lnk( psix_eiv, 'U', -1. ) 810 CALL lbc_lnk( psiy_eiv, 'V', -1. ) 811 812 813 END SUBROUTINE ldf_slp_grif 342 814 343 815
Note: See TracChangeset
for help on using the changeset viewer.