diff --git a/CHANGELOG.md b/CHANGELOG.md index e8e778e1..72fb8e21 100644 --- a/CHANGELOG.md +++ b/CHANGELOG.md @@ -2,6 +2,10 @@ ## Version 2.3 (in development) +### New features + +* [Issue 276](https://github.com/MassimoCimmino/pygfunction/issues/276) - Added functions to the `boreholes` module for the generation of rectangular fields in a staggered configuration. + ### Bug fixes * [Issue 274](https://github.com/MassimoCimmino/pygfunction/issues/274) - Fixed scalar assignment from ndim-1 array. It is deprecated as of `numpy` version `1.25`. Only ndim-0 arrays can be treated as scalars. diff --git a/examples/regular_bore_field.py b/examples/regular_bore_field.py index aea2e78c..604870fd 100644 --- a/examples/regular_bore_field.py +++ b/examples/regular_bore_field.py @@ -2,6 +2,8 @@ """ Example of definition of a bore field using pre-defined configurations. """ +import matplotlib.pyplot as plt + import pygfunction as gt @@ -29,6 +31,14 @@ def main(): # Rectangular field of 4 x 3 boreholes rectangularField = gt.boreholes.rectangle_field(N_1, N_2, B, B, H, D, r_b) + # Rectangular field triangular field of 4 x 3 borehole rows + staggeredRectangularField = gt.boreholes.staggered_rectangle_field( + N_1, N_2, B, B, H, D, r_b, False) + + # Dense field triangular field of 4 x 3 borehole rows + denseRectangularField = gt.boreholes.dense_rectangle_field( + N_1, N_2, B, H, D, r_b, False) + # Box-shaped field of 4 x 3 boreholes boxField = gt.boreholes.box_shaped_field(N_1, N_2, B, B, H, D, r_b) @@ -41,11 +51,17 @@ def main(): # Circular field of 8 boreholes circleField = gt.boreholes.circle_field(N_b, R, H, D, r_b) + # Filled circle field of 20 boreholes + filledCircleField = gt.boreholes.filled_circle_field(20, 10, 100, 1, 0.05) + # ------------------------------------------------------------------------- # Draw bore fields # ------------------------------------------------------------------------- - for field in [rectangularField, boxField, UField, LField, circleField]: + for field in [ + rectangularField, staggeredRectangularField, denseRectangularField, + boxField, UField, LField, circleField, filledCircleField]: gt.boreholes.visualize_field(field) + plt.show() return diff --git a/pygfunction/boreholes.py b/pygfunction/boreholes.py index 5da7bf23..f9e29c74 100644 --- a/pygfunction/boreholes.py +++ b/pygfunction/boreholes.py @@ -713,6 +713,176 @@ def rectangle_field(N_1, N_2, B_1, B_2, H, D, r_b, tilt=0., origin=None): return borefield +def staggered_rectangle_field( + N_1, N_2, B_1, B_2, H, D, r_b, include_last_borehole, tilt=0., + origin=None): + """ + Build a list of boreholes in a rectangular bore field configuration, with + boreholes placed in a staggered configuration. + + Parameters + ---------- + N_1 : int + Number of borehole in the x direction. + N_2 : int + Number of borehole in the y direction. + B_1 : float + Distance (in meters) between adjacent boreholes in the x direction. + B_2 : float + Distance (in meters) between adjacent boreholes in the y direction. + H : float + Borehole length (in meters). + D : float + Borehole buried depth (in meters). + r_b : float + Borehole radius (in meters). + include_last_borehole : bool + If True, then each row of boreholes has equal numbers of boreholes. + If False, then the staggered rows have one borehole less so they are + contained within the imaginary 'box' around the borefield. + tilt : float, optional + Angle (in radians) from vertical of the axis of the borehole. The + orientation of the tilt is orthogonal to the origin coordinate. + Default is 0. + origin : tuple, optional + A coordinate indicating the origin of reference for orientation of + boreholes. + Default is the center of the rectangle. + + Returns + ------- + boreField : list of Borehole objects + List of boreholes in the rectangular bore field. + + Notes + ----- + Boreholes located at the origin will remain vertical. + + Examples + -------- + >>> boreField = gt.boreholes.rectangle_field(N_1=3, N_2=2, B_1=5., B_2=5., + H=100., D=2.5, r_b=0.05) + + The bore field is constructed line by line. For N_1=3 and N_2=3, the bore + field layout is as follows, if `include_last_borehole` is True:: + + 6 7 8 + 3 4 5 + 0 1 2 + + and if `include_last_borehole` is False:: + + 5 6 7 + 3 4 + 0 1 2 + + """ + borefield = [] + + if N_1 == 1 or N_2 == 1: + return rectangle_field(N_1, N_2, B_1, B_2, H, D, r_b, tilt, origin) + + if origin is None: + # When no origin is supplied, compute the origin to be at the center of + # the rectangle + if include_last_borehole: + x0 = (N_1 - 1) / 2 * B_1 + y0 = (N_2 - 1) / 2 * B_2 + else: + x0 = (N_1 - 0.5) / 2 * B_1 + y0 = (N_2 - 0.5) / 2 * B_2 + else: + x0, y0 = origin + + for j in range(N_2): + for i in range(N_1): + x = i * B_1 + (B_1 / 2 if j % 2 == 1 else 0) + y = j * B_2 + # The borehole is inclined only if it does not lie on the origin + if np.sqrt((x - x0)**2 + (y - y0)**2) > r_b: + orientation = np.arctan2(y - y0, x - x0) + if i < (N_1 - 1) or include_last_borehole or (j % 2 == 0): + borefield.append( + Borehole( + H, D, r_b, x, y, tilt=tilt, orientation=orientation)) + else: + if i < (N_1 - 1) or include_last_borehole or (j % 2 == 0): + borefield.append(Borehole(H, D, r_b, x, y)) + + return borefield + + +def dense_rectangle_field( + N_1, N_2, B, H, D, r_b, include_last_borehole, tilt=0., origin=None): + """ + Build a list of boreholes in a rectangular bore field configuration, with + boreholes placed in a staggered configuration with uniform spacing between + boreholes. + + Parameters + ---------- + N_1 : int + Number of borehole in the x direction. + N_2 : int + Number of borehole in the y direction. + B : float + Distance (in meters) between adjacent boreholes. + H : float + Borehole length (in meters). + D : float + Borehole buried depth (in meters). + r_b : float + Borehole radius (in meters). + include_last_borehole : bool + If True, then each row of boreholes has equal numbers of boreholes. + If False, then the staggered rows have one borehole less so they are + contained within the imaginary 'box' around the borefield. + tilt : float, optional + Angle (in radians) from vertical of the axis of the borehole. The + orientation of the tilt is orthogonal to the origin coordinate. + Default is 0. + origin : tuple, optional + A coordinate indicating the origin of reference for orientation of + boreholes. + Default is the center of the rectangle. + + Returns + ------- + boreField : list of Borehole objects + List of boreholes in the rectangular bore field. + + Notes + ----- + Boreholes located at the origin will remain vertical. + + Examples + -------- + >>> boreField = gt.boreholes.rectangle_field( + N_1=3, N_2=2, B_1=5., B_2=5., H=100., D=2.5, r_b=0.05, + include_last_borehole=True) + + The bore field is constructed line by line. For N_1=3 and N_2=3, the bore + field layout is as follows, if `include_last_borehole` is True:: + + 6 7 8 + 3 4 5 + 0 1 2 + + and if `include_last_borehole` is False:: + + 5 6 7 + 3 4 + 0 1 2 + + """ + if N_1 == 1: + # line field + return rectangle_field(N_1, N_2, B, B, H, D, r_b, tilt, origin) + return staggered_rectangle_field( + N_1, N_2, B, np.sqrt(3)/2 * B, H, D, r_b, include_last_borehole, + tilt=tilt, origin=origin) + + def L_shaped_field(N_1, N_2, B_1, B_2, H, D, r_b, tilt=0., origin=None): """ Build a list of boreholes in a L-shaped bore field configuration. @@ -1078,7 +1248,6 @@ def circle_field(N, R, H, D, r_b, tilt=0., origin=None): for i in range(N): x = R * np.cos(2 * pi * i / N) y = R * np.sin(2 * pi * i / N) - orientation = np.arctan2(y - y0, x - x0) # The borehole is inclined only if it does not lie on the origin if np.sqrt((x - x0)**2 + (y - y0)**2) > r_b: orientation = np.arctan2(y - y0, x - x0) @@ -1091,6 +1260,70 @@ def circle_field(N, R, H, D, r_b, tilt=0., origin=None): return borefield +def filled_circle_field(N, B, H, D, r_b, tilt=0., origin=None): + """ + Build a list of boreholes in a filled circular configuration, starting from the origin + and spiralling outside. + + Parameters + ---------- + N : int + Number of boreholes in the bore field. + B : float + Distance between consecutive circles of boreholes (in meters). + H : float + Borehole length (in meters). + D : float + Borehole buried depth (in meters). + r_b : float + Borehole radius (in meters). + tilt : float, optional + Angle (in radians) from vertical of the axis of the borehole. The + orientation of the tilt is towards the exterior of the bore field. + Default is 0. + origin : tuple + A coordinate indicating the origin of reference for orientation of + boreholes. + Default is the origin (0, 0). + + Returns + ------- + boreField : list of Borehole objects + List of boreholes in the filled circle shaped bore field. + + Notes + ----- + Boreholes located at the origin will remain vertical. + """ + + if origin is None: + # When no origin is supplied, compute the origin to be at the center of + # the rectangle + x0 = 0. + y0 = 0. + else: + x0, y0 = origin + + field = [Borehole(H, D, r_b, 0, 0)] + a = 6 # number of boreholes in a circle, 6 for the inner circle + B2 = B # radius of the borehole + counter = 0 + + for i in range(N - 1): + i = i + 1 + counter = counter + 1 + if counter > a: + a = a + 6 + B2 = B2 + B + counter = 1 + # The borehole is inclined only if it does not lie on the origin, which is always the case + x, y = B2 * np.sin(2 * np.pi * counter / a), B2 * np.cos(2 * np.pi * counter / a) + orientation = np.arctan2(y - y0, x - x0) + field.append(Borehole(H, D, r_b, x, y, tilt=tilt, orientation=orientation)) + + return field + + def field_from_file(filename): """ Build a list of boreholes given coordinates and dimensions provided in a diff --git a/tests/boreholes_test.py b/tests/boreholes_test.py index ebc3c5e9..8bc6bc64 100644 --- a/tests/boreholes_test.py +++ b/tests/boreholes_test.py @@ -94,6 +94,96 @@ def test_rectangular_field(N_1, N_2, B_1, B_2): ]) +# Test staggered_rectangle_field +@pytest.mark.parametrize("N_1, N_2, B_1, B_2, include_last_element", [ + (1, 1, 5., 5., True), # 1 by 1 + (2, 1, 5., 5., True), # 2 by 1 + (1, 2, 5., 5., True), # 1 by 2 + (2, 2, 5., 7.5, True), # 2 by 2 (different x/y spacings) + (10, 9, 7.5, 5., True), # 10 by 9 (different x/y spacings), + (1, 1, 5., 5., False), # 1 by 1 + (2, 1, 5., 5., False), # 2 by 1 + (1, 2, 5., 5., False), # 1 by 2 + (2, 2, 5., 7.5, False), # 2 by 2 (different x/y spacings) + (10, 9, 7.5, 5., False), # 10 by 9 (different x/y spacings) +]) +def test_staggered_rectangular_field(N_1, N_2, B_1, B_2, include_last_element): + H = 150. # Borehole length [m] + D = 4. # Borehole buried depth [m] + r_b = 0.075 # Borehole radius [m] + # Generate the bore field + field = gt.boreholes.staggered_rectangle_field( + N_1, N_2, B_1, B_2, H, D, r_b, + include_last_borehole=include_last_element) + # Evaluate the borehole to borehole distances + x = np.array([b.x for b in field]) + y = np.array([b.y for b in field]) + dis = np.sqrt( + np.subtract.outer(x, x)**2 + np.subtract.outer(y, y)**2)[ + ~np.eye(len(field), dtype=bool)] + + if include_last_element or N_1 == 1 or N_2 == 1: + expected_nBoreholes = N_1 * N_2 + elif N_2 % 2 == 0: + expected_nBoreholes = N_2 * (2 * N_1 - 1) / 2 + else: + expected_nBoreholes = (N_2 - 1) * (2 * N_1 - 1) / 2 + N_1 + + assert np.all( + [len(field) == expected_nBoreholes, + np.allclose(H, [b.H for b in field]), + np.allclose(D, [b.D for b in field]), + np.allclose(r_b, [b.r_b for b in field]), + len(field) == 1 or np.isclose( + np.min(dis), min(B_1, np.sqrt(B_2**2 + 0.25 * B_1**2))), + ]) + + +# Test dense_rectangle_field +@pytest.mark.parametrize("N_1, N_2, B, include_last_element", [ + (1, 1, 5., True), # 1 by 1 + (2, 1, 5., True), # 2 by 1 + (1, 2, 5., True), # 1 by 2 + (2, 2, 5., True), # 2 by 2 + (10, 9, 7.5, True), # 10 by 9 + (10, 10, 7.5, True), # 10 by 10 + (1, 1, 5., False), # 1 by 1 + (2, 1, 5., False), # 2 by 1 + (1, 2, 5., False), # 1 by 2 + (2, 2, 5., False), # 2 by 2 + (10, 9, 7.5, False), # 10 by 9 + (10, 10, 7.5, False), # 10 by 10 +]) +def test_dense_rectangle_field(N_1, N_2, B, include_last_element): + H = 150. # Borehole length [m] + D = 4. # Borehole buried depth [m] + r_b = 0.075 # Borehole radius [m] + # Generate the bore field + field = gt.boreholes.dense_rectangle_field( + N_1, N_2, B, H, D, r_b, include_last_borehole=include_last_element) + # Evaluate the borehole to borehole distances + x = np.array([b.x for b in field]) + y = np.array([b.y for b in field]) + dis = np.sqrt( + np.subtract.outer(x, x)**2 + np.subtract.outer(y, y)**2)[ + ~np.eye(len(field), dtype=bool)] + + if include_last_element or N_1 == 1 or N_2 == 1: + expected_nBoreholes = N_1 * N_2 + elif N_2 % 2 == 0: + expected_nBoreholes = N_2 * (2 * N_1 - 1) / 2 + else: + expected_nBoreholes = (N_2 - 1) * (2 * N_1 - 1) / 2 + N_1 + + assert np.all( + [len(field) == expected_nBoreholes, + np.allclose(H, [b.H for b in field]), + np.allclose(D, [b.D for b in field]), + np.allclose(r_b, [b.r_b for b in field]), + len(field) == 1 or np.isclose(np.min(dis), B) + ]) + + # Test L_shaped_field @pytest.mark.parametrize("N_1, N_2, B_1, B_2", [ (1, 1, 5., 5.), # 1 by 1 @@ -217,3 +307,31 @@ def test_circle_field(N, R): len(field) == 1 or np.isclose(np.min(dis), B_min), len(field) == 1 or np.max(dis) <= (2 + 1e-6) * R, ]) + +# Test circle_field +@pytest.mark.parametrize("N, B", [ + (1, 5.), # 1 borehole + (2, 5.), # 2 boreholes + (3, 7.5), # 3 boreholes + (10, 9.), # 10 boreholes + ]) +def test_filled_circle_field(N, B): + H = 150. # Borehole length [m] + D = 4. # Borehole buried depth [m] + r_b = 0.075 # Borehole radius [m] + # Generate the bore field + field = gt.boreholes.filled_circle_field(N, B, H, D, r_b) + # Evaluate the borehole to borehole distances + x = np.array([b.x for b in field]) + y = np.array([b.y for b in field]) + dis = np.sqrt( + np.subtract.outer(x, x)**2 + np.subtract.outer(y, y)**2)[ + ~np.eye(len(field), dtype=bool)] + + assert np.all( + [len(field) == N, + np.allclose(H, [b.H for b in field]), + np.allclose(D, [b.D for b in field]), + np.allclose(r_b, [b.r_b for b in field]), + len(field) == 1 or np.isclose(np.min(dis), B), + ])