Calculation of g-functions with mixed parallel and series connectionsΒΆ

This example demonstrates the use of the g-function module and the networks module to calculate g-functions considering the piping connections between the boreholes, based on the method of Cimmino 1. For boreholes connected in series, it is considered the outlet fluid temperature of the upstream borehole is equal to the inlet fluid temperature of the downstream borehole. The total rate of heat extraction in the bore field is constant.

The following script generates the g-functions of a field of 5 equally spaced borehole on a straight line and connected in series. The boreholes have different lengths. The g-function considering piping conections is compared to the g-function obtained using a boundary condition of uniform borehole wall temperature.

The script is located in: pygfunction/examples/mixed_inlet_conditions.py

  1# -*- coding: utf-8 -*-
  2""" Example of calculation of g-functions using mixed inlet temperatures.
  3
  4    The g-functions of a field of 5 boreholes of different lengths connected
  5    in series are calculated for 2 boundary conditions: (a) uniform borehole
  6    wall temperature, and (b) series connections between boreholes. The
  7    g-function for case (b) is based on the effective borehole wall
  8    temperature, rather than the average borehole wall temperature.
  9
 10"""
 11import matplotlib.pyplot as plt
 12import numpy as np
 13from matplotlib.ticker import AutoMinorLocator
 14from scipy import pi
 15
 16import pygfunction as gt
 17
 18
 19def main():
 20    # -------------------------------------------------------------------------
 21    # Simulation parameters
 22    # -------------------------------------------------------------------------
 23
 24    # Borehole dimensions
 25    D = 4.0             # Borehole buried depth (m)
 26    # Borehole length (m)
 27    H_boreholes = np.array([75.0, 100.0, 125.0, 150.0, 75.0])
 28    H_mean = np.mean(H_boreholes)
 29    r_b = 0.075         # Borehole radius (m)
 30    B = 7.5             # Borehole spacing (m)
 31
 32    # Pipe dimensions
 33    r_out = 0.02        # Pipe outer radius (m)
 34    r_in = 0.015        # Pipe inner radius (m)
 35    D_s = 0.05          # Shank spacing (m)
 36    epsilon = 1.0e-6    # Pipe roughness (m)
 37
 38    # Pipe positions
 39    # Single U-tube [(x_in, y_in), (x_out, y_out)]
 40    pos_pipes = [(-D_s, 0.), (D_s, 0.)]
 41
 42    # Ground properties
 43    alpha = 1.0e-6      # Ground thermal diffusivity (m2/s)
 44    k_s = 2.0           # Ground thermal conductivity (W/m.K)
 45
 46    # Grout properties
 47    k_g = 1.0           # Grout thermal conductivity (W/m.K)
 48
 49    # Pipe properties
 50    k_p = 0.4           # Pipe thermal conductivity (W/m.K)
 51
 52    # Fluid properties
 53    m_flow_network = 0.25   # Total fluid mass flow rate in network (kg/s)
 54    # All boreholes are in series
 55    m_flow_borehole = m_flow_network
 56    # The fluid is propylene-glycol (20 %) at 20 degC
 57    fluid = gt.media.Fluid('MPG', 20.)
 58    cp_f = fluid.cp     # Fluid specific isobaric heat capacity (J/kg.K)
 59    rho_f = fluid.rho   # Fluid density (kg/m3)
 60    mu_f = fluid.mu     # Fluid dynamic viscosity (kg/m.s)
 61    k_f = fluid.k       # Fluid thermal conductivity (W/m.K)
 62
 63    # g-Function calculation options
 64    nSegments = 8
 65    options = {'nSegments': nSegments,
 66               'disp': True,
 67               'profiles': True}
 68    # The similarities method is used since the 'equivalent' method does not
 69    # apply if boreholes are connected in series
 70    method = 'similarities'
 71
 72    # Geometrically expanding time vector.
 73    dt = 100*3600.                  # Time step
 74    tmax = 3000. * 8760. * 3600.    # Maximum time
 75    Nt = 25                         # Number of time steps
 76    ts = H_mean**2/(9.*alpha)   # Bore field characteristic time
 77    time = gt.utilities.time_geometric(dt, tmax, Nt)
 78
 79    # -------------------------------------------------------------------------
 80    # Borehole field
 81    # -------------------------------------------------------------------------
 82
 83    boreField = []
 84    bore_connectivity = []
 85    for i, H in enumerate(H_boreholes):
 86        x = i*B
 87        borehole = gt.boreholes.Borehole(H, D, r_b, x, 0.)
 88        boreField.append(borehole)
 89        # Boreholes are connected in series: The index of the upstream
 90        # borehole is that of the previous borehole
 91        bore_connectivity.append(i - 1)
 92
 93    # -------------------------------------------------------------------------
 94    # Initialize pipe model
 95    # -------------------------------------------------------------------------
 96
 97    # Pipe thermal resistance
 98    R_p = gt.pipes.conduction_thermal_resistance_circular_pipe(
 99        r_in, r_out, k_p)
100    # Fluid to inner pipe wall thermal resistance (Single U-tube)
101    m_flow_pipe = m_flow_borehole
102    h_f = gt.pipes.convective_heat_transfer_coefficient_circular_pipe(
103        m_flow_pipe,  r_in, mu_f, rho_f, k_f, cp_f, epsilon)
104    R_f = 1.0/(h_f*2*pi*r_in)
105
106    # Single U-tube, same for all boreholes in the bore field
107    UTubes = []
108    for borehole in boreField:
109        SingleUTube = gt.pipes.SingleUTube(
110            pos_pipes, r_in, r_out, borehole, k_s, k_g, R_f + R_p)
111        UTubes.append(SingleUTube)
112    network = gt.networks.Network(
113        boreField, UTubes, bore_connectivity=bore_connectivity,
114        m_flow_network=m_flow_network, cp_f=cp_f, nSegments=nSegments)
115
116    # -------------------------------------------------------------------------
117    # Evaluate the g-functions for the borefield
118    # -------------------------------------------------------------------------
119
120    # Calculate the g-function for uniform temperature
121    gfunc_Tb = gt.gfunction.gFunction(
122        boreField, alpha, time=time, boundary_condition='UBWT',
123        options=options, method=method)
124
125    # Calculate the g-function for mixed inlet fluid conditions
126    gfunc_equal_Tf_mixed = gt.gfunction.gFunction(
127        network, alpha, time=time, boundary_condition='MIFT', options=options,
128        method=method)
129
130    # -------------------------------------------------------------------------
131    # Plot g-functions
132    # -------------------------------------------------------------------------
133
134    ax = gfunc_Tb.visualize_g_function().axes[0]
135    ax.plot(np.log(time/ts), gfunc_equal_Tf_mixed.gFunc, 'r-.')
136    ax.legend(['Uniform temperature', 'Mixed inlet temperature'])
137    plt.tight_layout()
138
139    # For the mixed inlet fluid temperature condition, draw the temperatures
140    # and heat extraction rates
141    gfunc_equal_Tf_mixed.visualize_temperatures()
142    gfunc_equal_Tf_mixed.visualize_temperature_profiles()
143    gfunc_equal_Tf_mixed.visualize_heat_extraction_rates()
144    gfunc_equal_Tf_mixed.visualize_heat_extraction_rate_profiles()
145
146    return
147
148
149# Main function
150if __name__ == '__main__':
151    main()

References

1

Cimmino, M. (2018). g-Functions for bore fields with mixed parallel and series connections considering the axial fluid temperature variations. Proceedings of the IGSHPA Sweden Research Track 2018. Stockholm, Sweden. pp. 262-270.