0

I am trying to write write an algorithm which maximizes the sum of Euclidian distance between 2D points. Based on (and thanks to) this discussion, I want to put the whole thing in for loop so that this can be done for multiple points.

Here is what I attempted to do:

from ortools.sat.python import cp_model
model = cp_model.CpModel()

distance_range = 20
number_of_points = 3

a = []    # list to collect all points
# create set of required (x, y) points 
for i in range(number_robots):
    x_i = model.NewIntVar(1, distance_range , 'x_%i' %i)
    y_i = model.NewIntVar(1, distance_range , 'y_%i' %i)
    a.append((x_i, y_i))

# find difference 
x_diff = model.NewIntVar(-distance_range , distance_range , 'x_diff')
y_diff = model.NewIntVar(-distance_range , distance_range , 'y_diff')
difference = []
for i in range(2):
    for j in range(2):
         model.AddAbsEquality(x_diff, a[i][0] - a[j][0])
         model.AddAbsEquality(y_diff, a[i][1] - a[j][1])
         difference.append((x_diff, y_diff))

# square the difference value
x_diff_sq = model.NewIntVar(0, distance_range^2 , '')
y_diff_sq = model.NewIntVar(0, distance_range^2 , '')
difference_sq = []

for i in range(2):
    model.AddMultiplicationEquality(x_diff_sq, a_diff[i][0], a_diff[i][0])
    model.AddMultiplicationEquality(y_diff_sq, a_diff[i][1], a_diff[i][1])
    difference_sq.append((x_diff_sq, y_diff_sq))

# sum of square distances, (x_i^2 + y_i^2)
distance_sq= model.NewIntVar(0, 2*distance_range , '')
for i in range(2):
    model.AddAbsEquality(distance_sq, difference_sq[i][0]+ difference_sq[i][1])

#Objective
model.Maximize(distance_sq)

# Creates a solver and solves the model.
solver = cp_model.CpSolver()
status = solver.Solve(model)


if status == cp_model.OPTIMAL or status == cp_model.FEASIBLE:
    for i in range(number_of_points):
        print("x_%i" %i, solver.Value(a[i][0]))
        print("y_%i" %i, solver.Value(a[i][1]))
else:
    print('No solution found.')

I tried executing the above given code but I never get the distance maximized. All the points have the same coordinate values. So I tried putting the below given constraint for 3 points

model.AllDifferent([a[0][0], a[1][0], a[2][0])    # all values of x should be different
model.AllDifferent([a[0][1], a[1][1], a[2][1])    # all values of y should be different

Then I am not getting any optimal solution. How can I tackle this?

I was able to use the above code to maximize/minimize the distance between a fixed point and multiple variables (decision variables) but few of the points are overlapping and that is why I am trying to enforce this maximization objective between the points.

Laurent Perron
  • 8,594
  • 1
  • 8
  • 22
Ken Adams
  • 25
  • 4
  • As CP-SAT only works with integers, it may be better to use the squared distance. – Erwin Kalvelagen Feb 10 '23 at 15:32
  • The sum of distances is unbounded. Can you explain the complete problem setting ? –  Feb 10 '23 at 15:34
  • @YvesDaoust I think the sum of distances is bounded, can be seen during its declaration. But I also tried to add an extra linear constraint defining its bounding range but I did not see any change any output. – Ken Adams Feb 10 '23 at 17:38
  • @YvesDaoust YES, I can describe you the problem setting. Its very easy actually. I want to maximize the distance between set of points but none of them have any initial value. The only constraint is its bounding range. Also we can add the maximum separation between 2 points. Do you think the issue is because of not having any initial values to the points? – Ken Adams Feb 10 '23 at 17:41
  • @ErwinKalvelagen Yes, CP-SAT only takes integer value. I am using the square distances right? But the format and the OR TOOLs syntax in general is getting me confused. I am very new linear programming. – Ken Adams Feb 10 '23 at 17:43
  • `distance_range^2` is probably something different than what you want. – Erwin Kalvelagen Feb 12 '23 at 10:20
  • Shouldn't ``x_diff`` and ``y_diff`` as well as ``x_diff_sq`` and ``y_diff_sq`` be arrays with an entry for each point? If not, at least you will need to create a new IntVar for each difference and square term in the respective inner loop. In the ``model.AddMultiplicationEquality(x_diff_sq, a_diff[i][0], a_diff[i][0])`` lines, don't you mean ``model.AddMultiplicationEquality(x_diff_sq, difference[i][0], difference[i][0])``? – Christopher Hamkins Feb 13 '23 at 12:41

1 Answers1

2

Here is a working version with the square distance:

from ortools.sat.python import cp_model

model = cp_model.CpModel()

distance_range = 20
number_of_points = 3

# create set of required (x, y) points
x = [
    model.NewIntVar(1, distance_range, f'x_{i}')
    for i in range(number_of_points)
]
y = [
    model.NewIntVar(1, distance_range, f'y_{i}')
    for i in range(number_of_points)
]

# Build intermediate variables and get the list of square differences on each dimension.
objective_terms = []
for i in range(number_of_points - 1):
    for j in range(i + 1, number_of_points):
        x_diff = model.NewIntVar(-distance_range, distance_range, f'x_diff{i}')
        y_diff = model.NewIntVar(-distance_range, distance_range, f'y_diff{i}')
        model.Add(x_diff == x[i] - x[j])
        model.Add(y_diff == y[i] - y[j])
        x_diff_sq = model.NewIntVar(0, distance_range ** 2, f'x_diff_sq{i}')
        y_diff_sq = model.NewIntVar(0, distance_range ** 2, f'y_diff_sq{i}')
        model.AddMultiplicationEquality(x_diff_sq, x_diff, x_diff)
        model.AddMultiplicationEquality(y_diff_sq, y_diff, y_diff)
        objective_terms.append(x_diff_sq)
        objective_terms.append(y_diff_sq)

#Objective
model.Maximize(sum(objective_terms))

# Creates a solver and solves the model.
solver = cp_model.CpSolver()
solver.parameters.log_search_progress = True
solver.parameters.num_workers = 8
status = solver.Solve(model)

# Prints the solution.
if status == cp_model.OPTIMAL or status == cp_model.FEASIBLE:
    for i in range(number_of_points):
        print(f'x_{i}={solver.Value(x[i])}')
        print(f'y_{i}={solver.Value(y[i])}')
else:
    print('No solution found.')

So many wrong things in your code:

  • it does not compile (number_robots vs number_of_points)
  • it is badly indented
  • you reuse the same variables (xdiff, y_diff, x_diff_sq, y_diff_sq) across multiple robots
  • you count distances twice
  • distance ^ 2 is not the square of distance, it is a boolean operator. Use **

So, because we are using square of distances, the solution is degenerate (2 points in a corner, the other in the opposite corner.

x_0=1
y_0=1
x_1=20
y_1=20
x_2=1
y_2=1

To fix the problem, let's introduce an integer variable that will be less that the square root of the distance between 2 robots. To increase precision, we will add a scaling factor.

here is the intermediate code:

objective_terms = []
scaling = 100
for i in range(number_of_points - 1):
    for j in range(i + 1, number_of_points):
        x_diff = model.NewIntVar(-distance_range, distance_range, f'x_diff{i}')
        y_diff = model.NewIntVar(-distance_range, distance_range, f'y_diff{i}')
        model.Add(x_diff == x[i] - x[j])
        model.Add(y_diff == y[i] - y[j])
        x_diff_sq = model.NewIntVar(0, distance_range ** 2, f'x_diff_sq{i}')
        y_diff_sq = model.NewIntVar(0, distance_range ** 2, f'y_diff_sq{i}')
        model.AddMultiplicationEquality(x_diff_sq, x_diff, x_diff)
        model.AddMultiplicationEquality(y_diff_sq, y_diff, y_diff)

        # we create a dist variable such that dist * dist < x_diff_sq + y_diff_sq.
        # To improve the precision, we scale x_diff_sq and y_diff_sq by a constant.
        scaled_distance = model.NewIntVar(0, distance_range * 2 * scaling, '')
        scaled_distance_square = model.NewIntVar(
            0, 2 * distance_range * distance_range * scaling * scaling, '')
        model.AddMultiplicationEquality(scaled_distance_square,
                                        scaled_distance, scaled_distance)
        model.Add(scaled_distance_square <= x_diff_sq * scaling +
                  y_diff_sq * scaling)
        objective_terms.append(scaled_distance)

which gives a reasonable answer:

x_0=1
y_0=1
x_1=20
y_1=1
x_2=1
y_2=20

I also rewrote the example to maximize the minimum euclidian distance between two robots. The code is available here.

Laurent Perron
  • 8,594
  • 1
  • 8
  • 22