Changeset 2163
- Timestamp:
- 2010-10-06T11:30:20+02:00 (14 years ago)
- Location:
- branches/DEV_r1924_nocs_latphys/NEMO/OPA_SRC
- Files:
-
- 4 edited
Legend:
- Unmodified
- Added
- Removed
-
branches/DEV_r1924_nocs_latphys/NEMO/OPA_SRC/LDF/ldfslp.F90
r2142 r2163 41 41 REAL(wp), PUBLIC, DIMENSION(jpi,jpj,jpk) :: alpha, beta !: alpha,beta at T points 42 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 43 45 44 46 REAL(wp), DIMENSION(jpi,jpj,jpk) :: omlmask ! mask of the surface mixed layer at T-pt … … 382 384 zdenr, zrhotmp, zdndt, zdddt, & ! " " 383 385 zdnds, zddds, znum, zden, & ! " " 384 zs xe, za_sxe, zslopec, zdsloper,& ! " "386 zslope, za_sxe, zslopec, zdsloper,& ! " " 385 387 zfact, zepsln, zatempw,zatempu,zatempv, & ! " " 386 ze1ur,ze2vr,ze3wr,zdxt,zdxs,zdyt,zdys,zdzt,zdzs,zvolf 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 387 394 REAL(wp), DIMENSION(jpi,jpj,jpk) :: & 388 395 zdkt,zdks 396 REAL(wp), DIMENSION(jpi,jpj) :: & 397 zr_ml_basew 389 398 !!---------------------------------------------------------------------- 390 399 391 400 ! 0. Local constant initialization 392 401 ! -------------------------------- 393 394 zslopec=0.004 395 zdsloper=1000.0 402 zr_slpmax = 1.0_wp/slpmax 403 404 ! zslopec=0.004 405 ! zdsloper=1000.0 396 406 zepsln=1e-25 397 407 … … 495 505 496 506 507 ! --------------------------------------- 508 ! 1. Horizontal tracer gradients at T-level jk 509 ! --------------------------------------- 497 510 DO jk = 1, jpkm1 498 511 DO jj = 1, jpjm1 … … 531 544 ENDIF 532 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, jpjm1 561 DO ji = 1, fs_jpim1 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 533 568 wslp2(:,:,:)=0.0 534 569 tfw(:,:,:) = 0.0 … … 542 577 ftvd(:,:,:) = 0.0 543 578 fsvd(:,:,:) = 0.0 544 545 DO jk = 2, jpk 546 ! 547 ! 1. Vertical tracer gradient at level jk 548 ! --------------------------------------- 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 ! interior (1=<jk=<jpk-1) 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) 557 587 DO ip=0,1 558 588 DO kp=0,1 559 589 560 590 DO jk = 1, jpkm1 561 562 591 DO jj = 1, jpjm1 563 564 592 DO ji = 1, fs_jpim1 ! vector opt. 565 593 … … 571 599 zdzt=zdkt(ji+ip,jj,jk+kp)*ze3wr 572 600 zdzs=zdks(ji+ip,jj,jk+kp)*ze3wr 573 ! Calculate the density gradient terms 574 zsxe=-(alpha(ji+ip,jj,jk)*zdxt+beta(ji+ip,jj,jk)*zdxs) / & 575 (-ABS(alpha(ji+ip,jj,jk)*zdzt+beta(ji+ip,jj,jk)*zdzs)-zepsln) 576 ! Check value and taper if necessary 577 za_sxe = MIN(1.0,(abs(zsxe)-zslopec)*zdsloper) 578 za_sxe = MAX(-1.0,za_sxe) 579 zfact = 0.5*(1.0-za_sxe)*umask(ji,jj,jk)*tmask(ji+ip,jj,jk+kp) 580 zsxe = zsxe*zfact 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,jj) 620 zsx_ml_base(ji,jj,ip,kp)=zsx(ji,jj,jk+1-kp,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 zfact = 1 - 1/(1 + (jk+1-kp)/nmln(ji,jj)) 637 zsx(ji,jj,jk,ip,kp) = zfact*zsx(ji,jj,jk,ip,kp) + & 638 & (1.0_wp-zfact)*(fsdepw(ji,jj,jk+1-kp)*zr_ml_basew(ji,jj))*zsx_ml_base(ji,jj,ip,kp) 639 END DO 640 641 END DO 642 643 END DO 644 END DO 645 END DO 646 647 ! Use zsx to calculate fluxes and save averaged slopes psix_eiv at psi-points 648 DO ip=0,1 649 DO kp=0,1 650 651 DO jk = 1, jpkm1 652 653 DO jj = 1, jpjm1 654 655 DO ji = 1, fs_jpim1 ! vector opt. 656 657 ze1ur=1.0/e1u(ji,jj) 658 zdxt=zdit(ji,jj,jk)*ze1ur 659 zdxs=zdis(ji,jj,jk)*ze1ur 660 661 ze3wr=1.0/fse3w(ji+ip,jj,jk+kp) 662 zdzt=zdkt(ji+ip,jj,jk+kp)*ze3wr 663 zdzs=zdks(ji+ip,jj,jk+kp)*ze3wr 664 zslope=zsx(ji+ip,jj,jk,1-ip,kp) 581 665 582 666 zvolf = 0.25_wp*e1u(ji,jj)*e2u(ji,jj)*fse3u(ji,jj,jk) 583 667 584 ftu(ji,jj,jk)= ftu(ji,jj,jk)+zs xe*zdzt*zvolf*ze1ur585 fsu(ji,jj,jk)= fsu(ji,jj,jk)+zs xe*zdzs*zvolf*ze1ur668 ftu(ji,jj,jk)= ftu(ji,jj,jk)+zslope*zdzt*zvolf*ze1ur 669 fsu(ji,jj,jk)= fsu(ji,jj,jk)+zslope*zdzs*zvolf*ze1ur 586 670 ftud(ji,jj,jk)=ftud(ji,jj,jk)+fsahtu(ji,jj,jk)*zdxt*zvolf*ze1ur 587 671 fsud(ji,jj,jk)=fsud(ji,jj,jk)+fsahtu(ji,jj,jk)*zdxs*zvolf*ze1ur 588 tfw(ji+ip,jj,jk+kp)=tfw(ji+ip,jj,jk+kp)+(zvolf*ze3wr)*zs xe*zdxt589 sfw(ji+ip,jj,jk+kp)=sfw(ji+ip,jj,jk+kp)+(zvolf*ze3wr)*zs xe*zdxs672 tfw(ji+ip,jj,jk+kp)=tfw(ji+ip,jj,jk+kp)+(zvolf*ze3wr)*zslope*zdxt 673 sfw(ji+ip,jj,jk+kp)=sfw(ji+ip,jj,jk+kp)+(zvolf*ze3wr)*zslope*zdxs 590 674 wslp2(ji+ip,jj,jk+kp)=wslp2(ji+ip,jj,jk+kp)+ & 591 & ((zvolf*ze3wr)/(e1t(ji+ip,jj)*e2t(ji+ip,jj)))*(zsxe)**2 675 & ((zvolf*ze3wr)/(e1t(ji+ip,jj)*e2t(ji+ip,jj)))*(zslope)**2 676 psix_eiv(ji+ip,jj,jk+kp) = psix_eiv(ji+ip,jj,jk+kp) + 0.25_wp*zslope 592 677 593 678 END DO … … 598 683 END DO 599 684 685 ! ---------------------------------------------------------------------- 686 ! y-z plane 687 ! ---------------------------------------------------------------------- 688 689 ! calculate limited triad y-slopes zsy in interior (1=<jk=<jpk-1) 600 690 DO jp=0,1 601 691 DO kp=0,1 … … 603 693 DO jk = 1, jpkm1 604 694 DO jj = 1, jpjm1 605 DO ji = fs_2, fs_jpim1 ! vector opt.695 DO ji = 1, fs_jpim1 ! vector opt. 606 696 607 697 ze2vr=1.0/e2v(ji,jj) … … 612 702 zdzt=zdkt(ji,jj+jp,jk+kp)*ze3wr 613 703 zdzs=zdks(ji,jj+jp,jk+kp)*ze3wr 614 ! Calculate the density gradient terms 615 zsxe=-(alpha(ji,jj+jp,jk)*zdyt+beta(ji,jj+jp,jk)*zdys) / & 616 (-ABS(alpha(ji,jj+jp,jk)*zdzt+beta(ji,jj+jp,jk)*zdzs)-zepsln) 617 ! Check value and taper if necessary 618 za_sxe = MIN(1.0,(abs(zsxe)-zslopec)*zdsloper) 619 za_sxe = MAX(-1.0,za_sxe) 620 zfact = 0.5*(1.0-za_sxe)*vmask(ji,jj,jk)*tmask(ji,jj+jp,jk+kp) 621 zsxe = zsxe*zfact 704 ! Calculate the horizontal and vertical density gradient 705 zdyrho = alpha(ji,jj+jp,jk)*zdyt+beta(ji,jj+jp,jk)*zdys 706 zabs_dzrho = ABS(alpha(ji,jj+jp,jk)*zdzt+beta(ji,jj+jp,jk)*zdzs)+zepsln 707 ! Limit by slpmax, and mask by psi-point 708 zsy(ji,jj+jp,jk,1-jp,kp) = vmask(ji,jj,jk+kp) & 709 & *zdyrho/MAX(zabs_dzrho,zr_slpmax*ABS(zdyrho)) 710 END DO 711 END DO 712 END DO 713 714 END DO 715 END DO 716 717 ! calculate slope of triad immediately below mixed-layer base 718 DO jp=0,1 719 DO kp=0,1 720 DO jj = 1, jpjm1 721 DO ji = 1, fs_jpim1 722 jk = nmln(ji,jj) 723 zsy_ml_base(ji,jj,jp,kp)=zsy(ji,jj,jk+1-kp,jp,kp) 724 END DO 725 END DO 726 END DO 727 END DO 728 729 ! Below ML use limited zsy as is 730 ! Inside ML replace by linearly reducing zsy_ml_base towards surface 731 DO jp=0,1 732 DO kp=0,1 733 734 DO jk = 1, jpkm1 735 736 DO jj = 1, jpjm1 737 738 DO ji = 1, fs_jpim1 ! vector opt. 739 zfact = 1 - 1/(1 + (jk+1-kp)/nmln(ji,jj)) 740 zsy(ji,jj,jk,jp,kp) = zfact*zsy(ji,jj,jk,jp,kp) +& 741 & (1.0_wp-zfact)*(fsdepw(ji,jj,jk+1-kp)*zr_ml_basew(ji,jj))*zsy_ml_base(ji,jj,jp,kp) 742 END DO 743 744 END DO 745 746 END DO 747 END DO 748 END DO 749 750 ! Use zsy to calculate fluxes and save averaged slopes psiy_eiv at psi-points 751 DO jp=0,1 752 DO kp=0,1 753 754 DO jk = 1, jpkm1 755 756 DO jj = 1, jpjm1 757 758 DO ji = 1, fs_jpim1 ! vector opt. 759 760 ze2vr=1.0/e2v(ji,jj) 761 zdyt=zdjt(ji,jj,jk)*ze2vr 762 zdys=zdjs(ji,jj,jk)*ze2vr 763 764 ze3wr=1.0/fse3w(ji,jj+jp,jk+kp) 765 zdzt=zdkt(ji,jj+jp,jk+kp)*ze3wr 766 zdzs=zdks(ji,jj+jp,jk+kp)*ze3wr 767 zslope=zsy(ji,jj+jp,jk,1-jp,kp) 622 768 623 769 zvolf = 0.25_wp*e1v(ji,jj)*e2v(ji,jj)*fse3v(ji,jj,jk) 624 770 625 ftv(ji,jj,jk)= ftv(ji,jj,jk)+zs xe*zdzt*zvolf*ze2vr626 fsv(ji,jj,jk)= fsv(ji,jj,jk)+zs xe*zdzs*zvolf*ze2vr771 ftv(ji,jj,jk)= ftv(ji,jj,jk)+zslope*zdzt*zvolf*ze2vr 772 fsv(ji,jj,jk)= fsv(ji,jj,jk)+zslope*zdzs*zvolf*ze2vr 627 773 ftvd(ji,jj,jk)=ftvd(ji,jj,jk)+fsahtv(ji,jj,jk)*zdyt*zvolf*ze2vr 628 774 fsvd(ji,jj,jk)=fsvd(ji,jj,jk)+fsahtv(ji,jj,jk)*zdys*zvolf*ze2vr 629 tfw(ji,jj+jp,jk+kp)=tfw(ji,jj+jp,jk+kp)+(zvolf*ze3wr)*zs xe*zdyt630 sfw(ji,jj+jp,jk+kp)=sfw(ji,jj+jp,jk+kp)+(zvolf*ze3wr)*zs xe*zdys775 tfw(ji,jj+jp,jk+kp)=tfw(ji,jj+jp,jk+kp)+(zvolf*ze3wr)*zslope*zdyt 776 sfw(ji,jj+jp,jk+kp)=sfw(ji,jj+jp,jk+kp)+(zvolf*ze3wr)*zslope*zdys 631 777 wslp2(ji,jj+jp,jk+kp)=wslp2(ji,jj+jp,jk+kp)+ & 632 & ((zvolf*ze3wr)/(e1t(ji,jj+jp)*e2t(ji,jj+jp)))*(zsxe)**2 778 & ((zvolf*ze3wr)/(e1t(ji,jj+jp)*e2t(ji,jj+jp)))*(zslope)**2 779 psiy_eiv(ji,jj+jp,jk+kp) = psiy_eiv(ji,jj+jp,jk+kp) + 0.25_wp*zslope 633 780 634 781 END DO … … 654 801 CALL lbc_lnk( ftvd, 'V', -1. ) 655 802 CALL lbc_lnk( fsvd, 'V', -1. ) 803 CALL lbc_lnk( psix_eiv, 'U', -1. ) 804 CALL lbc_lnk( psiy_eiv, 'V', -1. ) 656 805 657 806 -
branches/DEV_r1924_nocs_latphys/NEMO/OPA_SRC/LDF/ldftra.F90
r2079 r2163 68 68 & ln_traldf_level, ln_traldf_hor , ln_traldf_iso, & 69 69 & ln_traldf_vis, ln_traldf_grif , & 70 & rn_aht_0 , rn_ahtb_0 , rn_aeiv_0 70 & rn_aht_0 , rn_ahtb_0 , rn_aeiv_0, & 71 & rn_slpmax 71 72 !!---------------------------------------------------------------------- 72 73 … … 92 93 ENDIF 93 94 95 slpmax = rn_slpmax 94 96 ! ! convert DOCTOR namelist names into OLD names 95 97 aht0 = rn_aht_0 -
branches/DEV_r1924_nocs_latphys/NEMO/OPA_SRC/LDF/ldftra_oce.F90
r2061 r2163 25 25 REAL(wp), PUBLIC :: rn_ahtb_0 = 0._wp !: lateral background eddy diffusivity (m2/s) 26 26 REAL(wp), PUBLIC :: rn_aeiv_0 = 2000._wp !: eddy induced velocity coefficient (m2/s) 27 REAL(wp), PUBLIC :: rn_slpmax = 0.01_wp !: slope limit 27 28 28 29 REAL(wp), PUBLIC :: aht0, ahtb0, aeiv0 !!: OLD namelist names … … 43 44 !! 'key_traldf_eiv' eddy induced velocity 44 45 !!---------------------------------------------------------------------- 45 LOGICAL, PUBLIC, PARAMETER :: lk_traldf_eiv = .TRUE. !: eddy induced velocity flag 46 LOGICAL, PUBLIC, PARAMETER :: lk_traldf_eiv = .TRUE. !: eddy induced velocity flag 47 REAL(wp), PUBLIC :: slpmax !: slope limit 46 48 47 49 # if defined key_traldf_c3d -
branches/DEV_r1924_nocs_latphys/NEMO/OPA_SRC/TRA/traldf_iso_grif.F90
r2060 r2163 178 178 END DO 179 179 ! 180 psix_eiv(:,:,:) = psix_eiv(:,:,:)*fsaeiu(:,:,:) 181 psiy_eiv(:,:,:) = psiy_eiv(:,:,:)*fsaeiv(:,:,:) 182 180 183 END SUBROUTINE tra_ldf_iso_grif 181 184
Note: See TracChangeset
for help on using the changeset viewer.