v2.1.1

In GHEtool v2.1.1 there are two major code changes that reduce the computational time significantly. One has to do with the way the sizing methodology (L3/L4) is implemented and another with the new Gfunction class. Both improvements are explained below.

Improvement in sizing

Previously, in v2.1.0, for the sizing methodology L3 and L4 (i.e. monthly and hourly), the temperatures where calculated every time step. This however (especially for long simulation periods) requires a lot of time due to the convolution step, especially for an hourly sizing. In v2.1.1 this is changed so that only the first and last year are calculated, since only these years are relevant for the sizing 1. This means that however long the simulation period may be, only two years are calculated. In the table below, the time required for one sizing iteration (i.e. one g-value convolution) for an hourly sizing is shown in µs. The code to come up with these numbers is added below.

Speed improvement for g-value convolution in hourly sizing

Simulation period [years]	Required time old method [µs]	Required time new method [µs]
5 years	15625 µs	0 µs
15 years	15625 µs	0 µs
25 years	31250 µs	0 µs
35 years	46875 µs	0 µs
45 years	62500 µs	0 µs
55 years	78125 µs	0 µs
65 years	78125 µs	0 µs
75 years	78125 µs	0 µs
85 years	109375 µs	0 µs
95 years	125000 µs	0 µs

Gfunction class

Due to the implementation of the GFunction class in GHEtool, a substantial speed improvement is made w.r.t. GHEtool v2.1.0 for computationally expensive tasks. In the tables below, one can find this speed improvement for the different sizing methods and for several of the examples documents. The results can be recreated by running the code below. The computational times as shown in the table below, are an average of 5 runs.

Speed benchmark sizing methods

Sizing method	Time v2.1.0 [ms]	Time v2.1.1 [ms]	Improvement [%]
L2 (three pulse) sizing	1.84 ms	1.43 ms	28%
L3 (monthly) sizing	12.34 ms	5.8 ms	113%
L4 (hourly) sizing	4.21 ms	3.63 ms	16%

Speed benchmark examples

Example	Time v2.1.0 [ms]	Time v2.1.1 [ms]	Improvement [%]
Main functionalities	3.38 ms	2.57 ms	32%
Optimise load profile	15.07 ms	0.63 ms	2305%
Sizing with Rb calculation	9.97 ms	9.99 ms	0%
Effect borefield configuration	1.57 ms	1.5 ms	4%

import numpy as np
import pygfunction as gt
from GHEtool import GroundConstantTemperature, Borefield, HourlyGeothermalLoad, FOLDER
from scipy.signal import convolve
from math import pi
from time import process_time_ns


def test_new_calc_method(simulation_period: int):
    """
    Test the new calculation method which is just considering the first and last year.

    Parameters
    ----------
    simulation_period : int
        simulation period [years]

    Returns
    -------
        None

    Raises
    -------
        AssertionError
    """

    h = 110

    # initiate ground data
    data = GroundConstantTemperature(3, 10)

    # initiate pygfunction borefield model
    borefield_gt = gt.boreholes.rectangle_field(10, 10, 6, 6, 110, 1, 0.075)

    # initiate borefield
    borefield = Borefield(100)

    # set borehole thermal equivalent resistance
    borefield.Rb = 0.12

    # set ground data in borefield
    borefield.set_ground_parameters(data)

    # set pygfunction borefield model
    borefield.set_borefield(borefield_gt)


    # load the hourly profile
    load = HourlyGeothermalLoad(simulation_period=simulation_period)
    load.load_hourly_profile(f'hourly_profile.csv', header=True, separator=";")

    borefield.load = load

    # borefield.g-function is a function that uses the precalculated data to interpolate the correct values of the
    # g-function. This dataset is checked over and over again and is correct
    g_values = borefield.gfunction(borefield.time_L4, borefield.H)

    # get process time at start of new method
    dt1 = process_time_ns()
    # determine load
    loads_short = borefield.hourly_cooling_load - borefield.hourly_heating_load
    # reverse the load
    loads_short_rev = loads_short[::-1]
    # init results vector
    results = np.zeros(loads_short.size * 2)
    # calculation of needed differences of the g-function values. These are the weight factors in the calculation
    # of Tb.
    g_value_differences = np.diff(g_values, prepend=0)

    # convolution to get the results for the first year
    results[:8760] = convolve(loads_short * 1000, g_value_differences[:8760])[:8760]
    # sum up g_values until the pre last year
    g_sum_n1 = g_value_differences[:8760 * (borefield.simulation_period - 1)].reshape(borefield.simulation_period - 1, 8760).sum(axis=0)
    # add up last year
    g_sum = g_sum_n1 + g_value_differences[8760 * (borefield.simulation_period - 1):]
    # add zero at start and reverse the order
    g_sum_n2 = np.concatenate((np.array([0]), g_sum_n1[::-1]))[:-1]
    # determine results for the last year by the influence of the year (first term) and the previous years (last term)
    results[8760:] = convolve(loads_short * 1000, g_sum)[:8760] + convolve(loads_short_rev * 1000, g_sum_n2)[:8760][::-1]
    # calculation the borehole wall temperature for every month i
    t_b = results / (2 * pi * borefield.ground_data.k_s) / (borefield.H * borefield.number_of_boreholes) + borefield._Tg(borefield.H)

    # get process time
    dt2 = process_time_ns()
    # determine hourly load
    hourly_load = np.tile(borefield.hourly_cooling_load - borefield.hourly_heating_load, borefield.simulation_period)
    # calculation of needed differences of the g-function values. These are the weight factors in the calculation
    # of Tb.
    g_value_differences = np.diff(g_values, prepend=0)

    # convolution to get the monthly results
    results = convolve(hourly_load * 1000, g_value_differences)[:hourly_load.size]

    # calculation the borehole wall temperature for every month i
    t_b_new = results / (2 * pi * borefield.ground_data.k_s) / (h * borefield.number_of_boreholes) + borefield._Tg(h)

    # print time for the different methods
    print(f'simulation period: {simulation_period}; old: {(process_time_ns() - dt2)/1000:.0f} µs;'
          f'new: {(dt2 - dt1)/1000:.0f} µs')
    # compare results to ensure they are the same
    assert np.allclose(t_b[:8760], t_b_new[:8760])
    assert np.allclose(t_b[8760:], t_b_new[8760*(borefield.simulation_period - 1):])


if __name__ == "__main__":
    for sim_year in np.arange(5, 101, 10):
        test_new_calc_method(sim_year)

from GHEtool import *
import numpy as np
import pygfunction as gt
import time
import os, contextlib

from GHEtool.Examples.main_functionalities import main_functionalities
from GHEtool.Examples.sizing_with_Rb_calculation import sizing_with_Rb
from GHEtool.Examples.effect_of_borehole_configuration import effect_borefield_configuration
from GHEtool import HourlyGeothermalLoad

# disable the plot function by monkey patching over it
Borefield._plot_temperature_profile = lambda *args, **kwargs: None


# disable the print outputs
def supress_stdout(func):
    def wrapper(*a, **ka):
        with open(os.devnull, 'w') as devnull:
            with contextlib.redirect_stdout(devnull):
                return func(*a, **ka)
    return wrapper


@supress_stdout
def run_without_messages(callable) -> None:
    """
    This function runs the callable without messages.

    Parameters
    ----------
    callable : Callable
        Function to be called

    Returns
    -------
    None
    """
    callable()


def optimise_load_profile() -> None:
    """
    This is a benchmark for the optimise load profile method.

    Returns
    -------
    None
    """
    # initiate ground data
    data = GroundData(3, 10, 0.2)

    # initiate pygfunction borefield model
    borefield_gt = gt.boreholes.rectangle_field(10, 10, 6, 6, 110, 1, 0.075)

    # initiate borefield
    borefield = Borefield()

    # set ground data in borefield
    borefield.set_ground_parameters(data)

    # set pygfunction borefield
    borefield.set_borefield(borefield_gt)

    # load the hourly profile
    load = HourlyGeothermalLoad()
    load.load_hourly_profile("hourly_profile.csv", header=True, separator=";")
    borefield.load = load

    # optimise the load for a 10x10 field (see data above) and a fixed depth of 150m.
    borefield.optimise_load_profile(depth=150, print_results=False)


def size_L2() -> None:
    """
    This is a benchmark for the L2 sizing method.

    Returns
    -------
    None
    """
    number_of_iterations = 5
    max_value_cooling = 700
    max_value_heating = 800

    monthly_load_cooling_array = np.empty((number_of_iterations, 12))
    monthly_load_heating_array = np.empty((number_of_iterations, 12))
    peak_load_cooling_array = np.empty((number_of_iterations, 12))
    peak_load_heating_array = np.empty((number_of_iterations, 12))

    # populate arrays with random values
    for i in range(number_of_iterations):
        for j in range(12):
            monthly_load_cooling_array[i, j] = np.random.randint(0, max_value_cooling)
            monthly_load_heating_array[i, j] = np.random.randint(0, max_value_heating)
            peak_load_cooling_array[i, j] = np.random.randint(monthly_load_cooling_array[i, j], max_value_cooling)
            peak_load_heating_array[i, j] = np.random.randint(monthly_load_heating_array[i, j], max_value_heating)

    # initiate borefield model
    data = GroundData(3, 10, 0.2)
    borefield_gt = gt.boreholes.rectangle_field(10, 12, 6, 6, 110, 1, 0.075)

    # Monthly loading values
    peak_cooling = np.array([0., 0, 34., 69., 133., 187., 213., 240., 160., 37., 0., 0.])  # Peak cooling in kW
    peak_heating = np.array([160., 142, 102., 55., 0., 0., 0., 0., 40.4, 85., 119., 136.])  # Peak heating in kW

    # annual heating and cooling load
    annual_heating_load = 300 * 10 ** 3  # kWh
    annual_cooling_load = 160 * 10 ** 3  # kWh

    # percentage of annual load per month (15.5% for January ...)
    monthly_load_heating_percentage = np.array(
        [0.155, 0.148, 0.125, .099, .064, 0., 0., 0., 0.061, 0.087, 0.117, 0.144])
    monthly_load_cooling_percentage = np.array([0.025, 0.05, 0.05, .05, .075, .1, .2, .2, .1, .075, .05, .025])

    # resulting load per month
    monthly_load_heating = annual_heating_load * monthly_load_heating_percentage  # kWh
    monthly_load_cooling = annual_cooling_load * monthly_load_cooling_percentage  # kWh

    # create the borefield object
    borefield = Borefield(simulation_period=20,
                          peak_heating=peak_heating,
                          peak_cooling=peak_cooling,
                          baseload_heating=monthly_load_heating,
                          baseload_cooling=monthly_load_cooling)
    borefield.set_ground_parameters(data)
    borefield.set_borefield(borefield_gt)

    # set temperature boundaries
    borefield.set_max_avg_fluid_temperature(16)  # maximum temperature
    borefield.set_min_avg_fluid_temperature(0)  # minimum temperature

    # size according to L2 method
    for i in range(number_of_iterations):
        borefield.set_baseload_cooling(monthly_load_cooling_array[i])
        borefield.set_baseload_heating(monthly_load_heating_array[i])
        borefield.set_peak_cooling(peak_load_cooling_array[i])
        borefield.set_peak_heating(peak_load_heating_array[i])
        borefield.size(100, L2_sizing=True)


def size_L3() -> None:
    """
    This is a benchmark for the L3 sizing method.

    Returns
    -------
    None
    """
    number_of_iterations = 5
    max_value_cooling = 700
    max_value_heating = 800

    monthly_load_cooling_array = np.empty((number_of_iterations, 12))
    monthly_load_heating_array = np.empty((number_of_iterations, 12))
    peak_load_cooling_array = np.empty((number_of_iterations, 12))
    peak_load_heating_array = np.empty((number_of_iterations, 12))

    # populate arrays with random values
    for i in range(number_of_iterations):
        for j in range(12):
            monthly_load_cooling_array[i, j] = np.random.randint(0, max_value_cooling)
            monthly_load_heating_array[i, j] = np.random.randint(0, max_value_heating)
            peak_load_cooling_array[i, j] = np.random.randint(monthly_load_cooling_array[i, j], max_value_cooling)
            peak_load_heating_array[i, j] = np.random.randint(monthly_load_heating_array[i, j], max_value_heating)

    # initiate borefield model
    data = GroundData(3, 10, 0.2)
    borefield_gt = gt.boreholes.rectangle_field(10, 12, 6, 6, 110, 1, 0.075)

    # Monthly loading values
    peak_cooling = np.array([0., 0, 34., 69., 133., 187., 213., 240., 160., 37., 0., 0.])  # Peak cooling in kW
    peak_heating = np.array([160., 142, 102., 55., 0., 0., 0., 0., 40.4, 85., 119., 136.])  # Peak heating in kW

    # annual heating and cooling load
    annual_heating_load = 300 * 10 ** 3  # kWh
    annual_cooling_load = 160 * 10 ** 3  # kWh

    # percentage of annual load per month (15.5% for January ...)
    monthly_load_heating_percentage = np.array(
        [0.155, 0.148, 0.125, .099, .064, 0., 0., 0., 0.061, 0.087, 0.117, 0.144])
    monthly_load_cooling_percentage = np.array([0.025, 0.05, 0.05, .05, .075, .1, .2, .2, .1, .075, .05, .025])

    # resulting load per month
    monthly_load_heating = annual_heating_load * monthly_load_heating_percentage  # kWh
    monthly_load_cooling = annual_cooling_load * monthly_load_cooling_percentage  # kWh

    # create the borefield object
    borefield = Borefield(simulation_period=20,
                          peak_heating=peak_heating,
                          peak_cooling=peak_cooling,
                          baseload_heating=monthly_load_heating,
                          baseload_cooling=monthly_load_cooling)
    borefield.set_ground_parameters(data)
    borefield.set_borefield(borefield_gt)

    # set temperature boundaries
    borefield.set_max_avg_fluid_temperature(16)  # maximum temperature
    borefield.set_min_avg_fluid_temperature(0)  # minimum temperature

    # size according to L3 method
    for i in range(number_of_iterations):
        borefield.set_baseload_cooling(monthly_load_cooling_array[i])
        borefield.set_baseload_heating(monthly_load_heating_array[i])
        borefield.set_peak_cooling(peak_load_cooling_array[i])
        borefield.set_peak_heating(peak_load_heating_array[i])
        borefield.size(100, L3_sizing=True)


def size_L4() -> None:
    """
    This is a benchmark for the L4 sizing method.

    Returns
    -------
    None
    """
    # initiate ground data
    data = GroundData(3, 10, 0.2)

    # initiate pygfunction borefield model
    borefield_gt = gt.boreholes.rectangle_field(10, 10, 6, 6, 110, 1, 0.075)

    # initiate borefield
    borefield = Borefield()

    # set ground data in borefield
    borefield.set_ground_parameters(data)

    # set pygfunction borefield
    borefield.set_borefield(borefield_gt)

    # load the hourly profile
    load = HourlyGeothermalLoad()
    load.load_hourly_profile("hourly_profile.csv", header=True, separator=";")
    borefield.load = load

    borefield.size(L4_sizing=True)


def benchmark(callable, name: str) -> None:
    """
    This function calls the callable five times and outputs the time needed to run the callable.

    Parameters
    ----------
    callable : Callable
        Function to be called
    name : str
        Name of the function

    Returns
    -------
    None
    """
    diff, diff_without = 0., 0.

    for i in range(5):
        GFunction.DEFAULT_STORE_PREVIOUS_VALUES = True

        start_time = time.time()
        run_without_messages(callable)
        end_time = time.time()
        diff = diff * i/(i+1) + (end_time - start_time)/(i+1)

        GFunction.DEFAULT_STORE_PREVIOUS_VALUES = False

        start_time_without = time.time()
        run_without_messages(callable)
        end_time_without = time.time()
        diff_without = diff_without * i/(i+1) + (end_time_without - start_time_without)/(i+1)

    print(f'{name} took  {round(diff_without, 2)} ms in v2.1.0 and '
          f'{round(diff, 2)} ms in v2.1.1. This is an improvement of {round((diff_without-diff)/diff*100)}%.')


# run examples
benchmark(main_functionalities, "Main functionalities")
benchmark(optimise_load_profile, "Optimise load profile")
benchmark(sizing_with_Rb, "Sizing with Rb calculation")
benchmark(effect_borefield_configuration, "Effect borehole configuration")

# run benchmark sizing methods
benchmark(size_L2, "Sizing with L2 method")
benchmark(size_L3, "Sizing with L3 method")
benchmark(size_L4, "Sizing with L4 (hourly) method")

References

1

Peere, W., Picard, D., Cupeiro Figueroa, I., Boydens, W., and Helsen, L. (2021) Validated combined first and last year borefield sizing methodology. In Proceedings of International Building Simulation Conference 2021. Brugge (Belgium), 1-3 September 2021. https://doi.org/10.26868/25222708.2021.30180