*2.1. Data*

The hourly integrated land-atmosphere interaction observation dataset from the National Tibetan Plateau Science Data Center (2005–2016) (http://data.tpdc.ac.cn/zh-hans/ data/b9ab35b2-81fb-4330-925f-4d9860ac47c3, accessed on 15 March 2020) [27,28] is employed as observations in this study, including meteorological gradient data, radiation data, and soil and turbulent flow data. The directly observed SH in this study comes from the turbulent flow data during the period 2006–2016, sampled by the open-path Eddy Covariance (EC) turbulent flux measurement system consisting of an ultrasonic anemometer and an infrared gas analyzer [28]. This dataset contains 6 stations (Figure 1 and Table 1), which are the Muztagh Ata Westerly Observation and Research Station, the Chinese Academy of Sciences (MAWORS), the Ngari Desert Observation and Research Station (NADORS), the Nagqu Station of Plateau Climate and Environment (BJ), the Nam Co Monitoring and Research Station for Multisphere Interactions (NAMORS), the Qomolangma Atmospheric and Environmental Observation and Research Station (QOMS), and the Southeast Tibet Observation and Research Station for the Alpine Environment (SETORS). Additionally, wind speed, air temperature, and surface temperature, derived from the meteorological gradient data in this dataset at these six stations, are also used in this study.

**Figure 1.** Geographical distributions of the six meteorological stations over the TP and the time series of the annual mean SH (units: W m−2) at these six stations during the period 2006–2016. The blue dotted lines indicate the decreasing and increasing trends at MAWORS and SETORS, respectively.

**Table 1.** List of the stations in the hourly integrated observational dataset, including the station names, latitudes, longitudes, time periods, and the time difference between local time and Beijing time.


*2.2. Methods*

For the calculation of the SH, the bulk transfer equation has been widely used [22,29–32] with the formula as:

$$SH = \rho \mathbb{C}\_p \mathbb{C}\_{DH} V\_{10} (Ts - Ta) \tag{1}$$

where *ρ* is the density of air, taken as a constant value of 0.8 kg m−<sup>3</sup> [33]; *Cp* is the specific heat of dry air under constant pressure, which is 1005 J Kg−<sup>1</sup> K−1; *CDH* refers to the heat transfer coefficient, which is usually prescribed as a constant magnitude of 4 × <sup>10</sup><sup>−</sup>3; *<sup>V</sup>*<sup>10</sup> is the wind speed at the height of 10 m; *Ts* is the ground surface temperature; and *Ta* refers to the air temperature at the height of 1.5 m. In the following, only the calculated *SH* at QOMS and NAMORS can be obtained by using Equation (1) based on wind speed, ground temperature, and air temperature, respectively, owing to the limitation of data.

Additionally, unitary linear regression analysis is applied in this study to calculate the trend in SH:

$$\mathbf{x}\_{i} = a + bt\_{i} \ (i = 1, 2, \cdots, n) \tag{2}$$

where *xi* denotes the meteorological variable with the sample size of *n*, and *ti* is the time corresponding to *xi*, and the linear regression equation between *xi* and *ti* can be solved as follows:

$$\begin{cases} b = \frac{\sum\_{i=1}^{n} \mathbf{x}\_{i} t\_{i} - \frac{1}{n} \left( \sum\_{i=1}^{n} \mathbf{x}\_{i} \right) \left( \sum\_{i=1}^{n} t\_{i} \right)}{\sum\_{i=1}^{n} t\_{i}^{2} - \frac{1}{n} \left( \sum\_{i=1}^{n} t\_{i} \right)} \\ a = \overline{\mathbf{x}} - b\overline{\mathbf{f}} \end{cases} \tag{3}$$

where *a* is the regression constant and *b* is the regression coefficient. The *x* is the average value of the meteorological variable *xi*, and *t* is the average value of the time.

The root mean square error (*RMSE*) denotes the extent to which the data deviate from the true value and tends to be applied for assessing the data reliability, which is used in Section 5 to evaluate the accuracy of new *SH* calculated by using new *CDH*.

$$RMSE = \sqrt{\frac{\sum\_{i=1}^{N} (y\_i - x\_i)^2}{N}} \tag{4}$$

where *N* is the total number of samples, and *xi* and *yi* denote the observed samples and the calculated samples, respectively.

#### **3. Annual and Seasonal Mean of the SH Diurnal Variations over the TP**

#### *3.1. Annual Mean*

Firstly, the annual mean time series of the observed SH at six stations are shown in Figure 1. Although the data length is not uniform for each station (Table 1), we can still see their common features. SH shows obvious interannual variabilities, especially in NADORS and NAMORS, while SH in MAWORS and SETORS displays a clear decreasing and increasing trends with −0.94 and 1.67 W m−<sup>2</sup> yr<sup>−</sup>1, respectively. This may be the result of multiscale variabilities in SH, in which the effect of diurnal variation may be important.

The general characteristics of the diurnal variation in the observed SH over the TP can be understood from its climatic state (Figure 2). The magnitudes of annual mean SH at each station are basically negative at night with no clear fluctuations, while they are positive in the day with larger magnitudes and obvious variations, and their peaks appear in the afternoon, which confirms previous results [22,32,34]. Of course, the SH diurnal variations over the TP also show some discrepancies, with different amplitudes and peak timings at those six stations. This results presumably from differences in the underlying surface, altitudes, and climate conditions among these stations. QOMS, with sparse and short surface vegetation, has the largest diurnal amplitude of SH, reaching a peak of 184.35 W m−<sup>2</sup> at 15:00, while MAWORS, influenced by the westerly wind all year round and being surrounded by large-scale modern glaciers, shows the smallest SH amplitude, with the peak value only reaching 112.49 Wm−<sup>2</sup> at 15:00. Three other stations (BJ, NADORS, and NAMORS), with altitudes higher than 4000 m, have similar SH amplitudes to the average mean of those six stations, but the peak timing at BJ is earlier, at 14:00. Additionally, SETORS, which is located in a forested valley, close to the southeastern TP, deviates from the average mean, and it has a weak amplitude comparable to that of MAWORS, but its peak occurs much earlier (about 12:00) than other stations. Obviously, six stations distributed along the east–west direction are scattered sparsely across the entire TP, and the time difference of 2–3 h among them cannot be ignored (Table 1). If the Beijing time of all stations is changed to local time, then we will find that the peaks of most stations appear around 12:00 or 13:00 (local time), except SETORS, which is covered by dense vegetation (mainly temperate needleleaf trees and alpine meadows), whose peak timing is around 10:00 (local time) (The following times are Beijing time unless otherwise specified). The conditions of the underlying surface may change the upward radiation flux, affect the peak timing of the diurnal variation in surface net radiation, and then affect the peak timing of SH. Further examination (figure omitted) confirms that the peak timing of the net shortwave radiation flux in SETORS is indeed much earlier (about 13:00) than that

of other stations (about 14:00 and 15:00), which is consistent with the peak timing in SH diurnal variation.

**Figure 2.** Climatological annual mean of diurnal variations of the observed SH (W m<sup>−</sup>2) at six stations over the TP during the period 2006–2016. The colored lines represent SH at the six stations, respectively, and the black line indicates their average. The time in *x* axis is Beijing time.

### *3.2. Seasonal Mean*

Figure 3 shows the seasonal mean of SH diurnal variations at the six stations, which indicates that the diurnal variations in four seasons generally have similar peak timing but different diurnal amplitudes at each station. The peak timing of SH in four seasons is mostly at 15:00 in MAWORS, NADORS, QOMS, and NAMORS, and at 14:00 and 12:00 at BJ and SETORS, respectively, which is consistent with Section 3.1 above.

**Figure 3.** Climatological seasonal mean of the diurnal variation of the observed SH (solid lines, W m−2) at (**a**) BJ, (**b**) MAWORS, (**c**) NADORS, (**d**) SETORS, (**e**) QOMS, and (**f**) NAMORS, and the calculated SH (dashed lines, W m<sup>−</sup>2) at (**e**) QOMS and (**f**) NAMORS.

The diurnal amplitudes of SH at three stations (MAWORS, NADORS, and QOMS) are consistently largest in spring, followed by summer and autumn, and smallest in winter, which is consistent with previous studies [19,34]. In contrast, winter SH amplitudes at two other stations (BJ and SETORS) are exceptionally strong; the strongest being in SETORS. Both stations are located in the eastern TP and are well covered by vegetation, which allows the surface to capture more net radiation absorption in winter, resulting in stronger SH and latent heat [35]. In spring, the precipitation at SETORS accounts for 25.4% of the whole year, which is 14–18% higher than that of two other stations (QOMS and NAMORS), which measured precipitation by the same method. In a word, more precipitation in spring will weaken the SH to a large extent. Moreover, the amplitudes of BJ and SETORS in summer are the weakest throughout the year; the stronger relative humidity (up to 63.22% and 80.45%) at these two stations and the resulting larger evapotranspiration will probably weaken the SH to a great extent, resulting in an abnormal small amplitude in diurnal variation of the summer SH. For NAMORS, the amplitude of SH diurnal variation in winter is stronger than that in autumn, which may be related to the similar diurnal variations of the difference between surface temperature and air temperature in winter and autumn at this station.

In the past, due to the lack of flux observations over the TP region, most studies on SH were based on the calculated SH by using Equation (1) or reanalysis datasets [36,37], but uncertainties and biases always existed. The diurnal variations of the calculated SH and the observed SH at QOMS and NAMORS are also displayed in Figure 3e,f to detect the bias between them. We find that the diurnal variation of calculated SH shows similar seasonal distribution with the observed one, showing the greatest amplitude in spring and the smallest in winter at QOMS, and an abnormally stronger amplitude in winter at NAMORS. However, significant differences exist in the diurnal variation amplitude between calculated SH and observed SH; that is, the former is about 64–100% larger than the latter, which suggests that the calculated SH over the TP used in previous studies may be largely overestimated. Additionally, the peak timing of the calculated SH also shows different features compared with that of observed SH, mainly in that the latter is at 16:00, which is one hour later than the former at QOMS, suggesting that a significant phase shift occurs between the calculated SH and observed SH, while this phenomenon is not evident at NAMORS.

#### **4. Monthly Changes of the SH Diurnal Variation over the TP**

In order to further understand the SH diurnal variation over the TP in detail, we need to shorten the time scale to obtain the monthly mean of the SH diurnal variation.

#### *4.1. Monthly Changes of the Diurnal Variation in Observed SH*

Figure 4 shows the monthly mean of the diurnal variation of the observed SH. Obviously, it can be seen that the diurnal amplitudes of the SH among stations are significantly different, with smaller values at MAWORS and SETORS and larger values at the other four stations, which is also consistent with the results in Section 3. Furthermore, the diurnal variation characteristics of the SH differ significantly from month to month for each station. Generally, the duration of positive SH increases from January, with the longest duration in May or June about 13 h, from 8:00 to 21:00, at BJ, NADORS, and NAMORS, while it is about 10 h, from 10:00 to 22:00, at MAWORS, and it then decreases in wintertime. However, at QOMS, which has the higher altitude, the particularly long duration of positive SH in July or August is up to about 19 h from 2:00 to 21:00. The special feature of SETORS is that the changes of the positive SH duration are less distinctive, with more scattered positive values at night.

For the timing of peak values in SH diurnal variation (Figure 5), it also fluctuates month by month and varies from station to station, almost within a time range of 1 h. For example, at MAWORS (Figure 5a), the timing of the SH peak fluctuates between 15:00 and 16:00 for most months; at SETORS and NAMORS (Figure 5c), the peak timing fluctuates during 12:00–13:00 and 14:00–15:00, respectively. It should also be noted that although the peak timing of monthly SH fluctuates steadily within a certain range, there are still interannual variations at some stations. For example, the peak timing at BJ (Figure 5a) varies during 13:00–14:00 before 2013, but during 14:00–15:00 after 2013. Similarly, the peak timing at NADORS (Figure 5b) almost fluctuates between 14:00 and 15:00 before 2014, but between 15:00 and 16:00 after 2014. Especially in QOMS (Figure 5b), this phenomenon is more obvious; that is, before 2015, the peak timing fluctuates between 14:00 and 15:00 for most months in the eight years, while it fluctuates between 15:00 and 16:00 after 2015. Surface radiation flux tends to largely affect the peak timing changes of SH. During the period 2006–2016, the incoming radiation may have changed due to the variations in atmospheric conditions, such as aerosol, cloud cover, and water vapor, and the outgoing radiation may have changed due to the changes in underlying surface conditions, such as vegetation growth and soil moisture. These variations of surface radiation flux will be reflected in the land–air temperature difference, and then result in the SH changes.

**Figure 4.** Monthly mean of the diurnal variation in observed SH (shaded, W m−2) at (**a**) BJ, (**b**) MAWORS, (**c**) NADORS, (**d**) SETORS, (**e**) QOMS, and (**f**) NAMORS during the period 2006–2016.

**Figure 5.** Peak timing (hour) of monthly diurnal variation in SH at (**a**) BJ and MAWORS, (**b**) NADORS and QOMS, and (**c**) NAMORS and SETORS.

Another important characteristic (Figure 4) is that the SH diurnal variations at QOMS and BJ have two distinct centers per year, as detailed in Table 2. In fact, the other four stations (MAWORS, NADORS, SETORS, and NAMORS) have similar characteristics to the two centers, but the second center is inconspicuous, so only the information of the first center of these stations is presented in Table 2. We can clearly see that, in general, the first center of the SH diurnal variations at all stations appears in spring afternoon for most years, except for SETORS, which appears in January to March. However, the second center appears in different months at QOMS and BJ, occurring almost at 14:00–15:00 in September or October for QOMS but at 14:00–15:00 in October to December for BJ. One thing of note is that the first center value is about 27% stronger than the second one on average, which indicates the spring SH is dominant on multiple time scales. Another point worth noting is that the SH diurnal variation in SETORS shows an increasing trend during the period 2007–2016 (Figure 4d), which is obviously manifested by the increasing peak values during the day. This is consistent with the trend in annual mean SH (Figure 1), indicating the internal consistency of different timescales.


**Table 2.** The month and the specific timing of the occurrence of the two large centers in SH diurnal variation at six stations from 2006 to 2016.


**Table 2.** *Cont.*

Furthermore, the main factors leading to the obvious centers of SH are discussed by taking QOMS (Figure 4e) as an example. Considering Equation (1) and Figure 6, we can clearly see that two centers in land–air temperature difference are similar to those in SH compared to wind speed, which shows many irregular centers. For the diurnal variation in land–air temperature difference, the first centers always appear in March to May, and the second ones with a smaller value almost occur in September or October, which corresponds well to those in SH (Table 2) each year. This indicates that the land–air temperature difference is the main factor causing the diurnal variation features of the SH at QOMS.

**Figure 6.** Monthly mean of the diurnal variation in (**a**) land–air temperature difference (Ts-Ta, unit: ◦C) and (**b**) wind speed at the height of 10 m (V10, unit: m/s) at QOMS during the period 2006–2016.

### *4.2. Monthly Changes of the Diurnal Variation in Calculated SH*

From the above, we have understood that there is a certain difference in the seasonal mean diurnal variation between the calculated SH and the observed SH, so it is necessary to explore the condition for the monthly mean diurnal variation between them, which is of great significance to improve the SH calculation.

Figure 7 shows the monthly diurnal variation of the calculated SH at QOMS and NAMORS. Generally, comparing the diurnal variation characteristics of the observed (Figures 4 and 5) and calculated (Figure 7) SH, we can see that the calculated SH at these two stations have two large-value centers in each year, and the second center is more obvious at QOMS than at NAMORS, similar to the observed characteristics. However, some clear differences still exist between them. The diurnal amplitude in the calculated SH at QOMS and NAMORS can be up to 2.86 and 3.09 times stronger than those in observed SH, respectively. Moreover, the differences between calculated SH and observed SH also exist in the peak timing of the SH diurnal variation and the peak timing fluctuations range. For QOMS, the peak timing shows a clear delay in calculated SH. Before 2011, it is basically stable around 17:00, and then remains at around 16:00 after 2013, which is relatively delayed by 1~3 h compared with the observed one. For NAMORS, the range of the peak timing fluctuations in calculated SH is much greater than that in the observations. Especially before 2009, the peak timing of calculated SH changes during 14:00–17:00 in most months; From 2009 to 2013, it fluctuates within 2 h during 14:00–16:00, and then fluctuates within 1 h during 14:00–15:00 after 2013, but the peak timing of observed SH mainly varies within the range of about one hour.

**Figure 7.** Monthly mean of the diurnal variation in calculated SH (shaded, W m<sup>−</sup>2) at (**a**) QOMS and (**b**) NAMORS during the period 2006–2016. The curve lines in (**c**) represent the peak timing (hour) of monthly SH diurnal variation at each station, respectively.

Heat transfer coefficient CDH is widely used as a constant value of 4 × <sup>10</sup>−<sup>3</sup> in the SH calculations over the TP, but the above facts show that there is a significant deviation between the observed SH and the calculated SH by Equation (1) on the diurnal scale, not only in the diurnal amplitude, but also in the peak timing and its fluctuation range. From Equation (1), we can clearly see that SH is calculated by using the heat transfer coefficient (CDH), wind speed (V10), and difference between surface temperature and air temperature (Ts-Ta). Actually, V10, Ts, and Ta are originally from the same dataset and the same stations (QOMS and NAMORS), hence the CDH may be the key factor resulting in the deviation in diurnal variation between the observed and calculated SH at these two stations.
