Changeset 3274 for palm/trunk/SOURCE/init_3d_model.f90
- Timestamp:
- Sep 24, 2018 3:42:55 PM (6 years ago)
- File:
-
- 1 edited
-
palm/trunk/SOURCE/init_3d_model.f90 (modified) (20 diffs)
Legend:
- Unmodified
- Added
- Removed
-
palm/trunk/SOURCE/init_3d_model.f90
r3241 r3274 25 25 ! ----------------- 26 26 ! $Id$ 27 ! Modularization of all bulk cloud physics code components 28 ! 29 ! 3241 2018-09-12 15:02:00Z raasch 27 30 ! unused variables removed 28 31 ! … … 255 258 ! precipitation_amount, precipitation_rate, prr moved to arrays_3d. 256 259 ! Initialization of nc_1d, nr_1d, pt_1d, qc_1d, qr_1d, q_1d moved to 257 ! microphysics_init.260 ! bcm_init. 258 261 ! 259 262 ! 1845 2016-04-08 08:29:13Z raasch … … 498 501 USE arrays_3d 499 502 503 USE basic_constants_and_equations_mod, & 504 ONLY: c_p, g, l_v, pi, r_d, exner_function, exner_function_invers, & 505 ideal_gas_law_rho, ideal_gas_law_rho_pt, barometric_formula 506 507 USE bulk_cloud_model_mod, & 508 ONLY: bulk_cloud_model, bcm_init, bcm_init_arrays 509 500 510 USE chemistry_model_mod, & 501 511 ONLY: chem_emissions 502 512 503 USE cloud_parameters, &504 ONLY: cp, l_v, r_d505 506 USE constants, &507 ONLY: pi508 509 513 USE control_parameters 510 514 … … 530 534 USE lsf_nudging_mod, & 531 535 ONLY: lsf_init, ls_forcing_surf, nudge_init 532 533 USE microphysics_mod, &534 ONLY: microphysics_init535 536 536 537 USE model_1d_mod, & … … 738 739 #endif 739 740 740 IF ( cloud_physics ) THEN741 !742 !-- Liquid water content743 #if defined( __nopointer )744 ALLOCATE ( ql(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )745 #else746 ALLOCATE ( ql_1(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )747 #endif748 749 !750 !-- 3D-cloud water content751 IF ( .NOT. microphysics_morrison ) THEN752 #if defined( __nopointer )753 ALLOCATE( qc(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )754 #else755 ALLOCATE( qc_1(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )756 #endif757 ENDIF758 !759 !-- Precipitation amount and rate (only needed if output is switched)760 ALLOCATE( precipitation_amount(nysg:nyng,nxlg:nxrg) )761 762 !763 !-- 3d-precipitation rate764 ALLOCATE( prr(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )765 766 IF ( microphysics_morrison ) THEN767 !768 !-- 3D-cloud drop water content, cloud drop concentration arrays769 #if defined( __nopointer )770 ALLOCATE( nc(nzb:nzt+1,nysg:nyng,nxlg:nxrg), &771 nc_p(nzb:nzt+1,nysg:nyng,nxlg:nxrg), &772 qc(nzb:nzt+1,nysg:nyng,nxlg:nxrg), &773 qc_p(nzb:nzt+1,nysg:nyng,nxlg:nxrg), &774 tnc_m(nzb:nzt+1,nysg:nyng,nxlg:nxrg), &775 tqc_m(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )776 #else777 ALLOCATE( nc_1(nzb:nzt+1,nysg:nyng,nxlg:nxrg), &778 nc_2(nzb:nzt+1,nysg:nyng,nxlg:nxrg), &779 nc_3(nzb:nzt+1,nysg:nyng,nxlg:nxrg), &780 qc_1(nzb:nzt+1,nysg:nyng,nxlg:nxrg), &781 qc_2(nzb:nzt+1,nysg:nyng,nxlg:nxrg), &782 qc_3(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )783 #endif784 ENDIF785 786 IF ( microphysics_seifert ) THEN787 !788 !-- 3D-rain water content, rain drop concentration arrays789 #if defined( __nopointer )790 ALLOCATE( nr(nzb:nzt+1,nysg:nyng,nxlg:nxrg), &791 nr_p(nzb:nzt+1,nysg:nyng,nxlg:nxrg), &792 qr(nzb:nzt+1,nysg:nyng,nxlg:nxrg), &793 qr_p(nzb:nzt+1,nysg:nyng,nxlg:nxrg), &794 tnr_m(nzb:nzt+1,nysg:nyng,nxlg:nxrg), &795 tqr_m(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )796 #else797 ALLOCATE( nr_1(nzb:nzt+1,nysg:nyng,nxlg:nxrg), &798 nr_2(nzb:nzt+1,nysg:nyng,nxlg:nxrg), &799 nr_3(nzb:nzt+1,nysg:nyng,nxlg:nxrg), &800 qr_1(nzb:nzt+1,nysg:nyng,nxlg:nxrg), &801 qr_2(nzb:nzt+1,nysg:nyng,nxlg:nxrg), &802 qr_3(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )803 #endif804 ENDIF805 806 ENDIF807 808 741 IF ( cloud_droplets ) THEN 809 742 ! … … 868 801 ! 869 802 !-- Density profile calculation for anelastic approximation 870 t_surface = pt_surface * ( surface_pressure / 1000.0_wp )**( r_d / c p )803 t_surface = pt_surface * ( surface_pressure / 1000.0_wp )**( r_d / c_p ) 871 804 IF ( TRIM( approximation ) == 'anelastic' ) THEN 872 805 DO k = nzb, nzt+1 873 806 p_hydrostatic(k) = surface_pressure * 100.0_wp * & 874 ( 1 - ( g * zu(k) ) / ( c p * t_surface )&875 )**( c p / r_d )807 ( 1 - ( g * zu(k) ) / ( c_p * t_surface ) & 808 )**( c_p / r_d ) 876 809 rho_air(k) = ( p_hydrostatic(k) * & 877 810 ( 100000.0_wp / p_hydrostatic(k) & 878 )**( r_d / c p )&811 )**( r_d / c_p ) & 879 812 ) / ( r_d * pt_init(k) ) 880 813 ENDDO … … 887 820 DO k = nzb, nzt+1 888 821 p_hydrostatic(k) = surface_pressure * 100.0_wp * & 889 ( 1 - ( g * zu(nzb) ) / ( c p * t_surface )&890 )**( c p / r_d )822 ( 1 - ( g * zu(nzb) ) / ( c_p * t_surface ) & 823 )**( c_p / r_d ) 891 824 rho_air(k) = ( p_hydrostatic(k) * & 892 825 ( 100000.0_wp / p_hydrostatic(k) & 893 )**( r_d / c p )&826 )**( r_d / c_p ) & 894 827 ) / ( r_d * pt_init(nzb) ) 895 828 ENDDO … … 923 856 momentumflux_input_conversion(k) = rho_air_zw(k) 924 857 ELSEIF ( TRIM( flux_input_mode ) == 'dynamic' ) THEN 925 heatflux_input_conversion(k) = 1.0_wp / c p858 heatflux_input_conversion(k) = 1.0_wp / c_p 926 859 waterflux_input_conversion(k) = 1.0_wp / l_v 927 860 momentumflux_input_conversion(k) = 1.0_wp … … 933 866 momentumflux_output_conversion(k) = drho_air_zw(k) 934 867 ELSEIF ( TRIM( flux_output_mode ) == 'dynamic' ) THEN 935 heatflux_output_conversion(k) = c p868 heatflux_output_conversion(k) = c_p 936 869 waterflux_output_conversion(k) = l_v 937 870 momentumflux_output_conversion(k) = 1.0_wp … … 1065 998 IF ( humidity ) THEN 1066 999 q => q_1; q_p => q_2; tq_m => q_3 1067 IF ( humidity ) THEN 1068 vpt => vpt_1 1069 IF ( cloud_physics ) THEN 1070 ql => ql_1 1071 IF ( .NOT. microphysics_morrison ) THEN 1072 qc => qc_1 1073 ENDIF 1074 IF ( microphysics_morrison ) THEN 1075 qc => qc_1; qc_p => qc_2; tqc_m => qc_3 1076 nc => nc_1; nc_p => nc_2; tnc_m => nc_3 1077 ENDIF 1078 IF ( microphysics_seifert ) THEN 1079 qr => qr_1; qr_p => qr_2; tqr_m => qr_3 1080 nr => nr_1; nr_p => nr_2; tnr_m => nr_3 1081 ENDIF 1082 ENDIF 1083 ENDIF 1000 vpt => vpt_1 1084 1001 IF ( cloud_droplets ) THEN 1085 1002 ql => ql_1 … … 1102 1019 !-- Initialize surface arrays 1103 1020 CALL init_surface_arrays 1021 ! 1022 !-- Allocate microphysics module arrays 1023 IF ( bulk_cloud_model ) THEN 1024 CALL bcm_init_arrays 1025 ENDIF 1104 1026 ! 1105 1027 !-- Allocate land surface model arrays … … 1297 1219 !-- Set inital w to 0 1298 1220 w = 0.0_wp 1299 !1300 !-- Initialize the remaining quantities1301 IF ( humidity ) THEN1302 IF ( cloud_physics .AND. microphysics_morrison ) THEN1303 DO i = nxlg, nxrg1304 DO j = nysg, nyng1305 qc(:,j,i) = 0.0_wp1306 nc(:,j,i) = 0.0_wp1307 ENDDO1308 ENDDO1309 ENDIF1310 1311 IF ( cloud_physics .AND. microphysics_seifert ) THEN1312 DO i = nxlg, nxrg1313 DO j = nysg, nyng1314 qr(:,j,i) = 0.0_wp1315 nr(:,j,i) = 0.0_wp1316 ENDDO1317 ENDDO1318 ENDIF1319 1320 ENDIF1321 1221 1322 1222 IF ( passive_scalar ) THEN … … 1376 1276 ENDDO 1377 1277 ENDDO 1378 IF ( cloud_physics .AND. microphysics_morrison ) THEN1379 DO i = nxlg, nxrg1380 DO j = nysg, nyng1381 qc(:,j,i) = 0.0_wp1382 nc(:,j,i) = 0.0_wp1383 ENDDO1384 ENDDO1385 ENDIF1386 IF ( cloud_physics .AND. microphysics_seifert ) THEN1387 DO i = nxlg, nxrg1388 DO j = nysg, nyng1389 qr(:,j,i) = 0.0_wp1390 nr(:,j,i) = 0.0_wp1391 ENDDO1392 ENDDO1393 ENDIF1394 1278 ENDIF 1395 1279 … … 1485 1369 ENDDO 1486 1370 ENDDO 1487 IF ( cloud_physics .AND. microphysics_morrison ) THEN1488 DO i = nxlg, nxrg1489 DO j = nysg, nyng1490 qc(:,j,i) = 0.0_wp1491 nc(:,j,i) = 0.0_wp1492 ENDDO1493 ENDDO1494 ENDIF1495 1496 IF ( cloud_physics .AND. microphysics_seifert ) THEN1497 DO i = nxlg, nxrg1498 DO j = nysg, nyng1499 qr(:,j,i) = 0.0_wp1500 nr(:,j,i) = 0.0_wp1501 ENDDO1502 ENDDO1503 ENDIF1504 1505 1371 ENDIF 1506 1372 … … 1598 1464 !-- Store initial profile of mixing ratio and potential 1599 1465 !-- temperature 1600 IF ( cloud_physics.OR. cloud_droplets ) THEN1466 IF ( bulk_cloud_model .OR. cloud_droplets ) THEN 1601 1467 hom(:,1,27,:) = SPREAD( q(:,nys,nxl), 2, statistic_regions+1 ) 1602 1468 hom(:,1,28,:) = SPREAD( pt(:,nys,nxl), 2, statistic_regions+1 ) … … 1645 1511 !-- In case of iterative solvers, p must get an initial value 1646 1512 IF ( psolver(1:9) == 'multigrid' .OR. psolver == 'sor' ) p = 0.0_wp 1647 1648 !1649 !-- Treating cloud physics, liquid water content and precipitation amount1650 !-- are zero at beginning of the simulation1651 IF ( cloud_physics ) THEN1652 ql = 0.0_wp1653 qc = 0.0_wp1654 1655 precipitation_amount = 0.0_wp1656 ENDIF1657 1513 ! 1658 1514 !-- Impose vortex with vertical axis on the initial velocity profile … … 1693 1549 tq_m = 0.0_wp 1694 1550 q_p = q 1695 IF ( cloud_physics .AND. microphysics_morrison ) THEN1696 tqc_m = 0.0_wp1697 qc_p = qc1698 tnc_m = 0.0_wp1699 nc_p = nc1700 ENDIF1701 IF ( cloud_physics .AND. microphysics_seifert ) THEN1702 tqr_m = 0.0_wp1703 qr_p = qr1704 tnr_m = 0.0_wp1705 nr_p = nr1706 ENDIF1707 1551 ENDIF 1708 1552 … … 1963 1807 IF ( humidity ) THEN 1964 1808 q_p = q 1965 IF ( cloud_physics .AND. microphysics_morrison ) THEN1966 qc_p = qc1967 nc_p = nc1968 ENDIF1969 IF ( cloud_physics .AND. microphysics_seifert ) THEN1970 qr_p = qr1971 nr_p = nr1972 ENDIF1973 1809 ENDIF 1974 1810 IF ( passive_scalar ) s_p = s … … 1982 1818 IF ( humidity ) THEN 1983 1819 tq_m = 0.0_wp 1984 IF ( cloud_physics .AND. microphysics_morrison ) THEN1985 tqc_m = 0.0_wp1986 tnc_m = 0.0_wp1987 ENDIF1988 IF ( cloud_physics .AND. microphysics_seifert ) THEN1989 tqr_m = 0.0_wp1990 tnr_m = 0.0_wp1991 ENDIF1992 1820 ENDIF 1993 1821 IF ( passive_scalar ) ts_m = 0.0_wp … … 2417 2245 CALL init_ocean 2418 2246 2419 E LSE2247 ENDIF 2420 2248 ! 2421 2249 !-- Initialize quantities for handling cloud physics 2422 2250 !-- This routine must be called before lpm_init, because 2423 !-- otherwise, array pt_d_t, needed in data_output_dvrp (called by2251 !-- otherwise, array d_exner, needed in data_output_dvrp (called by 2424 2252 !-- lpm_init) is not defined. 2425 CALL init_cloud_physics 2426 ! 2427 !-- Initialize bulk cloud microphysics 2428 CALL microphysics_init 2253 IF ( .NOT. ocean ) THEN 2254 2255 ALLOCATE( hyp(nzb:nzt+1) ) 2256 ALLOCATE( d_exner(nzb:nzt+1) ) 2257 ALLOCATE( exner(nzb:nzt+1) ) 2258 ALLOCATE( hyrho(nzb:nzt+1) ) 2259 ! 2260 !-- Check temperature in case of too large domain height 2261 DO k = nzb, nzt+1 2262 IF ( ( pt_surface * exner_function(surface_pressure * 100.0_wp) - g/c_p * zu(k) ) < 0.0_wp ) THEN 2263 WRITE( message_string, * ) 'absolute temperature < 0.0 at zu(', k, & 2264 ') = ', zu(k) 2265 CALL message( 'init_bulk_cloud_model', 'PA0142', 1, 2, 0, 6, 0 ) 2266 ENDIF 2267 ENDDO 2268 2269 ! 2270 !-- Calculate vertical profile of the hydrostatic pressure (hyp) 2271 hyp = barometric_formula(zu, pt_surface * exner_function(surface_pressure * 100.0_wp), surface_pressure * 100.0_wp) 2272 d_exner = exner_function_invers(hyp) 2273 exner = 1.0_wp / exner_function_invers(hyp) 2274 hyrho = ideal_gas_law_rho_pt(hyp, pt_init) 2275 ! 2276 !-- Compute reference density 2277 rho_surface = ideal_gas_law_rho(surface_pressure * 100.0_wp, pt_surface * exner_function(surface_pressure * 100.0_wp)) 2278 2279 ENDIF 2280 ! 2281 !-- If required, initialize quantities needed for the microphysics module 2282 IF ( bulk_cloud_model ) THEN 2283 CALL bcm_init 2429 2284 ENDIF 2430 2285
Note: See TracChangeset
for help on using the changeset viewer.
