MassimoCimmino
commented
4 days ago A system of equations can be solved to obtain the heat extraction rate and the borehole wall temperature every time step:
g_dt = gFunc.gFunc[0]
a_in, a_b = network.coefficients_network_heat_extraction_rate(
m_flow_network, cp_f, 1)
a_in = a_in.item()
# The effective borehole wall temperature is applied to all boreholes
a_b = np.sum(a_b)
# Solve the system of equations:
# 1) Q_tot = a_in * T_f_in + a_b * T_b
# 2) T_b = T_b0 - Q_tot / (2*pi*k_s*H*nB) * g(dt)
denom = (2 * pi * k_s * H * nBoreholes)
T_b[i] = (T_b0 - a_in * T_f_in[i] * g_dt / denom) / (
1 + a_b * g_dt / denom)
Q_tot[i] = a_in * T_f_in[i] + a_b * T_b[i]Below is a modified version of the example fluid_temperature_multiple_boreholes with an inlet fluid temperature signal.
# -*- coding: utf-8 -*-
""" Example of simulation of a geothermal system with multiple boreholes.
The g-function of a bore field is calculated for boundary condition of
mixed inlet fluid temperature into the boreholes. Then, the borehole
wall temperature variations resulting from a time-varying inlet fluid
temperature profile are simulated using the aggregation method of Claesson
and Javed (2012). Predicted outlet fluid temperatures of double U-tube
borehole are evaluated.
"""
import matplotlib.pyplot as plt
import numpy as np
from scipy.constants import pi
import pygfunction as gt
def main():
# -------------------------------------------------------------------------
# Simulation parameters
# -------------------------------------------------------------------------
# Borehole dimensions
D = 4.0 # Borehole buried depth (m)
H = 150.0 # Borehole length (m)
r_b = 0.075 # Borehole radius (m)
# Bore field geometry (rectangular array)
N_1 = 6 # Number of boreholes in the x-direction (columns)
N_2 = 4 # Number of boreholes in the y-direction (rows)
B = 7.5 # Borehole spacing, in both directions (m)
# Pipe dimensions
r_out = 0.0211 # Pipe outer radius (m)
r_in = 0.0147 # Pipe inner radius (m)
D_s = 0.052 # Shank spacing (m)
epsilon = 1.0e-6 # Pipe roughness (m)
# Pipe positions
# Double U-tube [(x_in1, y_in1), (x_in2, y_in2),
# (x_out1, y_out1), (x_out2, y_out2)]
pos = [(-D_s, 0.), (0., -D_s), (D_s, 0.), (0., D_s)]
# Ground properties
alpha = 1.0e-6 # Ground thermal diffusivity (m2/s)
k_s = 2.0 # Ground thermal conductivity (W/m.K)
T_g = 10.0 # Undisturbed ground temperature (degC)
# Grout properties
k_g = 1.0 # Grout thermal conductivity (W/m.K)
# Pipe properties
k_p = 0.4 # Pipe thermal conductivity (W/m.K)
# Fluid properties
m_flow_borehole = 0.25 # Total fluid mass flow rate per borehole (kg/s)
m_flow_network = m_flow_borehole*N_1*N_2 # Total fluid mass flow rate (kg/s)
# The fluid is propylene-glycol (20 %) at 20 degC
fluid = gt.media.Fluid('MPG', 20.)
cp_f = fluid.cp # Fluid specific isobaric heat capacity (J/kg.K)
rho_f = fluid.rho # Fluid density (kg/m3)
mu_f = fluid.mu # Fluid dynamic viscosity (kg/m.s)
k_f = fluid.k # Fluid thermal conductivity (W/m.K)
# g-Function calculation options
options = {'nSegments': 8,
'disp': True}
# Simulation parameters
dt = 3600. # Time step (s)
tmax = 1.*8760. * 3600. # Maximum time (s)
Nt = int(np.ceil(tmax/dt)) # Number of time steps
time = dt * np.arange(1, Nt+1)
# Load aggregation scheme
LoadAgg = gt.load_aggregation.ClaessonJaved(dt, tmax)
# -------------------------------------------------------------------------
# Initialize bore field and pipe models
# -------------------------------------------------------------------------
# The field is a retangular array
boreField = gt.boreholes.rectangle_field(N_1, N_2, B, B, H, D, r_b)
nBoreholes = len(boreField)
# Pipe thermal resistance
R_p = gt.pipes.conduction_thermal_resistance_circular_pipe(
r_in, r_out, k_p)
# Fluid to inner pipe wall thermal resistance (Double U-tube in parallel)
m_flow_pipe = m_flow_borehole / 2
h_f = gt.pipes.convective_heat_transfer_coefficient_circular_pipe(
m_flow_pipe, r_in, mu_f, rho_f, k_f, cp_f, epsilon)
R_f = 1.0/(h_f*2*pi*r_in)
# Double U-tube (parallel), same for all boreholes in the bore field
UTubes = []
for borehole in boreField:
UTube = gt.pipes.MultipleUTube(
pos, r_in, r_out, borehole, k_s, k_g, R_f + R_p,
nPipes=2, config='parallel')
UTubes.append(UTube)
# Build a network object from the list of UTubes
network = gt.networks.Network(
boreField, UTubes, m_flow_network=m_flow_network, cp_f=cp_f)
# -------------------------------------------------------------------------
# Calculate g-function
# -------------------------------------------------------------------------
# Get time values needed for g-function evaluation
time_req = LoadAgg.get_times_for_simulation()
# Calculate g-function
gFunc = gt.gfunction.gFunction(
network, alpha, time=time_req, boundary_condition='MIFT',
options=options)
# Initialize load aggregation scheme
LoadAgg.initialize(gFunc.gFunc/(2*pi*k_s))
# -------------------------------------------------------------------------
# Simulation
# -------------------------------------------------------------------------
# Evaluate inlet fluid temperature
one_year = 8760 * 3600.
one_day = 24 * 3600.
T_f_in = 15. - 10 * np.cos(2 * np.pi * time / one_year) - 5 * np.sin(np.pi * time / one_day)**2
T_b = np.zeros(Nt)
T_f_out = np.zeros(Nt)
Q_tot = np.zeros(Nt)
g_dt = gFunc.gFunc[0]
a_in, a_b = network.coefficients_network_heat_extraction_rate(
m_flow_network, cp_f, 1)
a_in = a_in.item()
# The effective borehole wall temperature is applied to all boreholes
a_b = np.sum(a_b)
for i, t in enumerate(time):
# Increment time step by (1)
LoadAgg.next_time_step(t)
# Evaluate borehole wall temperature assuming no extraction
T_b0 = T_g - LoadAgg.temporal_superposition().item()
# Solve the system of equations:
# 1) Q_tot = a_in * T_f_in + a_b * T_b
# 2) T_b = T_b0 - Q_tot / (2*pi*k_s*H*nB) * g(dt)
denom = (2 * pi * k_s * H * nBoreholes)
T_b[i] = (T_b0 - a_in * T_f_in[i] * g_dt / denom) / (
1 + a_b * g_dt / denom)
Q_tot[i] = a_in * T_f_in[i] + a_b * T_b[i]
# Apply current load (in watts per meter of borehole)
LoadAgg.set_current_load(Q_tot[i] / (H * nBoreholes))
# Evaluate outlet fluid temperature
T_f_out[i] = network.get_network_outlet_temperature(
T_f_in[i], T_b[i], m_flow_network, cp_f, nSegments=1)
# -------------------------------------------------------------------------
# Plot hourly heat extraction rates and temperatures
# -------------------------------------------------------------------------
# Configure figure and axes
fig = gt.utilities._initialize_figure()
ax1 = fig.add_subplot(211)
# Axis labels
ax1.set_xlabel(r'Time [hours]')
ax1.set_ylabel(r'Total heat extraction rate [W/m]')
gt.utilities._format_axes(ax1)
# Plot heat extraction rates
hours = np.arange(1, Nt+1) * dt / 3600.
ax1.plot(hours, Q_tot / (nBoreholes * H))
ax2 = fig.add_subplot(212)
# Axis labels
ax2.set_xlabel(r'Time [hours]')
ax2.set_ylabel(r'Temperature [degC]')
gt.utilities._format_axes(ax2)
# Plot temperatures
ax2.plot(hours, T_b, label='Borehole wall')
ax2.plot(hours, T_f_in, '--',
label='Inlet, double U-tube (parallel)')
ax2.plot(hours, T_f_out, '-.',
label='Outlet, double U-tube (parallel)')
ax2.legend()
# Adjust to plot window
plt.tight_layout()
return
# Main function
if __name__ == '__main__':
main()Caution is needed when using monitored inlet fluid temperatures during simulations. The variation of the borehole wall temperature depends on the net heat extraction rate. This means that the error can accumulate during a simulation if, for example, the model underestimates the heat extraction rate for a given inlet fluid temperature. The can be a drift between the trends of the monitored and simulated temperatures in some cases.
Thanks for developing this great package!
Recently I want to use pygfunction to simulate a pipe network where the time series of inlet/outlet temperature are known.
I want to use the time series of inlet temperature to calculate the outlet temperature for validation.
However, I only see the simulation of fluid temperature where the rate of heat extraction is known as time series. Is there any way to use inlet temperature as input?