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