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

References