 Timestamp:
 Jul 11, 2019 1:57:56 PM (4 years ago)
 Location:
 palm/trunk/SOURCE
 Files:

 2 edited
Legend:
 Unmodified
 Added
 Removed

palm/trunk/SOURCE/basic_constants_and_equations_mod.f90
r4084 r4088 25 25 !  26 26 ! $Id$ 27 ! Comment of barometric formula improved, function for ideal gas law revised 28 ! 29 ! 4084 20190710 17:09:11Z knoop 27 30 ! Changed precomputed fractions to be variable based 28 31 ! … … 236 239 ! 237 240 ! compute density according to ideal gas law: 238 ideal_gas_law_rho_pt_0d = p / (r_d * (1.0_wp / exner_function_invers(p)) * t)241 ideal_gas_law_rho_pt_0d = p / (r_d * exner_function(p) * t) 239 242 240 243 END FUNCTION ideal_gas_law_rho_pt_0d … … 255 258 ! 256 259 ! compute density according to ideal gas law: 257 ideal_gas_law_rho_pt_1d = p / (r_d * (1.0_wp / exner_function_invers(p)) * t)260 ideal_gas_law_rho_pt_1d = p / (r_d * exner_function(p) * t) 258 261 259 262 END FUNCTION ideal_gas_law_rho_pt_1d … … 334 337 ! Description: 335 338 !  336 !> Compute the barometric formula for scalar arguments. 339 !> Compute the barometric formula for scalar arguments. The calculation is 340 !> based on the assumption of a polytropic atmosphere and neutral 341 !> stratification, where the temperature lapse rate is g/cp. 337 342 !! 338 343 FUNCTION barometric_formula_0d( z, t_0, p_0) … … 354 359 ! Description: 355 360 !  356 !> Compute the barometric formula for 1D array arguments. 361 !> Compute the barometric formula for 1D array arguments. The calculation is 362 !> based on the assumption of a polytropic atmosphere and neutral 363 !> stratification, where the temperature lapse rate is g/cp. 357 364 !! 358 365 FUNCTION barometric_formula_1d( z, t_0, p_0) 
palm/trunk/SOURCE/init_3d_model.f90
r4048 r4088 25 25 !  26 26 ! $Id$ 27 ! Pressure and density profile calculations revised using basic functions 28 ! 29 ! 4048 20190621 21:00:21Z knoop 27 30 ! Further modularization of particle code components 28 31 ! … … 845 848 ALLOCATE( drho_air_zw(nzb:nzt+1) ) 846 849 ! 847 ! Density profile calculation for anelastic approximation 848 t_surface = pt_surface * ( surface_pressure / 1000.0_wp )**( r_d / c_p ) 849 IF ( TRIM( approximation ) == 'anelastic' ) THEN 850 ! Density profile calculation for anelastic and Boussinesq approximation 851 ! In case of a Boussinesq approximation, a constant density is calculated 852 ! mainly for output purposes. This density do not need to be considered 853 ! in the model's system of equations. 854 ! t_surface = pt_surface * ( surface_pressure / 1000.0_wp )**( r_d / c_p ) 855 IF ( TRIM( approximation ) == 'anelastic' ) THEN 850 856 DO k = nzb, nzt+1 851 p_hydrostatic(k) = surface_pressure * 100.0_wp * & 852 ( 1  ( g * zu(k) ) / ( c_p * t_surface ) & 853 )**( c_p / r_d ) 854 rho_air(k) = ( p_hydrostatic(k) * & 855 ( 100000.0_wp / p_hydrostatic(k) & 856 )**( r_d / c_p ) & 857 ) / ( r_d * pt_init(k) ) 857 p_hydrostatic(k) = barometric_formula(zu(k), pt_surface * & 858 exner_function(surface_pressure * 100.0_wp), & 859 surface_pressure * 100.0_wp) 860 861 rho_air(k) = ideal_gas_law_rho_pt(p_hydrostatic(k), pt_init(k)) 858 862 ENDDO 863 859 864 DO k = nzb, nzt 860 865 rho_air_zw(k) = 0.5_wp * ( rho_air(k) + rho_air(k+1) ) 861 866 ENDDO 867 862 868 rho_air_zw(nzt+1) = rho_air_zw(nzt) & 863 869 + 2.0_wp * ( rho_air(nzt+1)  rho_air_zw(nzt) ) 870 864 871 ELSE 865 872 DO k = nzb, nzt+1 866 p_hydrostatic(k) = surface_pressure * 100.0_wp * & 867 ( 1  ( g * zu(nzb) ) / ( c_p * t_surface ) & 868 )**( c_p / r_d ) 869 rho_air(k) = ( p_hydrostatic(k) * & 870 ( 100000.0_wp / p_hydrostatic(k) & 871 )**( r_d / c_p ) & 872 ) / ( r_d * pt_init(nzb) ) 873 p_hydrostatic(k) = barometric_formula(zu(nzb), pt_surface * & 874 exner_function(surface_pressure * 100.0_wp), & 875 surface_pressure * 100.0_wp) 876 877 rho_air(k) = ideal_gas_law_rho_pt(p_hydrostatic(k), pt_init(nzb)) 873 878 ENDDO 879 874 880 DO k = nzb, nzt 875 881 rho_air_zw(k) = 0.5_wp * ( rho_air(k) + rho_air(k+1) ) 876 882 ENDDO 883 877 884 rho_air_zw(nzt+1) = rho_air_zw(nzt) & 878 885 + 2.0_wp * ( rho_air(nzt+1)  rho_air_zw(nzt) ) 886 879 887 ENDIF 880 888 !
Note: See TracChangeset
for help on using the changeset viewer.