Using Long Term Simulations to Understand Heat Transfer Processes during Steady Flow Conditions in Combined Sewers

This paper describes a new heat transfer parameterisation between wastewater and in-sewer air based on understanding the physical phenomena observed in free surface wastewater and in-sewer air. Long-term wastewater and in-sewer air temperature data were collected and studied to indicate the importance of considering the heat exchange with in-sewer air and the relevant seasonal changes. The new parameterisation was based on the physical flow condition variations. Accurate modelling of wastewater temperature in linked combined sewers is needed to assess the feasibility of in-sewer heat recovery. Historically, the heat transfer coefficient between wastewater and in-sewer air has been estimated using simple empirical relationships. The newly developed parameterisation was implemented and validated using independent long-term flow and temperature datasets. Predictive accuracy of wastewater temperatures was investigated using a Taylor diagram, where absolute errors and correlations between modelled and observed values were plotted for different site sizes and seasons. The newly developed coefficient improved wastewater temperature modelling accuracy, compared with the older empirical approaches, which resulted in predicting more potential for heat recovery from large sewer networks. For individual locations, the RMSE between observed and predicted temperatures ranged between 0.15 and 0.5 °C with an overall average of 0.27 °C. Previous studies showed higher RMSE ranges, e.g., between 0.12 and 7.8 °C, with overall averages of 0.35, 0.42 and 2 °C. The new coefficient has also provided stable values at various seasons and minimised the number of required model inputs.


Introduction
This section highlights the importance of modelling heat transfer in sewers, addresses the expected challenges with existing techniques and proposes a new parameterisation for the heat transfer between wastewater and in-sewer air. The section also introduces advanced methods for analysing the modelling accuracy using a Taylor diagram.

The Need for Modeling Heat Transfer in Sewers
The increasing demand for clean energy sources caused by more stringent carbon emission regulations is creating an urgent need for alternative low carbon energy sources in urban areas. Recovering heat from sewers may provide an attractive potential source of lowcarbon heat due to the relatively high wastewater temperature in sewers, the predictable daily dry weather flow patterns and the proximity of numerous potential heat users in many urban areas. Heat recovery from sewers has more significant potential than other forms of energy, such as chemical energy, e.g., Hao et al. (2019) [1] have estimated theoretical energy densities of 4.64 and 1.54 kWh/m 3 for heat and chemical energy potentials, respectively. These estimates neglect the efficiency factors when extracting the energy. Recovering heat from wastewater may also qualify for government subsidies, e.g., the renewable heat incentive (RHI) in the UK (Bulteau et al., 2019) [2]. Short-term measurement programmes of wastewater temperatures in combined sewers recorded temperatures reaching 27 • C (Schilperoort andClemens, 2009 andCipolla andMaglionico, 2014) [3,4]. A case study by Abdel-Aal et al. (2018) [5] indicated that 7 to 18% of the heat demand for a 79,500 PE Belgian city could be met by sewer heat recovery. This was based on assuming a 100% efficient heat exchange system.
Quantifying the potential heat recovery from sewer pipes requires a reliable prediction of the wastewater temperature variation to ensure the changing environmental conditions do not significantly vary the amount of potential heat recovery without compromising the operation of "end of system" wastewater treatment plants. It is also important to always maintain a wastewater temperature in a sewer pipe above the freezing point in order to avoid physical damage and blockage. The modelling of wastewater temperature requires knowledge of wastewater hydraulics, surrounding temperatures (soil and in-sewer air) and heat transfer coefficients. Durrenmatt (2006), Durrenmatt andWanner (2008 and [6][7][8] developed a deterministic model called TEMPEST that estimates wastewater temperature variations along with the longitudinal profile of a sewer pipe based on heat transfer laws and the principles of mass balance. Heat transfer parameters, including thermal conductivity and penetration depth of surrounding soil, were calibrated and validated in TEMPEST using data observed only over a few days in the same season [8]. Furthermore, TEMPEST can be highly parameterised since it requires fifteen input parameters, including specific variables, such as chemical oxygen demand and degradation rate. The TEMPEST model was tested onsite by utilising a limited dataset measured between 14 February and 22 March 2008. Abdel-Aal (2015) [9] developed a more efficient computational model that predicts the total amount of potential heat recovery from sewer networks and the consequent wastewater temperature variations along the pipes. The latter model, which will be referred to as the "2015 model", was calibrated and validated using long-duration datasets (February to July in 2012 and 2013) in a number of combined sewers.
Abdel-Aal (2015) [9] found that calibrated coefficients were inconsistent between the winter and summer data, which indicated that the heat transfer between wastewater and in-sewer air is not negligible. Elías-Maxil et al. (2017) [10] developed a parsimonious computational model to estimate water temperature using similar heat and mass balance principles to that of TEMPEST. The main difference between the TEMPEST model and that developed by Elías-Maxil et al. (2017) [10] was that the latter neglected the heat transfer between water and air. Elías-Maxil et al. (2017) [10] tested their model in the field on an out-of-service sewer with an internal diameter of 0.23 m using controlled inflows of cold and hot water. The current work, therefore, focuses on developing a new model for predicting sewage temperature changes, focusing on the interaction between in-sewer air temperature, air velocity and wastewater temperature to reflect the physical interaction of the in-sewer air/wastewater boundary layer.

Complexity of Modeling Heat Transfer in Sewer Pipes
Heat transfer in sewer pipes is usually a complex process since the pipe is partially filled, some of the parameters are difficult to quantify (e.g., in-sewer air velocity), and the relevant parameters are subject to seasonal variations. Previous studies concerned with modelling wastewater temperatures have not addressed in detail the physical interaction between wastewater and in-sewer air, which can improve the modelling accuracy and result in predicting more potential heat recovery, as will be shown in Section 3.4. The challenge of quantifying in-sewer air velocity can compromise the model accuracy since in-sewer air temperature is an influential parameter, as proven by Abdel-Aal et al. (2015) [9] using a predictive machine learning model. Therefore, this paper develops a new model for predicting sewage temperature changes, focusing on the physical interaction between in-sewer air and wastewater boundary layer while accounting for the seasonal variations presented in observed long-term field data.

Challenges Associated with Modeling In-Sewer Air Velocity
The heat transfer between wastewater and in-sewer air has previously been estimated using empirical methods as a function of in-sewer air velocity. The velocity of in-sewer air is not usually known and can also be difficult to predict accurately for large sewer networks. Pescod and Price (1982) [11] developed an empirical model to estimate in-sewer air velocity based on observations from a single laboratory flume. Other researchers have studied the airflow in sewers; for example, Edwini-Bonsu and Steffler (2006) [12] simulated the ventilation process in sewers by predicting the air movement and the pressure variation in the system. The latter authors estimated the in-sewer airflow rate as a function of acceleration due to gravity, available headspace, wastewater velocity and wastewater depth. However, the study was a general framework for drop structures, and the model was not validated using observed data. This may limit the model applicability in complex systems, e.g., large sewer networks. A conceptual model by Qian et al. (2018) [13] was developed to simulate the in-sewer air movement in a small prototype sewer network by considering air pressure, wastewater drag and friction forces. Their model assumed in-sewer air velocities to predict the airflow rate while considering a single direction for the airflow, which may not be the case in real sewer systems. Qian et al.'s (2018) [13] model has a number of assumptions that may result in discrepancies, especially if the model is applied to large sewer networks. These include the pipe connections having no further branches and the limited number of in-sewer air leakage sources. Some of the model discrepancy causes were found, by Qian et al. (2018) [13], to be mainly due to changes in drop diameters, breakup lengths and wastewater terminal velocities. The models by Qian et al. (2018) [13] and Edwini-Bonsu and Steffler (2006) [12] were aimed at predicting the in-sewer airflow by assuming in-sewer air velocity (u a ) mainly for controlling odours, and therefore, their focus was not to develop a methodology for estimating u a .
Durrenmatt and Wanner (2014) [8] incorporated u a in their model and assumed it is 50% of wastewater velocity. The latter authors considered u a to be a key parameter in their TEMPEST model since it was utilised to estimate the airflow, air transfer coefficient and heat transfer between wastewater and in-sewer air. Durrenmatt and Wanner (2014) [8] assumed that the airflow in sewers was based mainly on the pipe structure, which makes the insewer air velocity independent of external forces and variations of airflow streamwise. This assumption may accumulate simulation errors if the model to be deployed on large sewer networks. Abdel-Aal (2015) [9] calibrated u a, which was initially found to be unrealistically high (100 s m/s) before limiting its calibration range between −0.5 and 0.5 m/s, indicating that the way his model was formulated could be improved. Therefore, there is a need for developing a new parameterisation that can incorporate the effects of in-sewer air on wastewater temperature while relying on measurable parameters. This parameterisation needs to incorporate the widely varying temperatures of in-sewer air and wastewater over the year.

Introducing a New Parameterisation for the Heat Exchange between Wastewater and In-Sewer Air
The development of the new heat transfer coefficient (h Fr ) will enable the estimation of the heat transfer parameter at the air-wastewater surface by utilising physically based parameters that can reflect changes in the wastewater and in-sewer air dynamics at the interface. The 2015 model, which was based on the heat transfer coefficient being dependent on u a , as will be shown later by Equation (1), showed wide variability in the root-meansquared error (RMSE) and calibrated parameters, which indicates a weakness in the model structure. Therefore, utilising a relation that better simulates the physical variation of the relevant heat transfer parameters expressed by the dimensionless Fr number is expected to enhance the modelling accuracy.

Taylor Diagram
An advanced way to compare the accuracy between various models is to plot a Taylor diagram, which was developed by Taylor (2001) [14]. This can be achieved by graphically plotting various statistics that interpret the model accuracy in terms of correlation, RMSE and standard deviation in a single diagram for various models. A Taylor diagram can also be used to represent the normalised statistics and thus facilitates the fair comparison between different models. Therefore, this study will utilise the Taylor diagram to compare the performance of the model developed in this paper with previously developed models.

Novelty, Aim and Objectives
Novelty: Introduce and examine the reliability of a new heat transfer coefficient (h Fr ) between wastewater and in-sewer air.
Aim: Develop a new heat transfer coefficient between wastewater and in-sewer air that is based on physical variations in wastewater hydraulics. Objectives: 1.
Analyse the fluid dynamics at the interface between wastewater and in-sewer air to develop a dimensionless heat transfer parameter; 2.
Utilise long-term measured seasonal datasets for calibrating the new parameterisation for the heat transfer coefficient; 3.
Compare the accuracy of the new model developed in this paper and existing models by employing advanced analysis techniques, e.g., the Taylor diagram.

Methodology
This section describes the field sites and shows the collected seasonal field data, which include measured values for the wastewater, in-sewer air and soil temperatures, as well as in-sewer wastewater flow rates. The section explains the techniques used for modelling wastewater temperature and the process followed for calibrating and validating the key heat transfer parameters in sewers, including h Fr . Data analysis techniques used to examine the model performance over a long-term interval, including the probability density function (PDF) and Taylor diagram, are also described.

Data Collection
The data collection was performed in collaboration with Aquafin, a Belgian water company in Antwerp. The authors specified the measuring parameters based on the modelling requirements, defined the sensor specifications and processed the measured data. Four sewer pipes were utilised to measure the temperatures of wastewater (upstream and downstream), in-sewer air and soil, wastewater flow rate, velocity and depth. Since the underlying heat transfer mechanism varies depending on the sewer size (Abdel-Aal, 2015) [9] the four sites in this study were selected to represent the range of sewer pipes in a catchment area; small residential sewers and large-scale collector sewers. The site category depends on the dry-weather flow rate (DWF) and the position in the sewer network. The residential sewers have a higher diameter to flow ratios than larger collectors and are closer to residential homes, which results in higher wastewater temperatures. Collector sewers have no direct connections to residential homes and are mainly designed to transport wastewater from the downstream of the residential areas to the treatment plant. Table 1 describes the sites from which data were recorded between March 2012 and July 2012 for residential sewers and between February 2013 and June 2013 for collector sewers. Residential sewers were located 1.5 km apart in the same catchment, while the large collectors were 3 km apart and located 15 km from the residential sites. Temperatures of wastewater and in-sewer air were observed every 20 min by Tinytag (PBRF-5006-5m) sensors with ±0.06 • C accuracy and less than 0.05 • C resolution. The measurement frequency was found to be adequate to account for the temperature variation. The in-sewer air temperature sensors were installed 1 m below the manhole cover, while the wastewater temperature sensors were fixed just above the sewer pipe invert level and were always covered with wastewater. Thermocouples of the EJB378 K-type, calibrated by BERCU EJB, were used to measure soil temperatures 3 m underground with ±0.7 • C accuracy. Wastewater flow rate, velocity and depth were recorded every 2 min using Isco 2150 area-velocity type meters that were calibrated onsite by Studiebureau Patrick Casier, Belgium. The flow meters were visited weekly when the wastewater velocity and depth were checked, and sensors cleaned if necessary. Wastewater depth was observed using an Isco 2150 flow meters through submerged pressure transducers with ±0.003 m accuracy. Wastewater velocity was measured, with an accuracy of ±0.03 m/s, using Isco 2150 area velocity type meter employing the Doppler ultrasonic method with a 500 kHz frequency. All sensors were purchased from Antwerp, Belgium. More details on the measurement campaign are available in Abdel-Aal (2015) [9]. DWF for each site was estimated by observing the flow rate variations when there are no relatively high flow rates (spikes), such as the periods between 15 and 30 March at sites 1 and 2 (Figures 1 and 2). These flow values at dry periods were then averaged over their timespan to define the average DWF for the relevant site. The modelling, calibration and validation for the urban sites were based on DWF periods to avoid complications associated with sudden increase in flow rates and potentially high wastewater temperature changes caused by rainfall. The omission of rainfall events was achieved by removing (filtering out) flow rates and their corresponding parameters, e.g., observed temperatures, that are greater than the DWF. The collector sewers did not show significant differences in terms of calibration and validation results when DWF or the total flow was considered. Hence, there was no need to omit the rainfall events in the collector sewer data. The plotting of observed data was smoothened using a moving average of 15 points (i.e., over 5 h for temperatures and 30 min for flow rate) to avoid noisy data. Figures 1 and 2 show the variation of wastewater temperatures and flow rates of the smoothed data for sites 1 and 2, respectively.  Expected patterns of wastewater temperatures and flow rates were observed at sites 1 and 2, where the flow peaks reflect the infrequent rainfall events. The rainfall runoff peaks often show short-term dips in the sewer flow temperature, but not always. The occurrence and size of these dips would be dependent on the temperature of the rainfall runoff as well as the volume of rainfall runoff, but this has not been studied further. Wastewater temperatures at sites 1 and 2 generally increased before plateauing, e.g., on 10 April and 9 June in Figures 1 and 2. This is thought to be due to the change in in-sewer air temperature, which plateaus or dips around the same periods of time. Figure 3 shows the local soil temperature, measured 3 km from site 1 at 3.7 m depth below ground level, and the variations of in-sewer air at site 1. The soil temperature at the above depth is expected to be the same or very close to that in the vicinity of the four sites. In-sewer air temperatures for sites 1 and 2 were close in terms of patterns and values; hence only site 1 data were plotted in Figure 3. The pattern of observed in-sewer air temperature is very close to that of the wastewater temperatures at site 1 ( Figure 1). Unlike soil temperatures, in-sewer air temperatures vary significantly during the day and increase considerably from winter to summer (Figure 3), which suggests that the diurnal wastewater temperature variation was less influenced by soil temperatures. This indicates that the temperature of in-sewer air largely affects that of wastewater. Therefore, the heat transfer process between wastewater and in-sewer air is thought to be important when modelling wastewater temperature in sewers over long durations. Figures 4 and 5 show the observed wastewater temperatures and flow rates for sites 3 and 4, respectively. Some technical issues related to data monitoring occurred between 7 and 21 June 2013 in site 4, hence the gap in Figure 4.  A relatively steady drop in wastewater flow rate was observed, in Figures 4 and 5, during wintertime until around the start of April at sites 3 and 4. It is likely that from early April onwards, the groundwater infiltration was lower than that in other months. The lower temperature of the infiltrated water may cause less increase in wastewater temperature up until April, compared to sites 1 and 2. There was a sharp drop in the wastewater flow rate at site 3 between 1 February and 8 March, which is likely due to the earlier snow that melted over February to reflect this drop. This explains the sharp drop in wastewater temperature on around 20 March, as the melted snow was mixed with surface runoff and drained into the sewers at a lower temperature. The sudden flow rate drop also at site 3 on 13 March may be caused by a blocked sensor.  [6,9,15,16] utilised a relationship, developed by Flinspach (1973) [17], to estimate the heat transfer coefficient (h wa ) between wastewater and in-sewer air knowing the relative in-sewer air velocity as given by Equation (1).

Existing Technique of Modeling Heat Exchange between Wastewater and In-Sewer Air
where u wa is the relative in-sewer air velocity to that of wastewater (m/s). The derivation of Equation (1) is not fully explained by Flinspach (1973) [17]; hence the origin of the calibration coefficient (5.85) is unclear. In addition to the challenges associated with an accurate estimation of in-sewer air velocity (Section 1.2), a closer examination of Equation (1) suggests that it does not fully consider the heat transfer processes between the wastewater and the in-sewer air, as the velocity of the wastewater free surface relative to the in-sewer air in contact with the free surface can reach zero, as shown in Figure 7 (u w = u a ). The wastewater and in-sewer air velocities can be estimated based on the interface fluid dynamics. Abdel-Aal (2015) [9] described the heat transfer process in sewer pipes, as shown in Figure 6, and utilised Equation (1) to estimate h wa in his 2015 model. The challenge in estimating the in-sewer air velocity has led Abdel-Aal (2015) [9] to calibrate for the in-sewer air velocity. Figure 6. Illustration of heat transfer processes in a sewer pipe. q is the heat transfer rate (Watt) between wastewater (w) and in-sewer air (a) and soil (s).

New Parameterisation for the Heat Exchange between Wastewater and In-Sewer Air
Unlike Abdel-Aal (2015) and Abdel-Aal et al. (2014) [9,15,16] simulations, the model developed in this paper attempts to account for the impact of the in-sewer air velocity profile, close to the interface, on the heat transfer processes and develop a heat transfer parameterisation that accounts for the free surface air boundary layer. Convective heat transfer of a fluid moving on a surface relies on the velocity boundary layer caused by the shear stresses between the fluid and the surface (Schlichting andGersten, 2003 andIncropera et al., 2007) [18,19]. In this study, in-sewer air is assumed to be moving over the "rough" wastewater surface and therefore, it is necessary to study the boundary layer of the in-sewer air created by the moving wastewater surface. It is important to consider the effect of the "roughness" of the wastewater surface waves on the shearing of the adjacent airflow, as shown in Figure 7. Since the exact behaviour of wastewater interaction with in-sewer air is poorly understood, as a first approximation, it is assumed that the characteristics of wastewater surface movement in a sewer pipe have a similar effect on the airflow over a static rough surface with turbulent conditions. This assumption is reasonable to reflect the irregularity of the wastewater surface at the interface zone and develop a better understanding of the local air and wastewater velocity profiles. Based on the aforementioned assumptions, Figure 7 illustrates the shear stress and velocity boundary layers of in-sewer air and wastewater. For simplicity, the rough wastewater surface was assumed to have a fluid roughness that is related to the geometry of the surface waves. This is analogous to the fluid roughness of a static boundary in a channel or river bed, which can be represented by a layer of spheres packed together (Schlichting and Gersten, 2003) [18]. The wastewater shear stress (τ w ) below the interface zone of Figure 7 becomes zero when the velocity gradient of the wastewater along the vertical axis ( du w dy ) is zero. Einstein and El-Samni (1949) [20] found that, for a static bed, the fluid roughness is around 20% of the static roughness size (k W ), depicted as ∆ in Figure 7. The shear stresses and velocities of in-sewer air and wastewater close to those boundaries are controlled by the effective fluid roughness (∆), i.e., 20% of the wastewater wave amplitude. The scale of ∆ implies that the impact of in-sewer air velocity is less effective when analyzing the heat transfer at the interface zone. Therefore, and since wastewater hydraulic data are available in this work, the focus will be on investigating measurable parameters that can help describe the wastewater shear stress boundary layer, which influences the in-sewer air velocity gradient and hence impacts the physical heat transfer between wastewater and in-sewer air. Accounting for eddy tangential velocity and combining the shear stress and moment equations, the shear stress of wastewater in a boundary layer flow can be expressed, based on the Prandtl eddy model approach, by Equation (2) [18]: where τ w is the wastewater shear stress in the boundary layer (N/m 2 ), ρ w is the wastewater density (kg/m 3 ), K is Von Karman's constant, y is the relevant vertical distance from the zero-velocity plane (m), and u w is the wastewater velocity (m/s). Similarly, Equation (3) can be used to describe the shear stress of in-sewer air: where τ a is the in-sewer air shear stress in the boundary layer (N/m 2 ), ρ a is air density (kg/m 3 ), and u a is the in-sewer air velocity (m/s). The shear stresses of wastewater and in-sewer air are equal at the effective interface between the two fluids ( Figure 7). The density of in-sewer air (ρ a ) is significantly less than that of wastewater, and since both τ w and τ a are equal at the interface, the variation of wastewater velocity over the vertical axis ( du w dy ) is much less than that of in-sewer air ( du a dy ) as indicated by Equations (2) and (3). Assuming the in-sewer air velocity (u a ) follows a logarithmic profile, u a can be expressed as follows in terms of u * , which is the square root of the in-sewer air shear stress over its density, as shown by Equation (4): where u a is the in-sewer air velocity (m/s), K is Von Karman's constant, ∆ is the depth from the zero-velocity plane to the point where τ a reaches its maximum (m), and y is the vertical distance in the in-sewer air from the zero-velocity. Equation (4) demonstrates that the shearing of the in-sewer air velocity profile is influenced by the wastewater fluid roughness height, i.e., the wave amplitude. Since shear stresses of wastewater and in-sewer air are equal at the interface surface, Equations (2) and (3) show that the shear stress of in-sewer air is dependent on the wastewater surface wave amplitude and its velocity. Therefore, to enhance the accuracy of modelling the heat transfer between wastewater and in-sewer air, a parameterisation that includes the influence of the wastewater surface velocity and wave amplitude can be developed. Analysis of water surface elevation data collected by Romanova (2013) [21] indicates that for a variety of pipe roughness conditions at shallow depths (e.g., less than a third of the pipe diameter), the standard deviation of the surface wave height is related to Froude number. This limitation inflow depth corresponds to dry weather flow in combined sewers. This indication leads to the utilisation of the dimensionless Froude number (Fr), given by Equation (5), since it is related to wastewater surface velocity and surface wave height. Adding a dimensionless calibrating factor (f c ) to obtain an optimal value of the new heat transfer coefficient (h FR ), Equation (6) can be developed to express h Fr as a function of Fr in a partially filled sewer pipe. Although the wastewater hydraulic parameters were measured in this study, they can be retrieved from existing calibrated hydrodynamic models when deployed, e.g., in large sewer networks as performed in Section 2.8.
h Fr = f c × u w g·D w (6) h Fr is the new heat transfer coefficient between wastewater and in-sewer air as a function of Froude number (W/m 2 ·K), f c is the calibrating parameter for h Fr (W/m 2 ·K), u w is the wastewater velocity (m/s), g is the acceleration due to gravity (m/s 2 ), and D w is the hydraulic depth of wastewater in a sewer pipe (m).

Modelling
The model developed in this study is referred to as the "2020 model". The calibration of heat transfer parameters, using wastewater Froude numbers, is based on the principle of energy balance over a pipe length and incorporates heat transfer equations assuming slowly changing DWF conditions, that is, quasi-steady flow conditions relative to the timescale of the heat fluxes. The model assumes that wastewater temperature variation along a sewer pipe profile is due to heat being exchanged with the in-sewer air through the flow surface area and with the surrounding soil through the wetted pipe wall. This modelling technique is efficient since it only accounts for the most significant processes. The model efficiency is crucial for simulating wastewater temperature within large sewer networks to investigate optimum locations for potential heat recovery applications.
Referring to Figure 6, Equation (7) can be developed: . mc p,w dT = q wa + q ws where . m is mass flow rate (kg/s), c p,w is the specific heat capacity for wastewater (J/kg·K), dT is the temperature difference between upstream and downstream ends (K), q wa is heat transfer between wastewater and in-sewer air (Watt), and q ws is the heat transfer between wastewater and soil (Watt).
Equation (8) was originally formulated by Abdel-Aal (2015) [9] on the basis of Equation (7) to estimate wastewater temperature variation along the sewer pipe profiles and was also implemented in this study.
where T is the temperature (K), m is the nodal expression of the wastewater temperature location along the pipe length, R is thermal resistivity (m.K/W) between wastewater and in-sewer air (wa) and between wastewater and soil (ws), ∆x is the increment length streamwise (m), ρ w is the wastewater density (kg/m 3 ), Q is the wastewater volumetric flow rate (m 3 /s) and c p,w is the specific heat capacity for wastewater (J/kg·K).
Referring to Equation (8) and the newly proposed Equation (6), the expressions for the thermal resistivity between wastewater and in-sewer air ( R wa ) and between wastewater and soil (R ws ) are given by Equations (9) and Equation (10), respectively: where w, a, p and s denote for wastewater, air, pipe and soil, respectively, h Fr is the new convective heat transfer coefficient between wastewater and in-sewer air as a function of the Fr number (Equation (6)) (W/m 2 ·K), b is the surface width of wastewater in a sewer pipe (m), t p is the pipe wall thickness (m), d s is the soil penetration depth (m), wet.p is the sewer pipe wetted perimeter (m), and k is thermal conductivity (W/m·K).

Comparison with Existing Sewer Pipe Models
The model developed in this paper was compared with existing models developed by Abdel-Aal (2015) [9], Durrenmatt and Wanner (2014) [8] and Elías-Maxil et al. (2017) [10]. The same measured parameters, i.e., wastewater flow and temperatures of wastewater, in-sewer air and soil, were utilised for running the former models. The Durrenmatt and Wanner (2014) [8] model (known as TEMPEST) is given by Equation (11): Elías-Maxil et al.'s (2017) [10] model, given by Equation (12 ), was based on TEMPEST while omitting some parameters: where A W is the wetted cross-sectional area of wastewater (m 2 ), T W is the wastewater temperature (K), q is the heat flux between the pipe wall and the wastewater (pw) (W/m 2 ), between wastewater and in-sewer air (wa) (W/m 2 ), due to evaporation (ew) (W/m 2 ) and produced from biochemical reactions (COD) (W/m 3 ). b is the wastewater surface width (m).

Calibration
Calibration of the key heat transfer parameters was based on the sensitivity of the wastewater temperature model to the parameters, the availability of the parameter values in published literature and the ability to measure them onsite. The heat transfer coefficient between wastewater and in-sewer air plays a major role in simulating wastewater temperatures, as explained in Section 1, and hence it was important to calibrate the f c parameter (Equation (6)). Sensitivity analysis performed by Abdel-Aal (2015) [9] indicated that parameters related to soil, i.e., penetration depth (d s ) and soil thermal conductivity (k s ), were more important when compared to pipe thermal conductivity. Soil penetration depth is defined as the distance from the wetted area of the sewer pipe invert level at which soil temperature becomes independent of the heat transferred between wastewater and soil through the pipe wall. Determination of d s may prove impractical since it would usually require accurate measurement of the temperature profile of the surrounding soil. Soil thermal conductivity (k s ) is difficult to determine since it varies with the soil type and its saturation conditions (Mitchell and Soga, 1993) [22]. This led Abdel-Aal (2015) [9] to calibrate both d s and k s . However, Abdel-Aal (2015) [9] concluded that calibrating both d s and k s as separate parameters may compromise the accuracy of the calibrated parameter values, as they were not independent. Hence, the obtained values of calibrated d s and k s may not necessarily reflect their actual values. Therefore, it was decided to calibrate both d s and k s as a single parameter in the form of the d s /k s ratio. Hence, the calibration parameters in this work are the calibrating parameter (f c ) for heat transfer coefficient between wastewater and in-sewer air (h Fr ) and the ratio of soil penetration depth to soil thermal conductivity (d s /k s ).
Optimisation of constrained nonlinear function in MATLAB 2020a [23] software was employed to calibrate f c and d s /k s . An optimisation function with constraints was selected since the values of calibrated parameters were known to be positive and ought to be within a realistic range. The f c parameter was assumed to be varying between 0.01 and 100. This was based on observing the calibration results of a wider range (0.01-500), which revealed that f c would not exceed 100. The d s /k s range was estimated based on the available values in the published literature. Hence, k s was assumed to be between 0.1 and 2.5 W/m.K (Mitchell and Soga, 1993) [22]. The penetration depth is a function of the thermal diffusivity, which is related to the soil thermal conductivity and thermal capacity. A range between 0.001 to 5 m was assumed for d s to ensure full coverage of possible penetration depths, which is wider than that found in Durrenmatt and Wanner (2014) [8] (0.01 to 1 m). Considering the above d s and k s ranges and since low k s values lead to shorter d s , the calibrating parameter of d s /k s ranged from 0.01 to 2 m 2 ·K/W. The optimisation function operates by iterating values within the above pre-set ranges to find the minimum of a target function. The target, in this case, was set to minimise the difference between observed and modelled downstream wastewater temperatures as given by Equation (13). This target was interpreted by the norm of the difference between observed and modelled values as given by Equation (14), which was set to zero in an attempt to minimise it. Setting the target of Equation (13) to zero assumes that the distribution of the modelling error (Equation (15)) is close to normal. Abdel-Aal (2015) showed that the modelling error was approximately normally distributed in most cases, using a simpler form of sewer pipe modelling and a more limited dataset. Therefore, it is reasonable that Equation (13) can be used for the target function: where V is the square root of the sum of the squared error array, as shown by Equation (14): "error" is the modelling error as given by Equation (15): where N is a total number of time steps, n is the observed sample number, DS denotes downstream, WW denotes wastewater, O is observed temperature, and M is modelled (or predicted) temperature. The calibrating parameter values obtained by Abdel-Aal (2015) [9] were found to be varying with calendar months, and hence, the calibration process in this study was performed on a monthly basis to examine temporal variations. The calibration was performed using Equation (8) on residential and collector sewer data separately. Both sites 1 and 2 represent residential sewer characteristics, and site 1 data were arbitrarily selected for calibrating the residential sewer model, while site 2 data were used for validating the resulting model. For the collector sewer model, site 4 was used for calibration since it provided more data, while site 3 data were used for validation.

Summarising the Models' Performance
This section explains how the modelling accuracy was assessed using the probability density function (PDF) and Taylor diagram plots.

Probability Density Function (PDF)
An effective way to compare modelling results is to plot the PDF of the modelling error Equation (15). The PDF was plotted for the model results during both calibration and validation for each calendar month.

Taylor Diagrams
Taylor diagram [14] was used to present the correlation between observed and modelled wastewater temperature time-series, their normalised standard deviation and normalised centered root mean square errors. The formulae used to compute the above parameters are explained in this section. The correlation coefficient (R) between observed (O) and modelled (M) wastewater temperatures is given by Equation (16). The advantage of presenting correlation in this work is that the calibrated models (residential and collector) were based on time-series data, and hence plotting the correlation provides a clear representation of the overall relationship between modelled and observed values across the given calendar month, which is around 2000 time steps: 16) where N is the total number of time steps, n is the observed sample number, M and O are the modelled and observed wastewater temperatures, respectively ( • C), M and O are the mean values for modelled and observed wastewater temperatures, respectively ( • C), σ M and σ O are the standard deviations for modelled and observed wastewater temperatures, respectively ( • C). The root-mean-squared error (RMSE) was computed using Equation (17): In order to eliminate the effects of the overall bias (B), which is defined by Equation (19), the means of modelled ( M) and observed ( O) wastewater temperatures were subtracted to produce the centered root mean square error (CRMSE). Hence, and based on Equation (17), Equation (18) was formulated to show how CRMSE can be computed: CRMSE can also be expressed in terms of standard deviations and correlation coefficient, as shown by Equation (20): The standard deviation for modelled (σ M ) and observed (σ O ) temperatures can be expressed in a similar manner, as shown by Equations (21) and (22): Thus, CRMSE, σ M , σ O and R for all modelled values can be plotted on a single Taylor diagram. Since the wastewater temperatures were observed at different months and sites, where the behavior of heat transfer mechanism is likely to vary (Abdel-Aal, 2015) [9], and in order to compare the statistics of various models in one Taylor diagram for all sites, the model statistics were normalised. Normalisation can be achieved by dividing the statistics by the relevant observed standard deviation (σ O ). This will also account for the uncertainty on the observed values, which is quantified by σ O . Hence, normalised centered root mean square error (ĈRMSE), normalised standard deviations (σ M andσ O ) are introduced by dividing CRMSE, σ M and σ O by σ O as shown in Equations (23)-(25), respectively. The value ofσ O was used as a reference point in the Taylor diagram:

Estimating the Significance of the 2020 Model on a Sewer Network
The new heat transfer coefficient and the d s /k s ratio were utilised to predict the viable potential heat recovery from a large 3000 sewer pipe network that serves 79,500 PE, which is described by Abdel-Aal et al. (2018) [5] in more detail. In order to account for the seasonal variations, the months of January, March and May were considered when estimating the potential heat recovery. It was assumed that heat was being recovered over a working day (24 h) in each of these months to estimate the consequent wastewater temperature variations across the network and at the wastewater treatment plant (WwTP) influent. Annual potential heat recovery was then estimated based on the daily predicted amounts. It was assumed that a viable heat recovery would result in minimum wastewater temperatures of 9 • C at the WwTP influent and 5 • C at the sewer pipes. This was based on previous work done by the aforementioned authors. The 2015 and 2020 models used the previously developed and the new heat transfer parameters, respectively, to observe the differences in the predicted amounts of viable heat recovery.

Results
This section shows the values of the calibrated parameters obtained by the optimisation routine using the field data. In order to assess the accuracy of the model, PDFs of modelling errors were plotted using data from sites 1, 2, 3, and 4. A Taylor diagram was plotted to show R,ĈRMSE andσ M for all sites over the four seasons for previous models, including the 2015 model and the new 2020 model results.

Calibration Results
The calibrated parameters for the heat transfer coefficient between wastewater and insewer air (f c ) and the calibrated ratio of soil penetration depth to soil thermal conductivity (d s /k s ) are shown in Table 2 for both sewer categories. The new heat transfer coefficient between wastewater and in-sewer air (h Fr ) was computed, using Equation (6), by considering f c , monthly average u w and D w values and the results are presented in Table 2. Unlike the 2015 model, f c and consequently h Fr values for the 2020 model have close values in winter and then decrease in summer. This is likely to be due to the dependency of the heat transfer coefficient in the 2015 model (h wa ) on parameters that were calibrated either by Flinspach (1973) [17], i.e., the 5.85 value, or by Abdel-Aal (2015) [9], i.e., the in-sewer air velocity, rather than being influenced by the physical variation in wastewater hydraulics, represented by the Fr number. Table 2. Values of the calibrated ratio of soil penetration depth to soil thermal conductivity (d s /k s ), calibrated heat transfer coefficient factor ( f c ) for heat transfer between wastewater and in-sewer air, the heat transfer coefficient (h Fr ) and the thermal resistivity values between wastewater and in-sewer (R wa ) and soil (R ws ) for the 2020 model. The probability density function (PDF) was plotted for modelled results using calibration data from site 1, for the residential sewer model, and from site 4, for the collector sewer model, as shown in Figure 8. The PDF plots show the distribution of modelling errors (Equation (15)) to represent the relative probability of error ranges when modelling wastewater temperatures.

Month k s /d s (W/m 2 ·K) f c (W/m 2 ·K) h Fr (W/m 2 ·K) R wa (m·K/W) R ws (m·K/W)
The distribution of errors for site 1 was relatively close for the majority of the months. However, May has shown a slightly skewed PDF for site 1, which may be due to steeper wastewater temperature variations (Figure 1) where calibration performance may have been compromised by such variations. The summer months (May-July) have shown lower h Fr values, which is likely to be caused by the variation in air thermal conductivity. It is clear that the collector sewer shows a more normal PDF of errors compared to the residential sewer, where the symmetry is almost around zero modelling error for site 4.

Validation Results
Observed data from the real sewer pipes (Table 1) was utilised to validate the 2020 model. Figure 9 shows the validation results by plotting a PDF of modelling errors using data from residential (site 2) and collector (site 3) sewers. Similar to the calibration results, the residential sewer PDF from the validation data was less symmetric than that of collector sewers ( Figure 9). However, site 3 presented symmetry at around −0.1 • C modelling error, which occurred in all months.

Assessing Model Performance through the Taylor Diagram
The correlation coefficient (R), normalised centered root mean square error (ĈRMSE) and normalised standard deviation (σ M ), for 2020 (developed in this paper), 2015 and Elías-Maxil models were all plotted on a Taylor diagram using equations 16 to 22, as shown in Figure 10. For clarity, an additional Taylor diagram was plotted in Figure 11 to include the TEMPEST model results. For each model, the statistics are shown for sites 1, 2, 3 and 4 in all months. The horizontal and vertical axes representσ M and the blue dashed contour lines reflect theĈRMSE values. The arc between the horizontal and vertical axes shows the correlations between observed and modelled values. Figure 10 shows that a model may have a relatively low correlation yet produce lowσ M andĈRMSE. A highly accurate model would show statistic values near the reference point. Some plotted values in Figure 10 overlap and cannot appear, e.g., 3 June. The normalised standard deviation (σ M ) was mostly close to 1 for the 2020 model, where 60% are ofσ M values are greater than 0.9 and less than 1.1. April, particularly at site 1, which was used for validating the 2015 model, presented lower correlation and higherĈRMSE values than other months. April 1 in the 2015 model had a CRMSE value of 0.41 • C and a low σ O of 0.35 • C, which results in a relatively highĈRMSE. April can be a challenging month for calibration since the variation of observed wastewater temperature was high, as shown by having the second-lowest σ O . Furthermore, 1 April in the 2015 model presented a relatively high variation of the modelled results, expressed by σ M . The performance of a model, in terms ofĈRMSE is influenced by the difference between σ O and σ M , the higher the difference, the higher theĈRMSE (Equation (18)), which indicates a larger model discrepancy. March 1 in the 2015 model had similar characteristics to 1 April of the same model, yet σ O was slightly higher than that of April and hence the relatively lowerĈRMSE. The highĈRMSE values are due to CRMSE being close to or above the observed standard deviation (σ O ) (Equation (23)). The 2020 model has generally shown less discrepancy than the 2015 model, as demonstrated by the average RMSE of 0.27 • C, compared to 0.35 • C for the 2015 model. Furthermore, the ranges of some key statistics were considerably reduced in the 2020 model. For example, The RMSE values ranged between 0.15 and 0.50 • C for all months at all sites in the 2020 model, which is a tighter range than its 2015 model equivalent where RMSE varied from 0.14 to 0.83 • C. The 2020 model also has shown less variation in the bias range (−0.14 to 0.40 • C) than that of the 2015 model (−0.50 to 0.74 • C). The three models, presented in Figure 10, seem more capable of predicting sites 3 and 4 compared to sites 1 and 2. The Elías-Maxil model has presented higher modelling errors than the 2020 and 2015 models with maximumĈRMSE of around 1.5, which is 50% higher than its equivalent in the 2020 model. Nevertheless, the Elías-Maxil model showed a much higher correlation andĈRMSE values than that of TEMPEST, as shown in Figure 11. This suggests that omitting some parameters from TEMPEST may improve the modelling accuracy.  The RMSE for each month, computed using Equation (17), are plotted in Figure 12 for all models. TEMPEST has shown relatively large RMSEs, and hence it has a separate vertical axis. It is clear from Figure 12 that the 2020 model showed the lowest RMSE in most months or slightly higher than the 2015 models, e.g., in February and March at site 3.

The Significance of the 2020 Model in Estimating Potential Heat Recovery
Comparing the outputs of the 2015 and 2020 models using the 3000 sewer-pipe network (Section 2.8), viable annual heat recovery from each model is shown in Table 3.

Discussion
Observed wastewater temperatures showed an expected overall pattern of variation in residential sewers (sites 1 and 2) and collector sewers (sites 3 and 4), where wastewater temperature patterns followed those of in-sewer air temperature. The wastewater flow rates in the residential sewers generally presented a low, stable DWF with some flow spikes that reflect the rainfall events. The collector sewers have shown different wastewater flow variations, where sites 3 and 4 presented steady drops in flow rates from February to April. This is likely to be caused by higher infiltration rates during the winter months. This type of long-term flow data, combined with temperature data, may be utilised in the future to make more accurate predictions of groundwater infiltration.
The residential sewers showed similar values of heat transfer coefficient between wastewater and in-sewer air (h Fr ) during March and April, which was followed by a noticeable drop, from 43 W/m 2 ·K in April to 7 W/m 2 ·K in May, reflecting the higher thermal resistivity between wastewater and in-sewer air (R wa ). This suggests that the behaviour of heat transfer between wastewater and in-sewer air in residential sites is season-dependent. The average in-sewer air temperature in residential sewers was around 20 • C in summer (May to July), which is approximately twice that in winter (March to April). Assuming the moisture content in sewers is constant, the lower temperature in winter can increase the relative humidity, as can be observed from a hygrometric chart. For example, the in-sewer air relative humidity of 50% at 20 • C would increase to 90% if the temperature drops to 10 • C. The higher relative humidity can increase the air thermal conductivity, which may explain the larger winter values for the calibrated parameter (f c ) in residential sewers. Since Fr was almost constant in all months of each sewer category, the variation of heat transfer coefficient (h Fr ) is directly proportional to f c . Unlike the urban sewers, the collector sites have generally presented more consistent f c values in winter and summer. The average in-sewer air temperature in collector sewers was around 11 • C in summer, excluding June, and 9 • C in winter, which limits the variation in relative humidity and hence the closer f c values. However, June has shown a slightly lower f c value in collector sewers, which may be due to the relatively high in-sewer air temperature of around 15 • C. It is also worth noting that in-sewer air has less influence on the heat exchange in collector sewers than that in residential sewers as the ratio of wastewater depth to diameters in the former is around 40% compared to approximately 5% in residential sewers.
An additional process that may be necessary to consider when simulating wastewater temperatures is the heat transfer between the in-sewer air and the ambient air. This could particularly improve the modelling accuracy in April at residential sewers, which resulted in the largest modelling discrepancy in the 2015, 2020 and Elías-Maxil models. Preliminary calibration of the thermal resistivity between in-sewer air and ambient air (R aa ) in April, using the 2020 model at urban sewers, was performed to investigate the impact of the heat transfer within the air on the modelling performance. It was estimated that if R aa was around 0.001 m.K/W and the temperature difference between in-sewer air and ambient air was 3 • C, the 2020 model accuracy, represented by the bias, can be improved by around 50%. However, more details, e.g., local ambient temperatures and in-sewer air pressure, is needed to perform a more in-depth calibration of R aa . Unfortunately, this level of detail was not available and is beyond the scope of the paper. A simplified heat transfer model in a partially filled pipe can be created in the future, e.g., using a computational fluid dynamic (CFD) software package, to investigate impact of heat exchanged between in-sewer and ambient air.
The calibrated ratios of soil thermal conductivity to penetration depth (k s /d s ) were often larger in colder months in residential and collector sewers. The high k s /d s (or low d s /k s ) ratios result in low thermal resistivity between wastewater and soil (R ws ) (Equation (10)). Higher soil thermal resistivity was observed in residential sewers during warmer months (Table 2). This is likely since the temperature difference between wastewater and soil was generally higher, by 3 • C than that between wastewater and in-sewer air. Hence, the energy balance (Equation (7)) would imply higher soil thermal resistivity. Nevertheless, key parameters related to soil thermal resistivity, e.g., soil structure, type, void size and saturation level, were unavailable which made it difficult to confirm the exact physical causes of soil thermal resistivity variations. The high flow rate in collector sewers indicates high thermal energy, which means more heat needs to be dissipated, thus the lower thermal resistivity in collector sites.
The 2020 model has generally improved the accuracy of predicting wastewater temperatures, as proven by the majority of modelled values showing correlation coefficients higher than 0.9 ( Figure 10). TheĈRMSE values of the 2020 model were relatively low, with around 60% being equal to or less than 0.4, compared to 40% in the case of the Elías-Maxil model. AĈRMSE value of 0.4 is considered low for the purpose of predicting wastewater temperatures; for example, March at site 3 (March 3) has aĈRMSE of 0.36, and a corresponding CRMSE of 0.10 • C, which is close to the measurement accuracy. Therefore, reducing the model error would require more accurate sensors. The CRMSE has the advantage of representing the absolute error of the model, which is independent of the overall bias. Thus, when neglecting the bias between modelled and observed temperatures, e.g., in May at site 1 of the 2020 model, the CRMSE was relatively low (0.3 • C) and hence showed better accuracy than that interpreted by the PDF plot ( Figure 8).
The RMSE of the 2020 model ranged between 0.15 and 0.5 • C, which was relatively low in comparison to the 2015 model (0.12 to 0.87 • C), Elías-Maxil (0.16 to 0.8 • C) and TEMPEST (0.8 • C to 7.8 • C). Overall, the 2020 model average RMSE (among all months) was 0.27 • C, which is slightly lower than that of previous models, e.g., 0.35 • C in the 2015 model, 0.42 • C in Elías-Maxil. Nevertheless, the Elías-Maxil model presented higher absolute errors, as shown by the highĈRMSE, of up to 1.5 in comparison to 1 in the 2020 model, and low correlation values that reached 0.25 in comparison to 0.72 in the 2020 model. The TEMPEST model has presented the highest overall RMSE values with an overall average of 2.0 • C. It is worth noting that TEMPEST was validated on shorter datasets and limited time periods that excluded months like April, which particularly presented higher variation in observed temperatures.
Overall, the 2020 model has presented more accurate results that showed considerable differences when predicting the potential viable heat recovery in comparison to the 2015 model. Abdel-Aal (2015) [9] performed a sensitivity analysis on all model parameters for urban and large collector sewers. The sensitivity of the model to d s and k s is much less (e.g., by 50%) than that of the heat transfer coefficient between wastewater and in-sewer air. Therefore, and through close examination of the model behaviour, the reason behind the improved accuracy is the new relation developed in this paper.
The 2020 model requires few input parameters that are accessible, measurable or predictable through existing hydrodynamic models; these are mainly the upstream wastewater temperature, depth and velocity in addition to the in-sewer air temperature. This is important for the model's practicality and deployment in large sewer networks, which can save considerable computational time. Enhancing the model's practicality facilitates its implementation in existing hydrodynamic software packages.
Collector sewers have generally shown lower errors than residential sewers. This was expected because of the higher wastewater flow rates, around 35-fold, than that of residential sewers. The higher flow rates demand more energy to be dissipated from wastewater to drop its temperature by 1 • C. This is due to the heat capacity of wastewater, i.e., the multiplication product of specific heat capacity (c p,w ) and the wastewater mass flow rate. Therefore, the model is less sensitive to the calibrated parameters when collector sewers' data were utilised for calibration, using either model, in comparison to that of the residential sewers. This explains the lower collector sewerĈRMSE values and the more symmetric collector sewer PDF plots.

Conclusions
A novel relationship to estimate the heat transfer coefficient between wastewater and in-sewer air was developed to better predict the behaviour of heat exchange mechanisms at the wastewater-air boundary in combined sewers. Case study results predicting temperatures with the novel relationship indeed indicate the importance of including in-sewer air-water boundary heat transfer. The new relations enabled more accurate predictions of the potential heat recovery from sewers, which can result in estimating up to 25% additional viable heat recovery in winter.
The new heat transfer coefficient between wastewater and in-sewer air (h Fr ), which ranged between 5 and 58 W/m 2 ·K, was estimated and calibrated as a function of the wastewater Froude number using seasonal data from small residential and large collector sewers. The model was validated on larger independent datasets for both residential and collector sewers. The calibrated parameter for the heat transfer coefficient between wastewater and in-sewer air (f c ) had a stable value for periods of several months, indicating that the proposed parameterisation was an improvement on previous approaches and appeared to work well in both residential and collector sewers. Warmer weather presented lower but stable f c values to reflect the change in the heat transfer processes during warmer periods. The new parameterisation improved the modelling accuracy as presented by the absolute errors and is based on available hydrodynamic data that can facilitate the simulation process.
A normalised Taylor diagram was effective in summarising modelling errors graphically and comparing between various models, where the majority of the 2020 model results were close to the observed reference point, which demonstrates a generally high modelling accuracy. This accuracy is adequate for the purpose of the 2020 model, which is to assess the viability of heat recovery from individual sewers in large sewer networks. The PDF of the 2020 model results has shown over-prediction and underprediction of wastewater temperatures, which can overall minimise the modelling error when it is implemented to assess the viability of heat recovery from large networks.
Observed wastewater temperatures show seasonal variations in both residential and collector sewers and showed wastewater temperature patterns following air temperature patterns. Hence, considering seasonal data as well as in-sewer air-water boundary heat transfer is important when studying the process of heat transfer in sewers.
Although the modelling error in April at urban sewers was reduced by the 2020 model, it has particularly shown the largest discrepancy, and this requires further investigations. It is envisaged that additional processes, such as heat transfer between in-sewer air and ambient air, may need to be considered in order to further improve the modelling accuracy.  Institutional Review Board Statement: Ethical review and approval were waived for this study, due to studies not involving humans or animals.