Changeset 2252


Ignore:
Timestamp:
Jun 7, 2017 9:35:37 AM (4 years ago)
Author:
knoop
Message:

air density now depending on surface_pressure even in boussinesq mode

Location:
palm/trunk/SOURCE
Files:
2 edited

Legend:

Unmodified
Added
Removed
  • palm/trunk/SOURCE/flow_statistics.f90

    r2233 r2252  
    2525! -----------------
    2626! $Id$
     27! perturbation pressure now depending on flux_output_mode
     28!
     29! 2233 2017-05-30 18:08:54Z suehring
    2730!
    2831! 2232 2017-05-30 17:47:52Z suehring
     
    637640                sums_l(k,10,tn) = sums_l(k,10,tn) + kh(k,j,i) * rmask(j,i,sr)  &
    638641                                                              * flag
    639                 sums_l(k,40,tn) = sums_l(k,40,tn) + p(k,j,i)  * flag
     642                sums_l(k,40,tn) = sums_l(k,40,tn) + ( p(k,j,i)                 &
     643                                         / momentumflux_output_conversion(k) ) &
     644                                                              * flag
    640645
    641646                sums_l(k,33,tn) = sums_l(k,33,tn) + &
     
    14151420
    14161421                   sums_ll(k,2) = sums_ll(k,2) + 0.5_wp * w(k,j,i)             &
    1417                                                * ( p(k,j,i) + p(k+1,j,i) ) * flag
     1422                                       * ( ( p(k,j,i) + p(k+1,j,i) )           &
     1423                                         / momentumflux_output_conversion(k) ) &
     1424                                       * flag
    14181425
    14191426                ENDDO
     
    17501757       hom(:,1,70,sr) = sums(:,70)     ! q*2
    17511758       hom(:,1,71,sr) = sums(:,71)     ! prho
    1752        hom(:,1,72,sr) = hyp * 1E-4_wp  ! hyp in dbar
     1759       hom(:,1,72,sr) = hyp * 1E-2_wp  ! hyp in hPa
    17531760       hom(:,1,73,sr) = sums(:,73)     ! nr
    17541761       hom(:,1,74,sr) = sums(:,74)     ! qr
     
    19311938       IF ( hom(k_surface_level,1,18,sr) > 1.0E-8_wp  .AND.  z_i(1) /= 0.0_wp )&
    19321939       THEN
    1933           hom(nzb+8,1,pr_palm,sr) = &
     1940          hom(nzb+8,1,pr_palm,sr) =                                            &
    19341941             ( g / hom(k_surface_level+1,1,4,sr) *                             &
    1935              ( hom(k_surface_level,1,18,sr) / heatflux_output_conversion(nzb) )&
     1942             ( hom(k_surface_level,1,18,sr) /                                  &
     1943             ( heatflux_output_conversion(nzb) * rho_air(nzb) ) )              &
    19361944             * ABS( z_i(1) - mean_surface_level_height(sr) ) )**0.333333333_wp
    19371945       ELSE
  • palm/trunk/SOURCE/init_3d_model.f90

    r2233 r2252  
    2525! -----------------
    2626! $Id$
     27! rho_air now depending on surface_pressure even in Boussinesq mode
     28!
     29! 2233 2017-05-30 18:08:54Z suehring
    2730!
    2831! 2232 2017-05-30 17:47:52Z suehring
     
    657660!
    658661!-- Density profile calculation for anelastic approximation
     662    t_surface = pt_surface * ( surface_pressure / 1000.0_wp )**( r_d / cp )
    659663    IF ( TRIM( approximation ) == 'anelastic' ) THEN
    660        t_surface = pt_surface * ( surface_pressure / 1000.0_wp )**( r_d / cp )
    661664       DO  k = nzb, nzt+1
    662665          p_hydrostatic(k)    = surface_pressure * 100.0_wp *                  &
     
    674677                            + 2.0_wp * ( rho_air(nzt+1) - rho_air_zw(nzt)  )
    675678    ELSE
    676        rho_air     = 1.0_wp
    677        rho_air_zw  = 1.0_wp
     679       DO  k = nzb, nzt+1
     680          p_hydrostatic(k)    = surface_pressure * 100.0_wp *                  &
     681                                ( 1 - ( g * zu(nzb) ) / ( cp * t_surface )       &
     682                                )**( cp / r_d )
     683          rho_air(k)          = ( p_hydrostatic(k) *                           &
     684                                  ( 100000.0_wp / p_hydrostatic(k)             &
     685                                  )**( r_d / cp )                              &
     686                                ) / ( r_d * pt_init(nzb) )
     687       ENDDO
     688       DO  k = nzb, nzt
     689          rho_air_zw(k) = 0.5_wp * ( rho_air(k) + rho_air(k+1) )
     690       ENDDO
     691       rho_air_zw(nzt+1)  = rho_air_zw(nzt)                                    &
     692                            + 2.0_wp * ( rho_air(nzt+1) - rho_air_zw(nzt)  )
    678693    ENDIF
    679694
Note: See TracChangeset for help on using the changeset viewer.