Calculation of g-functions with equal inlet fluid temperatureΒΆ

This example demonstrates the use of the g-function module and the pipes module to calculate g-functions using a boundary condition of equal inlet fluid temperature into all boreholes, based on the method of Cimmino 1. The total rate of heat extraction in the bore field is constant.

The following script generates the g-functions of a rectangular field of 6 x 4 boreholes. The g-function using a boundary condition of equal inlet fluid temperature is compared to the g-functions obtained using boundary conditions of uniform heat extraction rate and of uniform borehole wall temperature.

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

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

References

1

Cimmino, M. (2015). The effects of borehole thermal resistances and fluid flow rate on the g-functions of geothermal bore fields. International Journal of Heat and Mass Transfer, 91, 1119-1127.