## 2 Analytical Methods
### Use of Symmetry: From Single to Multiple Boreholes  
Our developed solution for a single borehole can be extended to multiple boreholes by summing up each and every borehole (a.k.a., thermal superposition), as long as the borholes are placed in a symmetrical configuration, as seen in Figure 2. Note that \"_thermal superposition_\" should not be confused with the \"_method of superposition_\" in Section 2.3: the former is to add the thermal contribution of each borehole, and the latter is the summation of transient and steady-state parts of a PDE to remove its non-homogeneity. The thermal superposition is generally expressed as:  
\[T_{total}=\sum_{i=1}^{N_{total}}T_{i}=T_{1}+T_{2}+T_{3}+...+T_{N_{total}}, \tag{11}\]  
where \(T_{total}\) is the temperature distribution for a total of \(N_{total}\) boreholes. \(T_{i}\) represents the temperature of the \(i\)-th borehole; the determination of each borehole is based on both the general solution, Equation (10), and its location. In this way, the 1D radial temperature solution for a single borehole builds up a 2D solution for multiple boreholes.  
Recall that a few arrangements of boreholes are studied as shown in Figure 2: single borehole, double boreholes, and multiple boreholes (3-by-3 and N-by-N configurations. Setting 3-by-3 boreholes as an example here). In this study, we are particularly interested in calculating the temperature profile along the center line at \(y=0\), which is equivalent to the line at \(x=0\) due to its symmetry. That is, \(T_{total}(x,y=0,t)\). \(T_{total}\) is denoted as the temperature profile of all boreholes in the domain.  
When there is only one borehole in the system, the total temperature profile along the center line is mirrored at \(x=0\). Our previously calculated 1-D radial transient temperature solution can be readily implemented as follows:  
\[T_{total}(x,y=0,t)=T_{1}(r=x,t), \tag{12}\]  
where the subscript 1 indicates the first and the only borehole in this scenario.  
The case of double boreholes means that two boreholes are placed at the center of the domain with a distance of \(d\) apart. The total temperature profile along the center line can be calculated by summing the 1D radial transient temperature solution for each borehole:  
\[T_{total}(x,y=0,t)=T_{1}(r=x+\frac{d}{2},t)+T_{2}(r=x-\frac{d}{2},t). \tag{13}\]  
When more than two boreholes are present in the system (e.g., four or nine boreholes), there will be boreholes above and below the center line contributing \(T_{total}\). Consequently, the Pythagorean theorem needs to be applied for examining the correct radial distance from the upper and lower boreholes. In the case of nine boreholes, we find:\[T_{total}(x,y=0,t)\] \[= T_{1}(r=\sqrt{d^{2}+(x+d)^{2}},t)+T_{2}(r=\sqrt{d^{2}+x^{2}},t)+T_{3 }(r=\sqrt{d^{2}+(x-d)^{2}},t)\] \[+T_{4}(r=x+d,t)+T_{5}(r=x,t)+T_{6}(r=x-d,t)\] \[+T_{7}(r=\sqrt{d^{2}+(x+d)^{2}},t)+T_{8}(r=\sqrt{d^{2}+x^{2}},t)+T _{9}(r=\sqrt{d^{2}+(x-d)^{2}},t) \tag{14}\] \[= 2\times\left[T_{1}(r=\sqrt{d^{2}+(x+d)^{2}},t)+T_{2}(r=\sqrt{d^{ 2}+x^{2}},t)+T_{3}(r=\sqrt{d^{2}+(x-d)^{2}},t)\right]\] \[+T_{4}(r=x+d,t)+T_{5}(r=x,t)+T_{6}(r=x-d,t).\]  
The boreholes are labeled from the top left corner to the bottom right corner, one row after another. Since the top boreholes (No. 1, 2, 3) and bottom boreholes (No. 7, 8, 9) contribute to the center-line temperature in the same manner, we can double the top ones for simplicity, as expressed in Equation (14).  
It is rather clear to observe that the closed-form analytical solutions to single, double, and multiple GBS have simple algebraic expressions, without any requirement of temporal and spatial iterations that usually appear in numerical models. As seen in Equation (13) and Equation (14), increasing the number of boreholes only adds more terms in the thermal superposition, which can be rapidly computed and do not add much computational time or power.