Volume of liquid drained from pipe at different locations

Question

This might be more of an algorithm question, but I'm writing it in Python.

I have a set of data on a pipeline that gains and loses altitude as it progresses. My data is two columns, the measure along the pipeline, and the elevation at that measure. There are tens of thousands of rows in my dataset. (these would be columns instead of rows)

Measure: 1, 2, 3, 4, 5, 6, 7, 8, 9, 10

Elevation: 5, 7, 9, 15, 12, 13, 18, 14, 23, 9

In this script, the pipe would be assume to be capped at both ends. The goal is to calculate the total volume of liquid that would drain from a leak at any point in the pipe. Pressure/flow rate don't matter. The main part I'm trying to account for would be all of the catches/valleys (like in a bathroom sink) that liquid would remain in even while the rest of the pipe drains, like this:

https://www.youtube.com/watch?v=o82yNzLIKYo

The pipe radius, and location of the leak would be user-set parameters.

I'm really looking for a push in the right direction, I want to figure this out on my own as much as possible. I'm fine with the programming, but any advice on the actual logic would be helpful, thank you in advanced. enter image description here

Let's say in this graph that a leak appears at point 9 on the x-axis, and the pipe has a known radius r. I'm trying to figure out how to get my script to output the total quantity of liquid in terms of r will be emptied out, regardless of time. And if through damage a leak developed in a pipe, air would come in and water would come out, but not all of the water because of the various catches and different elevations of the pipeline.

Answer 1

If I understand the problem correctly, I think this can be achieved by traversing the pipe left and right from the point of the leak. At each point the current water level is compared with the pipe elevation, either resulting in a submersed point in which the water level remains the same, or a beach and a new dry peak. Interpolation is necessary to calculate the location of the beach.

An implementation is shown below. The bulk of the algorithm is in the traverse function. Hopefully the comments provide adequate description.

#!/usr/bin/python3

import numpy as np
import matplotlib.pyplot as pp

# Positions and elevations
n = 25
h = np.random.random((n, ))
x = np.linspace(0, 1, h.size)

# Index of leak
leak = np.random.randint(0, h.size)

# Traverse a pipe with positions (x) and elevations (h) from a leak index
# (leak) in a direction (step, +1 or -1). Return the positions of the changes
# in water level (y) the elevations at these changes (g) and the water level
# values (w).
def traverse(x, h, leak, step):
    # End of the index range for the traversal
    end = h.size if step == 1 else -1
    # Initialise 1-element output arrays with values at the leak
    y, g, w = [x[leak]], [h[leak]], [h[leak]]
    # Loop from the index adjacent to the leak
    for i in range(leak + step, end, step):
        if w[-1] > h[i]:
            # The new height is less than the old water level. Location i is
            # submerged. No new points are created and the water level stays
            # the same.
            y.append(x[i])
            g.append(h[i])
            w.append(w[-1])
        else:
            # The new height is greater than the old water level. We have a
            # "beach" and a "dry peak".
            # ...
            # Calculate the location of the beach as the position where the old
            # water level intersects the pipe section from [i-step] to [i].
            # This is added as a new point. The elevation and water level are
            # the same as the old water level.
            # ...
            # The if statement is not strictly necessary. It just prevents
            # duplicate points being generated.
            if w[-1] != h[i-step]:
                t = (w[-1] - h[i-step])/(h[i] - h[i-step])
                b = x[i-step] + (x[i] - x[i-step])*t
                y.append(b)
                g.append(w[-1])
                w.append(w[-1])
            # ...
            # Add the dry peak.
            y.append(x[i])
            g.append(h[i])
            w.append(h[i])
    # Convert from list to numpy array and return
    return np.array(y), np.array(g), np.array(w)

# Traverse left and right
yl, gl, wl = traverse(x, h, leak, -1)
yr, gr, wr = traverse(x, h, leak, 1)

# Combine, reversing the left arrays and deleting the repeated start point
y = np.append(yl[:0:-1], yr)
g = np.append(gl[:0:-1], gr)
w = np.append(wl[:0:-1], wr)

# Output the total volume of water by integrating water level minus elevation
print('Total volume =', np.trapz(w - g, y), 'm^3 per unit cross sectional area')

# Display
pp.plot(x, h, '.-', label='elevation')
pp.plot(y, w, '.-', label='water level')
pp.plot([x[leak]], [h[leak]], 'o', label='leak')
pp.legend()
pp.show()

Answer 2

For a constant radius pipe and with a radius much smaller than the variation of elevation, i.e. the section of the pipe is always filled by water. I think, in this case, it doesn't work if the pipe is capped at the ends, some air have to get in to let the water get out. The part of the pipe remaining filled by water is between the left free surface (green circle) and the right free surface (red square). For simplicity, it is assumed that the both ends of the pipe are the points of maximum elevation, otherwise the pipe will empty itself. The equilibrium may be unstable.

import numpy as np
import matplotlib.pyplot as plt
%matplotlib inline

def find_first_intersection(x, y, y_leak):
    for i in range(len(x)-1):
        dy_left = y[i] - y_leak
        dy_right = y[i+1] -  y_leak

        if  dy_left*dy_right < 0:
            x_free = x[i] + (y_leak - y[i])*(x[i+1] - x[i])/(y[i+1] - y[i])
            break

    return x_free

# Generate random data
x = np.linspace(0, 1, 10)
y = np.random.rand(*np.shape(x))
y[0], y[-1] = 1.1, 1.1
x_leak = np.random.rand(1)

# Look for the free surfaces
y_leak = np.interp(x_leak, x, y)

x_free_left = find_first_intersection(x, y, y_leak)
x_free_right =  find_first_intersection(x[::-1], y[::-1], y_leak)

# Plot
plt.plot(x, y, '-', label='pipe');
plt.plot(x_leak, y_leak, 'sk', label='leak')
plt.axhline(y=y_leak, linestyle=':', label='surface level');

plt.plot(x_free_left, y_leak, 'o', label="left free surface");
plt.plot(x_free_right, y_leak, 's', label="right free surface");

plt.legend(bbox_to_anchor=(1.5, 1.)); plt.xlabel('x'); plt.ylabel('y');

I added some annotations on the graph. I think it is confusing that the water will remain in the "confusing part" because I think this is valid only for very small diameter pipe. For a larger pipe the water here will flow through the leak and then estimating the remaining filled part of the pipe is more complicated...