Changeset 2037 for palm/trunk/SOURCE/flow_statistics.f90
- Timestamp:
- Oct 26, 2016 11:15:40 AM (8 years ago)
- File:
-
- 1 edited
Legend:
- Unmodified
- Added
- Removed
-
palm/trunk/SOURCE/flow_statistics.f90
r2032 r2037 20 20 ! Current revisions: 21 21 ! ----------------- 22 ! 22 ! Anelastic approximation implemented 23 23 ! 24 24 ! Former revisions: … … 214 214 215 215 USE arrays_3d, & 216 ONLY: ddzu, ddzw, e, hyp, km, kh, nr, ol, p, prho, prr, pt, q, qc, ql,& 217 qr, qs, qsws, qswst, rho_ocean, s, sa, ss, ssws, sswst, saswsb, & 218 saswst, shf, td_lsa_lpt, td_lsa_q, td_sub_lpt, td_sub_q, & 219 time_vert, ts, tswst, u, ug, us, usws, uswst, vsws, v, vg, vpt, & 220 vswst, w, w_subs, zw 216 ONLY: ddzu, ddzw, e, heatflux_output_conversion, hyp, km, kh, & 217 momentumflux_output_conversion, nr, ol, p, prho, prr, pt, q, & 218 qc, ql, qr, qs, qsws, qswst, rho_air, rho_air_zw, rho_ocean, s, & 219 sa, ss, ssws, sswst, saswsb, saswst, shf, td_lsa_lpt, td_lsa_q, & 220 td_sub_lpt, td_sub_q, time_vert, ts, tswst, u, ug, us, usws, & 221 uswst, vsws, v, vg, vpt, vswst, w, w_subs, & 222 waterflux_output_conversion, zw 221 223 222 224 USE cloud_parameters, & … … 326 328 sums_l(:,11,:) = sums_l_l(:,sr,:) ! mixing length from diffusivities 327 329 !-- WARNING: next line still has to be adjusted for OpenMP 328 sums_l(:,21,0) = sums_wsts_bc_l(:,sr) ! heat flux from advec_s_bc 330 sums_l(:,21,0) = sums_wsts_bc_l(:,sr) * & 331 heatflux_output_conversion ! heat flux from advec_s_bc 329 332 sums_l(nzb+9,pr_palm,0) = sums_divold_l(sr) ! old divergence from pres 330 333 sums_l(nzb+10,pr_palm,0) = sums_divnew_l(sr) ! new divergence from pres … … 358 361 ! 359 362 !-- Swap the turbulent quantities evaluated in advec_ws. 360 sums_l(:,13,i) = sums_wsus_ws_l(:,i) ! w*u* 361 sums_l(:,15,i) = sums_wsvs_ws_l(:,i) ! w*v* 363 sums_l(:,13,i) = sums_wsus_ws_l(:,i) & 364 * momentumflux_output_conversion ! w*u* 365 sums_l(:,15,i) = sums_wsvs_ws_l(:,i) & 366 * momentumflux_output_conversion ! w*v* 362 367 sums_l(:,30,i) = sums_us2_ws_l(:,i) ! u*2 363 368 sums_l(:,31,i) = sums_vs2_ws_l(:,i) ! v*2 … … 373 378 374 379 DO i = 0, threads_per_task-1 375 sums_l(:,17,i) = sums_wspts_ws_l(:,i) ! w*pt* 380 sums_l(:,17,i) = sums_wspts_ws_l(:,i) & 381 * heatflux_output_conversion ! w*pt* 376 382 IF ( ocean ) sums_l(:,66,i) = sums_wssas_ws_l(:,i) ! w*sa* 377 IF ( humidity ) sums_l(:,49,i) = sums_wsqs_ws_l(:,i) ! w*q* 383 IF ( humidity ) sums_l(:,49,i) = sums_wsqs_ws_l(:,i) & 384 * waterflux_output_conversion ! w*q* 378 385 IF ( passive_scalar ) sums_l(:,116,i) = sums_wsss_ws_l(:,i) ! w*s* 379 386 ENDDO … … 722 729 ( u(k+1,j,i) - u(k,j,i) ) * ddzu(k+1) & 723 730 + ( w(k,j,i) - w(k,j,i-1) ) * ddx & 724 ) * rmask(j,i,sr) 731 ) * rmask(j,i,sr) & 732 * rho_air_zw(k) & 733 * momentumflux_output_conversion(k) 725 734 ! 726 735 !-- Momentum flux w"v" … … 730 739 ( v(k+1,j,i) - v(k,j,i) ) * ddzu(k+1) & 731 740 + ( w(k,j,i) - w(k,j-1,i) ) * ddy & 732 ) * rmask(j,i,sr) 741 ) * rmask(j,i,sr) & 742 * rho_air_zw(k) & 743 * momentumflux_output_conversion(k) 733 744 ! 734 745 !-- Heat flux w"pt" … … 736 747 - 0.5_wp * ( kh(k,j,i) + kh(k+1,j,i) )& 737 748 * ( pt(k+1,j,i) - pt(k,j,i) ) & 749 * rho_air_zw(k) & 750 * heatflux_output_conversion(k) & 738 751 * ddzu(k+1) * rmask(j,i,sr) 739 752 … … 754 767 - 0.5_wp * ( kh(k,j,i) + kh(k+1,j,i) )& 755 768 * ( vpt(k+1,j,i) - vpt(k,j,i) ) & 769 * rho_air_zw(k) & 770 * heatflux_output_conversion(k) & 756 771 * ddzu(k+1) * rmask(j,i,sr) 757 772 sums_l(k,48,tn) = sums_l(k,48,tn) & 758 773 - 0.5_wp * ( kh(k,j,i) + kh(k+1,j,i) )& 759 774 * ( q(k+1,j,i) - q(k,j,i) ) & 775 * rho_air_zw(k) & 776 * waterflux_output_conversion(k)& 760 777 * ddzu(k+1) * rmask(j,i,sr) 761 778 … … 765 782 * ( ( q(k+1,j,i) - ql(k+1,j,i) )& 766 783 - ( q(k,j,i) - ql(k,j,i) ) ) & 784 * rho_air_zw(k) & 785 * waterflux_output_conversion(k)& 767 786 * ddzu(k+1) * rmask(j,i,sr) 768 787 ENDIF … … 784 803 IF ( use_surface_fluxes ) THEN 785 804 sums_l(nzb,12,tn) = sums_l(nzb,12,tn) + & 805 momentumflux_output_conversion(nzb) * & 786 806 usws(j,i) * rmask(j,i,sr) ! w"u" 787 807 sums_l(nzb,14,tn) = sums_l(nzb,14,tn) + & 808 momentumflux_output_conversion(nzb) * & 788 809 vsws(j,i) * rmask(j,i,sr) ! w"v" 789 810 sums_l(nzb,16,tn) = sums_l(nzb,16,tn) + & 811 heatflux_output_conversion(nzb) * & 790 812 shf(j,i) * rmask(j,i,sr) ! w"pt" 791 813 sums_l(nzb,58,tn) = sums_l(nzb,58,tn) + & … … 799 821 IF ( humidity ) THEN 800 822 sums_l(nzb,48,tn) = sums_l(nzb,48,tn) + & 823 waterflux_output_conversion(nzb) * & 801 824 qsws(j,i) * rmask(j,i,sr) ! w"q" (w"qv") 802 825 sums_l(nzb,45,tn) = sums_l(nzb,45,tn) + ( & 803 826 ( 1.0_wp + 0.61_wp * q(nzb,j,i) ) * & 804 827 shf(j,i) + 0.61_wp * pt(nzb,j,i) * & 805 qsws(j,i) ) 828 qsws(j,i) ) & 829 * heatflux_output_conversion(nzb) 806 830 IF ( cloud_droplets ) THEN 807 831 sums_l(nzb,45,tn) = sums_l(nzb,45,tn) + ( & 808 832 ( 1.0_wp + 0.61_wp * q(nzb,j,i) - & 809 833 ql(nzb,j,i) ) * shf(j,i) + & 810 0.61_wp * pt(nzb,j,i) * qsws(j,i) ) 834 0.61_wp * pt(nzb,j,i) * qsws(j,i) ) & 835 * heatflux_output_conversion(nzb) 811 836 ENDIF 812 837 IF ( cloud_physics ) THEN 813 838 ! 814 839 !-- Formula does not work if ql(nzb) /= 0.0 815 sums_l(nzb,51,tn) = sums_l(nzb,51,tn) + & 840 sums_l(nzb,51,tn) = sums_l(nzb,51,tn) + & 841 waterflux_output_conversion(nzb) * & 816 842 qsws(j,i) * rmask(j,i,sr) ! w"q" (w"qv") 817 843 ENDIF … … 860 886 IF ( use_top_fluxes ) THEN 861 887 sums_l(nzt:nzt+1,12,tn) = sums_l(nzt:nzt+1,12,tn) + & 888 momentumflux_output_conversion(nzt:nzt+1) * & 862 889 uswst(j,i) * rmask(j,i,sr) ! w"u" 863 890 sums_l(nzt:nzt+1,14,tn) = sums_l(nzt:nzt+1,14,tn) + & 891 momentumflux_output_conversion(nzt:nzt+1) * & 864 892 vswst(j,i) * rmask(j,i,sr) ! w"v" 865 893 sums_l(nzt:nzt+1,16,tn) = sums_l(nzt:nzt+1,16,tn) + & 894 heatflux_output_conversion(nzt:nzt+1) * & 866 895 tswst(j,i) * rmask(j,i,sr) ! w"pt" 867 896 sums_l(nzt:nzt+1,58,tn) = sums_l(nzt:nzt+1,58,tn) + & … … 876 905 IF ( humidity ) THEN 877 906 sums_l(nzt,48,tn) = sums_l(nzt,48,tn) + & 907 waterflux_output_conversion(nzt) * & 878 908 qswst(j,i) * rmask(j,i,sr) ! w"q" (w"qv") 879 909 sums_l(nzt,45,tn) = sums_l(nzt,45,tn) + ( & 880 910 ( 1.0_wp + 0.61_wp * q(nzt,j,i) ) * & 881 911 tswst(j,i) + 0.61_wp * pt(nzt,j,i) * & 882 qswst(j,i) ) 912 qswst(j,i) ) & 913 * heatflux_output_conversion(nzt) 883 914 IF ( cloud_droplets ) THEN 884 915 sums_l(nzt,45,tn) = sums_l(nzt,45,tn) + ( & 885 916 ( 1.0_wp + 0.61_wp * q(nzt,j,i) - & 886 917 ql(nzt,j,i) ) * tswst(j,i) + & 887 0.61_wp * pt(nzt,j,i) * qswst(j,i) ) 918 0.61_wp * pt(nzt,j,i) * qswst(j,i) )& 919 * heatflux_output_conversion(nzt) 888 920 ENDIF 889 921 IF ( cloud_physics ) THEN … … 891 923 !-- Formula does not work if ql(nzb) /= 0.0 892 924 sums_l(nzt,51,tn) = sums_l(nzt,51,tn) + & ! w"q" (w"qv") 925 waterflux_output_conversion(nzt) * & 893 926 qswst(j,i) * rmask(j,i,sr) 894 927 ENDIF … … 943 976 vpt(k+1,j,i) - hom(k+1,1,44,sr) ) 944 977 sums_l(k,46,tn) = sums_l(k,46,tn) + pts * w(k,j,i) * & 978 heatflux_output_conversion(k) * & 945 979 rmask(j,i,sr) 946 980 sums_l(k,54,tn) = sums_l(k,54,tn) + ql(k,j,i) * rmask(j,i,sr) … … 953 987 hom(k+1,1,42,sr) ) 954 988 sums_l(k,52,tn) = sums_l(k,52,tn) + pts * w(k,j,i) * & 989 waterflux_output_conversion(k) * & 955 990 rmask(j,i,sr) 956 991 sums_l(k,75,tn) = sums_l(k,75,tn) + qc(k,j,i) * & … … 971 1006 vpt(k+1,j,i) - hom(k+1,1,44,sr) ) 972 1007 sums_l(k,46,tn) = sums_l(k,46,tn) + pts * w(k,j,i) * & 1008 heatflux_output_conversion(k) * & 973 1009 rmask(j,i,sr) 974 1010 ELSE IF ( ws_scheme_sca .AND. sr == 0 ) THEN 975 sums_l(k,46,tn) = ( 1.0_wp + 0.61_wp * & 976 hom(k,1,41,sr) ) * & 977 sums_l(k,17,tn) + & 978 0.61_wp * hom(k,1,4,sr) * & 979 sums_l(k,49,tn) 1011 sums_l(k,46,tn) = ( ( 1.0_wp + 0.61_wp * & 1012 hom(k,1,41,sr) ) * & 1013 sums_l(k,17,tn) + & 1014 0.61_wp * hom(k,1,4,sr) * & 1015 sums_l(k,49,tn) & 1016 ) * heatflux_output_conversion(k) 980 1017 END IF 981 1018 END IF … … 996 1033 sums_l(k,37,tn) = sums_l(k,37,tn) + w(k,j,i) * 0.5_wp * & 997 1034 ( ust**2 + vst**2 + w(k,j,i)**2 ) & 1035 * momentumflux_output_conversion(k) & 998 1036 * rmask(j,i,sr) 999 1037 ENDDO … … 1017 1055 sums_l(k,13,tn) = sums_l(k,13,tn) + 0.5_wp * & 1018 1056 ( w(k,j,i-1) + w(k,j,i) ) & 1057 * momentumflux_output_conversion(k) & 1019 1058 * ust * rmask(j,i,sr) 1020 1059 ! … … 1022 1061 sums_l(k,15,tn) = sums_l(k,15,tn) + 0.5_wp * & 1023 1062 ( w(k,j-1,i) + w(k,j,i) ) & 1063 * momentumflux_output_conversion(k) & 1024 1064 * vst * rmask(j,i,sr) 1025 1065 ENDDO … … 1038 1078 ( pt(k,j,i) - hom(k,1,4,sr) + & 1039 1079 pt(k+1,j,i) - hom(k+1,1,4,sr) ) & 1080 * heatflux_output_conversion(k) & 1040 1081 * w(k,j,i) * rmask(j,i,sr) 1041 1082 IF ( humidity ) THEN … … 1043 1084 q(k+1,j,i) - hom(k+1,1,41,sr) ) 1044 1085 sums_l(k,49,tn) = sums_l(k,49,tn) + pts * w(k,j,i) * & 1086 waterflux_output_conversion(k) * & 1045 1087 rmask(j,i,sr) 1046 1088 ENDIF … … 1147 1189 ( kh(k,j,i) + kh(k,j,i-1) ) & 1148 1190 * ( pt(k,j,i-1) - pt(k,j,i) ) & 1191 * rho_air_zw(k) & 1192 * heatflux_output_conversion(k) & 1149 1193 * ddx * rmask(j,i,sr) 1150 1194 sums_l(k,61,tn) = sums_l(k,61,tn) - 0.5_wp * & 1151 1195 ( kh(k,j,i) + kh(k,j-1,i) ) & 1152 1196 * ( pt(k,j-1,i) - pt(k,j,i) ) & 1197 * rho_air_zw(k) & 1198 * heatflux_output_conversion(k) & 1153 1199 * ddy * rmask(j,i,sr) 1154 1200 ! … … 1157 1203 ( u(k,j,i) - hom(k,1,1,sr) ) & 1158 1204 * 0.5_wp * ( pt(k,j,i-1) - hom(k,1,4,sr) + & 1159 pt(k,j,i) - hom(k,1,4,sr) ) 1205 pt(k,j,i) - hom(k,1,4,sr) ) & 1206 * heatflux_output_conversion(k) 1160 1207 pts = 0.5_wp * ( pt(k,j-1,i) - hom(k,1,4,sr) + & 1161 1208 pt(k,j,i) - hom(k,1,4,sr) ) … … 1163 1210 ( v(k,j,i) - hom(k,1,2,sr) ) & 1164 1211 * 0.5_wp * ( pt(k,j-1,i) - hom(k,1,4,sr) + & 1165 pt(k,j,i) - hom(k,1,4,sr) ) 1212 pt(k,j,i) - hom(k,1,4,sr) ) & 1213 * heatflux_output_conversion(k) 1166 1214 ENDDO 1167 1215 ENDDO … … 1489 1537 ENDIF 1490 1538 1539 hom(:,1,121,sr) = rho_air ! rho_air in Kg/m^3 1540 hom(:,1,122,sr) = rho_air_zw ! rho_air_zw in Kg/m^3 1541 1491 1542 hom(:,1,pr_palm,sr) = sums(:,pr_palm) 1492 1543 ! u*, w'u', w'v', t* (in last profile) … … 1597 1648 THEN 1598 1649 hom(nzb+8,1,pr_palm,sr) = & 1599 ( g / hom(k_surface_level+1,1,4,sr) * hom(k_surface_level,1,18,sr)& 1650 ( g / hom(k_surface_level+1,1,4,sr) * & 1651 ( hom(k_surface_level,1,18,sr) / heatflux_output_conversion(nzb) )& 1600 1652 * ABS( z_i(1) - mean_surface_level_height(sr) ) )**0.333333333_wp 1601 1653 ELSE … … 1708 1760 1709 1761 USE arrays_3d, & 1710 ONLY: ddzu, ddzw, e, hyp, km, kh, nr, p, prho, pt, q, qc, ql, qr, qs, & 1711 qsws, qswst, rho_ocean, s, sa, saswsb, saswst, shf, ss, ssws, sswst, & 1712 td_lsa_lpt, td_lsa_q, td_sub_lpt, td_sub_q, time_vert, ts, & 1713 tswst, u, ug, us, usws, uswst, vsws, v, vg, vpt, vswst, w, & 1714 w_subs, zw 1762 ONLY: ddzu, ddzw, e, heatflux_output_conversion, hyp, km, kh, & 1763 momentumflux_output_conversion, nr, p, prho, pt, q, qc, ql, qr, & 1764 qs, qsws, qswst, rho_air, rho_air_zw, rho_ocean, s, sa, saswsb, & 1765 saswst, shf, ss, ssws, sswst, td_lsa_lpt, td_lsa_q, td_sub_lpt, & 1766 td_sub_q, time_vert, ts, tswst, u, ug, us, usws, uswst, vsws, & 1767 v, vg, vpt, vswst, w, w_subs, waterflux_output_conversion, zw 1715 1768 1716 1769 … … 1828 1881 sums_l(:,11,:) = sums_l_l(:,sr,:) ! mixing length from diffusivities 1829 1882 !-- WARNING: next line still has to be adjusted for OpenMP 1830 sums_l(:,21,0) = sums_wsts_bc_l(:,sr) ! heat flux from advec_s_bc 1883 sums_l(:,21,0) = sums_wsts_bc_l(:,sr) * & 1884 heatflux_output_conversion ! heat flux from advec_s_bc 1831 1885 sums_l(nzb+9,pr_palm,0) = sums_divold_l(sr) ! old divergence from pres 1832 1886 sums_l(nzb+10,pr_palm,0) = sums_divnew_l(sr) ! new divergence from pres … … 1860 1914 ! 1861 1915 !-- Swap the turbulent quantities evaluated in advec_ws. 1862 sums_l(:,13,i) = sums_wsus_ws_l(:,i) ! w*u* 1863 sums_l(:,15,i) = sums_wsvs_ws_l(:,i) ! w*v* 1916 sums_l(:,13,i) = sums_wsus_ws_l(:,i) & 1917 * momentumflux_output_conversion ! w*u* 1918 sums_l(:,15,i) = sums_wsvs_ws_l(:,i) & 1919 * momentumflux_output_conversion ! w*v* 1864 1920 sums_l(:,30,i) = sums_us2_ws_l(:,i) ! u*2 1865 1921 sums_l(:,31,i) = sums_vs2_ws_l(:,i) ! v*2 … … 1881 1937 1882 1938 DO i = 0, threads_per_task-1 1883 sums_l(:,17,i) = sums_wspts_ws_l(:,i) ! w*pt* from advec_s_ws 1939 sums_l(:,17,i) = sums_wspts_ws_l(:,i) & 1940 * heatflux_output_conversion ! w*pt* from advec_s_ws 1884 1941 IF ( ocean ) sums_l(:,66,i) = sums_wssas_ws_l(:,i) ! w*sa* 1885 IF ( humidity ) sums_l(:,49,i) = sums_wsqs_ws_l(:,i) !w*q* 1942 IF ( humidity ) sums_l(:,49,i) = sums_wsqs_ws_l(:,i) & 1943 * waterflux_output_conversion !w*q* 1886 1944 IF ( passive_scalar ) sums_l(:,116,i) = sums_wsss_ws_l(:,i) !w*s* 1887 1945 ENDDO … … 2378 2436 + ( w(k,j,i) - w(k,j,i-1) ) * ddx & 2379 2437 ) & 2380 * rmask(j,i,sr) * rflags_invers(j,i,k+1) 2438 * rmask(j,i,sr) * rflags_invers(j,i,k+1) & 2439 * rho_air_zw(k) & 2440 * momentumflux_output_conversion(k) 2381 2441 ! 2382 2442 !-- Momentum flux w"v" … … 2387 2447 + ( w(k,j,i) - w(k,j-1,i) ) * ddy & 2388 2448 ) & 2389 * rmask(j,i,sr) * rflags_invers(j,i,k+1) 2449 * rmask(j,i,sr) * rflags_invers(j,i,k+1) & 2450 * rho_air_zw(k) & 2451 * momentumflux_output_conversion(k) 2390 2452 ! 2391 2453 !-- Heat flux w"pt" 2392 2454 s3 = s3 - 0.5_wp * ( kh(k,j,i) + kh(k+1,j,i) ) & 2393 2455 * ( pt(k+1,j,i) - pt(k,j,i) ) & 2456 * rho_air_zw(k) & 2457 * heatflux_output_conversion(k) & 2394 2458 * ddzu(k+1) * rmask(j,i,sr) & 2395 2459 * rflags_invers(j,i,k+1) … … 2435 2499 s1 = s1 - 0.5_wp * ( kh(k,j,i) + kh(k+1,j,i) ) & 2436 2500 * ( vpt(k+1,j,i) - vpt(k,j,i) ) & 2501 * rho_air_zw(k) & 2502 * heatflux_output_conversion(k) & 2437 2503 * ddzu(k+1) * rmask(j,i,sr) & 2438 2504 * rflags_invers(j,i,k+1) 2439 2505 s2 = s2 - 0.5_wp * ( kh(k,j,i) + kh(k+1,j,i) ) & 2440 2506 * ( q(k+1,j,i) - q(k,j,i) ) & 2507 * rho_air_zw(k) & 2508 * waterflux_output_conversion(k) & 2441 2509 * ddzu(k+1) * rmask(j,i,sr) & 2442 2510 * rflags_invers(j,i,k+1) … … 2459 2527 * ( ( q(k+1,j,i) - ql(k+1,j,i) ) & 2460 2528 - ( q(k,j,i) - ql(k,j,i) ) ) & 2529 * rho_air_zw(k) & 2530 * waterflux_output_conversion(k) & 2461 2531 * ddzu(k+1) * rmask(j,i,sr) & 2462 2532 * rflags_invers(j,i,k+1) … … 2506 2576 ! 2507 2577 !-- Subgridscale fluxes in the Prandtl layer 2508 s1 = s1 + usws(j,i) * rmask(j,i,sr) ! w"u" 2509 s2 = s2 + vsws(j,i) * rmask(j,i,sr) ! w"v" 2510 s3 = s3 + shf(j,i) * rmask(j,i,sr) ! w"pt" 2578 s1 = s1 + usws(j,i) * momentumflux_output_conversion(nzb) & 2579 * rmask(j,i,sr) ! w"u" 2580 s2 = s2 + vsws(j,i) * momentumflux_output_conversion(nzb) & 2581 * rmask(j,i,sr) ! w"v" 2582 s3 = s3 + shf(j,i) * heatflux_output_conversion(nzb) & 2583 * rmask(j,i,sr) ! w"pt" 2511 2584 s4 = s4 + 0.0_wp * rmask(j,i,sr) ! u"pt" 2512 2585 s5 = s5 + 0.0_wp * rmask(j,i,sr) ! v"pt" … … 2545 2618 DO i = nxl, nxr 2546 2619 DO j = nys, nyn 2547 s1 = s1 + qsws(j,i) * rmask(j,i,sr) ! w"q" (w"qv") 2620 s1 = s1 + qsws(j,i) * waterflux_output_conversion(nzb) & 2621 * rmask(j,i,sr) ! w"q" (w"qv") 2548 2622 s2 = s2 + ( ( 1.0_wp + 0.61_wp * q(nzb,j,i) ) * shf(j,i) & 2549 + 0.61_wp * pt(nzb,j,i) * qsws(j,i) ) 2623 + 0.61_wp * pt(nzb,j,i) * qsws(j,i) ) & 2624 * heatflux_output_conversion(nzb) 2550 2625 ENDDO 2551 2626 ENDDO … … 2564 2639 s1 = s1 + ( ( 1.0_wp + & 2565 2640 0.61_wp * q(nzb,j,i) - ql(nzb,j,i) ) * & 2566 shf(j,i) + 0.61_wp * pt(nzb,j,i) * qsws(j,i) ) 2641 shf(j,i) + 0.61_wp * pt(nzb,j,i) * qsws(j,i) )& 2642 * heatflux_output_conversion(nzb) 2567 2643 ENDDO 2568 2644 ENDDO … … 2582 2658 ! 2583 2659 !-- Formula does not work if ql(nzb) /= 0.0 2584 s1 = s1 + qsws(j,i) * rmask(j,i,sr) ! w"q" (w"qv") 2660 s1 = s1 + qsws(j,i) * waterflux_output_conversion(nzb) & 2661 * rmask(j,i,sr) ! w"q" (w"qv") 2585 2662 ENDDO 2586 2663 ENDDO … … 2624 2701 DO i = nxl, nxr 2625 2702 DO j = nys, nyn 2626 s1 = s1 + uswst(j,i) * rmask(j,i,sr) ! w"u" 2627 s2 = s2 + vswst(j,i) * rmask(j,i,sr) ! w"v" 2628 s3 = s3 + tswst(j,i) * rmask(j,i,sr) ! w"pt" 2703 s1 = s1 + uswst(j,i) * momentumflux_output_conversion(nzt:nzt+1) & 2704 * rmask(j,i,sr) ! w"u" 2705 s2 = s2 + vswst(j,i) * momentumflux_output_conversion(nzt:nzt+1) & 2706 * rmask(j,i,sr) ! w"v" 2707 s3 = s3 + tswst(j,i) * heatflux_output_conversion(nzt:nzt+1) & 2708 * rmask(j,i,sr) ! w"pt" 2629 2709 s4 = s4 + 0.0_wp * rmask(j,i,sr) ! u"pt" 2630 2710 s5 = s5 + 0.0_wp * rmask(j,i,sr) ! v"pt" … … 2663 2743 DO i = nxl, nxr 2664 2744 DO j = nys, nyn 2665 s1 = s1 + qswst(j,i) * rmask(j,i,sr) ! w"q" (w"qv") 2745 s1 = s1 + qswst(j,i) * waterflux_output_conversion(nzt) & 2746 * rmask(j,i,sr) ! w"q" (w"qv") 2666 2747 s2 = s2 + ( ( 1.0_wp + 0.61_wp * q(nzt,j,i) ) * tswst(j,i) +& 2667 0.61_wp * pt(nzt,j,i) * qswst(j,i) ) 2748 0.61_wp * pt(nzt,j,i) * qswst(j,i) ) & 2749 * heatflux_output_conversion(nzt) 2668 2750 ENDDO 2669 2751 ENDDO … … 2683 2765 0.61_wp * q(nzt,j,i) - ql(nzt,j,i) ) * & 2684 2766 tswst(j,i) + & 2685 0.61_wp * pt(nzt,j,i) * qswst(j,i) ) 2767 0.61_wp * pt(nzt,j,i) * qswst(j,i) ) & 2768 * heatflux_output_conversion(nzt) 2686 2769 ENDDO 2687 2770 ENDDO … … 2701 2784 ! 2702 2785 !-- Formula does not work if ql(nzb) /= 0.0 2703 s1 = s1 + qswst(j,i) * rmask(j,i,sr) ! w"q" (w"qv") 2786 s1 = s1 + qswst(j,i) * waterflux_output_conversion(nzt) & 2787 * rmask(j,i,sr) ! w"q" (w"qv") 2704 2788 ENDDO 2705 2789 ENDDO … … 2755 2839 !-- Energy flux w*e* (has to be adjusted?) 2756 2840 s3 = s3 + w(k,j,i) * 0.5_wp * ( ust**2 + vst**2 + w(k,j,i)**2 )& 2757 * rmask(j,i,sr) * rflags_invers(j,i,k+1) 2841 * rmask(j,i,sr) * rflags_invers(j,i,k+1) & 2842 * momentumflux_output_conversion(k) 2758 2843 ENDDO 2759 2844 ENDDO … … 2822 2907 s1 = s1 + 0.5_wp * ( vpt(k,j,i) - hom(k,1,44,sr) + & 2823 2908 vpt(k+1,j,i) - hom(k+1,1,44,sr) ) * & 2909 heatflux_output_conversion(k) * & 2824 2910 w(k,j,i) * rmask(j,i,sr) * rflags_invers(j,i,k+1) 2825 2911 ENDDO … … 2839 2925 s1 = s1 + 0.5_wp * ( ( q(k,j,i) - ql(k,j,i) ) - hom(k,1,42,sr) + & 2840 2926 ( q(k+1,j,i) - ql(k+1,j,i) ) - hom(k+1,1,42,sr) ) & 2927 * waterflux_output_conversion(k) & 2841 2928 * w(k,j,i) * rmask(j,i,sr) * rflags_invers(j,i,k+1) 2842 2929 ENDDO … … 2928 3015 s1 = s1 + 0.5_wp * ( vpt(k,j,i) - hom(k,1,44,sr) + & 2929 3016 vpt(k+1,j,i) - hom(k+1,1,44,sr) ) & 3017 * heatflux_output_conversion(k) & 2930 3018 * w(k,j,i) * rmask(j,i,sr) * rflags_invers(j,i,k+1) 2931 3019 ENDDO … … 2939 3027 !$acc parallel loop present( hom, sums_l ) 2940 3028 DO k = nzb, nzt_diff 2941 sums_l(k,46,tn) = ( 1.0_wp + 0.61_wp * hom(k,1,41,sr) ) * sums_l(k,17,tn) + & 2942 0.61_wp * hom(k,1,4,sr) * sums_l(k,49,tn) 3029 sums_l(k,46,tn) = ( ( 1.0_wp + 0.61_wp * hom(k,1,41,sr) ) * & 3030 sums_l(k,17,tn) + 0.61_wp * & 3031 hom(k,1,4,sr) * sums_l(k,49,tn) & 3032 ) * heatflux_output_conversion(k) 2943 3033 ENDDO 2944 3034 !$acc end parallel loop … … 2993 3083 s1 = s1 + 0.5_wp * ( w(k,j,i-1) + w(k,j,i) ) & 2994 3084 * ust * rmask(j,i,sr) & 3085 * momentumflux_output_conversion(k) & 2995 3086 * rflags_invers(j,i,k+1) 2996 3087 ! … … 2998 3089 s2 = s2 + 0.5_wp * ( w(k,j-1,i) + w(k,j,i) ) & 2999 3090 * vst * rmask(j,i,sr) & 3091 * momentumflux_output_conversion(k) & 3000 3092 * rflags_invers(j,i,k+1) 3001 3093 ENDDO 3002 3094 ENDDO 3003 3095 sums_l(k,13,tn) = s1 3004 sums_l(k,15,tn) = s 13096 sums_l(k,15,tn) = s2 3005 3097 ENDDO 3006 3098 !$acc end parallel loop … … 3021 3113 s1 = s1 + 0.5_wp * ( pt(k,j,i) - hom(k,1,4,sr) + & 3022 3114 pt(k+1,j,i) - hom(k+1,1,4,sr) ) & 3115 * heatflux_output_conversion(k) & 3023 3116 * w(k,j,i) * rmask(j,i,sr) & 3024 3117 * rflags_invers(j,i,k+1) … … 3039 3132 s1 = s1 + 0.5_wp * ( q(k,j,i) - hom(k,1,41,sr) + & 3040 3133 q(k+1,j,i) - hom(k+1,1,41,sr) ) & 3134 * waterflux_output_conversion(k) & 3041 3135 * w(k,j,i) * rmask(j,i,sr) & 3042 3136 * rflags_invers(j,i,k+1) … … 3167 3261 ( kh(k,j,i) + kh(k,j,i-1) ) & 3168 3262 * ( pt(k,j,i-1) - pt(k,j,i) ) & 3263 * rho_air_zw(k) & 3264 * heatflux_output_conversion(k) & 3169 3265 * ddx * rmask(j,i,sr) 3170 3266 sums_l(k,61,tn) = sums_l(k,61,tn) - 0.5_wp * & 3171 3267 ( kh(k,j,i) + kh(k,j-1,i) ) & 3172 3268 * ( pt(k,j-1,i) - pt(k,j,i) ) & 3269 * rho_air_zw(k) & 3270 * heatflux_output_conversion(k) & 3173 3271 * ddy * rmask(j,i,sr) 3174 3272 ! … … 3177 3275 ( u(k,j,i) - hom(k,1,1,sr) ) * 0.5_wp * & 3178 3276 ( pt(k,j,i-1) - hom(k,1,4,sr) + & 3179 pt(k,j,i) - hom(k,1,4,sr) ) 3277 pt(k,j,i) - hom(k,1,4,sr) ) & 3278 * heatflux_output_conversion(k) 3180 3279 pts = 0.5_wp * ( pt(k,j-1,i) - hom(k,1,4,sr) + & 3181 3280 pt(k,j,i) - hom(k,1,4,sr) ) … … 3183 3282 ( v(k,j,i) - hom(k,1,2,sr) ) * 0.5_wp * & 3184 3283 ( pt(k,j-1,i) - hom(k,1,4,sr) + & 3185 pt(k,j,i) - hom(k,1,4,sr) ) 3284 pt(k,j,i) - hom(k,1,4,sr) ) & 3285 * heatflux_output_conversion(k) 3186 3286 ENDDO 3187 3287 ENDDO … … 3468 3568 END IF 3469 3569 3570 hom(:,1,121,sr) = rho_air ! rho_air in Kg/m^3 3571 hom(:,1,122,sr) = rho_air_zw ! rho_air_zw in Kg/m^3 3572 3470 3573 hom(:,1,pr_palm,sr) = sums(:,pr_palm) 3471 3574 ! u*, w'u', w'v', t* (in last profile) … … 3576 3679 IF ( hom(nzb,1,18,sr) > 1.0E-8_wp .AND. z_i(1) /= 0.0_wp ) THEN 3577 3680 hom(nzb+8,1,pr_palm,sr) = & 3578 ( g / hom(k_surface_level+1,1,4,sr) * hom(k_surface_level,1,18,sr)& 3681 ( g / hom(k_surface_level+1,1,4,sr) * & 3682 ( hom(k_surface_level,1,18,sr) / heatflux_output_conversion(nzb) )& 3579 3683 * ABS( z_i(1) - mean_surface_level_height(sr) ) )**0.333333333_wp 3580 3684 ELSE
Note: See TracChangeset
for help on using the changeset viewer.