Changeset 10475
- Timestamp:
- 2019-01-08T19:00:51+01:00 (6 years ago)
- File:
-
- 1 edited
-
NEMO/trunk/src/ICE/icedyn_adv_umx.F90 (modified) (30 diffs)
Legend:
- Unmodified
- Added
- Removed
-
NEMO/trunk/src/ICE/icedyn_adv_umx.F90
r10446 r10475 193 193 zamsk = 0._wp 194 194 ! 195 zhvar(:,:,:) = pv_ip(:,:,:) * z1_ai (:,:,:)195 zhvar(:,:,:) = pv_ip(:,:,:) * z1_aip(:,:,:) 196 196 CALL adv_umx( zamsk, kn_umx, jt, kt, zdt, zudy, zvdx, zua_ho, zva_ho, zcu_box, zcv_box, zhvar, pv_ip ) ! mp volume 197 197 ENDIF … … 271 271 CASE ( 20 ) !== centered second order ==! 272 272 ! 273 CALL cen2( pamsk, jt, kt, pdt, pt, pu, pv, puc, pvc, ptc,zt_ups, zfu_ups, zfv_ups, zfu_ho, zfv_ho )273 CALL cen2( pamsk, jt, kt, pdt, pt, pu, pv, zt_ups, zfu_ups, zfv_ups, zfu_ho, zfv_ho ) 274 274 ! 275 275 CASE ( 1:5 ) !== 1st to 5th order ULTIMATE-MACHO scheme ==! 276 276 ! 277 CALL macho( pamsk, kn_umx, jt, kt, pdt, pt, pu, pv, pu c, pvc, pubox, pvbox, ptc, zt_ups, zfu_ups, zfv_ups, zfu_ho, zfv_ho )277 CALL macho( pamsk, kn_umx, jt, kt, pdt, pt, pu, pv, pubox, pvbox, zt_ups, zfu_ups, zfv_ups, zfu_ho, zfv_ho ) 278 278 ! 279 279 END SELECT … … 282 282 ! new fluxes = u*H * u*a / u 283 283 ! ---------------------------- 284 IF( pamsk == 0. ) THEN284 IF( pamsk == 0._wp ) THEN 285 285 DO jl = 1, jpl 286 286 DO jj = 1, jpjm1 … … 440 440 & + pv (ji,jj ) - pv (ji ,jj-1 ) ) * pt(ji,jj,jl) * (1.-pamsk) 441 441 ! 442 pt_ups(ji,jj,jl) = ( pt 442 pt_ups(ji,jj,jl) = ( pt(ji,jj,jl) + ztra * pdt * r1_e1e2t(ji,jj) ) * tmask(ji,jj,1) 443 443 END DO 444 444 END DO … … 449 449 450 450 451 SUBROUTINE cen2( pamsk, jt, kt, pdt, pt, pu, pv, p uc, pvc, ptc, pt_ups, pfu_ups, pfv_ups, pfu_ho, pfv_ho )451 SUBROUTINE cen2( pamsk, jt, kt, pdt, pt, pu, pv, pt_ups, pfu_ups, pfv_ups, pfu_ho, pfv_ho ) 452 452 !!--------------------------------------------------------------------- 453 453 !! *** ROUTINE cen2 *** … … 462 462 REAL(wp), DIMENSION(:,:,:) , INTENT(in ) :: pt ! tracer fields 463 463 REAL(wp), DIMENSION(:,: ) , INTENT(in ) :: pu, pv ! 2 ice velocity components 464 REAL(wp), DIMENSION(:,:,:) , INTENT(in ) :: puc, pvc ! 2 ice velocity * A components465 REAL(wp), DIMENSION(:,:,:) , INTENT(in ) :: ptc ! tracer content at before time step466 464 REAL(wp), DIMENSION(:,:,:) , INTENT(in ) :: pt_ups ! upstream guess of tracer 467 465 REAL(wp), DIMENSION(:,:,:) , INTENT(in ) :: pfu_ups, pfv_ups ! upstream fluxes … … 478 476 DO jj = 1, jpjm1 479 477 DO ji = 1, fs_jpim1 480 pfu_ho(ji,jj,jl) = 0.5 * pu(ji,jj) * ( pt(ji,jj,jl) + pt(ji+1,jj,jl) ) 481 pfv_ho(ji,jj,jl) = 0.5 * pv(ji,jj) * ( pt(ji,jj,jl) + pt(ji,jj+1,jl) ) 482 END DO 483 END DO 484 END DO 478 pfu_ho(ji,jj,jl) = 0.5_wp * pu(ji,jj) * ( pt(ji,jj,jl) + pt(ji+1,jj ,jl) ) 479 pfv_ho(ji,jj,jl) = 0.5_wp * pv(ji,jj) * ( pt(ji,jj,jl) + pt(ji ,jj+1,jl) ) 480 END DO 481 END DO 482 END DO 483 ! 485 484 IF ( kn_limiter == 1 ) THEN 486 485 CALL nonosc( pamsk, pdt, pu, pv, pt, pt_ups, pfu_ups, pfv_ups, pfu_ho, pfv_ho ) … … 497 496 DO jj = 1, jpjm1 498 497 DO ji = 1, fs_jpim1 499 pfu_ho(ji,jj,jl) = 0.5 * pu(ji,jj) * ( pt(ji,jj,jl) + pt(ji+1,jj,jl) )498 pfu_ho(ji,jj,jl) = 0.5_wp * pu(ji,jj) * ( pt(ji,jj,jl) + pt(ji+1,jj,jl) ) 500 499 END DO 501 500 END DO … … 518 517 DO jj = 1, jpjm1 519 518 DO ji = 1, fs_jpim1 520 pfv_ho(ji,jj,jl) = 0.5 * pv(ji,jj) * ( zpt(ji,jj,jl) + zpt(ji,jj+1,jl) )519 pfv_ho(ji,jj,jl) = 0.5_wp * pv(ji,jj) * ( zpt(ji,jj,jl) + zpt(ji,jj+1,jl) ) 521 520 END DO 522 521 END DO … … 529 528 DO jj = 1, jpjm1 530 529 DO ji = 1, fs_jpim1 531 pfv_ho(ji,jj,jl) = 0.5 * pv(ji,jj) * ( pt(ji,jj,jl) + pt(ji,jj+1,jl) )530 pfv_ho(ji,jj,jl) = 0.5_wp * pv(ji,jj) * ( pt(ji,jj,jl) + pt(ji,jj+1,jl) ) 532 531 END DO 533 532 END DO … … 550 549 DO jj = 1, jpjm1 551 550 DO ji = 1, fs_jpim1 552 pfu_ho(ji,jj,jl) = 0.5 * pu(ji,jj) * ( zpt(ji,jj,jl) + zpt(ji+1,jj,jl) )551 pfu_ho(ji,jj,jl) = 0.5_wp * pu(ji,jj) * ( zpt(ji,jj,jl) + zpt(ji+1,jj,jl) ) 553 552 END DO 554 553 END DO … … 564 563 565 564 566 SUBROUTINE macho( pamsk, kn_umx, jt, kt, pdt, pt, pu, pv, pu c, pvc, pubox, pvbox, ptc, pt_ups, pfu_ups, pfv_ups, pfu_ho, pfv_ho )565 SUBROUTINE macho( pamsk, kn_umx, jt, kt, pdt, pt, pu, pv, pubox, pvbox, pt_ups, pfu_ups, pfv_ups, pfu_ho, pfv_ho ) 567 566 !!--------------------------------------------------------------------- 568 567 !! *** ROUTINE macho *** … … 581 580 REAL(wp), DIMENSION(:,:,:) , INTENT(in ) :: pt ! tracer fields 582 581 REAL(wp), DIMENSION(:,: ) , INTENT(in ) :: pu, pv ! 2 ice velocity components 583 REAL(wp), DIMENSION(:,:,:) , INTENT(in ) :: puc, pvc ! 2 ice velocity * A components584 582 REAL(wp), DIMENSION(:,: ) , INTENT(in ) :: pubox, pvbox ! upstream velocity 585 REAL(wp), DIMENSION(:,:,:) , INTENT(in ) :: ptc ! tracer content at before time step586 583 REAL(wp), DIMENSION(:,:,:) , INTENT(in ) :: pt_ups ! upstream guess of tracer 587 584 REAL(wp), DIMENSION(:,:,:) , INTENT(in ) :: pfu_ups, pfv_ups ! upstream fluxes … … 714 711 DO jj = 1, jpjm1 715 712 DO ji = 1, fs_jpim1 ! vector opt. 716 pt_u(ji,jj,jl) = 0.5_wp * umask(ji,jj,1) * ( pt(ji+1,jj,jl) + pt(ji,jj,jl) &717 & - SIGN( 1._wp, pu(ji,jj) ) * ( pt(ji+1,jj,jl) - pt(ji,jj,jl) ) )713 pt_u(ji,jj,jl) = 0.5_wp * umask(ji,jj,1) * ( pt(ji+1,jj,jl) + pt(ji,jj,jl) & 714 & - SIGN( 1._wp, pu(ji,jj) ) * ( pt(ji+1,jj,jl) - pt(ji,jj,jl) ) ) 718 715 END DO 719 716 END DO … … 727 724 zcu = pu(ji,jj) * r1_e2u(ji,jj) * pdt * r1_e1u(ji,jj) 728 725 pt_u(ji,jj,jl) = 0.5_wp * umask(ji,jj,1) * ( pt(ji+1,jj,jl) + pt(ji,jj,jl) & 729 & -zcu * ( pt(ji+1,jj,jl) - pt(ji,jj,jl) ) )726 & - zcu * ( pt(ji+1,jj,jl) - pt(ji,jj,jl) ) ) 730 727 END DO 731 728 END DO … … 741 738 !!rachid zdx2 = e1u(ji,jj) * e1t(ji,jj) 742 739 pt_u(ji,jj,jl) = 0.5_wp * umask(ji,jj,1) * ( ( pt (ji+1,jj,jl) + pt (ji,jj,jl) & 743 & -zcu * ( pt (ji+1,jj,jl) - pt (ji,jj,jl) ) ) &744 & + z1_6 * zdx2 * ( zcu*zcu - 1._wp ) * (ztu2(ji+1,jj,jl) + ztu2(ji,jj,jl) &740 & - zcu * ( pt (ji+1,jj,jl) - pt (ji,jj,jl) ) ) & 741 & + z1_6 * zdx2 * ( zcu*zcu - 1._wp ) * ( ztu2(ji+1,jj,jl) + ztu2(ji,jj,jl) & 745 742 & - SIGN( 1._wp, zcu ) * ( ztu2(ji+1,jj,jl) - ztu2(ji,jj,jl) ) ) ) 746 743 END DO … … 756 753 zdx2 = e1u(ji,jj) * e1u(ji,jj) 757 754 !!rachid zdx2 = e1u(ji,jj) * e1t(ji,jj) 758 pt_u(ji,jj,jl) = 0.5_wp * umask(ji,jj,1) * ( ( pt (ji+1,jj,jl) + pt (ji,jj,jl) &759 & - zcu* ( pt (ji+1,jj,jl) - pt (ji,jj,jl) ) ) &760 & + z1_6 * zdx2 * ( zcu*zcu - 1._wp ) * (ztu2(ji+1,jj,jl) + ztu2(ji,jj,jl) &761 & - 0.5_wp * zcu* ( ztu2(ji+1,jj,jl) - ztu2(ji,jj,jl) ) ) )755 pt_u(ji,jj,jl) = 0.5_wp * umask(ji,jj,1) * ( ( pt (ji+1,jj,jl) + pt (ji,jj,jl) & 756 & - zcu * ( pt (ji+1,jj,jl) - pt (ji,jj,jl) ) ) & 757 & + z1_6 * zdx2 * ( zcu*zcu - 1._wp ) * ( ztu2(ji+1,jj,jl) + ztu2(ji,jj,jl) & 758 & - 0.5_wp * zcu * ( ztu2(ji+1,jj,jl) - ztu2(ji,jj,jl) ) ) ) 762 759 END DO 763 760 END DO … … 773 770 !!rachid zdx2 = e1u(ji,jj) * e1t(ji,jj) 774 771 zdx4 = zdx2 * zdx2 775 pt_u(ji,jj,jl) = 0.5_wp * umask(ji,jj,1) * ( (pt (ji+1,jj,jl) + pt (ji,jj,jl) &776 & - zcu* ( pt (ji+1,jj,jl) - pt (ji,jj,jl) ) ) &777 & + z1_6 * zdx2 * ( zcu*zcu - 1._wp ) * (ztu2(ji+1,jj,jl) + ztu2(ji,jj,jl) &778 & - 0.5_wp * zcu* ( ztu2(ji+1,jj,jl) - ztu2(ji,jj,jl) ) ) &772 pt_u(ji,jj,jl) = 0.5_wp * umask(ji,jj,1) * ( ( pt (ji+1,jj,jl) + pt (ji,jj,jl) & 773 & - zcu * ( pt (ji+1,jj,jl) - pt (ji,jj,jl) ) ) & 774 & + z1_6 * zdx2 * ( zcu*zcu - 1._wp ) * ( ztu2(ji+1,jj,jl) + ztu2(ji,jj,jl) & 775 & - 0.5_wp * zcu * ( ztu2(ji+1,jj,jl) - ztu2(ji,jj,jl) ) ) & 779 776 & + z1_120 * zdx4 * ( zcu*zcu - 1._wp ) * ( zcu*zcu - 4._wp ) * ( ztu4(ji+1,jj,jl) + ztu4(ji,jj,jl) & 780 777 & - SIGN( 1._wp, zcu ) * ( ztu4(ji+1,jj,jl) - ztu4(ji,jj,jl) ) ) ) … … 793 790 DO ji = 1, fs_jpim1 794 791 IF( pt_u(ji,jj,jl) < 0._wp ) THEN 795 pt_u(ji,jj,jl) = 0.5_wp * umask(ji,jj,1) * ( pt(ji+1,jj,jl) + pt(ji,jj,jl) &796 & - SIGN( 1._wp, pu(ji,jj) ) * ( pt(ji+1,jj,jl) - pt(ji,jj,jl) ) )792 pt_u(ji,jj,jl) = 0.5_wp * umask(ji,jj,1) * ( pt(ji+1,jj,jl) + pt(ji,jj,jl) & 793 & - SIGN( 1._wp, pu(ji,jj) ) * ( pt(ji+1,jj,jl) - pt(ji,jj,jl) ) ) 797 794 ENDIF 798 795 END DO … … 871 868 DO jj = 1, jpjm1 872 869 DO ji = 1, fs_jpim1 873 pt_v(ji,jj,jl) = 0.5_wp * vmask(ji,jj,1) * ( ( pt(ji,jj+1,jl) + pt(ji,jj,jl) )&874 & - SIGN( 1._wp, pv(ji,jj) ) * ( pt(ji,jj+1,jl) - pt(ji,jj,jl) ) )870 pt_v(ji,jj,jl) = 0.5_wp * vmask(ji,jj,1) * ( pt(ji,jj+1,jl) + pt(ji,jj,jl) & 871 & - SIGN( 1._wp, pv(ji,jj) ) * ( pt(ji,jj+1,jl) - pt(ji,jj,jl) ) ) 875 872 END DO 876 873 END DO … … 882 879 DO ji = 1, fs_jpim1 883 880 zcv = pv(ji,jj) * r1_e1v(ji,jj) * pdt * r1_e2v(ji,jj) 884 pt_v(ji,jj,jl) = 0.5_wp * vmask(ji,jj,1) * ( ( pt(ji,jj+1,jl) + pt(ji,jj,jl) ) & 885 & - zcv * ( pt(ji,jj+1,jl) - pt(ji,jj,jl) ) ) 886 END DO 887 END DO 888 END DO 889 CALL lbc_lnk( 'icedyn_adv_umx', pt_v, 'V', 1. ) 881 pt_v(ji,jj,jl) = 0.5_wp * vmask(ji,jj,1) * ( pt(ji,jj+1,jl) + pt(ji,jj,jl) & 882 & - zcv * ( pt(ji,jj+1,jl) - pt(ji,jj,jl) ) ) 883 END DO 884 END DO 885 END DO 890 886 ! 891 887 CASE( 3 ) !== 3rd order central TIM ==! (Eq. 24) … … 896 892 zdy2 = e2v(ji,jj) * e2v(ji,jj) 897 893 !!rachid zdy2 = e2v(ji,jj) * e2t(ji,jj) 898 pt_v(ji,jj,jl) = 0.5_wp * vmask(ji,jj,1) * ( (pt (ji,jj+1,jl) + pt (ji,jj,jl) &899 & -zcv * ( pt (ji,jj+1,jl) - pt (ji,jj,jl) ) ) &894 pt_v(ji,jj,jl) = 0.5_wp * vmask(ji,jj,1) * ( ( pt (ji,jj+1,jl) + pt (ji,jj,jl) & 895 & - zcv * ( pt (ji,jj+1,jl) - pt (ji,jj,jl) ) ) & 900 896 & + z1_6 * zdy2 * ( zcv*zcv - 1._wp ) * ( ztv2(ji,jj+1,jl) + ztv2(ji,jj,jl) & 901 897 & - SIGN( 1._wp, zcv ) * ( ztv2(ji,jj+1,jl) - ztv2(ji,jj,jl) ) ) ) … … 911 907 zdy2 = e2v(ji,jj) * e2v(ji,jj) 912 908 !!rachid zdy2 = e2v(ji,jj) * e2t(ji,jj) 913 pt_v(ji,jj,jl) = 0.5_wp * vmask(ji,jj,1) * ( (pt (ji,jj+1,jl) + pt (ji,jj,jl) &914 & - zcv* ( pt (ji,jj+1,jl) - pt (ji,jj,jl) ) ) &915 & + z1_6 * zdy2 * ( zcv*zcv - 1._wp ) * ( ztv2(ji,jj+1,jl) + ztv2(ji,jj,jl) &916 & - 0.5_wp * zcv* ( ztv2(ji,jj+1,jl) - ztv2(ji,jj,jl) ) ) )909 pt_v(ji,jj,jl) = 0.5_wp * vmask(ji,jj,1) * ( ( pt (ji,jj+1,jl) + pt (ji,jj,jl) & 910 & - zcv * ( pt (ji,jj+1,jl) - pt (ji,jj,jl) ) ) & 911 & + z1_6 * zdy2 * ( zcv*zcv - 1._wp ) * ( ztv2(ji,jj+1,jl) + ztv2(ji,jj,jl) & 912 & - 0.5_wp * zcv * ( ztv2(ji,jj+1,jl) - ztv2(ji,jj,jl) ) ) ) 917 913 END DO 918 914 END DO … … 928 924 zdy4 = zdy2 * zdy2 929 925 pt_v(ji,jj,jl) = 0.5_wp * vmask(ji,jj,1) * ( ( pt (ji,jj+1,jl) + pt (ji,jj,jl) & 930 & - zcv* ( pt (ji,jj+1,jl) - pt (ji,jj,jl) ) ) &931 & + z1_6 * zdy2 * ( zcv*zcv - 1._wp ) * (ztv2(ji,jj+1,jl) + ztv2(ji,jj,jl) &932 & - 0.5_wp * zcv* ( ztv2(ji,jj+1,jl) - ztv2(ji,jj,jl) ) ) &926 & - zcv * ( pt (ji,jj+1,jl) - pt (ji,jj,jl) ) ) & 927 & + z1_6 * zdy2 * ( zcv*zcv - 1._wp ) * ( ztv2(ji,jj+1,jl) + ztv2(ji,jj,jl) & 928 & - 0.5_wp * zcv * ( ztv2(ji,jj+1,jl) - ztv2(ji,jj,jl) ) ) & 933 929 & + z1_120 * zdy4 * ( zcv*zcv - 1._wp ) * ( zcv*zcv - 4._wp ) * ( ztv4(ji,jj+1,jl) + ztv4(ji,jj,jl) & 934 930 & - SIGN( 1._wp, zcv ) * ( ztv4(ji,jj+1,jl) - ztv4(ji,jj,jl) ) ) ) … … 947 943 DO ji = 1, fs_jpim1 948 944 IF( pt_v(ji,jj,jl) < 0._wp ) THEN 949 pt_v(ji,jj,jl) = 0.5_wp * vmask(ji,jj,1) * ( ( pt(ji,jj+1,jl) + pt(ji,jj,jl) ) &950 & - SIGN( 1._wp, pv(ji,jj) ) * ( pt(ji,jj+1,jl) - pt(ji,jj,jl) ) )945 pt_v(ji,jj,jl) = 0.5_wp * vmask(ji,jj,1) * ( ( pt(ji,jj+1,jl) + pt(ji,jj,jl) ) & 946 & - SIGN( 1._wp, pv(ji,jj) ) * ( pt(ji,jj+1,jl) - pt(ji,jj,jl) ) ) 951 947 ENDIF 952 948 END DO … … 984 980 ! 985 981 INTEGER :: ji, jj, jl ! dummy loop indices 986 REAL(wp) :: zpos, zneg, zbig, z sml, zup, zdo, z1_dt ! local scalars987 REAL(wp) :: zau, zbu, zcu, zav, zbv, zcv, z sign, zcoef, zzt! - -982 REAL(wp) :: zpos, zneg, zbig, zup, zdo, z1_dt ! local scalars 983 REAL(wp) :: zau, zbu, zcu, zav, zbv, zcv, zcoef, zzt ! - - 988 984 REAL(wp), DIMENSION(jpi,jpj ) :: zbup, zbdo 989 985 REAL(wp), DIMENSION(jpi,jpj,jpl) :: zbetup, zbetdo, zti_ups, ztj_ups 990 986 !!---------------------------------------------------------------------- 991 987 zbig = 1.e+40_wp 992 zsml = epsi20993 988 994 989 ! antidiffusive flux : high order minus low order … … 1026 1021 DO jj = 2, jpjm1 1027 1022 DO ji = fs_2, fs_jpim1 1028 IF ( pfu_ho(ji,jj,jl) * ( pt_ups(ji+1,jj ,jl) - pt_ups(ji,jj,jl) ) <= 0..AND. &1029 & pfv_ho(ji,jj,jl) * ( pt_ups(ji ,jj+1,jl) - pt_ups(ji,jj,jl) ) <= 0.) THEN1023 IF ( pfu_ho(ji,jj,jl) * ( pt_ups(ji+1,jj ,jl) - pt_ups(ji,jj,jl) ) <= 0._wp .AND. & 1024 & pfv_ho(ji,jj,jl) * ( pt_ups(ji ,jj+1,jl) - pt_ups(ji,jj,jl) ) <= 0._wp ) THEN 1030 1025 ! 1031 IF( pfu_ho(ji,jj,jl) * ( zti_ups(ji+1,jj ,jl) - zti_ups(ji,jj,jl) ) <= 0..AND. &1032 & pfv_ho(ji,jj,jl) * ( ztj_ups(ji ,jj+1,jl) - ztj_ups(ji,jj,jl) ) <= 0.) THEN1033 pfu_ho(ji,jj,jl)=0. 1034 pfv_ho(ji,jj,jl)=0. 1026 IF( pfu_ho(ji,jj,jl) * ( zti_ups(ji+1,jj ,jl) - zti_ups(ji,jj,jl) ) <= 0._wp .AND. & 1027 & pfv_ho(ji,jj,jl) * ( ztj_ups(ji ,jj+1,jl) - ztj_ups(ji,jj,jl) ) <= 0._wp ) THEN 1028 pfu_ho(ji,jj,jl)=0._wp 1029 pfv_ho(ji,jj,jl)=0._wp 1035 1030 ENDIF 1036 1031 ! 1037 IF( pfu_ho(ji,jj,jl) * ( pt_ups(ji ,jj,jl) - pt_ups(ji-1,jj,jl) ) <= 0..AND. &1038 & pfv_ho(ji,jj,jl) * ( pt_ups(ji ,jj,jl) - pt_ups(ji,jj-1,jl) ) <= 0.) THEN1039 pfu_ho(ji,jj,jl)=0. 1040 pfv_ho(ji,jj,jl)=0. 1032 IF( pfu_ho(ji,jj,jl) * ( pt_ups(ji,jj,jl) - pt_ups(ji-1,jj ,jl) ) <= 0._wp .AND. & 1033 & pfv_ho(ji,jj,jl) * ( pt_ups(ji,jj,jl) - pt_ups(ji ,jj-1,jl) ) <= 0._wp ) THEN 1034 pfu_ho(ji,jj,jl)=0._wp 1035 pfv_ho(ji,jj,jl)=0._wp 1041 1036 ENDIF 1042 1037 ! … … 1076 1071 DO ji = fs_2, fs_jpim1 ! vector opt. 1077 1072 ! 1078 zup = MAX( zbup(ji,jj), zbup(ji-1,jj ), zbup(ji+1,jj ), zbup(ji ,jj-1), zbup(ji,jj+1) ) ! search max/min in neighbourhood1079 zdo = MIN( zbdo(ji,jj), zbdo(ji-1,jj ), zbdo(ji+1,jj ), zbdo(ji ,jj-1), zbdo(ji,jj+1) )1080 ! 1081 zpos = MAX( 0. , pfu_ho(ji-1,jj,jl) ) - MIN( 0., pfu_ho(ji ,jj,jl) ) & ! positive/negative part of the flux1082 & + MAX( 0. , pfv_ho(ji,jj-1,jl) ) - MIN( 0., pfv_ho(ji,jj ,jl) )1083 zneg = MAX( 0. , pfu_ho(ji ,jj,jl) ) - MIN( 0., pfu_ho(ji-1,jj,jl) ) &1084 & + MAX( 0. , pfv_ho(ji,jj ,jl) ) - MIN( 0., pfv_ho(ji,jj-1,jl) )1085 ! 1086 zpos = zpos - (pt(ji,jj,jl) * MIN( 0., pu(ji,jj) - pu(ji-1,jj) ) + pt(ji,jj,jl) * MIN( 0., pv(ji,jj) - pv(ji,jj-1) ) &1073 zup = MAX( zbup(ji,jj), zbup(ji-1,jj), zbup(ji+1,jj), zbup(ji,jj-1), zbup(ji,jj+1) ) ! search max/min in neighbourhood 1074 zdo = MIN( zbdo(ji,jj), zbdo(ji-1,jj), zbdo(ji+1,jj), zbdo(ji,jj-1), zbdo(ji,jj+1) ) 1075 ! 1076 zpos = MAX( 0._wp, pfu_ho(ji-1,jj ,jl) ) - MIN( 0._wp, pfu_ho(ji ,jj ,jl) ) & ! positive/negative part of the flux 1077 & + MAX( 0._wp, pfv_ho(ji ,jj-1,jl) ) - MIN( 0._wp, pfv_ho(ji ,jj ,jl) ) 1078 zneg = MAX( 0._wp, pfu_ho(ji ,jj ,jl) ) - MIN( 0._wp, pfu_ho(ji-1,jj ,jl) ) & 1079 & + MAX( 0._wp, pfv_ho(ji ,jj ,jl) ) - MIN( 0._wp, pfv_ho(ji ,jj-1,jl) ) 1080 ! 1081 zpos = zpos - (pt(ji,jj,jl) * MIN( 0., pu(ji,jj) - pu(ji-1,jj) ) + pt(ji,jj,jl) * MIN( 0., pv(ji,jj) - pv(ji,jj-1) ) & 1087 1082 & ) * ( 1. - pamsk ) 1088 zneg = zneg + (pt(ji,jj,jl) * MAX( 0., pu(ji,jj) - pu(ji-1,jj) ) + pt(ji,jj,jl) * MAX( 0., pv(ji,jj) - pv(ji,jj-1) ) &1083 zneg = zneg + (pt(ji,jj,jl) * MAX( 0., pu(ji,jj) - pu(ji-1,jj) ) + pt(ji,jj,jl) * MAX( 0., pv(ji,jj) - pv(ji,jj-1) ) & 1089 1084 & ) * ( 1. - pamsk ) 1090 1085 ! 1091 1086 ! ! up & down beta terms 1092 IF( zpos > 0. ) THEN ; zbetup(ji,jj,jl) = MAX( 0._wp, zup - pt_ups(ji,jj,jl) ) / zpos * e1e2t(ji,jj) * z1_dt1093 ELSE ; zbetup(ji,jj,jl) = 0.! zbig1087 IF( zpos > 0._wp ) THEN ; zbetup(ji,jj,jl) = MAX( 0._wp, zup - pt_ups(ji,jj,jl) ) / zpos * e1e2t(ji,jj) * z1_dt 1088 ELSE ; zbetup(ji,jj,jl) = 0._wp ! zbig 1094 1089 ENDIF 1095 1090 ! 1096 IF( zneg > 0. ) THEN ; zbetdo(ji,jj,jl) = MAX( 0._wp, pt_ups(ji,jj,jl) - zdo ) / zneg * e1e2t(ji,jj) * z1_dt1097 ELSE ; zbetdo(ji,jj,jl) = 0.! zbig1091 IF( zneg > 0._wp ) THEN ; zbetdo(ji,jj,jl) = MAX( 0._wp, pt_ups(ji,jj,jl) - zdo ) / zneg * e1e2t(ji,jj) * z1_dt 1092 ELSE ; zbetdo(ji,jj,jl) = 0._wp ! zbig 1098 1093 ENDIF 1099 1094 ! 1100 1095 ! if all the points are outside ice cover 1101 IF( zup == -zbig ) zbetup(ji,jj,jl) = 0. ! zbig1102 IF( zdo == zbig ) zbetdo(ji,jj,jl) = 0. ! zbig1096 IF( zup == -zbig ) zbetup(ji,jj,jl) = 0._wp ! zbig 1097 IF( zdo == zbig ) zbetdo(ji,jj,jl) = 0._wp ! zbig 1103 1098 ! 1104 1099 END DO … … 1115 1110 zau = MIN( 1._wp , zbetdo(ji,jj,jl) , zbetup(ji+1,jj,jl) ) 1116 1111 zbu = MIN( 1._wp , zbetup(ji,jj,jl) , zbetdo(ji+1,jj,jl) ) 1117 zcu = 0.5 + SIGN( 0.5, pfu_ho(ji,jj,jl) )1112 zcu = 0.5_wp + SIGN( 0.5_wp , pfu_ho(ji,jj,jl) ) 1118 1113 ! 1119 1114 zcoef = ( zcu * zau + ( 1._wp - zcu ) * zbu ) … … 1128 1123 zav = MIN( 1._wp , zbetdo(ji,jj,jl) , zbetup(ji,jj+1,jl) ) 1129 1124 zbv = MIN( 1._wp , zbetup(ji,jj,jl) , zbetdo(ji,jj+1,jl) ) 1130 zcv = 0.5 + SIGN( 0.5, pfv_ho(ji,jj,jl) )1125 zcv = 0.5_wp + SIGN( 0.5_wp , pfv_ho(ji,jj,jl) ) 1131 1126 ! 1132 1127 zcoef = ( zcv * zav + ( 1._wp - zcv ) * zbv )
Note: See TracChangeset
for help on using the changeset viewer.
