Connect series of bridges at midpoint by raising legs

Question

I have been trying for the past few weeks to solve a problem relating to a personal project I am working on. I have tried a few different approaches but none of them seem to give me the results I am looking for so I was hoping someone here could help at least point me in the right direction.

Situation: I want to build a series of bridges that will connect at the midpoint between each bridge. I am given the topography of the ground and the pillars that will serve as the legs of a bridge. (Each bridge has 2 legs to hold the road above it). The ground the legs stand on is not necessarily flat and so the angle between the legs of each bridge will vary. If I were to put a road on top of the legs of each bridge, since the angle of the bridge and the relative height of the legs are different, when the road is placed on the bridge, the road may not necessarily connect with the road that will be built on the following bridge at the midpoint between the 2 bridges. Because of this, I am allowed to raise any/all of the legs of the bridge within a set of bounds in order to make each of the roads placed on top of the legs of the bridges connect with the adjacent roads at the midpoints between the bridges. (The road extends beyond the legs of the bridge and terminates at the midpoint between any adjacent bridges)

Example: I have a list of points along the x-axis where the legs of the bridge will be placed: x = [0, 200, 300, 500, 600, 800] (The distance between each leg of a bridge is 200, the distance from the 2nd leg of a bridge to the first leg of an adjacent bridge is 100. Bridge 1 will rest on the legs that sit at 0 and 200 units, Bridge 2 will rest on the legs that sit at 300 and 500 units, etc... The midpoint where the roads should connect between bridges 1 and 2 lies at 250 units).

I also know that the base height of the legs of the bridge lie at these heights: y = [6, 1, 6, 9, 9, 12] above sea level (above 0).

Each of the legs can only be moved up within the bounds of bounds = (0, 10) (each leg can go up anywhere from 0 units to 10 units from its base height) for each leg in order to make the roads placed on top of the legs connect to adjacent roads at the midpoint between bridges.

I am trying to create an algorithm that returns a list of adjustments that need to be made to each leg in order to make sure each bridge connects to its neighbor bridges at their shared midpoints.

My current belief so far is that this is an optimization problem where I am trying to minimize the vertical distance between the ends of the roads of the bridges at the midpoints between the bridges. However, I have tried multiple different iterations and have had no luck so far. I am very new to optimization algorithms so this is probably the reason I have failed so far. (using scipy.optimize.minimize)

There is also a chance that this is not an optimization problem and it could be solved by other means and in that case I would love to listen to other approaches from anyone who could help point me in the right direction.

If anyone can help me solve this problem or at least help point me in the right direction I would really appreciate it.

Any feedback is greatly appreciated and I will be sure to respond quickly to any questions or comments you have. Thank you so much!

EDIT: For example: This is what our bridges would look like before any adjustments (Red is the legs of the bridge, green is the midpoints between the bridges, blue arrows are where the roads of the bridges would terminate, blue horizontal bars represent the base height of the bridge )

And this is what our bridge would look like with a series of successful adjustments: (Black arrows are where the new roads would rest on top of the legs)

As you can see, to make our bridges connect, one possible solution was to raise leg 2 and leg 5 (by roughly 4 units and 1 unit respectively)

Answer 1

By request, here's how to formulate this as a linear program and solve with OR-Tools.

For each leg, make a variable. The lower bound of this variable is the leg's current height. The upper bound of this variable is the leg's maximum height (e.g., current + 10).

For each adjacent pair of bridges, we constrain the midpoint between their neighboring legs to be at the same height for both bridges.

For an objective, I'll settle for minimizing the maximum difference in slopes. We use the slope instead of the angle because it's linear and differences in slopes decently approximates differences in angles. (It's definitely possible to optimize the difference in angles directly, but we'd have to switch to cvxpy or something that has support for non-linear objectives.)

from ortools.linear_solver import pywraplp

x = [0, 200, 300, 500, 600, 800]
y = [6, 1, 6, 9, 9, 120]
initial_bounds = [(0, 10), (0, 9), (0, 9.5), (0, 8), (0, 7), (0, 9)]

# Returns the y value at x2 of the line between (x0,y0) and (x1,y1).
def extrapolate(x0, y0, x1, y1, x2):
    return y0 + (x2 - x0) / (x1 - x0) * (y1 - y0)


def solve(is_phase_one, bounds):
    solver = pywraplp.Solver.CreateSolver("GLOP")
    if is_phase_one:
        objective = solver.NumVar(0, solver.infinity(), "objective")
        new_y = [
            solver.NumVar(-solver.infinity(), solver.infinity(), "y%d" % i)
            for i in range(len(y))
        ]
        for i in range(len(y)):
            solver.Add(new_y[i] >= y[i] + bounds[i][0] - objective)
            solver.Add(new_y[i] <= y[i] + bounds[i][1] + objective)
    else:
        objective = 0
        new_y = [
            solver.NumVar(y[i] + bounds[i][0], y[i] + bounds[i][1], "y%d" % i)
            for i in range(len(y))
        ]
    for i in range(3, len(y), 2):
        midpoint = (x[i - 2] + x[i - 1]) / 2
        solver.Add(
            extrapolate(x[i - 3], new_y[i - 3], x[i - 2], new_y[i - 2], midpoint)
            == extrapolate(x[i - 1], new_y[i - 1], x[i], new_y[i], midpoint)
        )
    if not is_phase_one:
        for i in range(1, len(y), 2):
            z = solver.NumVar(0, solver.infinity(), "")
            slope = (new_y[i] - new_y[i - 1]) / (x[i] - x[i - 1])
            solver.Add(z >= slope)
            solver.Add(z >= -slope)
            objective += z
    solver.Minimize(objective)
    status = solver.Solve()
    assert status == pywraplp.Solver.OPTIMAL
    if is_phase_one:
        return [
            (
                bounds[i][0] - objective.solution_value(),
                bounds[i][1] + objective.solution_value(),
            )
            for i in range(len(y))
        ]
    else:
        return [variable.solution_value() for variable in new_y]


relaxed_bounds = solve(True, initial_bounds)
print(relaxed_bounds)
heights = solve(False, relaxed_bounds)
print(heights)

Sample output:

[(-6.1999999999999975, 16.199999999999996), (-6.1999999999999975, 15.199999999999998), (-6.1999999999999975, 15.699999999999998), (-6.1999999999999975, 14.199999999999998), (-6.1999999999999975, 13.199999999999998), (-6.1999999999999975, 15.199999999999998)]
[-0.1999999999999975, 16.199999999999996, 16.800000000000033, 2.8000000000000025, 22.199999999999996, 113.8]