High Resolution modelling of the Soil Heat Flux with
the inclusion of a Vegetation Layer
Bachelor thesis by Casper Borgman Earth Sciences
Universiteit van Amsterdam
Supervisors: Willem Boutem & Judy Shamoun-Baranes Course Coordinator: Kenneth Rijsdijk
2
Abstract
A clear link can be drawn between atmospheric processes and bird migration. Therefore, it is useful to model these processes in order to understand more about the conditions that determine a birds' behaviour. This study is part of a larger project to model the sensible heat flux on high temporal and spatial resolution. This study presents a soil heat flux model with the inclusion of a vegetation layer on a resolution of 100 by 100 meters spatially and 1 hour temporally. The soil heat flux can influence the sensible heat flux, which is a proxy for the availability of thermals that birds use for soaring flight. Knowing as accurately as possible how large this flux is over a variety of conditions can thus be useful. The inclusion of a
vegetation layer and soil moisture content, amongst other parameters, add to the accuracy of the modelled flux because of their relatively strong influence on it. The simulated soil heat flux has realistic values when compared to other studies. There are still some issues that need to be addressed, such as a more detailed aerodynamic resistance, but the model is a proof of concept and works well with the limited input requirements.
3
Content
List of figures and tables ... 4
1 Introduction ... 5
2 Methods ... 7
2.1 Soil heat flux model ... 7
2.1.1 Theoretical background ... 7
2.1.2 Boundary conditions and vegetation layer ... 7
2.1.3 Model setup ... 8
2.2 Data preparation and parameterization ... 9
2.2.1 Data loading and adjustment ... 9
2.2.2 Parameterization ...10
2.3 Data analysis ...12
3 Results ...13
3.1 Model dynamics and output ...13
3.1.1 Influence on sensible heat flux ...14
3.2 Model sensitivity ...15
3.3 Comparison with non-mechanistic G estimates ...16
4 Discussion ...17
4.1 Discussion of results ...17
4.1.1 Research aim and questions ...17
4.2 Model validity ...18
5 Conclusion ...19
6 Acknowledgements and Evaluation ...20
6.1 Word of thanks ...20
6.2 Evaluation ...20
7 References ...21
4
List of figures and tables
Table 1: The data input and output of the soil heat flux model ...8
Table 2: Mapping extend of the data...9
Table 3: Data source, spatial resolution, temporal resolution, data format...9
Table 4: Texture data adjustments...10
Figure 1: Thermal conductivity based on bulk density and moisture content ...11
Figure 2: Ground heat flux maps of the Netherlands...13
Figure 3: Net radiation and soil heat flux...13
Figure 4: Min, mean and max layer temperatures of a grid cell...14
Figure 5: Average soil heat flux per hour for the entirety of the Netherlands...14
Figure 6: Box plots of the soil heat flux for different texture classes...15
Figure 7: Box plots of soil heat flux under different vegetation covers...16
Figure 8: Comparison of modelled soil heat flux against soil heat flux estimated from net radiation using LAI...16
5
1
Introduction
Many birds use soaring to save energy during flight. To soar, birds use thermal updrafts that have a larger upward speed than the sink rate of the bird (Pennycuick, 2008). They circle upwards in a thermal and then glide until another thermal updraft is encountered
(Hedenström, 1993).Birds that use thermal soaring sometimes prefer powered flapping over soaring depending on local weather conditions. When good convective conditions are present, birds are more likely to use soaring over flapping (Shamoun-Baranes et al., 2016, Hedenström 1993). The use of soaring thus heavily relies on the availability of thermals, which are very localized in the atmosphere, and their strength (Ákos et al., 2010,
Pennycuick, 2008). By tracking soaring birds during long-distance migration, researchers have shown large seasonal differences in the migratory performance of such birds. In order to relate variability in the behaviour of soaring birds to variability in thermal updrafts and atmospheric processes throughout their migratory journeys (Shamoun-Baranes et al., 2010), however researches are in need of models that accurately predict atmospheric conditions over large spatial domains.
Convection is particularly dependent on the energy fluxes in the atmosphere. There is a constant exchange of energy between the earth’s surface and the atmosphere. Equation 1 describes the energy balance of the atmosphere:
[W] (eq. 1) Where r is the albedo, which is used with incoming shortwave radiation to
calculate absorbed solar radiation. is the incoming longwave radiation, is the outgoing longwave radiation. Furthermore, SHF is the sensible heat flux, λE is latent heat and G is soil heat flux. Rising air occurs mainly when the sensible heat flux (SHF) is relatively large. The SHF can be seen as a proxy for determining the occurrence of thermal updrafts. On the contrary, when the latent heat flux (λE) is relatively large, there is less sensible heat. As can be derived from equation (1), the energy flow in and out of the soil, (G) is also of importance to the sensible heat flux. The soil heat flux, on which this study focuses, is dependent on a temperature gradient in the soil. The equation is as follows (Snyder and Paulo de Melo-Abreu, 2005):
[W m-2] (eq. 2) Wherein –Ks is the soil thermal conductivity and the thermal gradient, in which T is temperature and z depth.
Current large scale models often produce outputs on a spatial and temporal
resolution that is not sufficient for further insight into bird behaviour. In particular, the soil heat flux is often barely calculated per hour and not often on a high spatial resolution. The soil heat flux is assumed negligible on a daily basis (Seguin and Itier, 1983), partly because extensive hydrological data is needed to produce reliable values of this heat flux (Kustas and Daughtry,1989), and such data are often not available. However, in reality the soil heat flux will cause substantial variation in SHF and thermal availability for birds within days on a hour by hour basis as a large part of the net radiation can influence the soil (Kustas and Daughtry, 1989). During the morning the cool soil will mainly absorb energy, whereas at the end of the day the soil emits energy. This can shift the availability of thermals potentially by several hours. Furthermore, vegetation influences the soil heat flux as well. Under vegetative cover the soil heat flux can be nearly three times smaller than the flux of a bare soil (Monteith, 1973). Thus, the soil heat flux changes quite substantially during the day and over distance
6 (Choudhury et al., 1987). In addition, thermals are very localized and ephemeral, thus birds make in-flight decisions on a small timescale. Therefore, the modelling of soil heat flux on a good spatial and temporal resolution is important for understanding the behaviour of soaring birds.
This study focused on the creation of a soil heat flux model that uses existing hydrological and atmospheric data to get output on a better spatial and temporal resolution. The soil heat flux model is part of a sensible heat flux model that was created with
colleagues as a joint effort. The aim of the joint effort of this research is to create a model to produce sensible heat flux at high spatial resolution of 100 by 100 meters and temporal resolution of 1 hour. Accordingly, the soil heat flux model has the same resolutions. The main aim of this research was to create the model with the inclusion of a vegetation layer and to assess the variability of the soil heat flux and the sensitivity to different
circumstances. Furthermore, a brief look was taken at the potential influence of soil heat flux on the sensible heat flux. The main questions of this research were:
1. Can existing data from multiple sources be integrated to create a reliable soil heat flux model at one hour resolution on a spatial scale of 100 by 100 meters?
2. How does the soil heat flux vary as the parameters it is dependent on change under different circumstances?
3. Does the addition of a vegetation layer allow for more accurate modelling of the soil heat flux?
This report is built up as follows. First the methods of this research are given. Herein the model construction, data processing and analysis are discussed. Thereafter, the results are given followed by a discussion. Lastly, the conclusions of this research are presented.
7
2
Methods
2.1 Soil heat flux model
The model is built in Matlab according to the model setup recommended by Bouten (2015). Therefore, the model consists out of an initialization, a dynamic part and data storage and visualization. The model can be called as a function of the sensible heat flux model. This segment will describe how the model was constructed for this research.
2.1.1 Theoretical background
The primary functioning of the model is based on equation (2), which was given in the introduction. The soil heat flux is a flux of sensible heat and moves by conduction. It is dependent on soil moisture content, land cover and solar radiation, amongst others (Snyder and Pawu, 2011). The equation for heat flow through the soil is (Bouten, 2016):
[J m
-3 day-1] (eq. 3)
In which is the volumetric heat capacity, z is height in the soil column, and t is time. The soil can be divided into layers for the purpose of modelling. Each layer of the soil has a heat content, the amount of energy stored in a specified partition. The heat content per layer at a given time and layer is calculated as follows (Bouten, 2016):
[J m-3] (eq. 4) Wherein is the heat content. Temperature changes according to changes in the heat content. The energy flow in between layers results into a removal or addition of energy to a layer or the surface. This energy is then added or subtracted from the heat content. This results in a new temperature for the layer.
2.1.2 Boundary conditions and vegetation layer
The model requires boundary conditions as input to run. Net radiation is taken for the warming and cooling of the soil. Without this boundary condition the soil temperature would be in equilibrium, because no temperature influence would occur. However, not all radiation reaches the soil when vegetation is present. To account for that a fraction is introduced based on the Leaf Area Index, which is determined during the model initialisation. The formula used is (Yang et al., 1999):
[ ] (eq. 5) Where cf is the cover fraction, meaning the fraction that determines how much
radiation reaches the soil through vegetation. is the Leaf Area Index, and Zenith is the average solar zenith angle. The LAI has been calculated in the parameterization section of the model. The average solar zenith angle has been taken as 0.6 from Yang et al. (1999), because this allows for diffuse and direct solar radiation.
A layer of vegetation can also influence energy transfer. The temperature of a vegetation layer is different than the air temperature above it. This is due to aerodynamic resistance to heat transfer. The effect of aerodynamic resistance can be made clear by
8 introducing the vegetation layer as boundary condition. The temperature of the vegetation layer influences the soil temperature. For this purpose the following equation is used in the model (Zuang et al., 2015):
[W m
-2] (eq. 6)
Where is the mean air density at constant pressure (kg m3), Cp is the specific heat capacity of air, is the temperature of air, is the temperature of the surface and is the aerodynamic resistance.
Finally, it is assumed that energy flow does not continue endlessly into the soil below the lowest layer. The flow between the lowest layer and all soil beneath it is set to zero. This is the boundary at the bottom of the model.
2.1.3 Model setup
The main section of the model is dynamic, meaning the variables used in this section change over time. The formulae and processes described in section 2.1.1 and 2.1.2 are part of the dynamic part. The conductivity and heat capacity, which are dependent on texture and water content, are included in the dynamic section as well, but will be described in further detail in section 2.2.2
The model is three-dimensional; calculations are made for grid cells with coordinates x and y, over soil layers with depth z. The
model has six soil layers, because the recommended amount of soil layers while calculating the soil heat flux in the energy balance of the atmosphere is six (Benoit, 1976). The layers are 30 cm thick to prevent model oscillation. Oscillation occurs when changes are large, causing the flow to oscillate in each model iteration. The vegetation layer and surface layer calculations are done at the layer boundary, but the soil layer calculations
are performed for the whole layer. The model iterations are done with discrete steps of 1 hour, which will result in a 1 hour resolution.
Table 1 shows the data input needed to run the model and the output data of the model.
Data Input Model Output
Texture Soil heat flux
Net Radiation Temperature per layer
Leaf Area Index Flux per layer
Time Saturation
Saturation Deficiency Rootzone depth
Land cover
9
2.2 Data preparation and parameterization
The input data needs to be adjusted to the right extent and size before the model can run. This section describes which data and how data were processed for application in the model. The model for this study is run for the Netherlands, due to data availability. The model is run for the month of May in 2012.
2.2.1 Data loading and adjustment
The ESDAC texture map and Corine Land Cover map were retrieved and loaded into ArcGIS. Both maps were then clipped for the right extent of the Netherlands, as specified in table 2, for each wind direction. Both the texture and land cover map use the ETRS89 LAEA projection. The land cover and texture rasters were
transformed with ArcGIS to ascii text files for easy handling within Matlab. Table 3 shows the data used for the model and the data sources as well as resolution.
Direction Extent
Top 54° N
Bottom 50° N
Left 3° E
Right 8° E
Source Data provided Spatial
resolution Temporal resolution Data format ERA-interim (reanalysis) by ECMWF -Temperature 0.0125 by 0.0125 (degrees) 3-hourly .nc (netcdf) NHI (model) by Deltares - Soil moisture content -Rootzone depth -Soil moisture deficiency
1.0 by 1.0 (km) Daily .asc (ASCII)
Copernicus land monitoring services - Land-cover 100 by 100 (m) - .csv The Satellite Application Facility on Climate Monitoring (CM SAF) -Solar radiation 0.05 x 0.05 (degrees) 1-hourly .nc (netcdf)
European Soil Data Centre
-Soil texture 1000x1000 (m) Published
2006
.adf
Table 2. The clipping extend for the Netherlands
10 The data is then loaded into the Matlab workspace in an initialisation script. This script adjusts the data to the preferred size and adjusts data to remove false values and errors. In addition the net radiation is altered to J m-2 h-1 instead of W m-2 for calculations over time. Thereafter, a reference matrix is chosen to adjust the data to the same size. The
reference dataset is the loaded Corine Land Cover matrix. This dataset already has the target extent as well as spatial resolution. Much of the data has different sizes that need to be homogenized, because it needs to be the same size before any calculations in the model can happen. The texture data, temperature at 2 meters, root zone, root zone saturation, saturation deficiency, aerodynamic resistance and net radiation are then adjusted with interpolation to the same matrix size. The interpolation method used is the k-nearest neighbour algorithm. This method assigns the new object the integer value of the nearest neighbouring cell. Nearest neighbour interpolation prevents the assignment of incorrect values, such as unwanted float values for data classes.
The texture data needs to be adjusted further before it is ready to use. The data has six different classes (table 4). These classes are somewhat imprecise and no organic matter description is given, but the data is available for the whole extent of Europe. Based on the values given by the ESDAC new texture fraction values were assigned for use in the model, including soil organic matter, because organic matter lowers the conductivity of a soil (Abu-Hamdeh and Reeder, 2000). Furthermore , the 'no information' class is assumed to be completely covered by urban areas, otherwise there are large gaps in the data and no heat flux could be calculated.
Class Description Texture Assigned fraction
-9999
Water None 0
0 No information & Urban
area
None 1
1 Coarse <18% clay, >65%
sand
18% mineral, 77% sand, 5% organic matter
2 Medium 18-35% clay, >15%
sand
35% mineral, 59% sand, 6% organic matter
3 Medium fine <35% clay, <15%
sand
64% mineral, 30% sand, 6% organic matter
4 Fine 35% < clay < 60% 79% mineral, 14% sand, & organic
matter
5 Very fine Clay > 60% 90% mineral,
9 Peat None 45% mineral, 10% sand, 40% organic
matter
2.2.2 Parameterization
The thermal conductivity (Ks) and volumetric heat capacity of a soil are dependent on the materials the soil consists out of. Quartz, mineral soil and organic matter all conduct thermal energy differently. The pores in the soil itself are negligible, as air does not conduct very well at standard atmospheric pressure of 1 bar. However, water has a large influence on the thermal conductivity of a soil. The water fills the pores and increases the amount of contact in the soil layer, increasing conductivity (Hillel, 1980).
11 The thermal conductivity of the soil has been determined using the Johansen method (Johansen, 1977) as described in Balland and Arp (2005). The method will be described here briefly, but all formulae used are available in
Appendix A. Figure 1 shows the conductivity similar to the conductivity produced using this method. Balland and Arp (2005) propose a new method, however this method is more elaborate and would require more data and calculations and thus more processing time for a slight increase in accuracy. Therefore, to limit process time, the Johansen method for unfrozen soils was taken.
The Johansen method determines the conductivity of soils based on the saturation of the soil. The primary equation of this method is as follows (Balland and Arp, 2008):
[J m-1°C-1 day-1] (eq.7) Where is the thermal conductivity at full soil saturation, the conductivity of a dry soil, and the Kersten number. The Kersten number describes how much the
conductivity changes based on the saturation. The saturation is data loaded in the model and changes each model iteration, thus the conductivity also changes. The Kersten number is calculated using equation:
[ ] (eq.8)
Where is the degree of saturation given by
[ ] (eq.9)
In which is the volume of water and the volume of pores.
The volumetric heat capacity is a linear function of the water content (Appendix A). The heat capacity increases when more water is present as air, which has a low heat capacity, is replaced by water with a high volumetric heat capacity (Koorevaar et al., 1983).
The conductivities and volumetric heat capacity per assigned texture fraction were taken as recommended by Balland & Arp (2008) from Hillel (1982), De Vries (1963) and Alter (1969). The urban areas have been assigned a conductivity and heat capacity of concrete and asphalt. Open water was not assigned any value of conductivity.
The aerodynamic resistance was taken from Goudriaan (1977) with values of 50,100 and 500 respectively for grass, crops and forests. The parameterization has been done based on the Corine Land Cover map classification. The grassland aerodynamic resistance is assumed equal to the bare soil aerodynamic resistance to heat transfer for modelling purposes.
Figure 1. Soil thermal conductivity based on water content and bulk density. Source: Balland and Arp(2005)
12
2.3 Data analysis
The primary data to be analysed are the model outputs. As table 1 shows, these are the soil heat flux and the layer temperatures. The soil heat flux is positive during the daytime, meaning that energy goes into the soil. A negative soil heat flux means a release of energy, which happens when solar radiation is low or absent. The flux is visualised for a time during the day and a time during the night. It is also plotted along with the net radiation to show the relationship between both. In addition, it is possible to look at the hourly average values of the flux to derive a possible influence on the sensible heat flux. This is done by averaging the flux per hour for the Netherlands during the month of May 2012.
The heat flux and temperature can be compared to a study done by Jacobs et al. (2011) on a grassland area in the Netherlands. For this purpose a grid cell with a grass layer was taken. Jacobs et al. (2011) give the average soil heat flux per month as well as mean, minimum and maximum soil temperatures averaged over a year. Therefore, the mean minimum and maximum temperatures are visualised and the mean flux per month is calculated to assess similarity.
In addition, it is interesting to visualise and store the changes of the soil heat flux and soil temperatures over a range of conditions, such as different aerodynamic resistances and wet and dry days. Due to the nature of the model nearly all of the parameters should show a clear relation with the soil heat flux. The input parameters all influence the soil heat flux one way or another because the flux is directly dependent on the data inputs. The impact of different parameters is visualised using box plots. These show the soil heat flux according to assigned values of parameters, such as the texture classes. The median, quartile range and whiskers of the plot will show details about the sensitivity to parameters.
There are no data sets available against which the model can be compared. Testing the validity of the model can therefore be difficult. Another way to compare the results with reality is testing them against a non-mechanistic predictive formulae. Idso et al. (1975) found that the soil heat flux can be estimated as a ratio of the net radiation, with a ratio of 0.5 for dry soils and 0.3 for wet soils. Monteith (1973) suggests that the ratio is between 0.05 and 0.1 for fully vegetation covered soils. Another predictive formula is based on the Leaf Area Index, this formula is proposed by Choudhury et al. (1987). This formula is created using a linear regression between and LAI. The resulting formula is:
[W m-2] (eq.10) Where Rn is the net radiation. Thus, this formula is based solely on leaf area index and the relationship between soil heat flux and net radiation. In the case of Choudhury et al. (1987), the goodness of fit between measured and predicted soil heat flux is quite
reasonable (r=0.93). The modelled soil heat flux is compared to this formula. Furthermore, a comparison is made with the ratios formulated by Idso et al. (1975) and Monteith (1973).
13
3
Results
3.1 Model dynamics and output
The soil heat flux is the main output of the model. Figure 2 shows the soil heat flux map modelled for the Netherlands on the 2nd of May, 2012. This figure shows the difference between day and night. There is a clear distinction between the times of the day. In addition, the difference between vegetation is visible. The darker shaded grid cells on the left and brighter coloured grid cells on the right are the areas which experience the largest soil heat flux and the largest diurnal differences in the soil heat flux. These are areas with bare soil and grasslands. Furthermore, the Veluwe area is a distinct mostly forested area and is clearly visible on both maps in the middle surrounding cells x650,y450.
However, zooming in on a grid cell can explain model behaviour better than looking at an overview map, because the properties of grid cells are known. A visualisation of the soil heat flux versus the net radiation for the 1st through 10th of May in a random grid cell with crops can be seen in figure 3. This figure confirms that the soil heat flux follows the net radiation throughout the day. Between 80 and 150 hours, the soil heat flux forms a large part of the net radiation. During some hours, such as around hours 215-220, the soil heat flux is larger than the net radiation. This can
be due to modelling inaccuracies and boundary conditions.
The temperature mean, max and min for the month of may in 2012 is shown in figure 4. The temperature gradient of the soil layer temperature decreases with depth. Therefore there are also smaller fluxes between the lower layers. During the month of May the soil becomes warmer than at the start of the month, resulting in the shape of the mean in the figure.
Figure 2. Ground heat flux maps for the Netherlands. 2nd of may 2012. Left: 1:00 AM Right: 1 PM
14 3.1.1 Influence on sensible heat flux
Figure 5 shows the modelled mean heat flux per hour for all grid cells of the Netherlands combined for the month of May 2012. From this figure the average influence of soil heat flux on the strength of the sensible heat flux can be derived. During the morning, on average, the soil releases energy until roughly 8:00 am. The soil gains energy during the day until 8:00 pm, where after it starts releasing energy again during the night. In addition the means of the flux can be calculated to get an indication whether the soil is gaining net energy or losing net energy during a month. The mean of a grassland area is 2.71 W m-2 for the month of May. The mean for that same period for cropland is 1.48 W m-2. For forests the mean is -0.08 W m-2. This means that during this time the forest has a slight loss of energy.
However, as this is the average it does not represent a real location. In reality there will be differences in vegetation coverage, soil texture and water content per grid cell. The next section describes the model sensitivity to these differences.
Figure 4. Min, mean and max temperatures for the soil, with depth.
Figure 5. Average soil heat flux per hour for the entirety of the Netherlands
15
3.2 Model sensitivity
The soil moisture content has been brought up as potentially important to the soil heat flux. In figure 3, the fifth and sixth peaks of the soil heat flux are during relatively cloudy days. The last three peaks are during rainy days. During these last days, the thermal conductivity is higher than on dry days, because the soil moisture content increases after a rainfall event. The soil heat flux is notably lower on the wet days compared to the cloudy days. This is especially notable since the net radiation is larger on the rainy days compared to the cloudy days. This can be attributed to the greater thermal conductivity due to increased moisture content.
The thermal conductivity also depends on the chemical composition of the soil. The influence of each of the classes from table 4 on the soil heat flux can be seen in figure 6. The 'very fine' class is underrepresented and therefore excluded. This box plot shows that texture of a soil can incur slight changes in the soil heat flux. The coarse, medium and medium fine classes show similar median, quartile range and whisker range. There fine and coarse class show slight differences however. Yet, the soil texture does not necessarily cause the
difference in the soil heat flux. Different soils support different types of vegetation. In the Netherlands many grasslands are for example on clayey soils, which might explain some of the difference shown in this box plot. Therefore, the role of vegetation will be examined next.
The model sensitivity to vegetation can be looked at through the aerodynamic
resistance classes. The LAI, thus vegetation cover fraction, is automatically accounted for, as cells with high aerodynamic resistance also have a high LAI. Figure 7 shows box plots of the soil heat flux for different aerodynamic resistances from Goudriaan (1977). There is a large difference in soil heat visible, meaning that the vegetation layer has a large influence on the soil heat flux. Therefore, the soil heat flux is mainly dependent on the type of vegetation cover. The forest layer thus noticeably stops a fraction of radiation reaching the soil and dampens the temperature gradient that causes energy flow.
16
3.3 Comparison with non-mechanistic G estimates
In the methods section the formula to estimate G based on the LAI was described. Figure 8 shows the modelled soil heat flux against the estimated soil heat flux for a crop grid cell. The modelled grid cell shows a flux in the same order of magnitude as the estimated flux. There are some large deviations, however this is not unexpected. The formula used by Choudhury et al. (1987) only uses the LAI parameter, whereas the model uses multiple parameters to calculate G. The estimated and modelled G for all grid cells also shows the correct order of magnitude for multiple grid cells.
From figure 3, figure 5 and a comparison calculation can be derived that the flux under crop cover is roughly 15-20% of the net radiation. This relates to the relationship between vegetation and soil heat flux described by Monteith (1973), who argues that the soil heat flux under full vegetation cover is 5 to 10% of net radiation, and by Idso et al. (1975), who argue that it is 30 to 50% of net radiation on bare soil. The flux under forest cover, as can be verified by looking at figure 5, does not often exceed 40 W m-2, whereas the net radiation is often around 300 W m-2. The order of magnitude of the fluxes thus seem to align with their predicted fluxes.
Figure 7. Box plots of soil heat flux under different vegetation covers. Rh 50= grass/bare soil. Rh 100= crop cover. Rh 500 = forest cover
Figure 8. Comparison of modelled soil heat flux against soil heat flux estimated from net radiation using LAI.
17
4
Discussion
4.1 Discussion of results
The results show that the modelled soil heat flux is in the right order of magnitude and most likely quite accurate, and the model dynamics behave as expected. The modelled map of the Netherlands shows the spatial differences of the soil heat flux. The relationship between the net radiation and soil heat flux is mostly logical, however it should be remarked that the inclusion of a vegetation boundary layer changes the direct dependency of soil heat flux on the net radiation. The mean flux of particularly the grassland is of interest for comparison with the study done by Jacobs et al. (2011). When compared, the mean flux is relatively the same as what Jacobs found for a grassland in the Netherland, namely roughly 5 W m-2 over 27 years, whereas the modelled mean for May 2012 was 2.71 W m-2. However, this is still quite a large difference. This difference can arise due to different weather conditions, but also due to parameterization. Their study found deviations of roughly 2 W m-2 per month between their modelled mean soil heat flux and measured heat flux.
The comparison with the non-mechanistic estimates show that the modelled flux agrees with relationships described in the literature. There are some slight deviations, but that is to be expected as there are differences in parameters per grid cell. For example: the soil saturation, texture and vegetation can differ over distance. When using the formula by Choudhury et al. (1987) the model simulation also shows realistic values, even though the LAI is different for each type of land use, while the aerodynamic resistance is only different per groups of land use.
4.1.1 Research aim and questions
The aim of this research was to model the soil heat flux accurately with a high spatial and temporal resolution. The first question pertained the construction of a reliable soil heat flux model that includes a vegetation layer. The results prove that reliable calculation of the soil heat flux is possible through this model. When compared to non-mechanical estimation the modelled flux is in the right order of magnitude. In addition, comparison of mean flux with a study and comparison of the flux to the net radiation shows reasonable values. The model can be used to determine the influence of soil heat flux on the sensible heat flux by looking at hourly values of grid cells, which each have different characteristics.
The second research question pertained the sensitivity of the model to these different characteristics. The vegetation seems to have the largest influence on the soil heat flux. The vegetation blocks direct sunlight to the soil and the vegetation layer has different temperature dynamics than the air above it. The soil moisture has the second largest influence on the soil heat flux. The third important parameter, texture, does not matter so much as the vegetation that is above it.
Therefore, the answer to the third research question becomes clear. The addition of a vegetation layer allows for more accurate and realistic modelling of the soil heat flux
compared to modelling the flux without it. However, there are some remarks to be made about the current implementation of the vegetation layer, which will be discussed in the following section along with other remarks.
18
4.2 Model validity
The aerodynamic resistance Rh that is used to calculate the heat flux in the vegetation layer has only been divided into three broad categories taken from Goudriaan(1977). The
assigned values of Rh for grass, crops and forests will differ per vegetation type. The
equation for Rh in Appendix A shows the many parameters on which it is dependent(Yang et al., 2015). Due to time restraints and the scope of this research, a dynamic aerodynamic resistance is not implemented. Since the results show a large influence of the aerodynamic resistance, this will need to be resolved for more accurate results. However, the results show that the order of magnitude of the flux can be considered valid.
Secondly, the boundary condition of net radiation is missing the latent heat in the soil. A part of the net radiation will be used for latent heat in the soil. Especially when soil is wet this is expected to have implications. Currently, the conductivity and heat capacity change when water is added, but the latent heat flux in the soil is unaccounted for. The errors that arise can be quite large, up to 28 W m-2 (Mayoocchi and Bristow, 1995). Surfaces prone to high evaporation rates can even show errors up to 100 W m-2 (Buchan, 1989).
The volumetric heat capacity currently changes per model iteration based on how much water is present. However, due to this, energy can get lost or be gained when the heat content of the soil is updated. The implications of this are small, because the soil saturation changes daily. Moreover, the thermal conductivity right now is only correct for unfrozen soils. Balland and Arp (2005) describe a method for frozen soils, however this is not implemented in the model.
Another remark to be made is the effect of data input resolution. For example, the texture data has a resolution of one by one kilometer. Adjusting this data to 100 by 100 meter does only artificially increase the resolution, but not the real resolution. The same applies to data of temporal resolution of more than one hour. When this data is interpolated to a one hour timescale, the resolution is only artificially increased. Table 3 shows the datasets which experience this problem with spatial and temporal resolution. Land use, arguably the most important parameter of the model, is however on the wanted 100 by 100 meter resolution. The solar radiation, which is the main driving force behind the soil heat flux, is also available at the best resolution needed for this model. Therefore, the model resolutions of 1 hour and 100 by 100 meters still have added value over larger resolutions.
In addition, the ESDAC soil texture data map only has seven classes. Although water saturation has a larger influence than the textural composition of a soil, the data could still be more precise. The ESDAC data also has a 1km by 1km resolution, which means small pastures that could be potentially of importance to the sensible heat flux might not be recognized. Urban areas are also misrepresented in the model. There could be no correct aerodynamic resistance allocated nor does the thermal conductivity assume correct values. Urban areas are complex to model and could not be modelled effectively in this research.
The last remark is that the data used is from multiple sources. Therefore slight deviations can occur. Each provider uses different methods to get their data. The data they provide can. For example, for this research boundary conditions of net radiation and
temperature at 2 meters were used. The temperature at 2 meters is in part related to the net radiation. However, both are from different sources with a different resolution. Therefore, there might be small errors that arise. Furthermore, the edges of land close to bodies of water show slight errors, where at places data can be missing due to different mapping techniques and points.
19
5
Conclusion
In conclusion, the produced soil heat flux model works reasonably well and it is quite accurate using the available data. The inclusion of the vegetation layer seems to be useful but imprecise as of right now. The order of magnitude of the flux for vegetated surfaces seems to be correct however. The flux is sensitive to wet and dry circumstances as well. Not only through reduced net radiation on rainy days but also through thermal conductivity. It is not so much sensitive to textural changes in the soil.
The model is a proof of concept and can be expanded to be used for soil heat flux modelling for Europe, given that datasets are available. The model minimally requires the data as said in the methods and should work for a wide range of conditions, except for frozen soils. The ability of the model to produce accurate results can be enhanced by addressing the points brought forward in the discussion. It is imperative to expand on the aerodynamic resistance parameter.
All in all, the model achieves to calculate the soil heat flux on the resolution that was aimed for during this research. The 100 by 100 meter resolution on 1 hour timescale does add increased detail over larger scale models. The model is suitable to be used within the sensible heat flux model which in turn can be used for the analysis of soaring bird behaviour.
20
6
Acknowledgements and Evaluation
6.1 Word of thanks
I want to thank Maarten Mol, Bart Sweerts, and Daniel Kooij for the excellent teamwork on the creation of the sensible heat flux model. They have provided great feedback and work on this project. I also want to thank prof. dr. ir. Willem Bouten for the guidance and assistance during this project. In addition, I want to thank dr. Judy Shamoun-Baranes, Wouter
Vansteelant and dr. Adriaan Dokter for their enthusiasm and interest in this project.
Furthermore, thanks goes out to Hidde Leijnse from the KNMI, Rens van Beek from Utrecht University who gave access PCR-GLOBWB model data. Without the data the sensible heat flux model could not have been tested nor run at all. Furthermore, all organisations and institutions not yet named that have assisted with this project one way or another have my thanks.
6.2 Evaluation
The project involved quite a lot of teamwork, especially at the beginning but throughout the project as well. I had frequent meeting with Maarten, Bart and Daniel to discuss the project and how we would tackle problems such as data acquisition and model construction. In addition, as group we had on average one meeting per week with Willem Bouten, our supervisor. These meetings were often very useful, mostly because Willem urged us to prepare ahead of them. All problems that I faced could be addressed and discussed. The teamwork was pleasant as well. Often times we would work individually in the same working space at the university, so that whenever one of us had problems, we could deal with it as a team.
One of the main problems with programming and modelling are bugs in the code. These bugs can take many hours to solve at times, especially if there is no error message in the programming software. For example, I made a mistake with matrix multiplication.
Therefore, the model did not do the calculation I wanted correctly. However, this does not show up as an error and requires evaluation of every step in the model to find the mistake.
The workload was manageable, however due to another course the first month and a half of the project were sometimes quite stressful. All in all, I didn’t meet the expectations of the time planning of my proposal. Nearly everything took more time than expected. However, this was always seen as a possibility and therefore no major setback of this project was experienced.
What I have learnt is to never underestimate the writing of a program or model. Nearly every line of code has a thought out purpose and has a functioning in the greater whole of the model. One has to be mindful of this when writing code and has to be wary of bugs and minor faults. Another lesson learnt is that problems and issues do occur on projects like these, and that along the way your perception of what you are doing might change from your initial thoughts.
21
7
References
Abu-Hamdeh, N. H., & Reeder, R. C. (2000). Soil thermal conductivity effects of density, moisture, salt concentration, and organic matter. Soil Science Society of America Journal, 64(4), 1285-1290.
Ákos, Z., Nagy, M., Leven, S., Viscek, T. (2010). Thermal soaring flight of birds and unmanned aerial vehicles. Bioinspiration & Biomimetics, 5(4).
Alter, A. J. (1969). Water supply in cold regions. Cold Regions Research Engineering Laboratory Monograph, III-C5a.
ArcGIS Desktop (2012). Release 10.1. ESRI Redlands, CA: Environmental Systems Research.
Balland, V., & Arp, P. A. (2005). Modeling soil thermal conductivities over a wide range of conditions. Journal of Environmental Engineering and Science,4(6), 549-558.
Benoit, R. (1976). A Comprehensive Parameterization of the Atmospheric Boundary Layer for General Circulation Models. Dissertation Abstracts International, 37(11).
Bouten, W. (2016) Modelling and Simulating Spatial Processes in Landscape Ecosystems. Course manual for Modelling and Simulation, BSc major Earth Sciences, Universiteit van Amsterdam. Buchan, G. D. (1989). Soil heat flux and soil surface energy balance: a clarification of concepts. Proceedings of the Fourth Australasian Conference on Heat and Mass Transfer. University of Caterbury, Christchurch, New Zealand.
European Centre for Medium-Range Weather Forecasts (2016). ERA-Interim global reanalysis model. European Soil Data Centre (2004). European Soil Database v2.0 (vector and attribute data). European
Comission – DG JRC.
Choudhury, B. J., Idso, S. B., & Reginato, R. J. (1987). Analysis of an empirical model for soil heat flux under a growing wheat crop for estimating evaporation by an infrared-temperature based energy balance equation. Agricultural and Forest Meteorology, 39(4), 283-297.
Copernicus Land Monitoring Services (2012). Corine Land Cover vector data. Copernicus Programme,
European Commission.
Deltares (2016). Het Nederlands Hydrologisch Instrumentarium Landelijk Hydrologisch Model. Goudriaan, J. (1977). Crop micrometeorology: a simulation study (No. 551.5 G688). Centre for Agricultural Publishing and Documentation, Wageningen.
Hedenström, A. (1993). Migration by soaring or flapping flight in birds: the relative importance of energy cost and speed. Philosophical Transactions of the Royal Society of London B: Biological
Sciences, 342(1302), 353-361.
Hillel, D. (1980) Fundamentals of Soil Physics. Academic press.
Idso, S. B., Aase, J. K., & Jackson, R. D. (1975). Net radiation—soil heat flux relations as influenced by soil water content variations. Boundary-Layer Meteorology, 9(1), 113-122.
Jacobs, A. F.G., Heusinkveld, B.G., Holtslag, A.A.M. (2011) Long-term record and analysis of soil temperatures and soil heat fluxes in a grassland area, The Netherlands. Agricultural and Forest Meteorology, 151(7), 774-780.
Johansen, O. (1977). Thermal conductivity of soils (No. CRREL-TL-637). Cold Regions Research and Engineering Lab Hanover NH.
22 Koorevaar, P., Menelik, G., & Dirksen, C. (1983). Elements of soil physics (Vol. 13). Amsterdam: Elsevier.
Kustas, W. P., & Daughtry, C. S. (1990). Estimation of the soil heat flux/net radiation ratio from spectral data. Agricultural and Forest Meteorology, 49(3), 205-223.
MATLAB (2014). version 7.10.0 (R2014a). Natick, Massachusetts: The MathWorks Inc.
Mayocchi, C. L., & Bristow, K. L. (1995). Soil surface heat flux: some general questions and comments on measurements. Agricultural and forest meteorology, 75(1), 43-50.
Monteith, J. L., 1973. Principles of Environmental Physics. New York: Elsevier. 241 pp. Pennycuick, C. J. (2008). Modelling the flying bird (Vol. 5). Pennycuick: Elsevier.
Seguin, B., & Itier, B. (1983). Using midday surface temperature to estimate daily evaporation from satellite thermal IR data. International Journal of Remote Sensing, 4(2), 371-383.
Shamoun-Baranes, J., Bouten, W., & van Loon, E. E. (2010). Integrating meteorology into research on migration. Integrative and Comparative Biology, icq011.
Shamoun-Baranes, J., Bouten, W., van Loon, E., Meijer, C., Camphuysen, C. (2016). Flap or soar? How a flight generalist responds to its aerial environment. Philosophical Transactions of the Royal
Society B. * Not published
Snyder, R.L., Paulo de Melo-Abreu, J. (2005). Frost Protection: fundamentals, practise, and economics. Food and Agriculture Organization of the United Nations, Rome.
Snyder, R.L., Pawu, K.T. (2001) Soil Heat Flow and Temperature. University of California. The Satellite Application Facility on Climate Monitoring (2016) Solar radiation. EUMETSAT. de Vries, D. A. (1963). Thermal properties of soils. Physics of plant environment.
Van Wijk, W. R. (1963). Physics of plant environment. Physics of Plant Environment.
Yang, Z-L., Dai, Y., Dickinson, R.E., Shuttleworth, W. J. (1999) Sensitivity of ground heat flux to vegetation cover fraction and leaf area index. Journal of Geophysical Research. 104, 505-514. Zuang, Q., Wu, B., Yan, N., Zhu, W., Xing, Q. (2015) A method for sensible heat flux model parameterization based on radiametric surface temperature and environmental factors without involving the parameter KB-1. International Journal of Applied Earth Observation and Geoinformation, 42, 50-69.
23
8
Appendix A: Formulae
Energy balance of the earth:
[W] = Incoming shortwave radiation
= Incoming longwave radiation = Outgoing longwave radiation SHF= Sensible heat flux
= Latent heat = Soil heat flux
Soil heat flux(Snyder and Paulo de Melo-Abreu, 2005): [W m-2]
= Thermal conductivity T= Temperature
z = depth
Heat exchange in the soil(Bouten, 2015):
[J m-3 day-1] C= Volumetric heat capacity
Heat content of a soil layer(Bouten, 2015):
[J m-3]
Vegetation cover fraction(Yang et al. 1999):
[ ] cf = cover fraction
LAI= Leaf Area Index
24 Sensible heat flux of vegetation layer(Zuang et al. 2015):
[W m
-2]
= mean air density at constant pressure(kg m3 ), Cp =specific heat capacity of air,
= the temperature of air, = the temperature of the surface
Aerodynamic resistance(Zuang et al. 2015): [s -1 m]
= roughness length for momentum transport = reference height
d = zero plane displacement height = roughness length for heat k= von Karman’s constant u= friction velocity
(h)= stability correction functions for heat
Complete Johansen(1977) method for thermal conductivity Balland and Arp(2005): [J m-1°C-1 day-1] [] [] [J m -1°C-1 day-1] [J m-1°C-1 day-1] [J m-1°C-1 day-1]
Ke= Kersten number = saturation fraction = bulk density
volumetric heat capacity(Hillel, 1982):
[J °C-1 m-3]
