Changeset 10281
- Timestamp:
- 2018-11-07T12:22:39+01:00 (6 years ago)
- File:
-
- 1 edited
Legend:
- Unmodified
- Added
- Removed
-
NEMO/branches/2018/dev_r9947_SI3_advection/tests/ICEADV/MY_SRC/icedyn_adv_umx.F90
r10277 r10281 35 35 REAL(wp) :: z1_120 = 1._wp / 120._wp ! =1/120 36 36 37 REAL(wp), DIMENSION(:,:), ALLOCATABLE :: z1_a 38 39 ! alternate directions for upstream 40 ! clem: it gives results for Lipscomb test that are the same as "ll_upsxy=false" 41 LOGICAL :: ll_upsxy = .FALSE. 42 43 ! prelimiter 44 LOGICAL :: ll_prelimiter = .TRUE. 45 LOGICAL :: ll_prelimiter_zalesak = .TRUE. ! from: Zalesak(1979) eq. 14 => better than Devore 46 LOGICAL :: ll_prelimiter_devore = .FALSE. ! from: Devore eq. 11 47 48 ! iterate on the limiter (only nonosc) 49 LOGICAL :: ll_limiter_it2 = .TRUE. 50 51 37 52 !! * Substitutions 38 53 # include "vectopt_loop_substitute.h90" … … 77 92 ! ! must be >= 0.01 and the best seems to be 0.1 78 93 REAL(wp), DIMENSION(jpi,jpj) :: zudy, zvdx, zcu_box, zcv_box, zfu, zfv 79 REAL(wp), DIMENSION(jpi,jpj) :: z 1_a, zh_i, zh_s, zs_i, zo_i, ze_i, ze_s, z1_ap, zh_ip94 REAL(wp), DIMENSION(jpi,jpj) :: zh_i, zh_s, zs_i, zo_i, ze_i, ze_s, zh_ip 80 95 !!---------------------------------------------------------------------- 81 96 ! … … 112 127 END DO 113 128 129 IF(.NOT. ALLOCATED(z1_a )) ALLOCATE(z1_a (jpi,jpj)) 114 130 !---------------! 115 131 !== advection ==! … … 149 165 END DO 150 166 ! 151 IF ( ln_pnd_H12 ) THEN ! melt ponds 167 IF ( ln_pnd_H12 ) THEN ! melt ponds (must be the last ones to be advected because of z1_a...) 152 168 ! to avoid a problem with the limiter nonosc when A gets close to 0 153 169 pa_ip(:,:,jl) = pa_ip(:,:,jl) + zeps * tmask(:,:,1) 154 170 ! 155 WHERE( pa_ip(:,:,jl) > epsi20 ) ; z1_a p(:,:) = 1._wp / pa_ip(:,:,jl)156 ELSEWHERE ; z1_a p(:,:) = 0.171 WHERE( pa_ip(:,:,jl) > epsi20 ) ; z1_a(:,:) = 1._wp / pa_ip(:,:,jl) 172 ELSEWHERE ; z1_a(:,:) = 0. 157 173 END WHERE 158 174 ! … … 213 229 INTEGER :: ji, jj ! dummy loop indices 214 230 REAL(wp) :: ztra ! local scalar 215 !!clem REAL(wp) :: zeps = 1.e-02 ! local scalar216 231 INTEGER :: kn_limiter = 1 ! 1=nonosc ; 2=superbee ; 3=h3(rachid) 217 REAL(wp), DIMENSION(jpi,jpj) :: zfu_ho , zfv_ho , zt_u, zt_v, zpt 232 REAL(wp), DIMENSION(jpi,jpj) :: zfu_ho , zfv_ho , zt_u, zt_v, zpt, zz1_a 218 233 REAL(wp), DIMENSION(jpi,jpj) :: zfu_ups, zfv_ups, zt_ups ! only for nonosc 219 234 !!---------------------------------------------------------------------- 220 235 ! 221 ! add a constant value to avoid problems with zeros222 DO jj = 1, jpj223 DO ji = 1, jpi224 zpt(ji,jj) = pt(ji,jj) !!clem + zeps * tmask(ji,jj,1)225 END DO226 END DO227 228 236 ! upstream (_ups) advection with initial mass fluxes 229 237 ! --------------------------------------------------- 230 ! fluxes 231 DO jj = 1, jpjm1 232 DO ji = 1, fs_jpim1 233 zfu_ups(ji,jj) = MAX( puc(ji,jj), 0._wp ) * zpt(ji,jj) + MIN( puc(ji,jj), 0._wp ) * zpt(ji+1,jj) 234 zfv_ups(ji,jj) = MAX( pvc(ji,jj), 0._wp ) * zpt(ji,jj) + MIN( pvc(ji,jj), 0._wp ) * zpt(ji,jj+1) 235 END DO 236 END DO 238 IF( .NOT. ll_upsxy ) THEN 239 240 ! fluxes 241 DO jj = 1, jpjm1 242 DO ji = 1, fs_jpim1 243 zfu_ups(ji,jj) = MAX( puc(ji,jj), 0._wp ) * pt(ji,jj) + MIN( puc(ji,jj), 0._wp ) * pt(ji+1,jj) 244 zfv_ups(ji,jj) = MAX( pvc(ji,jj), 0._wp ) * pt(ji,jj) + MIN( pvc(ji,jj), 0._wp ) * pt(ji,jj+1) 245 END DO 246 END DO 247 248 ELSE 249 ! 250 IF ( pamsk == 1. ) THEN ; zz1_a(:,:) = 1._wp ; ! 1 if advection of A 251 ELSEIF( pamsk == 0. ) THEN ; zz1_a(:,:) = z1_a(:,:) ; ENDIF ! 1/A if advection of V 252 ! 253 IF( MOD( (kt - 1) / nn_fsbc , 2 ) == MOD( (jt - 1) , 2 ) ) THEN !== odd ice time step: adv_x then adv_y ==! 254 ! flux in x-direction 255 DO jj = 1, jpjm1 256 DO ji = 1, fs_jpim1 257 zfu_ups(ji,jj) = MAX( puc(ji,jj), 0._wp ) * pt(ji,jj) + MIN( puc(ji,jj), 0._wp ) * pt(ji+1,jj) 258 END DO 259 END DO 260 ! first guess of tracer content from u-flux 261 DO jj = 2, jpjm1 262 DO ji = fs_2, fs_jpim1 ! vector opt. 263 zpt(ji,jj) = ( ptc(ji,jj) - ( zfu_ups(ji,jj) - zfu_ups(ji-1,jj) ) * pdt * r1_e1e2t(ji,jj) ) * tmask(ji,jj,1) * zz1_a(ji,jj) 264 END DO 265 END DO 266 CALL lbc_lnk( zpt, 'T', 1. ) 267 ! 268 ! flux in y-direction 269 DO jj = 1, jpjm1 270 DO ji = 1, fs_jpim1 271 zfv_ups(ji,jj) = MAX( pvc(ji,jj), 0._wp ) * zpt(ji,jj) + MIN( pvc(ji,jj), 0._wp ) * zpt(ji,jj+1) 272 END DO 273 END DO 274 ! 275 ELSE !== even ice time step: adv_y then adv_x ==! 276 ! flux in y-direction 277 DO jj = 1, jpjm1 278 DO ji = 1, fs_jpim1 279 zfv_ups(ji,jj) = MAX( pvc(ji,jj), 0._wp ) * pt(ji,jj) + MIN( pvc(ji,jj), 0._wp ) * pt(ji,jj+1) 280 END DO 281 END DO 282 ! first guess of tracer content from v-flux 283 DO jj = 2, jpjm1 284 DO ji = fs_2, fs_jpim1 ! vector opt. 285 zpt(ji,jj) = ( ptc(ji,jj) - ( zfv_ups(ji,jj) - zfv_ups(ji,jj-1) ) * pdt * r1_e1e2t(ji,jj) ) * tmask(ji,jj,1) * zz1_a(ji,jj) 286 END DO 287 END DO 288 CALL lbc_lnk( zpt, 'T', 1. ) 289 ! 290 ! flux in x-direction 291 DO jj = 1, jpjm1 292 DO ji = 1, fs_jpim1 293 zfu_ups(ji,jj) = MAX( puc(ji,jj), 0._wp ) * zpt(ji,jj) + MIN( puc(ji,jj), 0._wp ) * zpt(ji+1,jj) 294 END DO 295 END DO 296 ! 297 ENDIF 298 299 ENDIF 237 300 ! guess after content field with upstream scheme 238 301 DO jj = 2, jpjm1 … … 251 314 CASE ( 20 ) !== centered second order ==! 252 315 ! 253 CALL cen2( kn_limiter, jt, kt, pdt, zpt, pu, pv, puc, pvc, ptc, zfu_ho, zfv_ho, &316 CALL cen2( kn_limiter, jt, kt, pdt, pt, pu, pv, puc, pvc, ptc, zfu_ho, zfv_ho, & 254 317 & zt_ups, zfu_ups, zfv_ups ) 255 318 ! 256 319 CASE ( 1:5 ) !== 1st to 5th order ULTIMATE-MACHO scheme ==! 257 320 ! 258 CALL macho( pamsk, kn_limiter, kn_umx, jt, kt, pdt, zpt, pu, pv, puc, pvc, pubox, pvbox, ptc, zt_u, zt_v, zfu_ho, zfv_ho, &321 CALL macho( pamsk, kn_limiter, kn_umx, jt, kt, pdt, pt, pu, pv, puc, pvc, pubox, pvbox, ptc, zt_u, zt_v, zfu_ho, zfv_ho, & 259 322 & zt_ups, zfu_ups, zfv_ups ) 260 323 ! … … 266 329 DO jj = 1, jpjm1 267 330 DO ji = 1, fs_jpim1 268 pfu(ji,jj) = zfu_ho(ji,jj) !!clem - zeps * puc(ji,jj)269 pfv(ji,jj) = zfv_ho(ji,jj) !!clem - zeps * pvc(ji,jj)331 pfu(ji,jj) = zfu_ho(ji,jj) 332 pfv(ji,jj) = zfv_ho(ji,jj) 270 333 END DO 271 334 END DO … … 278 341 DO jj = 2, jpjm1 279 342 DO ji = fs_2, fs_jpim1 280 ztra = ( - ( zfu_ho(ji,jj) - zfu_ho(ji-1,jj) + zfv_ho(ji,jj) - zfv_ho(ji,jj-1) ) & ! Div(uaH) or Div(ua) 281 !!clem & + ( puc (ji,jj) - puc (ji-1,jj) + pvc (ji,jj) - pvc (ji,jj-1) ) * zeps & ! epsi * Div(ua) or epsi * Div(u) 343 ztra = ( - ( zfu_ho(ji,jj) - zfu_ho(ji-1,jj) + zfv_ho(ji,jj) - zfv_ho(ji,jj-1) ) & ! Div(uaH) or Div(ua) 282 344 & ) * r1_e1e2t(ji,jj) 283 345 ptc(ji,jj) = ( ptc(ji,jj) + pdt * ztra ) * tmask(ji,jj,1) … … 430 492 ! 431 493 INTEGER :: ji, jj ! dummy loop indices 432 REAL(wp), DIMENSION(jpi,jpj) :: zzt 494 REAL(wp) :: ztra 495 REAL(wp), DIMENSION(jpi,jpj) :: zzt, zzfu_ho, zzfv_ho 433 496 !!---------------------------------------------------------------------- 434 497 ! … … 478 541 ! 479 542 ENDIF 480 IF( kn_limiter == 1 ) CALL nonosc_2d ( pdt, ptc, pt_ups, pfu_ups, pfv_ups, pfu_ho, pfv_ho ) 543 IF( kn_limiter == 1 ) THEN 544 IF( .NOT. ll_limiter_it2 ) THEN 545 CALL nonosc_2d ( pdt, ptc, pt_ups, pfu_ups, pfv_ups, pfu_ho, pfv_ho ) 546 ELSE 547 zzfu_ho(:,:) = pfu_ho(:,:) 548 zzfv_ho(:,:) = pfv_ho(:,:) 549 ! 1st iteration of nonosc (limit the flux with the upstream solution) 550 CALL nonosc_2d ( pdt, ptc, pt_ups, pfu_ups, pfv_ups, zzfu_ho, zzfv_ho ) 551 ! guess after content field with high order 552 DO jj = 2, jpjm1 553 DO ji = fs_2, fs_jpim1 554 ztra = - ( zzfu_ho(ji,jj) - zzfu_ho(ji-1,jj) + zzfv_ho(ji,jj) - zzfv_ho(ji,jj-1) ) * r1_e1e2t(ji,jj) 555 zzt(ji,jj) = ( ptc(ji,jj) + pdt * ztra ) * tmask(ji,jj,1) 556 END DO 557 END DO 558 CALL lbc_lnk( zzt, 'T', 1. ) 559 ! 2nd iteration of nonosc (limit the flux with the limited high order solution) 560 CALL nonosc_2d ( pdt, ptc, zzt, zzfu_ho, zzfv_ho, pfu_ho, pfv_ho ) 561 ENDIF 562 ENDIF 481 563 ! 482 564 END SUBROUTINE macho … … 733 815 734 816 735 SUBROUTINE nonosc_2d( pdt, ptc, pt_ ups, pfu_ups, pfv_ups, pfu_ho, pfv_ho )817 SUBROUTINE nonosc_2d( pdt, ptc, pt_low, pfu_low, pfv_low, pfu_ho, pfv_ho ) 736 818 !!--------------------------------------------------------------------- 737 819 !! *** ROUTINE nonosc *** 738 820 !! 739 !! ** Purpose : compute monotonic tracer fluxes from the upstream821 !! ** Purpose : compute monotonic tracer fluxes from the upstream 740 822 !! scheme and the before field by a nonoscillatory algorithm 741 823 !! 742 824 !! ** Method : ... ??? 743 !! warning : ptc and pt_ upsmust be masked, but the boundaries825 !! warning : ptc and pt_low must be masked, but the boundaries 744 826 !! conditions on the fluxes are not necessary zalezak (1979) 745 827 !! drange (1995) multi-dimensional forward-in-time and upstream- … … 747 829 !!---------------------------------------------------------------------- 748 830 REAL(wp) , INTENT(in ) :: pdt ! tracer time-step 749 REAL(wp), DIMENSION (jpi,jpj), INTENT(in ) :: ptc, pt_ ups! before field & upstream guess of after field750 REAL(wp), DIMENSION (jpi,jpj), INTENT(in ) :: pfv_ ups, pfu_ups! upstream flux831 REAL(wp), DIMENSION (jpi,jpj), INTENT(in ) :: ptc, pt_low ! before field & upstream guess of after field 832 REAL(wp), DIMENSION (jpi,jpj), INTENT(in ) :: pfv_low, pfu_low ! upstream flux 751 833 REAL(wp), DIMENSION (jpi,jpj), INTENT(inout) :: pfv_ho, pfu_ho ! monotonic flux 752 834 ! 753 835 INTEGER :: ji, jj ! dummy loop indices 754 836 REAL(wp) :: zpos, zneg, zbig, zsml, z1_dt ! local scalars 755 REAL(wp) :: zau, zbu, zcu, zav, zbv, zcv, zup, zdo 756 REAL(wp), DIMENSION(jpi,jpj) :: zbetup, zbetdo, zbup, zbdo, zdiv 837 REAL(wp) :: zau, zbu, zcu, zav, zbv, zcv, zup, zdo, zsign ! - - 838 REAL(wp), DIMENSION(jpi,jpj) :: zbetup, zbetdo, zbup, zbdo, zdiv, zti_low, ztj_low 757 839 !!---------------------------------------------------------------------- 758 840 zbig = 1.e+40_wp … … 771 853 DO jj = 1, jpjm1 772 854 DO ji = 1, fs_jpim1 ! vector opt. 773 pfu_ho(ji,jj) = pfu_ho(ji,jj) - pfu_ups(ji,jj) 774 pfv_ho(ji,jj) = pfv_ho(ji,jj) - pfv_ups(ji,jj) 775 END DO 776 END DO 777 855 pfu_ho(ji,jj) = pfu_ho(ji,jj) - pfu_low(ji,jj) 856 pfv_ho(ji,jj) = pfv_ho(ji,jj) - pfv_low(ji,jj) 857 END DO 858 END DO 859 860 ! extreme case where pfu_ho has to be zero 861 ! ---------------------------------------- 862 ! pfu_ho 863 ! * ---> 864 ! | | * | | 865 ! | | | * | 866 ! | | | | * 867 ! t_low : i-1 i i+1 i+2 868 IF( ll_prelimiter ) THEN 869 870 DO jj = 2, jpjm1 871 DO ji = fs_2, fs_jpim1 872 zti_low(ji,jj)= pt_low(ji+1,jj ) 873 ztj_low(ji,jj)= pt_low(ji ,jj+1) 874 END DO 875 END DO 876 CALL lbc_lnk_multi( zti_low, 'T', 1., ztj_low, 'T', 1. ) 877 878 IF( ll_prelimiter_zalesak ) THEN 879 DO jj = 2, jpjm1 880 DO ji = fs_2, fs_jpim1 881 IF( SIGN( 1., pfu_ho(ji,jj) ) /= SIGN( 1., pt_low (ji+1,jj) - pt_low (ji ,jj) ) ) THEN 882 IF( SIGN( 1., pfu_ho(ji,jj) ) /= SIGN( 1., zti_low(ji+1,jj) - zti_low(ji ,jj) ) ) pfu_ho(ji,jj) = 0. 883 IF( SIGN( 1., pfu_ho(ji,jj) ) /= SIGN( 1., pt_low (ji ,jj) - pt_low (ji-1,jj) ) ) pfu_ho(ji,jj) = 0. 884 ENDIF 885 IF( SIGN( 1., pfv_ho(ji,jj) ) /= SIGN( 1., pt_low (ji,jj+1) - pt_low (ji,jj ) ) ) THEN 886 IF( SIGN( 1., pfv_ho(ji,jj) ) /= SIGN( 1., ztj_low(ji,jj+1) - ztj_low(ji,jj ) ) ) pfv_ho(ji,jj) = 0. 887 IF( SIGN( 1., pfv_ho(ji,jj) ) /= SIGN( 1., pt_low (ji,jj ) - pt_low (ji,jj-1) ) ) pfv_ho(ji,jj) = 0. 888 ENDIF 889 END DO 890 END DO 891 CALL lbc_lnk_multi( pfu_ho, 'U', -1., pfv_ho, 'V', -1. ) ! lateral boundary cond. 892 893 ELSEIF( ll_prelimiter_devore ) THEN 894 z1_dt = 1._wp / pdt 895 DO jj = 2, jpjm1 896 DO ji = fs_2, fs_jpim1 897 zsign = SIGN( 1., pt_low(ji+1,jj) - pt_low(ji,jj) ) 898 pfu_ho(ji,jj) = zsign * MAX( 0. , MIN( ABS(pfu_ho(ji,jj)) , & 899 & zsign * ( pt_low (ji ,jj) - pt_low (ji-1,jj) ) * e1e2t(ji ,jj) * z1_dt , & 900 & zsign * ( zti_low(ji+1,jj) - zti_low(ji ,jj) ) * e1e2t(ji+1,jj) * z1_dt ) ) 901 902 zsign = SIGN( 1., pt_low(ji,jj+1) - pt_low(ji,jj) ) 903 pfv_ho(ji,jj) = zsign * MAX( 0. , MIN( ABS(pfv_ho(ji,jj)) , & 904 & zsign * ( pt_low (ji,jj ) - pt_low (ji,jj-1) ) * e1e2t(ji,jj ) * z1_dt , & 905 & zsign * ( ztj_low(ji,jj+1) - ztj_low(ji,jj ) ) * e1e2t(ji,jj+1) * z1_dt ) ) 906 END DO 907 END DO 908 CALL lbc_lnk_multi( pfu_ho, 'U', -1., pfv_ho, 'V', -1. ) ! lateral boundary cond. 909 910 ENDIF 911 912 ENDIF 913 778 914 ! Search local extrema 779 915 ! -------------------- 780 ! max/min of ptc & pt_ upswith large negative/positive value (-/+zbig) outside ice cover916 ! max/min of ptc & pt_low with large negative/positive value (-/+zbig) outside ice cover 781 917 DO jj = 1, jpj 782 DO ji = fs_2, fs_jpim1783 IF( ptc(ji,jj) == 0._wp .AND. pt_ ups(ji,jj) == 0._wp .AND. zdiv(ji,jj) == 0._wp ) THEN918 DO ji = 1, jpi 919 IF( ptc(ji,jj) == 0._wp .AND. pt_low(ji,jj) == 0._wp .AND. zdiv(ji,jj) == 0._wp ) THEN 784 920 zbup(ji,jj) = -zbig 785 921 zbdo(ji,jj) = zbig 786 922 ELSE 787 zbup(ji,jj) = MAX( ptc(ji,jj) , pt_ ups(ji,jj) )788 zbdo(ji,jj) = MIN( ptc(ji,jj) , pt_ ups(ji,jj) )923 zbup(ji,jj) = MAX( ptc(ji,jj) , pt_low(ji,jj) ) 924 zbdo(ji,jj) = MIN( ptc(ji,jj) , pt_low(ji,jj) ) 789 925 ENDIF 790 926 END DO … … 795 931 DO ji = fs_2, fs_jpim1 ! vector opt. 796 932 ! 797 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 798 zdo = MIN( zbdo(ji,jj), zbdo(ji-1,jj), zbdo(ji+1,jj), zbdo(ji,jj-1), zbdo(ji,jj+1) ) 799 ! 800 zpos = MAX( 0., pfu_ho(ji-1,jj) ) - MIN( 0., pfu_ho(ji ,jj) ) & 801 & + MAX( 0., pfv_ho(ji,jj-1) ) - MIN( 0., pfv_ho(ji,jj ) ) ! positive/negative part of the flux 933 !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 934 !zdo = MIN( zbdo(ji,jj), zbdo(ji-1,jj ), zbdo(ji+1,jj ), zbdo(ji ,jj-1), zbdo(ji ,jj+1) ) 935 ! 936 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 937 & zbup(ji-1,jj-1), zbup(ji+1,jj+1), zbup(ji+1,jj-1), zbup(ji-1,jj+1) ) 938 zdo = MIN( zbdo(ji,jj), zbdo(ji-1,jj ), zbdo(ji+1,jj ), zbdo(ji ,jj-1), zbdo(ji ,jj+1), & 939 & zbdo(ji-1,jj-1), zbdo(ji+1,jj+1), zbdo(ji+1,jj-1), zbdo(ji-1,jj+1) ) 940 ! 941 zpos = MAX( 0., pfu_ho(ji-1,jj) ) - MIN( 0., pfu_ho(ji ,jj) ) & ! positive/negative part of the flux 942 & + MAX( 0., pfv_ho(ji,jj-1) ) - MIN( 0., pfv_ho(ji,jj ) ) 802 943 zneg = MAX( 0., pfu_ho(ji ,jj) ) - MIN( 0., pfu_ho(ji-1,jj) ) & 803 944 & + MAX( 0., pfv_ho(ji,jj ) ) - MIN( 0., pfv_ho(ji,jj-1) ) 804 945 ! 805 946 ! ! up & down beta terms 806 !!clem zbetup(ji,jj) = ( zup - pt_ups(ji,jj) ) / ( zpos + zsml ) * e1e2t(ji,jj) * z1_dt 807 !!clem zbetdo(ji,jj) = ( pt_ups(ji,jj) - zdo ) / ( zneg + zsml ) * e1e2t(ji,jj) * z1_dt 808 IF( zpos >= epsi20 ) THEN 809 zbetup(ji,jj) = ( zup - pt_ups(ji,jj) ) / zpos * e1e2t(ji,jj) * z1_dt 810 ELSE 811 zbetup(ji,jj) = zbig 947 !!clem zbetup(ji,jj) = ( zup - pt_low(ji,jj) ) / ( zpos + zsml ) * e1e2t(ji,jj) * z1_dt 948 !!clem zbetdo(ji,jj) = ( pt_low(ji,jj) - zdo ) / ( zneg + zsml ) * e1e2t(ji,jj) * z1_dt 949 IF( zpos >= epsi20 ) THEN ; zbetup(ji,jj) = ( zup - pt_low(ji,jj) ) / zpos * e1e2t(ji,jj) * z1_dt 950 ELSE ; zbetup(ji,jj) = zbig 812 951 ENDIF 813 952 ! 814 IF( zneg >= epsi20 ) THEN 815 zbetdo(ji,jj) = ( pt_ups(ji,jj) - zdo ) / zneg * e1e2t(ji,jj) * z1_dt 816 ELSE 817 zbetdo(ji,jj) = zbig 953 IF( zneg >= epsi20 ) THEN ; zbetdo(ji,jj) = ( pt_low(ji,jj) - zdo ) / zneg * e1e2t(ji,jj) * z1_dt 954 ELSE ; zbetdo(ji,jj) = zbig 818 955 ENDIF 819 956 ! … … 830 967 zcu = 0.5 + SIGN( 0.5 , pfu_ho(ji,jj) ) 831 968 ! 832 pfu_ho(ji,jj) = pfu_ho(ji,jj) * ( zcu * zau + ( 1._wp - zcu ) * zbu ) + pfu_ ups(ji,jj)969 pfu_ho(ji,jj) = pfu_ho(ji,jj) * ( zcu * zau + ( 1._wp - zcu ) * zbu ) + pfu_low(ji,jj) 833 970 END DO 834 971 END DO … … 840 977 zcv = 0.5 + SIGN( 0.5 , pfv_ho(ji,jj) ) 841 978 ! 842 pfv_ho(ji,jj) = pfv_ho(ji,jj) * ( zcv * zav + ( 1._wp - zcv ) * zbv ) + pfv_ ups(ji,jj)979 pfv_ho(ji,jj) = pfv_ho(ji,jj) * ( zcv * zav + ( 1._wp - zcv ) * zbv ) + pfv_low(ji,jj) 843 980 END DO 844 981 END DO
Note: See TracChangeset
for help on using the changeset viewer.