Skip to content

Commit

Permalink
Add dtype option to gFunction solvers
Browse files Browse the repository at this point in the history 
This is for #85
  • Loading branch information
@MassimoCimmino
MassimoCimmino committed Apr 10, 2021
1 parent 8b5bc09 commit a0f78b1
Showing 1 changed file with 23 additions and 19 deletions.
42 changes: 23 additions & 19 deletions pygfunction/gfunction.py
View file
Original file line number Diff line number Diff line change
Expand Up @@ -75,6 +75,9 @@ class gFunction(object):
response factors. See documentation for
scipy.interpolate.interp1d.
Default is linear.
dtype : numpy dtype, optional
numpy data type used for matrices and vectors.
Default is numpy.double.
The ''similarities'' solver accepts the following method-specific
options:
disTol : float, optional
Expand Down Expand Up @@ -768,7 +771,7 @@ class _BaseSolver(object):
"""
def __init__(self, boreholes, network, time, boundary_condition,
nSegments=12, processes=None, disp=False, profiles=False,
kind='linear', **other_options):
kind='linear', dtype=np.double, **other_options):
self.boreholes = boreholes
self.network = network
# Convert time to a 1d array
Expand All @@ -779,6 +782,7 @@ def __init__(self, boreholes, network, time, boundary_condition,
self.disp = disp
self.profiles = profiles
self.kind = kind
self.dtype = dtype
# Check the validity of inputs
self._check_inputs()
# Initialize the solver with solver-specific options
Expand Down Expand Up @@ -835,11 +839,11 @@ def solve(self, time, alpha):
if self.boundary_condition == 'UHTR':
Q = 1
else:
Q = np.zeros((self.nSources, nt))
Q = np.zeros((self.nSources, nt), dtype=self.dtype)
if self.boundary_condition == 'UBWT':
Tb = np.zeros(nt)
Tb = np.zeros(nt, dtype=self.dtype)
else:
Tb = np.zeros((self.nSources, nt))
Tb = np.zeros((self.nSources, nt), dtype=self.dtype)
# Calculate segment to segment thermal response factors
h_ij = self.thermal_response_factors(time, alpha, kind=self.kind)
# Segment lengths
Expand Down Expand Up @@ -893,7 +897,7 @@ def solve(self, time, alpha):
# Spatial superposition: [Tb] = [Tb0] + [h_ij_dt]*[Qb]
# Energy conservation: sum([Q*Hb]) = sum([Hb])
# ---------------------------------------------------------
A = np.block([[h_dt, -np.ones((self.nSources, 1))],
A = np.block([[h_dt, -np.ones((self.nSources, 1), dtype=self.dtype)],
[Hb, 0.]])
B = np.hstack((-Tb_0, Htot))
# Solve the system of equations
Expand Down Expand Up @@ -922,10 +926,10 @@ def solve(self, time, alpha):
a_in, a_b = self.network.coefficients_borehole_heat_extraction_rate(
self.network.m_flow, self.network.cp, self.nSegments)
k_s = self.network.p[0].k_s
A = np.block([[h_dt, -np.eye(self.nSources), np.zeros((self.nSources, 1))],
[np.eye(self.nSources), a_b/(2.0*pi*k_s*np.atleast_2d(Hb).T), a_in/(2.0*pi*k_s*np.atleast_2d(Hb).T)],
[Hb, np.zeros(self.nSources + 1)]])
B = np.hstack((-Tb_0, np.zeros(self.nSources), Htot))
A = np.block([[h_dt, -np.eye(self.nSources, dtype=self.dtype), np.zeros((self.nSources, 1), dtype=self.dtype)],
[np.eye(self.nSources, dtype=self.dtype), a_b/(2.0*pi*k_s*np.atleast_2d(Hb).T), a_in/(2.0*pi*k_s*np.atleast_2d(Hb).T)],
[Hb, np.zeros(self.nSources + 1, dtype=self.dtype)]])
B = np.hstack((-Tb_0, np.zeros(self.nSources, dtype=self.dtype), Htot))
# Solve the system of equations
X = np.linalg.solve(A, B)
# Store calculated heat extraction rates
Expand Down Expand Up @@ -967,7 +971,7 @@ def segment_lengths(self):
"""
# Borehole lengths
H = np.array([b.H for b in self.boreSegments])
H = np.array([b.H for b in self.boreSegments], dtype=self.dtype)
return H

def borehole_segments(self):
Expand Down Expand Up @@ -1058,7 +1062,7 @@ def load_history_reconstruction(self, time, Q):
# Reconstructed time vector
t_reconstructed = np.hstack((0., np.cumsum(dt_reconstructed)))
# Accumulated heat extracted
f = np.hstack((np.zeros((nSources, 1)), np.cumsum(Q*dt, axis=1)))
f = np.hstack((np.zeros((nSources, 1), dtype=self.dtype), np.cumsum(Q*dt, axis=1)))
f = np.hstack((f, f[:,-1:]))
# Create interpolation object for accumulated heat extracted
sf = interp1d(t, f, kind='linear', axis=1)
Expand Down Expand Up @@ -1211,7 +1215,7 @@ def thermal_response_factors(self, time, alpha, kind='linear'):
# Initialize chrono
tic = tim.time()
# Initialize segment-to-segment response factors
h_ij = np.zeros((self.nSources, self.nSources, nt))
h_ij = np.zeros((self.nSources, self.nSources, nt), dtype=self.dtype)

for i in range(self.nSources):
# Segment to same-segment thermal response factor
Expand All @@ -1220,7 +1224,7 @@ def thermal_response_factors(self, time, alpha, kind='linear'):
func = partial(finite_line_source,
alpha=alpha, borehole1=b2, borehole2=b2)
# Evaluate the FLS solution at all times in parallel
h = np.array(pool.map(func, time))
h = np.array(pool.map(func, time), dtype=self.dtype)
h_ij[i, i, :] = h

# Segment to other segments thermal response factor
Expand All @@ -1229,7 +1233,7 @@ def thermal_response_factors(self, time, alpha, kind='linear'):
# Evaluate the FLS solution at all times in parallel
func = partial(finite_line_source,
alpha=alpha, borehole1=b1, borehole2=b2)
h = np.array(pool.map(func, time))
h = np.array(pool.map(func, time), dtype=self.dtype)
h_ij[i, j, :] = h
h_ij[j, i, :] = b2.H / b1.H * h_ij[i, j, :]

Expand All @@ -1244,7 +1248,7 @@ def thermal_response_factors(self, time, alpha, kind='linear'):
# Interp1d object for thermal response factors
h_ij = interp1d(
np.hstack((0., time)),
np.dstack((np.zeros((self.nSources,self.nSources)), h_ij)),
np.dstack((np.zeros((self.nSources,self.nSources), dtype=self.dtype), h_ij)),
kind=kind, copy=True, axis=2)
toc = tim.time()
if self.disp: print('{} sec'.format(toc - tic))
Expand Down Expand Up @@ -1382,7 +1386,7 @@ def thermal_response_factors(self, time, alpha, kind='linear'):
# Initialize chrono
tic = tim.time()
# Initialize segment-to-segment response factors
h_ij = np.zeros((self.nSources, self.nSources, nt))
h_ij = np.zeros((self.nSources, self.nSources, nt), dtype=self.dtype)

# Similarities for real sources
for s in range(self.nSimPos):
Expand All @@ -1401,7 +1405,7 @@ def thermal_response_factors(self, time, alpha, kind='linear'):
alpha=alpha, borehole1=b1, borehole2=b2,
reaSource=True, imgSource=True)
# Evaluate the FLS solution at all times in parallel
hPos = np.array(pool.map(func, np.atleast_1d(time)))
hPos = np.array(pool.map(func, np.atleast_1d(time)), dtype=self.dtype)
# Assign thermal response factors to similar segment pairs
for (i, j) in self.simPos[s]:
h_ij[j, i, :] = hPos
Expand All @@ -1419,7 +1423,7 @@ def thermal_response_factors(self, time, alpha, kind='linear'):
alpha=alpha, borehole1=b1, borehole2=b2,
reaSource=False, imgSource=True)
# Evaluate the FLS solution at all times in parallel
hNeg = np.array(pool.map(func, time))
hNeg = np.array(pool.map(func, time), dtype=self.dtype)
# Assign thermal response factors to similar segment pairs
for (i, j) in self.simNeg[s]:
h_ij[j, i, :] = h_ij[j, i, :] + hNeg
Expand All @@ -1435,7 +1439,7 @@ def thermal_response_factors(self, time, alpha, kind='linear'):

# Interp1d object for thermal response factors
h_ij = interp1d(np.hstack((0., time)),
np.dstack((np.zeros((self.nSources,self.nSources)), h_ij)),
np.dstack((np.zeros((self.nSources,self.nSources), dtype=self.dtype), h_ij)),
kind=kind, copy=True, axis=2)
toc = tim.time()
if self.disp: print('{} sec'.format(toc - tic))
Expand Down

0 comments on commit a0f78b1

Please sign in to comment.