Changeset 14072 for NEMO/trunk/src/ICE/icedyn_adv_umx.F90
- Timestamp:
- 2020-12-04T08:48:38+01:00 (4 years ago)
- File:
-
- 1 edited
-
NEMO/trunk/src/ICE/icedyn_adv_umx.F90 (modified) (52 diffs)
Legend:
- Unmodified
- Added
- Removed
-
NEMO/trunk/src/ICE/icedyn_adv_umx.F90
r14005 r14072 14 14 !! ultimate_x(_y) : compute a tracer value at velocity points using ULTIMATE scheme at various orders 15 15 !! macho : compute the fluxes 16 !! nonosc_ice : limit the fluxes using a non-oscillatory algorithm 16 !! nonosc_ice : limit the fluxes using a non-oscillatory algorithm 17 17 !!---------------------------------------------------------------------- 18 18 USE phycst ! physical constant … … 63 63 !!---------------------------------------------------------------------- 64 64 !! *** ROUTINE ice_dyn_adv_umx *** 65 !! 66 !! ** Purpose : Compute the now trend due to total advection of 65 !! 66 !! ** Purpose : Compute the now trend due to total advection of 67 67 !! tracers and add it to the general trend of tracer equations 68 68 !! using an "Ultimate-Macho" scheme 69 69 !! 70 !! Reference : Leonard, B.P., 1991, Comput. Methods Appl. Mech. Eng., 88, 17-74. 70 !! Reference : Leonard, B.P., 1991, Comput. Methods Appl. Mech. Eng., 88, 17-74. 71 71 !!---------------------------------------------------------------------- 72 72 INTEGER , INTENT(in ) :: kn_umx ! order of the scheme (1-5=UM or 20=CEN2) … … 103 103 REAL(wp), DIMENSION(jpi,jpj,nlay_s,jpl) :: ze_s, zes_max 104 104 ! 105 REAL(wp), ALLOCATABLE, DIMENSION(:,:,:) :: zuv_ho, zvv_ho, zuv_ups, zvv_ups, z1_vi, z1_vs 105 REAL(wp), ALLOCATABLE, DIMENSION(:,:,:) :: zuv_ho, zvv_ho, zuv_ups, zvv_ups, z1_vi, z1_vs 106 106 !! diagnostics 107 REAL(wp), DIMENSION(jpi,jpj) :: zdiag_adv_mass, zdiag_adv_salt, zdiag_adv_heat 107 REAL(wp), DIMENSION(jpi,jpj) :: zdiag_adv_mass, zdiag_adv_salt, zdiag_adv_heat 108 108 !!---------------------------------------------------------------------- 109 109 ! … … 131 131 ELSEWHERE ; ze_s(:,:,jk,:) = 0._wp 132 132 END WHERE 133 END DO 133 END DO 134 134 CALL icemax4D( ze_i , zei_max ) 135 135 CALL icemax4D( ze_s , zes_max ) … … 143 143 zcflnow(1) = MAXVAL( ABS( pu_ice(:,:) ) * rDt_ice * r1_e1u(:,:) ) 144 144 zcflnow(1) = MAX( zcflnow(1), MAXVAL( ABS( pv_ice(:,:) ) * rDt_ice * r1_e2v(:,:) ) ) 145 145 146 146 ! non-blocking global communication send zcflnow and receive zcflprv 147 147 CALL mpp_delay_max( 'icedyn_adv_umx', 'cflice', zcflnow(:), zcflprv(:), kt == nitend - nn_fsbc + 1 ) … … 157 157 zvdx(:,:) = pv_ice(:,:) * e1v(:,:) 158 158 ! 159 ! setup transport for each ice cat 159 ! setup transport for each ice cat 160 160 DO jl = 1, jpl 161 161 zu_cat(:,:,jl) = zudy(:,:) … … 190 190 ! record at_i before advection (for open water) 191 191 zati1(:,:) = SUM( pa_i(:,:,:), dim=3 ) 192 192 193 193 ! inverse of A and Ap 194 194 WHERE( pa_i(:,:,:) >= epsi20 ) ; z1_ai(:,:,:) = 1._wp / pa_i(:,:,:) … … 201 201 ! setup a mask where advection will be upstream 202 202 IF( ll_neg ) THEN 203 IF( .NOT. ALLOCATED(imsk_small) ) ALLOCATE( imsk_small(jpi,jpj,jpl) ) 204 IF( .NOT. ALLOCATED(jmsk_small) ) ALLOCATE( jmsk_small(jpi,jpj,jpl) ) 203 IF( .NOT. ALLOCATED(imsk_small) ) ALLOCATE( imsk_small(jpi,jpj,jpl) ) 204 IF( .NOT. ALLOCATED(jmsk_small) ) ALLOCATE( jmsk_small(jpi,jpj,jpl) ) 205 205 DO jl = 1, jpl 206 206 DO_2D( 1, 0, 1, 0 ) … … 232 232 CALL adv_umx( zamsk, kn_umx, jt, kt, zdt, zudy , zvdx, zua_ho , zva_ho , zcu_box, zcv_box, & 233 233 & zhvar, pv_i, zua_ups, zva_ups ) 234 !== Snw volume ==! 234 !== Snw volume ==! 235 235 zhvar(:,:,:) = pv_s(:,:,:) * z1_ai(:,:,:) 236 236 CALL adv_umx( zamsk, kn_umx, jt, kt, zdt, zudy , zvdx, zua_ho , zva_ho , zcu_box, zcv_box, & … … 260 260 CALL adv_umx( zamsk, kn_umx, jt, kt, zdt, zudy , zvdx, zua_ho , zva_ho , zcu_box, zcv_box, & 261 261 & zhvar, pv_i, zua_ups, zva_ups ) 262 !== Snw volume ==! 262 !== Snw volume ==! 263 263 zhvar(:,:,:) = pv_s(:,:,:) * z1_ai(:,:,:) 264 264 CALL adv_umx( zamsk, kn_umx, jt, kt, zdt, zudy , zvdx, zua_ho , zva_ho , zcu_box, zcv_box, & … … 316 316 & zhvar, pe_i(:,:,jk,:), zuv_ups, zvv_ups ) 317 317 END DO 318 !== Snow volume ==! 318 !== Snow volume ==! 319 319 zuv_ups = zua_ups 320 320 zvv_ups = zva_ups … … 374 374 zati2(:,:) = SUM( pa_i(:,:,:), dim=3 ) 375 375 DO_2D( 0, 0, 0, 0 ) 376 pato_i(ji,jj) = pato_i(ji,jj) - ( zati2(ji,jj) - zati1(ji,jj) ) & 376 pato_i(ji,jj) = pato_i(ji,jj) - ( zati2(ji,jj) - zati1(ji,jj) ) & 377 377 & - ( zudy(ji,jj) - zudy(ji-1,jj) + zvdx(ji,jj) - zvdx(ji,jj-1) ) * r1_e1e2t(ji,jj) * zdt 378 378 END_2D … … 406 406 END SUBROUTINE ice_dyn_adv_umx 407 407 408 408 409 409 SUBROUTINE adv_umx( pamsk, kn_umx, jt, kt, pdt, pu, pv, puc, pvc, pubox, pvbox, & 410 410 & pt, ptc, pua_ups, pva_ups, pua_ho, pva_ho ) 411 411 !!---------------------------------------------------------------------- 412 412 !! *** ROUTINE adv_umx *** 413 !! 414 !! ** Purpose : Compute the now trend due to total advection of 413 !! 414 !! ** Purpose : Compute the now trend due to total advection of 415 415 !! tracers and add it to the general trend of tracer equations 416 416 !! … … 434 434 !! 435 435 !! in eq. c), one can solve the equation for S (ln_advS=T), then dVS/dt = -div(uV * uS / u) 436 !! or for HS (ln_advS=F), then dVS/dt = -div(uA * uHS / u) 436 !! or for HS (ln_advS=F), then dVS/dt = -div(uA * uHS / u) 437 437 !! 438 438 !! ** Note : - this method can lead to tiny negative V (-1.e-20) => set it to 0 while conserving mass etc. … … 462 462 REAL(wp), DIMENSION(jpi,jpj,jpl), INTENT( out), OPTIONAL :: pua_ho, pva_ho ! high order u*a fluxes 463 463 ! 464 INTEGER :: ji, jj, jl ! dummy loop indices 464 INTEGER :: ji, jj, jl ! dummy loop indices 465 465 REAL(wp) :: ztra ! local scalar 466 466 REAL(wp), DIMENSION(jpi,jpj,jpl) :: zfu_ho , zfv_ho , zpt … … 468 468 !!---------------------------------------------------------------------- 469 469 ! 470 ! Upstream (_ups) fluxes 470 ! Upstream (_ups) fluxes 471 471 ! ----------------------- 472 472 CALL upstream( pamsk, jt, kt, pdt, pt, pu, pv, zt_ups, zfu_ups, zfv_ups ) 473 474 ! High order (_ho) fluxes 473 474 ! High order (_ho) fluxes 475 475 ! ----------------------- 476 476 SELECT CASE( kn_umx ) … … 506 506 zfv_ups(ji,jj,jl) = zfv_ups(ji,jj,jl) * pva_ups(ji,jj,jl) / pv(ji,jj) 507 507 ELSE 508 zfv_ho (ji,jj,jl) = 0._wp 509 zfv_ups(ji,jj,jl) = 0._wp 508 zfv_ho (ji,jj,jl) = 0._wp 509 zfv_ups(ji,jj,jl) = 0._wp 510 510 ENDIF 511 511 END_2D … … 551 551 DO jl = 1, jpl 552 552 DO_2D( 0, 0, 0, 0 ) 553 ztra = - ( zfu_ho(ji,jj,jl) - zfu_ho(ji-1,jj,jl) + zfv_ho(ji,jj,jl) - zfv_ho(ji,jj-1,jl) ) 554 ! 555 ptc(ji,jj,jl) = ( ptc(ji,jj,jl) + ztra * r1_e1e2t(ji,jj) * pdt ) * tmask(ji,jj,1) 553 ztra = - ( zfu_ho(ji,jj,jl) - zfu_ho(ji-1,jj,jl) + zfv_ho(ji,jj,jl) - zfv_ho(ji,jj-1,jl) ) 554 ! 555 ptc(ji,jj,jl) = ( ptc(ji,jj,jl) + ztra * r1_e1e2t(ji,jj) * pdt ) * tmask(ji,jj,1) 556 556 END_2D 557 557 END DO … … 563 563 !!--------------------------------------------------------------------- 564 564 !! *** ROUTINE upstream *** 565 !! 565 !! 566 566 !! ** Purpose : compute the upstream fluxes and upstream guess of tracer 567 567 !!---------------------------------------------------------------------- … … 572 572 REAL(wp), DIMENSION(:,:,:) , INTENT(in ) :: pt ! tracer fields 573 573 REAL(wp), DIMENSION(:,: ) , INTENT(in ) :: pu, pv ! 2 ice velocity components 574 REAL(wp), DIMENSION(jpi,jpj,jpl), INTENT( out) :: pt_ups ! upstream guess of tracer 575 REAL(wp), DIMENSION(jpi,jpj,jpl), INTENT( out) :: pfu_ups, pfv_ups ! upstream fluxes 574 REAL(wp), DIMENSION(jpi,jpj,jpl), INTENT( out) :: pt_ups ! upstream guess of tracer 575 REAL(wp), DIMENSION(jpi,jpj,jpl), INTENT( out) :: pfu_ups, pfv_ups ! upstream fluxes 576 576 ! 577 577 INTEGER :: ji, jj, jl ! dummy loop indices … … 638 638 ! 639 639 ENDIF 640 640 641 641 ENDIF 642 642 ! … … 655 655 END SUBROUTINE upstream 656 656 657 657 658 658 SUBROUTINE cen2( pamsk, jt, kt, pdt, pt, pu, pv, pt_ups, pfu_ups, pfv_ups, pfu_ho, pfv_ho ) 659 659 !!--------------------------------------------------------------------- 660 660 !! *** ROUTINE cen2 *** 661 !! 661 !! 662 662 !! ** Purpose : compute the high order fluxes using a centered 663 !! second order scheme 663 !! second order scheme 664 664 !!---------------------------------------------------------------------- 665 665 REAL(wp) , INTENT(in ) :: pamsk ! advection of concentration (1) or other tracers (0) … … 669 669 REAL(wp), DIMENSION(:,:,:) , INTENT(in ) :: pt ! tracer fields 670 670 REAL(wp), DIMENSION(:,: ) , INTENT(in ) :: pu, pv ! 2 ice velocity components 671 REAL(wp), DIMENSION(:,:,:) , INTENT(in ) :: pt_ups ! upstream guess of tracer 672 REAL(wp), DIMENSION(:,:,:) , INTENT(in ) :: pfu_ups, pfv_ups ! upstream fluxes 673 REAL(wp), DIMENSION(jpi,jpj,jpl), INTENT( out) :: pfu_ho, pfv_ho ! high order fluxes 671 REAL(wp), DIMENSION(:,:,:) , INTENT(in ) :: pt_ups ! upstream guess of tracer 672 REAL(wp), DIMENSION(:,:,:) , INTENT(in ) :: pfu_ups, pfv_ups ! upstream fluxes 673 REAL(wp), DIMENSION(jpi,jpj,jpl), INTENT( out) :: pfu_ho, pfv_ho ! high order fluxes 674 674 ! 675 675 INTEGER :: ji, jj, jl ! dummy loop indices … … 750 750 ENDIF 751 751 IF( np_limiter == 1 ) CALL nonosc_ice( pamsk, pdt, pu, pv, pt, pt_ups, pfu_ups, pfv_ups, pfu_ho, pfv_ho ) 752 752 753 753 ENDIF 754 754 755 755 END SUBROUTINE cen2 756 756 757 757 758 758 SUBROUTINE macho( pamsk, kn_umx, jt, kt, pdt, pt, pu, pv, pubox, pvbox, pt_ups, pfu_ups, pfv_ups, pfu_ho, pfv_ho ) 759 759 !!--------------------------------------------------------------------- 760 760 !! *** ROUTINE macho *** 761 !! 762 !! ** Purpose : compute the high order fluxes using Ultimate-Macho scheme 761 !! 762 !! ** Purpose : compute the high order fluxes using Ultimate-Macho scheme 763 763 !! 764 764 !! ** Method : ... 765 765 !! 766 !! Reference : Leonard, B.P., 1991, Comput. Methods Appl. Mech. Eng., 88, 17-74. 766 !! Reference : Leonard, B.P., 1991, Comput. Methods Appl. Mech. Eng., 88, 17-74. 767 767 !!---------------------------------------------------------------------- 768 768 REAL(wp) , INTENT(in ) :: pamsk ! advection of concentration (1) or other tracers (0) … … 774 774 REAL(wp), DIMENSION(:,: ) , INTENT(in ) :: pu, pv ! 2 ice velocity components 775 775 REAL(wp), DIMENSION(:,: ) , INTENT(in ) :: pubox, pvbox ! upstream velocity 776 REAL(wp), DIMENSION(:,:,:) , INTENT(in ) :: pt_ups ! upstream guess of tracer 777 REAL(wp), DIMENSION(:,:,:) , INTENT(in ) :: pfu_ups, pfv_ups ! upstream fluxes 778 REAL(wp), DIMENSION(jpi,jpj,jpl), INTENT( out) :: pfu_ho, pfv_ho ! high order fluxes 776 REAL(wp), DIMENSION(:,:,:) , INTENT(in ) :: pt_ups ! upstream guess of tracer 777 REAL(wp), DIMENSION(:,:,:) , INTENT(in ) :: pfu_ups, pfv_ups ! upstream fluxes 778 REAL(wp), DIMENSION(jpi,jpj,jpl), INTENT( out) :: pfu_ho, pfv_ho ! high order fluxes 779 779 ! 780 780 INTEGER :: ji, jj, jl ! dummy loop indices … … 807 807 ! !-- limiter in y --! 808 808 IF( np_limiter == 2 .OR. np_limiter == 3 ) CALL limiter_y( pdt, pv, pt, pfv_ups, pfv_ho ) 809 ! 809 ! 810 810 ! 811 811 ELSE !== even ice time step: adv_y then adv_x ==! … … 821 821 & + pt (ji,jj,jl) * ( pv (ji,jj ) - pv (ji,jj-1 ) ) * r1_e1e2t(ji,jj) & 822 822 & * pamsk & 823 & ) * pdt ) * tmask(ji,jj,1) 823 & ) * pdt ) * tmask(ji,jj,1) 824 824 END_2D 825 825 END DO … … 845 845 !!--------------------------------------------------------------------- 846 846 !! *** ROUTINE ultimate_x *** 847 !! 848 !! ** Purpose : compute tracer at u-points 847 !! 848 !! ** Purpose : compute tracer at u-points 849 849 !! 850 850 !! ** Method : ... 851 851 !! 852 !! Reference : Leonard, B.P., 1991, Comput. Methods Appl. Mech. Eng., 88, 17-74. 852 !! Reference : Leonard, B.P., 1991, Comput. Methods Appl. Mech. Eng., 88, 17-74. 853 853 !!---------------------------------------------------------------------- 854 854 REAL(wp) , INTENT(in ) :: pamsk ! advection of concentration (1) or other tracers (0) … … 857 857 REAL(wp), DIMENSION(:,: ) , INTENT(in ) :: pu ! ice i-velocity component 858 858 REAL(wp), DIMENSION(:,:,:) , INTENT(in ) :: pt ! tracer fields 859 REAL(wp), DIMENSION(jpi,jpj,jpl), INTENT( out) :: pt_u ! tracer at u-point 860 REAL(wp), DIMENSION(jpi,jpj,jpl), INTENT( out) :: pfu_ho ! high order flux 859 REAL(wp), DIMENSION(jpi,jpj,jpl), INTENT( out) :: pt_u ! tracer at u-point 860 REAL(wp), DIMENSION(jpi,jpj,jpl), INTENT( out) :: pfu_ho ! high order flux 861 861 ! 862 862 INTEGER :: ji, jj, jl ! dummy loop indices … … 897 897 ! 898 898 CASE( 1 ) !== 1st order central TIM ==! (Eq. 21) 899 ! 899 ! 900 900 DO jl = 1, jpl 901 901 DO_2D( 0, 0, 1, 0 ) … … 911 911 zcu = pu(ji,jj) * r1_e2u(ji,jj) * pdt * r1_e1u(ji,jj) 912 912 pt_u(ji,jj,jl) = 0.5_wp * umask(ji,jj,1) * ( pt(ji+1,jj,jl) + pt(ji,jj,jl) & 913 & - zcu * ( pt(ji+1,jj,jl) - pt(ji,jj,jl) ) ) 914 END_2D 915 END DO 916 ! 913 & - zcu * ( pt(ji+1,jj,jl) - pt(ji,jj,jl) ) ) 914 END_2D 915 END DO 916 ! 917 917 CASE( 3 ) !== 3rd order central TIM ==! (Eq. 24) 918 918 ! … … 983 983 ! 984 984 END SUBROUTINE ultimate_x 985 986 985 986 987 987 SUBROUTINE ultimate_y( pamsk, kn_umx, pdt, pt, pv, pt_v, pfv_ho ) 988 988 !!--------------------------------------------------------------------- 989 989 !! *** ROUTINE ultimate_y *** 990 !! 991 !! ** Purpose : compute tracer at v-points 990 !! 991 !! ** Purpose : compute tracer at v-points 992 992 !! 993 993 !! ** Method : ... 994 994 !! 995 !! Reference : Leonard, B.P., 1991, Comput. Methods Appl. Mech. Eng., 88, 17-74. 995 !! Reference : Leonard, B.P., 1991, Comput. Methods Appl. Mech. Eng., 88, 17-74. 996 996 !!---------------------------------------------------------------------- 997 997 REAL(wp) , INTENT(in ) :: pamsk ! advection of concentration (1) or other tracers (0) … … 1000 1000 REAL(wp), DIMENSION(:,: ) , INTENT(in ) :: pv ! ice j-velocity component 1001 1001 REAL(wp), DIMENSION(:,:,:) , INTENT(in ) :: pt ! tracer fields 1002 REAL(wp), DIMENSION(jpi,jpj,jpl), INTENT( out) :: pt_v ! tracer at v-point 1003 REAL(wp), DIMENSION(jpi,jpj,jpl), INTENT( out) :: pfv_ho ! high order flux 1002 REAL(wp), DIMENSION(jpi,jpj,jpl), INTENT( out) :: pt_v ! tracer at v-point 1003 REAL(wp), DIMENSION(jpi,jpj,jpl), INTENT( out) :: pfv_ho ! high order flux 1004 1004 ! 1005 1005 INTEGER :: ji, jj, jl ! dummy loop indices … … 1115 1115 ! 1116 1116 END SUBROUTINE ultimate_y 1117 1117 1118 1118 1119 1119 SUBROUTINE nonosc_ice( pamsk, pdt, pu, pv, pt, pt_ups, pfu_ups, pfv_ups, pfu_ho, pfv_ho ) 1120 1120 !!--------------------------------------------------------------------- 1121 1121 !! *** ROUTINE nonosc_ice *** 1122 !! 1123 !! ** Purpose : compute monotonic tracer fluxes from the upstream 1124 !! scheme and the before field by a non-oscillatory algorithm 1122 !! 1123 !! ** Purpose : compute monotonic tracer fluxes from the upstream 1124 !! scheme and the before field by a non-oscillatory algorithm 1125 1125 !! 1126 1126 !! ** Method : ... … … 1141 1141 !!---------------------------------------------------------------------- 1142 1142 zbig = 1.e+40_wp 1143 1143 1144 1144 ! antidiffusive flux : high order minus low order 1145 1145 ! -------------------------------------------------- … … 1157 1157 ! pfu_ho 1158 1158 ! * ---> 1159 ! | | * | | 1160 ! | | | * | 1159 ! | | * | | 1160 ! | | | * | 1161 1161 ! | | | | * 1162 ! t_ups : i-1 i i+1 i+2 1162 ! t_ups : i-1 i i+1 i+2 1163 1163 IF( ll_prelim ) THEN 1164 1164 1165 1165 DO jl = 1, jpl 1166 1166 DO_2D( 0, 0, 0, 0 ) … … 1200 1200 z1_dt = 1._wp / pdt 1201 1201 DO jl = 1, jpl 1202 1202 1203 1203 DO_2D( 1, 1, 1, 1 ) 1204 1204 IF ( pt(ji,jj,jl) <= 0._wp .AND. pt_ups(ji,jj,jl) <= 0._wp ) THEN … … 1244 1244 ! if all the points are outside ice cover 1245 1245 IF( zup == -zbig ) zbetup(ji,jj,jl) = 0._wp ! zbig 1246 IF( zdo == zbig ) zbetdo(ji,jj,jl) = 0._wp ! zbig 1246 IF( zdo == zbig ) zbetdo(ji,jj,jl) = 0._wp ! zbig 1247 1247 ! 1248 1248 END_2D … … 1250 1250 CALL lbc_lnk_multi( 'icedyn_adv_umx', zbetup, 'T', 1.0_wp, zbetdo, 'T', 1.0_wp ) ! lateral boundary cond. (unchanged sign) 1251 1251 1252 1252 1253 1253 ! monotonic flux in the y direction 1254 1254 ! --------------------------------- … … 1280 1280 END SUBROUTINE nonosc_ice 1281 1281 1282 1282 1283 1283 SUBROUTINE limiter_x( pdt, pu, pt, pfu_ups, pfu_ho ) 1284 1284 !!--------------------------------------------------------------------- 1285 1285 !! *** ROUTINE limiter_x *** 1286 !! 1287 !! ** Purpose : compute flux limiter 1286 !! 1287 !! ** Purpose : compute flux limiter 1288 1288 !!---------------------------------------------------------------------- 1289 1289 REAL(wp) , INTENT(in ) :: pdt ! tracer time-step … … 1295 1295 REAL(wp) :: Cr, Rjm, Rj, Rjp, uCFL, zpsi, zh3, zlimiter, Rr 1296 1296 INTEGER :: ji, jj, jl ! dummy loop indices 1297 REAL(wp), DIMENSION (jpi,jpj,jpl) :: zslpx ! tracer slopes 1297 REAL(wp), DIMENSION (jpi,jpj,jpl) :: zslpx ! tracer slopes 1298 1298 !!---------------------------------------------------------------------- 1299 1299 ! … … 1304 1304 END DO 1305 1305 CALL lbc_lnk( 'icedyn_adv_umx', zslpx, 'U', -1.0_wp) ! lateral boundary cond. 1306 1306 1307 1307 DO jl = 1, jpl 1308 1308 DO_2D( 0, 0, 0, 0 ) 1309 1309 uCFL = pdt * ABS( pu(ji,jj) ) * r1_e1e2t(ji,jj) 1310 1310 1311 1311 Rjm = zslpx(ji-1,jj,jl) 1312 1312 Rj = zslpx(ji ,jj,jl) … … 1319 1319 ENDIF 1320 1320 1321 zh3 = pfu_ho(ji,jj,jl) - pfu_ups(ji,jj,jl) 1321 zh3 = pfu_ho(ji,jj,jl) - pfu_ups(ji,jj,jl) 1322 1322 IF( Rj > 0. ) THEN 1323 1323 zlimiter = MAX( 0., MIN( zh3, MAX(-Rr * 0.5 * ABS(pu(ji,jj)), & … … 1371 1371 END SUBROUTINE limiter_x 1372 1372 1373 1373 1374 1374 SUBROUTINE limiter_y( pdt, pv, pt, pfv_ups, pfv_ho ) 1375 1375 !!--------------------------------------------------------------------- 1376 1376 !! *** ROUTINE limiter_y *** 1377 !! 1378 !! ** Purpose : compute flux limiter 1377 !! 1378 !! ** Purpose : compute flux limiter 1379 1379 !!---------------------------------------------------------------------- 1380 1380 REAL(wp) , INTENT(in ) :: pdt ! tracer time-step … … 1386 1386 REAL(wp) :: Cr, Rjm, Rj, Rjp, vCFL, zpsi, zh3, zlimiter, Rr 1387 1387 INTEGER :: ji, jj, jl ! dummy loop indices 1388 REAL(wp), DIMENSION (jpi,jpj,jpl) :: zslpy ! tracer slopes 1388 REAL(wp), DIMENSION (jpi,jpj,jpl) :: zslpy ! tracer slopes 1389 1389 !!---------------------------------------------------------------------- 1390 1390 ! … … 1410 1410 ENDIF 1411 1411 1412 zh3 = pfv_ho(ji,jj,jl) - pfv_ups(ji,jj,jl) 1412 zh3 = pfv_ho(ji,jj,jl) - pfv_ups(ji,jj,jl) 1413 1413 IF( Rj > 0. ) THEN 1414 1414 zlimiter = MAX( 0., MIN( zh3, MAX(-Rr * 0.5 * ABS(pv(ji,jj)), & … … 1524 1524 pe_s(ji,jj,1:nlay_s,jl) = pe_s(ji,jj,1:nlay_s,jl) * zfra 1525 1525 pv_s(ji,jj,jl) = pa_i(ji,jj,jl) * phs_max(ji,jj,jl) 1526 ENDIF 1527 ! 1526 ENDIF 1527 ! 1528 1528 ! ! -- check s_i -- ! 1529 1529 ! if s_i is larger than the surrounding 9 pts => put salt excess in the ocean … … 1537 1537 ENDIF 1538 1538 END_2D 1539 END DO 1539 END DO 1540 1540 ! 1541 1541 ! ! -- check e_i/v_i -- ! … … 1624 1624 SUBROUTINE icemax3D( pice , pmax ) 1625 1625 !!--------------------------------------------------------------------- 1626 !! *** ROUTINE icemax3D *** 1626 !! *** ROUTINE icemax3D *** 1627 1627 !! ** Purpose : compute the max of the 9 points around 1628 1628 !!---------------------------------------------------------------------- … … 1633 1633 !!---------------------------------------------------------------------- 1634 1634 DO jl = 1, jpl 1635 DO jj = Njs0-1, Nje0+1 1635 DO jj = Njs0-1, Nje0+1 1636 1636 DO ji = Nis0, Nie0 1637 1637 zmax(ji,jj) = MAX( epsi20, pice(ji,jj,jl), pice(ji-1,jj,jl), pice(ji+1,jj,jl) ) 1638 1638 END DO 1639 1639 END DO 1640 DO jj = Njs0, Nje0 1640 DO jj = Njs0, Nje0 1641 1641 DO ji = Nis0, Nie0 1642 1642 pmax(ji,jj,jl) = MAX( epsi20, zmax(ji,jj), zmax(ji,jj-1), zmax(ji,jj+1) ) … … 1648 1648 SUBROUTINE icemax4D( pice , pmax ) 1649 1649 !!--------------------------------------------------------------------- 1650 !! *** ROUTINE icemax4D *** 1650 !! *** ROUTINE icemax4D *** 1651 1651 !! ** Purpose : compute the max of the 9 points around 1652 1652 !!---------------------------------------------------------------------- … … 1659 1659 DO jl = 1, jpl 1660 1660 DO jk = 1, jlay 1661 DO jj = Njs0-1, Nje0+1 1661 DO jj = Njs0-1, Nje0+1 1662 1662 DO ji = Nis0, Nie0 1663 1663 zmax(ji,jj) = MAX( epsi20, pice(ji,jj,jk,jl), pice(ji-1,jj,jk,jl), pice(ji+1,jj,jk,jl) ) 1664 1664 END DO 1665 1665 END DO 1666 DO jj = Njs0, Nje0 1666 DO jj = Njs0, Nje0 1667 1667 DO ji = Nis0, Nie0 1668 1668 pmax(ji,jj,jk,jl) = MAX( epsi20, zmax(ji,jj), zmax(ji,jj-1), zmax(ji,jj+1) )
Note: See TracChangeset
for help on using the changeset viewer.
