Changeset 11536 for NEMO/trunk/src/ICE/icevar.F90
- Timestamp:
- 2019-09-11T15:54:18+02:00 (5 years ago)
- File:
-
- 1 edited
Legend:
- Unmodified
- Added
- Removed
-
NEMO/trunk/src/ICE/icevar.F90
r11229 r11536 32 32 !! - vt_s(jpi,jpj) 33 33 !! - at_i(jpi,jpj) 34 !! - st_i(jpi,jpj) 34 35 !! - et_s(jpi,jpj) total snow heat content 35 36 !! - et_i(jpi,jpj) total ice thermal content … … 46 47 !! ice_var_zapneg : remove negative ice fields 47 48 !! ice_var_roundoff : remove negative values arising from roundoff erros 48 !! ice_var_itd : convert 1-cat to jpl-cat49 !! ice_var_itd2 : convert N-cat to jpl-cat50 49 !! ice_var_bv : brine volume 51 50 !! ice_var_enthalpy : compute ice and snow enthalpies from temperature 52 51 !! ice_var_sshdyn : compute equivalent ssh in lead 52 !! ice_var_itd : convert N-cat to M-cat 53 53 !!---------------------------------------------------------------------- 54 54 USE dom_oce ! ocean space and time domain … … 104 104 ! 105 105 ! ! integrated values 106 vt_i(:,:) = SUM( v_i(:,:,:) , dim=3 ) 107 vt_s(:,:) = SUM( v_s(:,:,:) , dim=3 ) 108 at_i(:,:) = SUM( a_i(:,:,:) , dim=3 ) 109 et_s(:,:) = SUM( SUM( e_s(:,:,:,:), dim=4 ), dim=3 ) 110 et_i(:,:) = SUM( SUM( e_i(:,:,:,:), dim=4 ), dim=3 ) 106 vt_i(:,:) = SUM( v_i (:,:,:) , dim=3 ) 107 vt_s(:,:) = SUM( v_s (:,:,:) , dim=3 ) 108 st_i(:,:) = SUM( sv_i(:,:,:) , dim=3 ) 109 at_i(:,:) = SUM( a_i (:,:,:) , dim=3 ) 110 et_s(:,:) = SUM( SUM( e_s (:,:,:,:), dim=4 ), dim=3 ) 111 et_i(:,:) = SUM( SUM( e_i (:,:,:,:), dim=4 ), dim=3 ) 111 112 ! 112 113 at_ip(:,:) = SUM( a_ip(:,:,:), dim=3 ) ! melt ponds … … 138 139 tm_si(:,:) = SUM( t_si(:,:,:) * a_i(:,:,:) , dim=3 ) * z1_at_i(:,:) 139 140 om_i (:,:) = SUM( oa_i(:,:,:) , dim=3 ) * z1_at_i(:,:) 140 sm_i (:,:) = SUM( sv_i(:,:,:) , dim=3 )* z1_vt_i(:,:)141 sm_i (:,:) = st_i(:,:) * z1_vt_i(:,:) 141 142 ! 142 143 tm_i(:,:) = 0._wp … … 158 159 tm_s (:,:) = rt0 159 160 END WHERE 160 161 ! 162 ! ! mean melt pond depth 163 WHERE( at_ip(:,:) > epsi20 ) ; hm_ip(:,:) = vt_ip(:,:) / at_ip(:,:) 164 ELSEWHERE ; hm_ip(:,:) = 0._wp 165 END WHERE 166 ! 161 167 DEALLOCATE( z1_at_i , z1_vt_i , z1_vt_s ) 168 ! 162 169 ENDIF 163 170 ! … … 263 270 ! 264 271 ! integrated values 265 vt_i (:,:) = SUM( v_i , dim=3 )266 vt_s (:,:) = SUM( v_s , dim=3 )267 at_i (:,:) = SUM( a_i , dim=3 )272 vt_i (:,:) = SUM( v_i , dim=3 ) 273 vt_s (:,:) = SUM( v_s , dim=3 ) 274 at_i (:,:) = SUM( a_i , dim=3 ) 268 275 ! 269 276 END SUBROUTINE ice_var_glo2eqv … … 533 540 534 541 ! to be sure that at_i is the sum of a_i(jl) 535 at_i (:,:) = SUM( a_i(:,:,:), dim=3 ) 536 vt_i (:,:) = SUM( v_i(:,:,:), dim=3 ) 542 at_i (:,:) = SUM( a_i (:,:,:), dim=3 ) 543 vt_i (:,:) = SUM( v_i (:,:,:), dim=3 ) 544 !!clem add? 545 ! vt_s (:,:) = SUM( v_s (:,:,:), dim=3 ) 546 ! st_i (:,:) = SUM( sv_i(:,:,:), dim=3 ) 547 ! et_s(:,:) = SUM( SUM( e_s (:,:,:,:), dim=4 ), dim=3 ) 548 ! et_i(:,:) = SUM( SUM( e_i (:,:,:,:), dim=4 ), dim=3 ) 549 !!clem 537 550 538 551 ! open water = 1 if at_i=0 … … 652 665 WHERE( pe_i (1:npti,:,:) < 0._wp .AND. pe_i (1:npti,:,:) > -epsi06 ) pe_i (1:npti,:,:) = 0._wp ! e_i must be >= 0 653 666 WHERE( pe_s (1:npti,:,:) < 0._wp .AND. pe_s (1:npti,:,:) > -epsi06 ) pe_s (1:npti,:,:) = 0._wp ! e_s must be >= 0 654 IF 667 IF( ln_pnd_H12 ) THEN 655 668 WHERE( pa_ip(1:npti,:) < 0._wp .AND. pa_ip(1:npti,:) > -epsi10 ) pa_ip(1:npti,:) = 0._wp ! a_ip must be >= 0 656 669 WHERE( pv_ip(1:npti,:) < 0._wp .AND. pv_ip(1:npti,:) > -epsi10 ) pv_ip(1:npti,:) = 0._wp ! v_ip must be >= 0 … … 773 786 !! ** Purpose : converting N-cat ice to jpl ice categories 774 787 !!------------------------------------------------------------------- 775 SUBROUTINE ice_var_itd_1c1c( zhti, zhts, zati, zh_i, zh_s, za_i ) 788 SUBROUTINE ice_var_itd_1c1c( phti, phts, pati , ph_i, ph_s, pa_i, & 789 & ptmi, ptms, ptmsu, psmi, patip, phtip, pt_i, pt_s, pt_su, ps_i, pa_ip, ph_ip ) 776 790 !!------------------------------------------------------------------- 777 791 !! ** Purpose : converting 1-cat ice to 1 ice category 778 792 !!------------------------------------------------------------------- 779 REAL(wp), DIMENSION(:), INTENT(in) :: zhti, zhts, zati ! input ice/snow variables 780 REAL(wp), DIMENSION(:), INTENT(inout) :: zh_i, zh_s, za_i ! output ice/snow variables 781 !!------------------------------------------------------------------- 782 zh_i(:) = zhti(:) 783 zh_s(:) = zhts(:) 784 za_i(:) = zati(:) 793 REAL(wp), DIMENSION(:), INTENT(in) :: phti, phts, pati ! input ice/snow variables 794 REAL(wp), DIMENSION(:), INTENT(inout) :: ph_i, ph_s, pa_i ! output ice/snow variables 795 REAL(wp), DIMENSION(:), INTENT(in) :: ptmi, ptms, ptmsu, psmi, patip, phtip ! input ice/snow temp & sal & ponds 796 REAL(wp), DIMENSION(:), INTENT(inout) :: pt_i, pt_s, pt_su, ps_i, pa_ip, ph_ip ! output ice/snow temp & sal & ponds 797 !!------------------------------------------------------------------- 798 ! == thickness and concentration == ! 799 ph_i(:) = phti(:) 800 ph_s(:) = phts(:) 801 pa_i(:) = pati(:) 802 ! 803 ! == temperature and salinity and ponds == ! 804 pt_i (:) = ptmi (:) 805 pt_s (:) = ptms (:) 806 pt_su(:) = ptmsu(:) 807 ps_i (:) = psmi (:) 808 pa_ip(:) = patip(:) 809 ph_ip(:) = phtip(:) 810 785 811 END SUBROUTINE ice_var_itd_1c1c 786 812 787 SUBROUTINE ice_var_itd_Nc1c( zhti, zhts, zati, zh_i, zh_s, za_i ) 813 SUBROUTINE ice_var_itd_Nc1c( phti, phts, pati , ph_i, ph_s, pa_i, & 814 & ptmi, ptms, ptmsu, psmi, patip, phtip, pt_i, pt_s, pt_su, ps_i, pa_ip, ph_ip ) 788 815 !!------------------------------------------------------------------- 789 816 !! ** Purpose : converting N-cat ice to 1 ice category 790 817 !!------------------------------------------------------------------- 791 REAL(wp), DIMENSION(:,:), INTENT(in) :: zhti, zhts, zati ! input ice/snow variables 792 REAL(wp), DIMENSION(:) , INTENT(inout) :: zh_i, zh_s, za_i ! output ice/snow variables 793 !!------------------------------------------------------------------- 794 ! 795 za_i(:) = SUM( zati(:,:), dim=2 ) 796 ! 797 WHERE( za_i(:) /= 0._wp ) 798 zh_i(:) = SUM( zhti(:,:) * zati(:,:), dim=2 ) / za_i(:) 799 zh_s(:) = SUM( zhts(:,:) * zati(:,:), dim=2 ) / za_i(:) 800 ELSEWHERE 801 zh_i(:) = 0._wp 802 zh_s(:) = 0._wp 818 REAL(wp), DIMENSION(:,:), INTENT(in) :: phti, phts, pati ! input ice/snow variables 819 REAL(wp), DIMENSION(:) , INTENT(inout) :: ph_i, ph_s, pa_i ! output ice/snow variables 820 REAL(wp), DIMENSION(:,:), INTENT(in) :: ptmi, ptms, ptmsu, psmi, patip, phtip ! input ice/snow temp & sal & ponds 821 REAL(wp), DIMENSION(:) , INTENT(inout) :: pt_i, pt_s, pt_su, ps_i, pa_ip, ph_ip ! output ice/snow temp & sal & ponds 822 ! 823 REAL(wp), ALLOCATABLE, DIMENSION(:) :: z1_ai, z1_vi, z1_vs 824 ! 825 INTEGER :: idim 826 !!------------------------------------------------------------------- 827 ! 828 idim = SIZE( phti, 1 ) 829 ! 830 ! == thickness and concentration == ! 831 ALLOCATE( z1_ai(idim), z1_vi(idim), z1_vs(idim) ) 832 ! 833 pa_i(:) = SUM( pati(:,:), dim=2 ) 834 835 WHERE( ( pa_i(:) ) /= 0._wp ) ; z1_ai(:) = 1._wp / pa_i(:) 836 ELSEWHERE ; z1_ai(:) = 0._wp 803 837 END WHERE 838 839 ph_i(:) = SUM( phti(:,:) * pati(:,:), dim=2 ) * z1_ai(:) 840 ph_s(:) = SUM( phts(:,:) * pati(:,:), dim=2 ) * z1_ai(:) 841 ! 842 ! == temperature and salinity == ! 843 WHERE( ( pa_i(:) * ph_i(:) ) /= 0._wp ) ; z1_vi(:) = 1._wp / ( pa_i(:) * ph_i(:) ) 844 ELSEWHERE ; z1_vi(:) = 0._wp 845 END WHERE 846 WHERE( ( pa_i(:) * ph_s(:) ) /= 0._wp ) ; z1_vs(:) = 1._wp / ( pa_i(:) * ph_s(:) ) 847 ELSEWHERE ; z1_vs(:) = 0._wp 848 END WHERE 849 pt_i (:) = SUM( ptmi (:,:) * pati(:,:) * phti(:,:), dim=2 ) * z1_vi(:) 850 pt_s (:) = SUM( ptms (:,:) * pati(:,:) * phts(:,:), dim=2 ) * z1_vs(:) 851 pt_su(:) = SUM( ptmsu(:,:) * pati(:,:) , dim=2 ) * z1_ai(:) 852 ps_i (:) = SUM( psmi (:,:) * pati(:,:) * phti(:,:), dim=2 ) * z1_vi(:) 853 854 ! == ponds == ! 855 pa_ip(:) = SUM( patip(:,:), dim=2 ) 856 WHERE( pa_ip(:) /= 0._wp ) ; ph_ip(:) = SUM( phtip(:,:) * patip(:,:), dim=2 ) / pa_ip(:) 857 ELSEWHERE ; ph_ip(:) = 0._wp 858 END WHERE 859 ! 860 DEALLOCATE( z1_ai, z1_vi, z1_vs ) 804 861 ! 805 862 END SUBROUTINE ice_var_itd_Nc1c 806 863 807 SUBROUTINE ice_var_itd_1cMc( zhti, zhts, zati, zh_i, zh_s, za_i ) 864 SUBROUTINE ice_var_itd_1cMc( phti, phts, pati , ph_i, ph_s, pa_i, & 865 & ptmi, ptms, ptmsu, psmi, patip, phtip, pt_i, pt_s, pt_su, ps_i, pa_ip, ph_ip ) 808 866 !!------------------------------------------------------------------- 809 867 !! … … 826 884 !! 4) Iterate until ok (SUM(itest(:) = 4) 827 885 !! 828 !! ** Arguments : zhti: 1-cat ice thickness829 !! zhts: 1-cat snow depth830 !! zati: 1-cat ice concentration886 !! ** Arguments : phti: 1-cat ice thickness 887 !! phts: 1-cat snow depth 888 !! pati: 1-cat ice concentration 831 889 !! 832 890 !! ** Output : jpl-cat … … 834 892 !! (Example of application: BDY forcings when input are cell averaged) 835 893 !!------------------------------------------------------------------- 836 INTEGER :: ji, jk, jl ! dummy loop indices 837 INTEGER :: idim, i_fill, jl0 838 REAL(wp) :: zarg, zV, zconv, zdh, zdv 839 REAL(wp), DIMENSION(:), INTENT(in) :: zhti, zhts, zati ! input ice/snow variables 840 REAL(wp), DIMENSION(:,:), INTENT(inout) :: zh_i, zh_s, za_i ! output ice/snow variables 841 INTEGER , DIMENSION(4) :: itest 842 !!------------------------------------------------------------------- 843 ! 844 ! ---------------------------------------- 894 REAL(wp), DIMENSION(:), INTENT(in) :: phti, phts, pati ! input ice/snow variables 895 REAL(wp), DIMENSION(:,:), INTENT(inout) :: ph_i, ph_s, pa_i ! output ice/snow variables 896 REAL(wp), DIMENSION(:) , INTENT(in) :: ptmi, ptms, ptmsu, psmi, patip, phtip ! input ice/snow temp & sal & ponds 897 REAL(wp), DIMENSION(:,:), INTENT(inout) :: pt_i, pt_s, pt_su, ps_i, pa_ip, ph_ip ! output ice/snow temp & sal & ponds 898 ! 899 INTEGER , DIMENSION(4) :: itest 900 REAL(wp), ALLOCATABLE, DIMENSION(:) :: zfra 901 INTEGER :: ji, jk, jl 902 INTEGER :: idim, i_fill, jl0 903 REAL(wp) :: zarg, zV, zconv, zdh, zdv 904 !!------------------------------------------------------------------- 905 ! 906 ! == thickness and concentration == ! 845 907 ! distribution over the jpl ice categories 846 ! ----------------------------------------847 ! a gaussian distribution for ice concentration is used848 ! then we check whether the distribution fullfills849 ! volume and area conservation, positivity and ice categories bounds850 idim = SIZE( zhti , 1 )851 zh_i(1:idim,1:jpl) = 0._wp852 zh_s(1:idim,1:jpl) = 0._wp853 za_i(1:idim,1:jpl) = 0._wp908 ! a gaussian distribution for ice concentration is used 909 ! then we check whether the distribution fullfills 910 ! volume and area conservation, positivity and ice categories bounds 911 idim = SIZE( phti , 1 ) 912 ! 913 ph_i(1:idim,1:jpl) = 0._wp 914 ph_s(1:idim,1:jpl) = 0._wp 915 pa_i(1:idim,1:jpl) = 0._wp 854 916 ! 855 917 DO ji = 1, idim 856 918 ! 857 IF( zhti(ji) > 0._wp ) THEN919 IF( phti(ji) > 0._wp ) THEN 858 920 ! 859 921 ! find which category (jl0) the input ice thickness falls into 860 922 jl0 = jpl 861 923 DO jl = 1, jpl 862 IF ( ( zhti(ji) >= hi_max(jl-1) ) .AND. ( zhti(ji) < hi_max(jl) ) ) THEN924 IF ( ( phti(ji) >= hi_max(jl-1) ) .AND. ( phti(ji) < hi_max(jl) ) ) THEN 863 925 jl0 = jl 864 926 CYCLE … … 872 934 i_fill = i_fill - 1 873 935 ! 874 zh_i(ji,1:jpl) = 0._wp875 za_i(ji,1:jpl) = 0._wp936 ph_i(ji,1:jpl) = 0._wp 937 pa_i(ji,1:jpl) = 0._wp 876 938 itest(:) = 0 877 939 ! 878 940 IF ( i_fill == 1 ) THEN !-- case very thin ice: fill only category 1 879 zh_i(ji,1) = zhti(ji)880 za_i (ji,1) = zati (ji)941 ph_i(ji,1) = phti(ji) 942 pa_i(ji,1) = pati (ji) 881 943 ELSE !-- case ice is thicker: fill categories >1 882 944 ! thickness 883 945 DO jl = 1, i_fill - 1 884 zh_i(ji,jl) = hi_mean(jl)946 ph_i(ji,jl) = hi_mean(jl) 885 947 END DO 886 948 ! 887 949 ! concentration 888 za_i(ji,jl0) = zati(ji) / SQRT(REAL(jpl))950 pa_i(ji,jl0) = pati(ji) / SQRT(REAL(jpl)) 889 951 DO jl = 1, i_fill - 1 890 952 IF ( jl /= jl0 ) THEN 891 zarg = ( zh_i(ji,jl) - zhti(ji) ) / ( zhti(ji) * 0.5_wp )892 za_i(ji,jl) = za_i (ji,jl0) * EXP(-zarg**2)953 zarg = ( ph_i(ji,jl) - phti(ji) ) / ( phti(ji) * 0.5_wp ) 954 pa_i(ji,jl) = pa_i (ji,jl0) * EXP(-zarg**2) 893 955 ENDIF 894 956 END DO 895 957 ! 896 958 ! last category 897 za_i(ji,i_fill) = zati(ji) - SUM( za_i(ji,1:i_fill-1) )898 zV = SUM( za_i(ji,1:i_fill-1) * zh_i(ji,1:i_fill-1) )899 zh_i(ji,i_fill) = ( zhti(ji) * zati(ji) - zV ) / MAX( za_i(ji,i_fill), epsi10 )959 pa_i(ji,i_fill) = pati(ji) - SUM( pa_i(ji,1:i_fill-1) ) 960 zV = SUM( pa_i(ji,1:i_fill-1) * ph_i(ji,1:i_fill-1) ) 961 ph_i(ji,i_fill) = ( phti(ji) * pati(ji) - zV ) / MAX( pa_i(ji,i_fill), epsi10 ) 900 962 ! 901 963 ! correction if concentration of upper cat is greater than lower cat … … 903 965 IF ( jl0 /= jpl ) THEN 904 966 DO jl = jpl, jl0+1, -1 905 IF ( za_i(ji,jl) > za_i(ji,jl-1) ) THEN906 zdv = zh_i(ji,jl) * za_i(ji,jl)907 zh_i(ji,jl ) = 0._wp908 za_i (ji,jl ) = 0._wp909 za_i (ji,1:jl-1) = za_i(ji,1:jl-1) + zdv / MAX( REAL(jl-1) * zhti(ji), epsi10 )967 IF ( pa_i(ji,jl) > pa_i(ji,jl-1) ) THEN 968 zdv = ph_i(ji,jl) * pa_i(ji,jl) 969 ph_i(ji,jl ) = 0._wp 970 pa_i (ji,jl ) = 0._wp 971 pa_i (ji,1:jl-1) = pa_i(ji,1:jl-1) + zdv / MAX( REAL(jl-1) * phti(ji), epsi10 ) 910 972 END IF 911 973 END DO … … 915 977 ! 916 978 ! Compatibility tests 917 zconv = ABS( zati(ji) - SUM( za_i(ji,1:jpl) ) )979 zconv = ABS( pati(ji) - SUM( pa_i(ji,1:jpl) ) ) 918 980 IF ( zconv < epsi06 ) itest(1) = 1 ! Test 1: area conservation 919 981 ! 920 zconv = ABS( zhti(ji)*zati(ji) - SUM( za_i(ji,1:jpl)*zh_i(ji,1:jpl) ) )982 zconv = ABS( phti(ji)*pati(ji) - SUM( pa_i(ji,1:jpl)*ph_i(ji,1:jpl) ) ) 921 983 IF ( zconv < epsi06 ) itest(2) = 1 ! Test 2: volume conservation 922 984 ! 923 IF ( zh_i(ji,i_fill) >= hi_max(i_fill-1) ) itest(3) = 1! Test 3: thickness of the last category is in-bounds ?985 IF ( ph_i(ji,i_fill) >= hi_max(i_fill-1) ) itest(3) = 1 ! Test 3: thickness of the last category is in-bounds ? 924 986 ! 925 987 itest(4) = 1 926 988 DO jl = 1, i_fill 927 IF ( za_i(ji,jl) < 0._wp ) itest(4) = 0! Test 4: positivity of ice concentrations989 IF ( pa_i(ji,jl) < 0._wp ) itest(4) = 0 ! Test 4: positivity of ice concentrations 928 990 END DO 929 991 ! !---------------------------- … … 933 995 END DO 934 996 935 ! Add Snow in each category where za_i is not 0997 ! Add Snow in each category where pa_i is not 0 936 998 DO jl = 1, jpl 937 999 DO ji = 1, idim 938 IF( za_i(ji,jl) > 0._wp ) THEN939 zh_s(ji,jl) = zh_i(ji,jl) * ( zhts(ji) / zhti(ji) )1000 IF( pa_i(ji,jl) > 0._wp ) THEN 1001 ph_s(ji,jl) = ph_i(ji,jl) * ( phts(ji) / phti(ji) ) 940 1002 ! In case snow load is in excess that would lead to transformation from snow to ice 941 1003 ! Then, transfer the snow excess into the ice (different from icethd_dh) 942 zdh = MAX( 0._wp, ( rhos * zh_s(ji,jl) + ( rhoi - rau0 ) * zh_i(ji,jl) ) * r1_rau0 )1004 zdh = MAX( 0._wp, ( rhos * ph_s(ji,jl) + ( rhoi - rau0 ) * ph_i(ji,jl) ) * r1_rau0 ) 943 1005 ! recompute h_i, h_s avoiding out of bounds values 944 zh_i(ji,jl) = MIN( hi_max(jl), zh_i(ji,jl) + zdh )945 zh_s(ji,jl) = MAX( 0._wp, zh_s(ji,jl) - zdh * rhoi * r1_rhos )1006 ph_i(ji,jl) = MIN( hi_max(jl), ph_i(ji,jl) + zdh ) 1007 ph_s(ji,jl) = MAX( 0._wp, ph_s(ji,jl) - zdh * rhoi * r1_rhos ) 946 1008 ENDIF 947 1009 END DO 948 1010 END DO 949 1011 ! 1012 ! == temperature and salinity == ! 1013 DO jl = 1, jpl 1014 pt_i (:,jl) = ptmi (:) 1015 pt_s (:,jl) = ptms (:) 1016 pt_su(:,jl) = ptmsu(:) 1017 ps_i (:,jl) = psmi (:) 1018 ps_i (:,jl) = psmi (:) 1019 END DO 1020 ! 1021 ! == ponds == ! 1022 ALLOCATE( zfra(idim) ) 1023 ! keep the same pond fraction atip/ati for each category 1024 WHERE( pati(:) /= 0._wp ) ; zfra(:) = patip(:) / pati(:) 1025 ELSEWHERE ; zfra(:) = 0._wp 1026 END WHERE 1027 DO jl = 1, jpl 1028 pa_ip(:,jl) = zfra(:) * pa_i(:,jl) 1029 END DO 1030 ! keep the same v_ip/v_i ratio for each category 1031 WHERE( ( phti(:) * pati(:) ) /= 0._wp ) ; zfra(:) = ( phtip(:) * patip(:) ) / ( phti(:) * pati(:) ) 1032 ELSEWHERE ; zfra(:) = 0._wp 1033 END WHERE 1034 DO jl = 1, jpl 1035 WHERE( pa_ip(:,jl) /= 0._wp ) ; ph_ip(:,jl) = zfra(:) * ( ph_i(:,jl) * pa_i(:,jl) ) / pa_ip(:,jl) 1036 ELSEWHERE ; ph_ip(:,jl) = 0._wp 1037 END WHERE 1038 END DO 1039 DEALLOCATE( zfra ) 1040 ! 950 1041 END SUBROUTINE ice_var_itd_1cMc 951 1042 952 SUBROUTINE ice_var_itd_NcMc( zhti, zhts, zati, zh_i, zh_s, za_i ) 1043 SUBROUTINE ice_var_itd_NcMc( phti, phts, pati , ph_i, ph_s, pa_i, & 1044 & ptmi, ptms, ptmsu, psmi, patip, phtip, pt_i, pt_s, pt_su, ps_i, pa_ip, ph_ip ) 953 1045 !!------------------------------------------------------------------- 954 1046 !! … … 971 1063 !! b) removing 25% ice area from the higher cat (descendant loop jlmax=>jlmin) 972 1064 !! 973 !! ** Arguments : zhti: N-cat ice thickness974 !! zhts: N-cat snow depth975 !! zati: N-cat ice concentration1065 !! ** Arguments : phti: N-cat ice thickness 1066 !! phts: N-cat snow depth 1067 !! pati: N-cat ice concentration 976 1068 !! 977 1069 !! ** Output : jpl-cat … … 979 1071 !! (Example of application: BDY forcings when inputs have N-cat /= jpl) 980 1072 !!------------------------------------------------------------------- 981 INTEGER :: ji, jl, jl1, jl2 ! dummy loop indices 1073 REAL(wp), DIMENSION(:,:), INTENT(in) :: phti, phts, pati ! input ice/snow variables 1074 REAL(wp), DIMENSION(:,:), INTENT(inout) :: ph_i, ph_s, pa_i ! output ice/snow variables 1075 REAL(wp), DIMENSION(:,:), INTENT(in) :: ptmi, ptms, ptmsu, psmi, patip, phtip ! input ice/snow temp & sal & ponds 1076 REAL(wp), DIMENSION(:,:), INTENT(inout) :: pt_i, pt_s, pt_su, ps_i, pa_ip, ph_ip ! output ice/snow temp & sal & ponds 1077 ! 1078 INTEGER , ALLOCATABLE, DIMENSION(:,:) :: jlfil, jlfil2 1079 INTEGER , ALLOCATABLE, DIMENSION(:) :: jlmax, jlmin 1080 REAL(wp), ALLOCATABLE, DIMENSION(:) :: z1_ai, z1_vi, z1_vs, ztmp, zfra 1081 ! 1082 REAL(wp), PARAMETER :: ztrans = 0.25_wp 1083 INTEGER :: ji, jl, jl1, jl2 982 1084 INTEGER :: idim, icat 983 REAL(wp), PARAMETER :: ztrans = 0.25_wp 984 REAL(wp), DIMENSION(:,:), INTENT(in) :: zhti, zhts, zati ! input ice/snow variables 985 REAL(wp), DIMENSION(:,:), INTENT(inout) :: zh_i, zh_s, za_i ! output ice/snow variables 986 INTEGER , DIMENSION(:,:), ALLOCATABLE :: jlfil, jlfil2 987 INTEGER , DIMENSION(:) , ALLOCATABLE :: jlmax, jlmin 988 !!------------------------------------------------------------------- 989 ! 990 idim = SIZE( zhti, 1 ) 991 icat = SIZE( zhti, 2 ) 1085 !!------------------------------------------------------------------- 1086 ! 1087 idim = SIZE( phti, 1 ) 1088 icat = SIZE( phti, 2 ) 1089 ! 1090 ! == thickness and concentration == ! 992 1091 ! ! ---------------------- ! 993 1092 IF( icat == jpl ) THEN ! input cat = output cat ! 994 1093 ! ! ---------------------- ! 995 zh_i(:,:) = zhti(:,:) 996 zh_s(:,:) = zhts(:,:) 997 za_i(:,:) = zati(:,:) 1094 ph_i(:,:) = phti(:,:) 1095 ph_s(:,:) = phts(:,:) 1096 pa_i(:,:) = pati(:,:) 1097 ! 1098 ! == temperature and salinity and ponds == ! 1099 pt_i (:,:) = ptmi (:,:) 1100 pt_s (:,:) = ptms (:,:) 1101 pt_su(:,:) = ptmsu(:,:) 1102 ps_i (:,:) = psmi (:,:) 1103 pa_ip(:,:) = patip(:,:) 1104 ph_ip(:,:) = phtip(:,:) 998 1105 ! ! ---------------------- ! 999 1106 ELSEIF( icat == 1 ) THEN ! input cat = 1 ! 1000 1107 ! ! ---------------------- ! 1001 CALL ice_var_itd_1cMc( zhti(:,1), zhts(:,1), zati(:,1), zh_i(:,:), zh_s(:,:), za_i(:,:) ) 1108 CALL ice_var_itd_1cMc( phti(:,1), phts(:,1), pati (:,1), & 1109 & ph_i(:,:), ph_s(:,:), pa_i (:,:), & 1110 & ptmi(:,1), ptms(:,1), ptmsu(:,1), psmi(:,1), patip(:,1), phtip(:,1), & 1111 & pt_i(:,:), pt_s(:,:), pt_su(:,:), ps_i(:,:), pa_ip(:,:), ph_ip(:,:) ) 1002 1112 ! ! ---------------------- ! 1003 1113 ELSEIF( jpl == 1 ) THEN ! output cat = 1 ! 1004 1114 ! ! ---------------------- ! 1005 CALL ice_var_itd_Nc1c( zhti(:,:), zhts(:,:), zati(:,:), zh_i(:,1), zh_s(:,1), za_i(:,1) ) 1115 CALL ice_var_itd_Nc1c( phti(:,:), phts(:,:), pati (:,:), & 1116 & ph_i(:,1), ph_s(:,1), pa_i (:,1), & 1117 & ptmi(:,:), ptms(:,:), ptmsu(:,:), psmi(:,:), patip(:,:), phtip(:,:), & 1118 & pt_i(:,1), pt_s(:,1), pt_su(:,1), ps_i(:,1), pa_ip(:,1), ph_ip(:,1) ) 1006 1119 ! ! ----------------------- ! 1007 1120 ELSE ! input cat /= output cat ! … … 1012 1125 1013 1126 ! --- initialize output fields to 0 --- ! 1014 zh_i(1:idim,1:jpl) = 0._wp1015 zh_s(1:idim,1:jpl) = 0._wp1016 za_i(1:idim,1:jpl) = 0._wp1127 ph_i(1:idim,1:jpl) = 0._wp 1128 ph_s(1:idim,1:jpl) = 0._wp 1129 pa_i(1:idim,1:jpl) = 0._wp 1017 1130 ! 1018 1131 ! --- fill the categories --- ! … … 1024 1137 DO jl2 = 1, icat 1025 1138 DO ji = 1, idim 1026 IF( hi_max(jl1-1) <= zhti(ji,jl2) .AND. hi_max(jl1) > zhti(ji,jl2) ) THEN1139 IF( hi_max(jl1-1) <= phti(ji,jl2) .AND. hi_max(jl1) > phti(ji,jl2) ) THEN 1027 1140 ! fill the right category 1028 zh_i(ji,jl1) = zhti(ji,jl2)1029 zh_s(ji,jl1) = zhts(ji,jl2)1030 za_i(ji,jl1) = zati(ji,jl2)1141 ph_i(ji,jl1) = phti(ji,jl2) 1142 ph_s(ji,jl1) = phts(ji,jl2) 1143 pa_i(ji,jl1) = pati(ji,jl2) 1031 1144 ! record categories that are filled 1032 1145 jlmax(ji) = MAX( jlmax(ji), jl1 ) … … 1045 1158 IF( jl1 > 1 ) THEN 1046 1159 ! fill the lower cat (jl1-1) 1047 za_i(ji,jl1-1) = ztrans * za_i(ji,jl1)1048 zh_i(ji,jl1-1) = hi_mean(jl1-1)1160 pa_i(ji,jl1-1) = ztrans * pa_i(ji,jl1) 1161 ph_i(ji,jl1-1) = hi_mean(jl1-1) 1049 1162 ! remove from cat jl1 1050 za_i(ji,jl1 ) = ( 1._wp - ztrans ) * za_i(ji,jl1)1163 pa_i(ji,jl1 ) = ( 1._wp - ztrans ) * pa_i(ji,jl1) 1051 1164 ENDIF 1052 1165 IF( jl2 < jpl ) THEN 1053 1166 ! fill the upper cat (jl2+1) 1054 za_i(ji,jl2+1) = ztrans * za_i(ji,jl2)1055 zh_i(ji,jl2+1) = hi_mean(jl2+1)1167 pa_i(ji,jl2+1) = ztrans * pa_i(ji,jl2) 1168 ph_i(ji,jl2+1) = hi_mean(jl2+1) 1056 1169 ! remove from cat jl2 1057 za_i(ji,jl2 ) = ( 1._wp - ztrans ) * za_i(ji,jl2)1170 pa_i(ji,jl2 ) = ( 1._wp - ztrans ) * pa_i(ji,jl2) 1058 1171 ENDIF 1059 1172 END DO … … 1065 1178 IF( jlfil(ji,jl-1) /= 0 .AND. jlfil(ji,jl) == 0 ) THEN 1066 1179 ! fill high 1067 za_i(ji,jl) = ztrans * za_i(ji,jl-1)1068 zh_i(ji,jl) = hi_mean(jl)1180 pa_i(ji,jl) = ztrans * pa_i(ji,jl-1) 1181 ph_i(ji,jl) = hi_mean(jl) 1069 1182 jlfil(ji,jl) = jl 1070 1183 ! remove low 1071 za_i(ji,jl-1) = ( 1._wp - ztrans ) * za_i(ji,jl-1)1184 pa_i(ji,jl-1) = ( 1._wp - ztrans ) * pa_i(ji,jl-1) 1072 1185 ENDIF 1073 1186 END DO … … 1079 1192 IF( jlfil2(ji,jl+1) /= 0 .AND. jlfil2(ji,jl) == 0 ) THEN 1080 1193 ! fill low 1081 za_i(ji,jl) = za_i(ji,jl) + ztrans * za_i(ji,jl+1)1082 zh_i(ji,jl) = hi_mean(jl)1194 pa_i(ji,jl) = pa_i(ji,jl) + ztrans * pa_i(ji,jl+1) 1195 ph_i(ji,jl) = hi_mean(jl) 1083 1196 jlfil2(ji,jl) = jl 1084 1197 ! remove high 1085 za_i(ji,jl+1) = ( 1._wp - ztrans ) * za_i(ji,jl+1)1198 pa_i(ji,jl+1) = ( 1._wp - ztrans ) * pa_i(ji,jl+1) 1086 1199 ENDIF 1087 1200 END DO … … 1090 1203 DEALLOCATE( jlfil, jlfil2 ) ! deallocate arrays 1091 1204 DEALLOCATE( jlmin, jlmax ) 1205 ! 1206 ! == temperature and salinity == ! 1207 ! 1208 ALLOCATE( z1_ai(idim), z1_vi(idim), z1_vs(idim), ztmp(idim) ) 1209 ! 1210 WHERE( SUM( pa_i(:,:), dim=2 ) /= 0._wp ) ; z1_ai(:) = 1._wp / SUM( pa_i(:,:), dim=2 ) 1211 ELSEWHERE ; z1_ai(:) = 0._wp 1212 END WHERE 1213 WHERE( SUM( pa_i(:,:) * ph_i(:,:), dim=2 ) /= 0._wp ) ; z1_vi(:) = 1._wp / SUM( pa_i(:,:) * ph_i(:,:), dim=2 ) 1214 ELSEWHERE ; z1_vi(:) = 0._wp 1215 END WHERE 1216 WHERE( SUM( pa_i(:,:) * ph_s(:,:), dim=2 ) /= 0._wp ) ; z1_vs(:) = 1._wp / SUM( pa_i(:,:) * ph_s(:,:), dim=2 ) 1217 ELSEWHERE ; z1_vs(:) = 0._wp 1218 END WHERE 1219 ! 1220 ! fill all the categories with the same value 1221 ztmp(:) = SUM( ptmi (:,:) * pati(:,:) * phti(:,:), dim=2 ) * z1_vi(:) 1222 DO jl = 1, jpl 1223 pt_i (:,jl) = ztmp(:) 1224 END DO 1225 ztmp(:) = SUM( ptms (:,:) * pati(:,:) * phts(:,:), dim=2 ) * z1_vs(:) 1226 DO jl = 1, jpl 1227 pt_s (:,jl) = ztmp(:) 1228 END DO 1229 ztmp(:) = SUM( ptmsu(:,:) * pati(:,:) , dim=2 ) * z1_ai(:) 1230 DO jl = 1, jpl 1231 pt_su(:,jl) = ztmp(:) 1232 END DO 1233 ztmp(:) = SUM( psmi (:,:) * pati(:,:) * phti(:,:), dim=2 ) * z1_vi(:) 1234 DO jl = 1, jpl 1235 ps_i (:,jl) = ztmp(:) 1236 END DO 1237 ! 1238 DEALLOCATE( z1_ai, z1_vi, z1_vs, ztmp ) 1239 ! 1240 ! == ponds == ! 1241 ALLOCATE( zfra(idim) ) 1242 ! keep the same pond fraction atip/ati for each category 1243 WHERE( SUM( pati(:,:), dim=2 ) /= 0._wp ) ; zfra(:) = SUM( patip(:,:), dim=2 ) / SUM( pati(:,:), dim=2 ) 1244 ELSEWHERE ; zfra(:) = 0._wp 1245 END WHERE 1246 DO jl = 1, jpl 1247 pa_ip(:,jl) = zfra(:) * pa_i(:,jl) 1248 END DO 1249 ! keep the same v_ip/v_i ratio for each category 1250 WHERE( SUM( phti(:,:) * pati(:,:), dim=2 ) /= 0._wp ) 1251 zfra(:) = SUM( phtip(:,:) * patip(:,:), dim=2 ) / SUM( phti(:,:) * pati(:,:), dim=2 ) 1252 ELSEWHERE 1253 zfra(:) = 0._wp 1254 END WHERE 1255 DO jl = 1, jpl 1256 WHERE( pa_ip(:,jl) /= 0._wp ) ; ph_ip(:,jl) = zfra(:) * ( ph_i(:,jl) * pa_i(:,jl) ) / pa_ip(:,jl) 1257 ELSEWHERE ; ph_ip(:,jl) = 0._wp 1258 END WHERE 1259 END DO 1260 DEALLOCATE( zfra ) 1092 1261 ! 1093 1262 ENDIF
Note: See TracChangeset
for help on using the changeset viewer.