Changeset 2328 for palm/trunk/SOURCE/radiation_model_mod.f90

Timestamp:

Aug 3, 2017 12:34:22 PM (7 years ago)

Author:

maronga

Message:

some improvements in land surface model

File:

• palm/trunk/SOURCE/radiation_model_mod.f90 (modified) (15 diffs)

palm/trunk/SOURCE/radiation_model_mod.f90

@@ -25 +25 @@
 ! -----------------
 ! $Id$
+! Emissivity can now be set individually for each pixel.
+! Albedo type can be inferred from land surface model.
+! Added default albedo type for bare soil
+! 
+! 2318 2017-07-20 17:27:44Z suehring
 ! Get topography top index via Function call 
 ! 
@@ -211 +216 @@
 !
 !-- Predefined Land surface classes (albedo_type) after Briegleb (1992)
-    CHARACTER(37), DIMENSION(0:17), PARAMETER :: albedo_type_name = (/      &
+    CHARACTER(37), DIMENSION(0:18), PARAMETER :: albedo_type_name = (/      &
                                    'user defined                         ', & !  0 
                                    'ocean                                ', & !  1
@@ -229 +234 @@
                                    'sea ice                              ', & ! 15 
                                    'snow                                 ', & ! 16
-                                   'pavement/roads                       '  & ! 17
+                                   'pavement/roads                       ', & ! 17
+                                   'bare soil                            '  & ! 18
                                                          /)
 
-    INTEGER(iwp) :: albedo_type  = 5,    & !< Albedo surface type (default: short grassland)
-                    day,                 & !< current day of the year
-                    day_init     = 172,  & !< day of the year at model start (21/06)
-                    dots_rad     = 0       !< starting index for timeseries output 
+    INTEGER(iwp) :: albedo_type  = 9999999, & !< Albedo surface type
+                    day,                    & !< current day of the year
+                    day_init     = 172,     & !< day of the year at model start (21/06)
+                    dots_rad     = 0          !< starting index for timeseries output 
 
     LOGICAL ::  unscheduled_radiation_calls = .TRUE., & !< flag parameter indicating whether additional calls of the radiation code are allowed
@@ -261 +267 @@
                 decl_3,                          & !< declination coef. 3
                 dt_radiation = 0.0_wp,           & !< radiation model timestep
-                emissivity = 0.98_wp,            & !< NAMELIST surface emissivity
+                emissivity = 9999999.9_wp,       & !< NAMELIST surface emissivity
                 lambda = 0.0_wp,                 & !< longitude in degrees
                 lon = 0.0_wp,                    & !< longitude in radians
@@ -278 +284 @@
     REAL(wp), DIMENSION(:,:), ALLOCATABLE :: &
                 alpha,                       & !< surface broadband albedo (used for clear-sky scheme)
+                emis,                        & !< surface broadband emissitivity
                 rad_lw_out_change_0,         & !< change in LW out due to change in surface temperature
                 rad_net,                     & !< net radiation at the surface
@@ -285 +292 @@
 !-- Land surface albedos for solar zenith angle of 60° after Briegleb (1992)     
 !-- (shortwave, longwave, broadband):   sw,      lw,      bb,
-    REAL(wp), DIMENSION(0:2,1:17), PARAMETER :: albedo_pars = RESHAPE( (/& 
+    REAL(wp), DIMENSION(0:2,1:18), PARAMETER :: albedo_pars = RESHAPE( (/& 
                                    0.06_wp, 0.06_wp, 0.06_wp,            & !  1
                                    0.09_wp, 0.28_wp, 0.19_wp,            & !  2
@@ -302 +309 @@
                                    0.90_wp, 0.65_wp, 0.77_wp,            & ! 15
                                    0.95_wp, 0.70_wp, 0.82_wp,            & ! 16
-                                   0.08_wp, 0.08_wp, 0.08_wp             & ! 17
-                                 /), (/ 3, 17 /) )
+                                   0.08_wp, 0.08_wp, 0.08_wp,            & ! 17
+                                   0.17_wp, 0.17_wp, 0.17_wp             & ! 18
+                                 /), (/ 3, 18 /) )
 
     REAL(wp), DIMENSION(:,:,:), ALLOCATABLE, TARGET :: &
@@ -510 +518 @@
 !
 !-- Public variables and constants / NEEDS SORTING
-    PUBLIC decl_1, decl_2, decl_3, dots_rad, dt_radiation, force_radiation_call,&
+    PUBLIC albedo, albedo_type, decl_1, decl_2, decl_3, dots_rad, dt_radiation, emissivity, force_radiation_call,&
            lat, lon, rad_net, rad_net_av, radiation, radiation_scheme, rad_lw_in,        &
            rad_lw_in_av, rad_lw_out, rad_lw_out_av, rad_lw_out_change_0,       &
@@ -888 +896 @@
 
 !
+!--    Allocate array for storing emissivity
+       IF ( .NOT. ALLOCATED ( emis ) )  THEN
+          ALLOCATE ( emis(nysg:nyng,nxlg:nxrg) )
+          emis = emissivity
+       ENDIF
+
+!
 !--    Allocate array for storing the surface net radiation
        IF ( .NOT. ALLOCATED ( rad_net ) )  THEN
@@ -956 +971 @@
 !
 !--       Overwrite albedo if manually set in parameter file
-          IF ( albedo_type /= 0 .AND. albedo == 9999999.9_wp )  THEN
+          IF ( albedo_type /= 0  .AND.  albedo_type /= 9999999 .AND. albedo == 9999999.9_wp )  THEN
              albedo = albedo_pars(2,albedo_type)
           ENDIF
-    
+!
+!--       Write albedo to 2d array alpha to allow surface heterogeneities    
           alpha = albedo
   
@@ -1194 +1210 @@
              rad_sw_in(0,j,i)  = solar_constant * sky_trans * zenith(0)
              rad_sw_out(0,j,i) = alpha(j,i) * rad_sw_in(0,j,i)
-             rad_lw_out(0,j,i) = emissivity * sigma_sb * (pt(k,j,i) * exn)**4
+             rad_lw_out(0,j,i) = emis(j,i) * sigma_sb * (pt(k,j,i) * exn)**4
 
              IF ( cloud_physics )  THEN
@@ -1207 +1223 @@
 
 
-             rad_lw_out_change_0(j,i) = 3.0_wp * sigma_sb * emissivity         &
+             rad_lw_out_change_0(j,i) = 3.0_wp * sigma_sb * emis(j,i)         &
                                         * (pt(k,j,i) * exn) ** 3
 
@@ -1254 +1270 @@
              ENDIF
 
-             rad_lw_out(0,j,i) = emissivity * sigma_sb * (pt(k,j,i) * exn)**4
+             rad_lw_out(0,j,i) = emis(j,i) * sigma_sb * (pt(k,j,i) * exn)**4
 
              rad_sw_in(0,j,i) = ( rad_net(j,i) - rad_lw_in(0,j,i)              &
@@ -1533 +1549 @@
              rrtm_tsfc = pt(nzb,j,i) * (surface_pressure / 1000.0_wp )**0.286_wp
 
+!
+!--          Set surface emissivity
+             rrtm_emis = emis(j,i)
+
              IF ( lw_radiation )  THEN
                CALL rrtmg_lw( 1, nzt_rad      , rrtm_icld    , rrtm_idrv      ,&
@@ -1723 +1743 @@
                  rrtm_aldir(0,:,:) = aldif
                  rrtm_asdir(0,:,:) = asdif
+
+!
+!--        Bare soil
+           ELSEIF ( albedo_type == 18 )  THEN
+                 rrtm_aldir(0,:,:) = aldif
+                 rrtm_asdir(0,:,:) = asdif
+
 !
 !--        Land surfaces