Simplified two-step, cross-sectional approach for coupled heat and moisture transfer modeling of shallow, horizontal, ground-based heat Exchangers

This study evaluates a two-step, cross-sectional approach for designing shallow, unsaturated, horizontal, ground-based heat exchangers (GHXs). Numerical modeling was conducted for coupled heat and moisture transfer around GHXs under transient climatological conditions through a cooling season. Soil samples were collected in Texas and Wisconsin to measure thermal conductivity dry-out curves and soil-water characteristic curves for use in the modeling. Average daily meteorological data from central Texas was applied for the top boundary condition. Heat generation from the GHX was calculated by total condenser heat rejection from an environmental control unit based on ambient temperatures and unit specifications. In the first step of the modeling, results showed that the fluid temperature through the GHX loop was approximately 43 ºC, and rapid heat and moisture fluxes were observed around the GHX loops. High moisture flux along the upper surface was also observed due to high ambient temperatures that occur during the summer season. Using these results, exiting temperature of the GHX was estimated for the second cross-sectional modeling step. This two-step, cross-sectional modeling approach provides a systematic analysis of coupled heat and moisture transfer around shallow, horizontal, unsaturated GHX loops, thus simplifying high computational effort needed for full three-dimensional modeling of shallow GHX systems.


Introduction
As ground has a relatively constant temperature with an increase in depth from the surface [1][2][3], the ground is used as heat source and/or sink [4][5][6][7]. Ground-based heat exchanger (GHX) in a ground-source heat pump (GSHP) system has been widely used for cooling and/or heating residential and small commercial buildings using the constant ground temperature benefiting various aspects (e.g., reliability, environmental-friendly, cost-effective). To obtain stored geothermal energy in the earth during heating season, fluid (e.g., water, anti-freeze solution) is circulated through GHX because the ground temperature during winter is relatively higher than the GHX fluid temperature. On the other hand, high GHX fluid temperature during the summer cooling season is transferred to the ground, which has a relatively lower temperature. The fluid is injected from a GSHP to the GHX in closed loop system, while GHX in open loop system directly uses groundwater for heat exchange.
In a closed-loop GHX system, there are two typical GHX configurations: vertical and horizontal. Vertical GHX is installed in a deep borehole (typically 20 m to 300 m in depth) where the constant and higher geothermal energy are available to use, and horizontal GHX is buried in a shallow trench (typically 1 m to 5 m in depth) where the ground temperature fluctuates depending on environmental and seasonal conditions. Despite the potential for environmental and seasonal impacts on performance and efficiency of the horizontal GHX, the horizontal GHX has been widely used due to its easier installation process and lower trench excavation cost than the borehole drilling of the vertical GHX.
To experimentally estimate the thermal performance of a GHX system, thermal response tests (TRTs) are often conducted at a field site [8][9][10][11][12]. During TRTs, heat transfer in soil is quantified based on a heating load, exchanger flow rate, and fluid temperatures at the inlet and outlet, and the resulting values are indicative of the bulk thermal performance of geothermal exchange well system. Laboratory tests more directly and discretely quantify soil thermal properties [13][14][15]. Laboratory testing results indicate that soil thermal properties are dependent on the physical properties of the soil or rock specimen in question (e.g., moisture content, density, mineralogy) and broadly indicate that GHX performance can be improved with thermally enhanced and/or engineered backfill. More recently, [16,17] have reported that thermal conductivity and degree of saturation of the soil surrounding GHX decreases significantly or even permanently (i.e., not recovering during winter) when the GHX operates, due to coupled heat and moisture transfer from the GHX, and this potentially indicates that the performance can be fluctuated during the operation. Many numerical and analytical modeling studies have also indicated that performance and efficiency of horizontal GHX systems are affected by variable conditions, such as environmental and seasonal effects on the top surface of the porous medium and specifications and configurations of the GHX, as well as soil thermal properties [6,7,9,[18][19][20][21][22][23][24][25].
Although three-dimensional (3D) numerical modeling [6,7,9,19,20,[22][23][24] provides more realistic and accurate model solutions than two-dimensional (2D) modeling, computational effort involved in 3D modeling is significantly higher than 2D modeling due to increased complexity and interrogated domain. According to Wu et al. (2015), 2D modeling provides reasonable results with a benefit of reduced modeling time as 2D model domain represents the mid-section of the GHX system. This study proposes a simplified two-step, cross-sectional approach to simplify complexity of 3D model solutions based on 2D modeling of coupled heat and moisture transfer around horizontal GHXs along two distinct cross sections. For the first step, coupled heat and moisture transfer at the first cross section is modeled using the excess quantity of heat rejected from the operation of an environmental control unit. Then, the heat and moisture transfer at a second cross section is evaluated based on GHX temperatures estimated using modeling results at the previous section in addition to the fluid properties.

Materials
Two soil samples were collected from sites in central Texas and Wisconsin and classified as silty sand (SM) and poorly graded sand (SP), respectively, by ASTM D422 [26] and ASTM D2487 [27]. Fig. 1 presents the grain-size distribution curves for the two samples. Index properties, soil particle diameters corresponding to 10% (D10) and 50% finer (D50), coefficient of uniformity (Cu), and coefficient of curvature (Cc), are tabulated in Table 1. Maximum dry unit weight (γdmax) of the two samples were also measured by standard Proctor tests (ASTM D698). The physical properties of SM and SP were used to represent the ground and backfill model domains, respectively.  The soil-water characteristic curve (SWCC) describes the relationship between matric suction and moisture content, and the thermal conductivity dry-out curve (TCDC) defines the relationship between soil thermal conductivity (λ) and moisture content (commonly represented as volumetric or gravimetric water content, or as degree of saturation). As the SWCC and TCDC provide required information for modeling coupled heat and moisture transfer in unsaturated soil [24,[28][29][30], these relationships for SP and SM samples were measured using multiple specimen and modified Tempe cell tests, which have been widely used by previous researchers [14][15]. In the multiple specimen tests, the soil thermal conductivity was repeatedly measured at various moisture contents and soil densities (e.g., 3% moisture content compacted to 90% γdmax). The modified Tempe cell testing was used to continuously measure SWCC and TCDC from fully saturated condition to dry condition. As residual moisture is primarily retained on the particle surfaces as water of hydration [29], very large increase in matric suction is thus required to remove further moisture. Accordingly, once volumetric water content (VWC) reached around 0.01, the testing was stopped, and thermal conductivity at fully dried condition (oven dry at 105 °C for 24 h) was separately measured and added as an ending point. The measured SWCCs and TCDCs were then calibrated and validated with SWCCs and TCDCs previously reported in the literature, such as [14,15,28], to represent typical hydraulic and thermal properties of unsaturated soils at a field site. Fig. 2 represents the approximated SWCCs and TCDCs used as material properties in the modeling. When degree of saturation decreases most notably below critical degree of saturation where liquid bridges among particles break down, matric suction increases and thermal conductivity decreases.
These results indicate that performance and efficiency of the GHX may be significantly affected by the rapid changes in backfill hydraulic and thermal properties.

Governing Equations
SVHeat and SVFlux in the SVOffice™ finite-element software package (SoilVision Systems Ltd.) were used to numerically model coupled heat and moisture transfer around GHXs. Equation 1 demonstrates the governing equation for moisture and vapor flow driven by water pressure and temperature gradients in terms of conductivity terms. In equation 2, heat transfer is governed by conduction, convection of pore fluids, and latent heat associated with phase change.
where is hydraulic conductivity, is vapor conductivity by diffusion, is molecular weight of vapor (0.018016 kg·mol -1 ), is pore water pressure due to hydraulic head, R is universal gas constant (9.314 J·molK -1 ), * is vapor diffusivity, T is temperature, is soil-water vapor pressure, is unit weight of water, is volumetric water content, is volumetric ice content, is water density, is ice density, y is elevation, t is time, is thermal conductivity, is volumetric latent heat of water vaporization or condensation, is volumetric heat capacity of water, is water flow flux, C is volumetric heat capacity of soil, is volumetric latent heat of fusion of water, and is slope of the soilfreezing characteristic curve, SFCC, which is a relationship between freezing soil temperature (or suction) and unfrozen water content.

Model Geometry and Boundary Conditions
Dimensions of the model domain were 12-m width and 12-m depth to sufficiently observe coupled heat and moisture transfer around the GHXs (Fig. 3). A shallow trench (width = 3.048-m, depth = 3.658-m) was located at center of the model domain and included engineered backfill (3.048-m width and 0.66-m height) surrounding five GHXs (diameter = 15 cm, spacing = 75 cm) on the bottom of the trench at the center. Lengths of the GHXs were assumed as 50 m. Physical properties of SM and SP were applied for the ground and backfill, respectively. Initial gravimetric water contents of the ground and backfill were 27% and 22%, respectively. Air temperature, solar radiation, wind speed, relative humidity, and precipitation data in central Texas were collected from regional airports, and the average daily data was assigned for meteorological conditions on the top surface (Fig. 4). An energy balance was applied for the thermal flux balance at the ground surface (i.e., Qboundary = QN -QH -QL, where QN = net radiation flux, QH = sensible heat flux, QL = latent thermal flux). To simulate symmetrical heat and moisture transfer, zeroflux boundary conditions were applied at the left and right sides of the domain. Initial ground temperature at any depth, y, was also simulated using a pure harmonic function as follows [1]: where Ta is average ground temperature, A0 is annual amplitude of the ground temperature, y is depth, α is soil thermal diffusivity, ω is annual angular frequency, t is time, and t0 is time lag from an arbitrary starting date (usually January 1) to the occurrence of the minimum temperature in year.

Heat Generation in GHXs
An environmental control unit (ECU) rejects excess heat from the space being conditioned to the fluid in the GHX loop via a condenser (i.e., cooling unit). For heat generation conditions used in the five GHX model domains, total condenser heat rejection from a 5-ton ECU was evaluated based on ambient temperature, solarconducted load, compressor heat, and supplemental heat [ Fig. 5(a)]. Average load data represented as a solid line in Fig. 5(a) was combined with the air temperature in central Texas. When the ambient temperature exceeds 26.67 °C (80 °F), the ECU is operated thus resulting in substantial fuel savings, particularly during the peak summer cooling season [ Fig. 5(b)].
In this study, it is assumed that the five GHXs have the same heat generation conditions. However, with future work, these heat generation conditions will be improved with additional modeling based on computational fluid dynamics (CFD) that investigates how much heat in each GHX dissipates through the backfill (i.e., significant information on fluid temperatures in each GHX depending on the heat dissipation through backfill).

Two-Step Cross-Sectional Approach
This study applied a two-step, cross-sectional approach using two-dimensional (2D) modeling to simplify computational effort. For the first step, a mid-section of a three-dimensional GHX system was modeled in 2D. The cross section was assumed as an inlet of 1-m 3 cube. Second, the fluid temperature at an outlet of the cube was estimated based on a model equation. Using the estimated outlet temperature (i.e., inlet temperature of the next cross section), the second cross section was modeled. Fig. 6 shows modeling results at the 200 th day, which represents one of the worst-case scenarios (i.e., high ambient temperature and heat generation). As the fluid temperature through the GHX piping increased with ambient temperature (about 43 °C fluid temperature at the 200 th day), coupled heat and moisture transfer rapidly occurred near the ground loop piping [see vectors in Fig.  6(a) and (b)]. As a result, VWC around the GHXs decreased below 0.05 m 3 /m 3 that soil thermal conductivity significantly decreased (see TCDCs in Fig.2). In other words, a dry zone formed near the GHX pipes and performed as an insulator due to the decreased thermal conductivity, resulting in a rapid temperature change [see the rapid color change from light pink color to orange color in Fig 6(a)]. In addition, the high ambient temperature (about 30.2 °C on the 200 th day) led to coupled heat and moisture transfer near the ground surface. For example, red vectors in Fig. 6(a) indicate rapid evaporative moisture migration toward the ground surface. The decrease in moisture content near the GHX piping and moisture evaporation near the ground surface, especially during summer, imply that performance and efficiency of the ECU system will decrease due to the decrease in soil thermal conductivity.

Coupled Heat and Moisture Transfer at 2 nd Cross-Section
Wu et al. (2015) evaluated the length of horizontal GHXs using simplified heat transfer model based on heat pump capacity, heat pump efficiency, fluid flow rate, entering and leaving water temperatures, and thermal resistances (e.g., pipe, soil). In this study, one of the equations used in [24] was adopted to estimate fluid temperature at the outlet. Equation 4 describes a relationship between the heat exchanger load and fluid properties.
where is heat exchanger load, is fluid flow rate, is specific heat capacity the fluid at constant pressure, and ∆ is temperature difference of the fluid (i.e., entering water temperature -leaving water temperature).
As described in Fig. 6 above, the GHX fluid temperature on the 200 th day was 43 °C under of 151,900 kJ/day·m 3 . That is, if the 43 °C is assumed as entering water temperature at the inlet of a 1-m 3 cube, leaving water temperature at the outlet of a 1-m 3 cube can be derived using equation 4 and known and parameters. Table 2 summarizes the estimated entering (Tin) and leaving (Tout) fluid temperatures at 10 crosssections when = 0.584 kg/s and = 3.923 kJ/kg·°C obtained from [24]. Table 2. Entering (Tin) and leaving (Tout) fluid temperatures at 10 cross-sections. Due to heat dissipation from the GHX piping to the backfill, the fluid temperature at the outlet (i.e., Tout) decreased with an increase in the sections (about 2% every section). However, these estimations of Tout using equation 4 include uncertainties. For example, the equation does not consider heat dissipation through the backfill based on soil thermal conductivity. Further study can be extended with CFD to improve the simplified equation for more accurate and mechanistically sound Tout estimations at the n th cross section, especially considering soil thermal properties.  transfer was slightly less than the heat and moisture transfer at the first cross section (e.g., less evaporative moisture migration described as red vectors). Despite of the uncertainties of equation 4, this two-step crosssectional approach can be used to predict coupled heat and moisture transfer at the n th cross section.

Summary and Conclusions
This study proposes a two-step, cross-sectional approach for 2D coupled heat and moisture transfer modeling of a five-pipe horizontal GHX loop at two different crosssections.
Laboratory tests with literature support/validation were conducted to develop the physical, hydraulic, and thermal properties of the ground and backfill model domains. Environmental and meteorological data in central Texas were collected for application of the surface boundary condition. The collected ambient temperature was also used to evaluate total condenser heat rejection used for heat generation conditions in the GHX domains. In the first step model result, coupled heat and moisture transfer was significant around the GHX pipes carrying the 43 °C fluid. Based on the modeling results and fluid properties, the leaving water temperature at the outlet of a 1-m 3 cube was evaluated using a simplified heat transfer model equation. This two-step, cross-sectional approach provides significant information to predict how coupled heat and moisture transfer occurs at n th cross-sections, as well as simplification of high computational effort associated with full 3D modeling. However, the simplified model equation to estimate the outlet fluid temperature has limitations, such as a lack of consideration of heat dissipation through the soil. Experiments and/or 3D CFD modeling are recommended to continue exploration and validation of this simplified two-step methodology and to improve the model equation, potentially including additional terms that describe the specific effect of soil thermal properties on heat dissipation.