 Timestamp:
 Jul 11, 2018 9:59:11 AM (4 years ago)
 File:

 1 edited
Legend:
 Unmodified
 Added
 Removed

palm/trunk/SOURCE/radiation_model_mod.f90
r3116 r3117 28 28 !  29 29 ! $Id$ 30 ! Bugfix: water vapor was not transfered to RRTMG when cloud_physics = .F. 31 ! Bugfix: changed the calculation of RRTMG boundary conditions (Mohamed Salim) 32 ! Bugfix: dry residual atmosphere is replaced by standard RRTMG atmosphere 33 ! 34 ! 3116 20180710 14:31:58Z suehring 30 35 ! Output of long/shortwave radiation at surface 31 36 ! … … 370 375 USE control_parameters, & 371 376 ONLY: cloud_droplets, cloud_physics, coupling_char, dz, g, & 377 humidity, & 372 378 initializing_actions, io_blocks, io_group, & 373 379 latitude, longitude, large_scale_forcing, lsf_surf, & … … 618 624 619 625 REAL(wp), DIMENSION(:), ALLOCATABLE :: hyp_snd, & !< hypostatic pressure from sounding data (hPa) 620 q_snd, & !< mixing ratio from sounding data (kg/kg)  dummy at the moment621 626 rrtm_tsfc, & !< dummy array for storing surface temperature 622 627 t_snd !< actual temperature from sounding data (hPa) … … 829 834 REAL(wp), DIMENSION(:), ALLOCATABLE :: surfoutlw !< array of total lw radiation outgoing from nonvirtual surfaces surfaces after all reflection 830 835 REAL(wp), DIMENSION(:), ALLOCATABLE :: surfhf !< array of total radiation flux incoming to minus outgoing from local surface 836 REAL(wp), DIMENSION(:), ALLOCATABLE :: surfemitlwl !< array of emitted lw radiation for local surface used to calculate effective surface temperature for radiation model 831 837 832 838 ! block variables needed for calculation of the plant canopy model inside the urban surface model … … 2850 2856 CALL calc_mean_profile( pt, 4 ) 2851 2857 pt_av = hom(:, 1, 4, 0) 2852 2858 2859 IF ( humidity ) THEN 2860 CALL calc_mean_profile( q, 41 ) 2861 q_av = hom(:, 1, 41, 0) 2862 ENDIF 2853 2863 ! 2854 2864 ! Prepare profiles of temperature and H2O volume mixing ratio … … 2856 2866 2857 2867 IF ( cloud_physics ) THEN 2858 CALL calc_mean_profile( q, 41 ) 2859 ! average q is now in hom(:, 1, 41, 0) 2860 q_av = hom(:, 1, 41, 0) 2868 2861 2869 CALL calc_mean_profile( ql, 54 ) 2862 2870 ! average ql is now in hom(:, 1, 54, 0) … … 2872 2880 rrtm_tlay(0,k) = pt_av(k) * ( (hyp(k) ) / 100000._wp & 2873 2881 )**.286_wp 2874 rrtm_h2ovmr(0,k) = 0._wp 2875 ENDDO 2882 ENDDO 2883 2884 IF ( humidity ) THEN 2885 DO k = nzb+1, nzt+1 2886 rrtm_h2ovmr(0,k) = mol_mass_air_d_wv * q_av(k) 2887 ENDDO 2888 ELSE 2889 rrtm_h2ovmr(0,nzb+1:nzt+1) = 0.0_wp 2890 ENDIF 2876 2891 ENDIF 2877 2892 … … 2943 2958 ! Set surface temperature 2944 2959 rrtm_tsfc = t_rad_urb 2945 2960 2946 2961 IF ( lw_radiation ) THEN 2947 2962 CALL rrtmg_lw( 1, nzt_rad , rrtm_icld , rrtm_idrv ,& … … 3081 3096 rrtm_tlay(0,k) = pt(k,j,i) * ( (hyp(k) ) / 100000.0_wp & 3082 3097 )**0.286_wp 3083 rrtm_h2ovmr(0,k) = 0.0_wp3084 3098 ENDDO 3099 3100 IF ( humidity ) THEN 3101 DO k = nzb+1, nzt+1 3102 rrtm_h2ovmr(0,k) = mol_mass_air_d_wv * q(k,j,i) 3103 ENDDO 3104 ELSE 3105 rrtm_h2ovmr(0,nzb+1:nzt+1) = 0.0_wp 3106 ENDIF 3085 3107 ENDIF 3086 3108 … … 3742 3764 DEALLOCATE( hyp_snd ) 3743 3765 DEALLOCATE( t_snd ) 3744 DEALLOCATE( q_snd )3745 3766 DEALLOCATE ( rrtm_play ) 3746 3767 DEALLOCATE ( rrtm_plev ) … … 3748 3769 DEALLOCATE ( rrtm_tlev ) 3749 3770 3750 DEALLOCATE ( rrtm_h2ovmr )3751 3771 DEALLOCATE ( rrtm_cicewp ) 3752 3772 DEALLOCATE ( rrtm_cldfr ) … … 3844 3864 3845 3865 ! 3846 ! Save data above LES domain in hyp_snd, t_snd and q_snd 3847 ! Note: q_snd_tmp is not calculated at the moment (dry residual atmosphere) 3866 ! Save data above LES domain in hyp_snd, t_snd 3848 3867 ALLOCATE( hyp_snd(nzb+1:nzt_rad) ) 3849 3868 ALLOCATE( t_snd(nzb+1:nzt_rad) ) 3850 ALLOCATE( q_snd(nzb+1:nzt_rad) )3851 3869 hyp_snd = 0.0_wp 3852 3870 t_snd = 0.0_wp 3853 q_snd = 0.0_wp3854 3871 3855 3872 hyp_snd(nzt+2:nzt_rad) = hyp_snd_tmp(nz_snd_start+1:nz_snd_end) … … 3891 3908 ALLOCATE ( rrtm_tlay(0:0,nzb+1:nzt_rad+1) ) 3892 3909 ALLOCATE ( rrtm_tlev(0:0,nzb+1:nzt_rad+2) ) 3893 ALLOCATE ( rrtm_h2ovmr(0:0,nzb+1:nzt_rad+1) )3894 3910 3895 3911 DO k = nzt+8, nzt_rad 3896 3912 rrtm_tlay(0,k) = t_snd(k) 3897 rrtm_h2ovmr(0,k) = q_snd(k)3898 3913 ENDDO 3899 3914 rrtm_tlay(0,nzt_rad+1) = 2.0_wp * rrtm_tlay(0,nzt_rad) & … … 3905 3920 * ( rrtm_plev(0,k)  rrtm_play(0,k1) ) 3906 3921 ENDDO 3907 rrtm_h2ovmr(0,nzt_rad+1) = rrtm_h2ovmr(0,nzt_rad)3908 3922 3909 3923 rrtm_tlev(0,nzt_rad+2) = 2.0_wp * rrtm_tlay(0,nzt_rad+1) & … … 3994 4008 IMPLICIT NONE 3995 4009 3996 INTEGER(iwp), PARAMETER :: num_trace_gases = 9!< number of trace gases (absorbers)4010 INTEGER(iwp), PARAMETER :: num_trace_gases = 10 !< number of trace gases (absorbers) 3997 4011 3998 4012 CHARACTER(LEN=5), DIMENSION(num_trace_gases), PARAMETER :: & !< trace gas names 3999 4013 trace_names = (/'O3 ', 'CO2 ', 'CH4 ', 'N2O ', 'O2 ', & 4000 'CFC11', 'CFC12', 'CFC22', 'CCL4 ' /)4014 'CFC11', 'CFC12', 'CFC22', 'CCL4 ', 'H2O'/) 4001 4015 4002 4016 INTEGER(iwp) :: id, & !< NetCDF id … … 4036 4050 DEALLOCATE ( rrtm_cfc22vmr ) 4037 4051 DEALLOCATE ( rrtm_ccl4vmr ) 4052 DEALLOCATE ( rrtm_h2ovmr ) 4038 4053 ENDIF 4039 4054 … … 4049 4064 ALLOCATE ( rrtm_cfc22vmr(0:0,1:nzt_rad+1) ) 4050 4065 ALLOCATE ( rrtm_ccl4vmr(0:0,1:nzt_rad+1) ) 4066 ALLOCATE ( rrtm_h2ovmr(0:0,1:nzt_rad+1) ) 4051 4067 4052 4068 ! … … 4226 4242 rrtm_ccl4vmr(0,:) = trace_path_tmp(:) 4227 4243 4244 CASE ( 'H2O' ) 4245 4246 rrtm_h2ovmr(0,:) = trace_path_tmp(:) 4247 4228 4248 CASE DEFAULT 4229 4249 … … 4388 4408 REAL(wp) :: area_surfl !< total area of surfaces in local processor 4389 4409 REAL(wp) :: area_surf !< total area of surfaces in all processor 4410 REAL(wp) :: area_horl !< total horizontal area of domain in local processor 4411 REAL(wp) :: area_hor !< total horizontal area of domain in all processor 4390 4412 4391 4413 … … 4605 4627 surfoutsw = 0.0_wp 4606 4628 surfoutlw = surfoutll 4629 surfemitlwl = surfoutll 4607 4630 ! surfhf = surfinsw + surfinlw  surfoutsw  surfoutlw 4608 4631 … … 4826 4849 IF ( kdir(d) == 0 ) facearea(d) = facearea(d) * dz(1) 4827 4850 ENDDO 4851 ! calculate horizontal area 4852 ! !!! ATTENTION!!! uniform grid is assumed here 4853 area_hor = (nx+1) * (ny+1) * dx * dy 4828 4854 ! 4829 4855 ! absorbed/received SW & LW and emitted LW energy of all physical … … 4839 4865 d = surfl(id, i) 4840 4866 ! received SW & LW 4841 pinswl = pinswl + surfinsw(i) * facearea(d)4842 pinlwl = pinlwl + surfinlw (i) * facearea(d)4867 pinswl = pinswl + (surfinswdir(i) + surfinswdif(i)) * facearea(d) 4868 pinlwl = pinlwl + surfinlwdif(i) * facearea(d) 4843 4869 ! absorbed SW & LW 4844 4870 pabsswl = pabsswl + (1._wp  albedo_surf(i)) * & … … 4846 4872 pabslwl = pabslwl + emiss_surf(i) * surfinlw(i) * facearea(d) 4847 4873 ! emitted LW 4848 pemitlwl = pemitlwl + surf outlw(i) * facearea(d)4874 pemitlwl = pemitlwl + surfemitlwl(i) * facearea(d) 4849 4875 ! emissivity and area sum 4850 4876 emiss_sum_surfl = emiss_sum_surfl + emiss_surf(i) * facearea(d) … … 4856 4882 pabsswl = pabsswl + SUM(pcbinsw) 4857 4883 pabslwl = pabslwl + SUM(pcbinlw) 4884 pinswl = pinswl + SUM(pcbinswdir) + SUM(pcbinswdif) 4858 4885 ENDIF 4859 4886 ! … … 4871 4898 pinlw = pinlwl 4872 4899 pabssw = pabsswl 4873 pabslw l = pabslw4874 pemitlw l = pemitlw4900 pabslw = pabslwl 4901 pemitlw = pemitlwl 4875 4902 emiss_sum_surf = emiss_sum_surfl 4876 4903 area_surf = area_surfl … … 4878 4905 4879 4906 ! (1) albedo 4880 IF ( pinsw /= 0.0_wp ) albedo_urb = 1._wp  pabssw / pinsw4881 4907 IF ( pinsw /= 0.0_wp ) & 4908 albedo_urb = (pinsw  pabssw) / pinsw 4882 4909 ! (2) average emmsivity 4883 IF ( area_surf /= 0.0_wp ) emissivity_urb = emiss_sum_surf / area_surf4884 4910 IF ( area_surf /= 0.0_wp ) & 4911 emissivity_urb = emiss_sum_surf / area_surf 4885 4912 ! (3) temperature 4886 t_rad_urb = ((pemitlw  pabslw + emissivity_urb*pinlw)/(emissivity_urb*sigma_sb*area_surf))**0.25_wp 4913 t_rad_urb = ( (pemitlw  pabslw + emissivity_urb*pinlw) / & 4914 (emissivity_urb*sigma_sb * area_hor) )**0.25_wp 4887 4915 4888 4889 4916 CONTAINS 4890 4917 … … 5356 5383 ALLOCATE( skyvf(nsurfl) ) 5357 5384 ALLOCATE( skyvft(nsurfl) ) 5385 ALLOCATE( surfemitlwl(nsurfl) ) 5358 5386 5359 5387 ! … … 5362 5390 ! For now set an arbitrary initial value. 5363 5391 IF ( average_radiation ) THEN 5364 albedo_urb = 0. 5_wp5365 emissivity_urb = 0. 5_wp5392 albedo_urb = 0.1_wp 5393 emissivity_urb = 0.9_wp 5366 5394 t_rad_urb = pt_surface 5367 5395 ENDIF
Note: See TracChangeset
for help on using the changeset viewer.