Or-tools: CP-Sat: Python 3.7.2, ortools 7.0 wrong optimization result

[Phibedy]

My test-cases fail on 3.7.2 with version 7.0.6546. It says that it found an optimal solution of 32, but it should be 16. Using python 3.6.7 with version 7.0.6546 it works like expected.

I also get some deprecation warnings:

/usr/local/lib/python3.7/site-packages/google/protobuf/descriptor.py:47
  /usr/local/lib/python3.7/site-packages/google/protobuf/descriptor.py:47: DeprecationWarning: Using or importing the ABCs from 'collections' instead of from 'collections.abc' is deprecated, and in 3.8 it will stop working
    from google.protobuf.pyext import _message

/usr/local/lib/python3.7/site-packages/google/protobuf/internal/well_known_types.py:788
  /usr/local/lib/python3.7/site-packages/google/protobuf/internal/well_known_types.py:788: DeprecationWarning: Using or importing the ABCs from 'collections' instead of from 'collections.abc' is deprecated, and in 3.8 it will stop working
    collections.MutableMapping.register(Struct)

-- Docs: https://docs.pytest.org/en/latest/warnings.html

I added the test cases after upgrading ortools from 6.10 to 7.0. So I don't know if it worked before.
Does someone have experienced this behavior?

Python 3.6.8 says that the optimum is 30, it differs but still wrong.

What I do

Put testcase in different docker images (in brackets the calculated optimum):

• python:3.7-slim fails (32) is python 3.7.1
• python:3.6-slim fails (30) is python 3.6.8
• python:3.6.7 fails (32)
  Those are all debian packages but I am on ubuntu and the testcase works.

Weird that the result differs from 3.6.7, 3.6.8, 3.7.1. Maybe some other dependency differs, I will try to find a working docker image and a minimal working example. At the moment I can't share the code.

[Phibedy]

My problem:
I am calculating two cost functions after another. After I solved the the first one I call:

self.model._CpModel__model.solution_hint.Clear()
for i,field in enumerate(self.model._CpModel__model.variables):
    self.model._CpModel__model.solution_hint.vars.extend([i])
    self.model._CpModel__model.solution_hint.values.extend([self.solver._CpSolver__solution.solution[i]])

to speed up the second optimization.
On some systems the output of the second optimization will be the same as the given solution_hint even if it is not even close to the optimal solution.

If I don't supply an solution_hint it works on all systems.
The weird part is, that it does not work on another ubuntu 18.04 computer with the same virtual env. On some arch linux with python 3.7.1 it worked. So the question is, what could go wrong?

[Phibedy]

Unluckily I can't come up with small working example which shows the error atm.

What I found out, that if I force model.Add(costfunction==16) it returns 16 on all machines even though they would return 32 if I just call Minimize()

[Phibedy]

It was quite painful to write a "small" working example but here it is: If you comment out the solution hint it works.
What it does:
Create n Handlers, each handler has 2 slots (start,end) and contains a list which states if the handler is currently used. Now the time each handler can work is restricted and the minimal change in the is_used list is searched. Slots don't overlap and there is always one 0 element between them.
Therefore each slot adds 2 to the cost function. The minimal optimal solution is therefore n*2 as only one slot is used per handler.

the _f_ test cases show, that the optimum is valid and can be reached. The code is not nice at all but it reproduces the error on all python docker images I tried.

import multiprocessing
from ortools.sat.python import cp_model

COUNTS = 10
MIN_VAL = 2
MAX_VAL = 4

class Slot:
    def __init__(self,model):
        self.model = model
        self.start = model.NewIntVar(-1,COUNTS-1,name="")#inclusive
        self.duration = model.NewIntVar(0,COUNTS,"")
        self.end = model.NewIntVar(-1,COUNTS,name="") #exclusive
        self.is_used = model.NewBoolVar("")
        self.is_used_table = []
        for t in range(0,COUNTS):
            self.is_used_table.append(model.NewBoolVar(""))

        #add constraints
        #only want to asign a working duration if a shift start is set
        model.Add((self.start==-1)).OnlyEnforceIf(self.is_used.Not())
        model.Add((self.start!=-1)).OnlyEnforceIf(self.is_used)
        model.Add((self.duration==0)).OnlyEnforceIf(self.is_used.Not())
        model.Add((self.duration>0)).OnlyEnforceIf(self.is_used)
        model.Add((self.start<self.end)).OnlyEnforceIf(self.is_used)
        model.Add(self.end == self.start+self.duration)

        #nothing before start
        for t in range(0,COUNTS):
            if t > 0:
                model.Add(self.start==t).OnlyEnforceIf([k.Not() for k in self.is_used_table[:t]]+[self.is_used_table[t]])
        #nothing after end
        for t in range(1,COUNTS):
            model.Add(self.end==t).OnlyEnforceIf([k.Not() for k in self.is_used_table[t:]]+[self.is_used_table[t-1]])

        #if we work at the shift there has to be at least one slot, otherwise it could let working_table be only zeroes 0
        model.Add(sum(self.is_used_table)>=self.is_used)
        #number of slots marked as 1 and duration has to match
        model.Add(sum(self.is_used_table)==self.duration)

        #fix working table to start and End
        for t in range(0,COUNTS):
            inersect = self.is_used_table[t]#
            model.Add(t >= self.start).OnlyEnforceIf([inersect,self.is_used])
            model.Add(t < self.end).OnlyEnforceIf([inersect,self.is_used])


class Handler:

    def __init__(self,model,list_count):
        self.model=model
        self.is_used_table = [model.NewBoolVar(name="") for _ in range(COUNTS)]
        self.used_changes = []
        self.possible_slots = []

        for s in range(list_count):
            self.possible_slots.append(Slot(model))
            """
            Shifts are sorted from the first to the last one. Obsolete shifts are at the end
            """
            if s > 0:
                model.Add(self.possible_slots[s].is_used==0).OnlyEnforceIf(self.possible_slots[s-1].is_used.Not())
                #sort shifts and intervals shouldn't overlap but only if they are working at that shift
                #> because otherwise shifts
                model.Add(self.possible_slots[s].start >
                          self.possible_slots[s-1].end).OnlyEnforceIf([self.possible_slots[s-1].is_used,self.possible_slots[s].is_used])
                model.Add(self.possible_slots[s].start >
                          self.possible_slots[s-1].start).OnlyEnforceIf([self.possible_slots[s-1].is_used,self.possible_slots[s].is_used])

        for t in range(0,COUNTS):
            #Combine slots with is_used_table
            overlap_list = []
            for s in self.possible_slots:
                overlap_list.append(s.is_used_table[t])

            #We just want exactly one overlap or none
            model.Add(sum(overlap_list) == self.is_used_table[t])


    def __count_used_changes(self):
        if len(self.used_changes) != 0:
            raise ValueError("count_used_changes can only be called once")
        rt = self.is_used_table
        for i,r in enumerate(rt):
            if i > 0:
                is_diff = self.model.NewBoolVar("")
                self.model.Add(r!=rt[i-1]).OnlyEnforceIf(is_diff)
                self.model.Add(r==rt[i-1]).OnlyEnforceIf(is_diff.Not())
                self.used_changes.append(is_diff)

    def get_used_changes(self):
        if len(self.used_changes) == 0:
            self.__count_used_changes()
        return self.used_changes

def refine_x_days(n_d,force_solution=False,print_to_file=None):
    model = cp_model.CpModel()
    solver = cp_model.CpSolver()
    solver.parameters.max_time_in_seconds = 30
    solver.parameters.num_search_workers = max(1,multiprocessing.cpu_count()-1)
    handlers = []
    for d in range(n_d):
        wd = Handler(model,2)#,lambda t: has_role_restriction(d,t))
        handlers.append(wd)

    total_changes = []
    for wd in handlers:
        total_changes.extend(wd.get_used_changes())
        #dont allow ones at the start/end -> handler adds at least a cost of 2
        model.Add(wd.is_used_table[0] == 0)
        model.Add(wd.is_used_table[-1] == 0)
        model.AddSumConstraint([ws.duration for ws in wd.possible_slots],MIN_VAL,MAX_VAL)

    model.Minimize(sum([]))#basically 0
    if print_to_file != None:
        with open(print_to_file+"_before.txt","a") as f:
            f.write("######################################## \n")
            f.write(str(model.Proto()))
            f.write("######################################## \n")
    result = solver.Solve(model)
    assert result == cp_model.OPTIMAL
    model._CpModel__model.solution_hint.Clear()
    for i,field in enumerate(model._CpModel__model.variables):
        model._CpModel__model.solution_hint.vars.extend([i])
        model._CpModel__model.solution_hint.values.extend([solver._CpSolver__solution.solution[i]])

    if force_solution:
        model.Add(sum(total_changes) == n_d*2)
    model.Minimize(sum(total_changes))
    if print_to_file != None:
        with open(print_to_file+"_after.txt","a") as f:
            f.write("######################################## \n")
            f.write(str(model.Proto()))
            f.write("######################################## \n")
    result = solver.Solve(model)
    assert result == cp_model.FEASIBLE or result == cp_model.OPTIMAL
    if result == cp_model.OPTIMAL:
        assert solver.Value(sum(total_changes)) == n_d*2
    else:
        assert solver.Value(sum(total_changes)) <= n_d*2*2

#force right solution
def test_f_one():
    refine_x_days(1,True)
def test_f_two():
    refine_x_days(2,True)
def test_f_three():
    refine_x_days(3,True)
def test_f_four():
    refine_x_days(4,True)
def test_f_five():
    refine_x_days(5,True)
def test_f_six():
    refine_x_days(6,True)
def test_f_seven():
    refine_x_days(7,True)

FILE = "protos"
#don'tforce right solution
def test_one():
    refine_x_days(1,print_to_file=FILE)
def test_two():
    refine_x_days(2,print_to_file=FILE)
def test_three():
    refine_x_days(3,print_to_file=FILE)
def test_four():
    refine_x_days(4,print_to_file=FILE)
def test_five():
    refine_x_days(5,print_to_file=FILE)
def test_six():
    refine_x_days(6,print_to_file=FILE)
def test_seven():
    refine_x_days(7,print_to_file=FILE)

On my system ubuntu with the latest ortools and protobuf 3.7.1 it works. On another ubuntu 18.04 system with same virtual-env it does not work.

The output is not even deterministic, earch run the result looks different, only the forced ones are always passed. Example output:

minimal_example.py::test_f_one PASSED                                                                                                                            [  7%]
minimal_example.py::test_f_two PASSED                                                                                                                            [ 14%]
minimal_example.py::test_f_three PASSED                                                                                                                          [ 21%]
minimal_example.py::test_f_four PASSED                                                                                                                           [ 28%]
minimal_example.py::test_f_five PASSED                                                                                                                           [ 35%]
minimal_example.py::test_f_six PASSED                                                                                                                            [ 42%]
minimal_example.py::test_f_seven PASSED                                                                                                                          [ 50%]
minimal_example.py::test_one PASSED                                                                                                                              [ 57%]
minimal_example.py::test_two FAILED                                                                                                                              [ 64%]
minimal_example.py::test_three FAILED                                                                                                                            [ 71%]
minimal_example.py::test_four FAILED                                                                                                                             [ 78%]
minimal_example.py::test_five PASSED                                                                                                                             [ 85%]
minimal_example.py::test_six PASSED                                                                                                                              [ 92%]
minimal_example.py::test_seven FAILED                                                                                                                            [100%]

        result = solver.Solve(model)
        assert result == cp_model.FEASIBLE or result == cp_model.OPTIMAL
        if result == cp_model.OPTIMAL:
>           assert solver.Value(sum(total_changes)) == n_d*2
E           assert 28 == 14
E             -28
E             +14

Any help is appreciated :)

[lperron]

before doing that, can you print the generated model (str(model.Proto())) and compare ?

[Phibedy]

I updated the example, the model before setting the hint is the same on all machines (diff server_protos_before.txt client_protos_before.txt). After adding the hints not. But only the hints differ not the models themselves. The output of the first testcase and the output of the model with set hints:
https://pastebin.com/ifL4Rw5g

[lperron]

can you try with the following hinting code:

model.Proto().solution_hint.Clear()
for i in range(len(model.Proto().variables):
    model.Proto().solution_hint.vars.append(i)
    model.Proto().solution_hint.values.append(solver.ResponseProto().solution[i])

[lperron]

The following code works for me with python 3.7 on my mac

    from ortools.sat.python import cp_model

    COUNTS = 10
    MIN_VAL = 2
    MAX_VAL = 4

    class Slot:
        def __init__(self,model):
            self.model = model
            self.start = model.NewIntVar(-1,COUNTS-1,name="")#inclusive
            self.duration = model.NewIntVar(0,COUNTS,"")
            self.end = model.NewIntVar(-1,COUNTS,name="") #exclusive
            self.is_used = model.NewBoolVar("")
            self.is_used_table = []
            for t in range(0,COUNTS):
                self.is_used_table.append(model.NewBoolVar(""))

            #add constraints
            #only want to asign a working duration if a shift start is set
            model.Add((self.start==-1)).OnlyEnforceIf(self.is_used.Not())
            model.Add((self.start!=-1)).OnlyEnforceIf(self.is_used)
            model.Add((self.duration==0)).OnlyEnforceIf(self.is_used.Not())
            model.Add((self.duration>0)).OnlyEnforceIf(self.is_used)
            model.Add((self.start<self.end)).OnlyEnforceIf(self.is_used)
            model.Add(self.end == self.start+self.duration)

            #nothing before start
            for t in range(0,COUNTS):
                if t > 0:
                    model.Add(self.start==t).OnlyEnforceIf([k.Not() for k in self.is_used_table[:t]]+[self.is_used_table[t]])
            #nothing after end
            for t in range(1,COUNTS):
                model.Add(self.end==t).OnlyEnforceIf([k.Not() for k in self.is_used_table[t:]]+[self.is_used_table[t-1]])

            #if we work at the shift there has to be at least one slot, otherwise it could let working_table be only zeroes 0
            model.Add(sum(self.is_used_table)>=self.is_used)
            #number of slots marked as 1 and duration has to match
            model.Add(sum(self.is_used_table)==self.duration)

            #fix working table to start and End
            for t in range(0,COUNTS):
                inersect = self.is_used_table[t]#
                model.Add(t >= self.start).OnlyEnforceIf([inersect,self.is_used])
                model.Add(t < self.end).OnlyEnforceIf([inersect,self.is_used])


    class Handler:

        def __init__(self,model,list_count):
            self.model=model
            self.is_used_table = [model.NewBoolVar(name="") for _ in range(COUNTS)]
            self.used_changes = []
            self.possible_slots = []

            for s in range(list_count):
                self.possible_slots.append(Slot(model))
                """
                Shifts are sorted from the first to the last one. Obsolete shifts are at the end
                """
                if s > 0:
                    model.Add(self.possible_slots[s].is_used==0).OnlyEnforceIf(self.possible_slots[s-1].is_used.Not())
                    #sort shifts and intervals shouldn't overlap but only if they are working at that shift
                    #> because otherwise shifts
                    model.Add(self.possible_slots[s].start >
                              self.possible_slots[s-1].end).OnlyEnforceIf([self.possible_slots[s-1].is_used,self.possible_slots[s].is_used])
                    model.Add(self.possible_slots[s].start >
                              self.possible_slots[s-1].start).OnlyEnforceIf([self.possible_slots[s-1].is_used,self.possible_slots[s].is_used])

            for t in range(0,COUNTS):
                #Combine slots with is_used_table
                overlap_list = []
                for s in self.possible_slots:
                    overlap_list.append(s.is_used_table[t])

                #We just want exactly one overlap or none
                model.Add(sum(overlap_list) == self.is_used_table[t])


        def __count_used_changes(self):
            if len(self.used_changes) != 0:
                raise ValueError("count_used_changes can only be called once")
            rt = self.is_used_table
            for i,r in enumerate(rt):
                if i > 0:
                    is_diff = self.model.NewBoolVar("")
                    self.model.Add(r!=rt[i-1]).OnlyEnforceIf(is_diff)
                    self.model.Add(r==rt[i-1]).OnlyEnforceIf(is_diff.Not())
                    self.used_changes.append(is_diff)

        def get_used_changes(self):
            if len(self.used_changes) == 0:
                self.__count_used_changes()
            return self.used_changes

    def refine_x_days(n_d,force_solution=False,print_to_file=None):
        model = cp_model.CpModel()
        handlers = []
        for d in range(n_d):
            wd = Handler(model,2)#,lambda t: has_role_restriction(d,t))
            handlers.append(wd)

        total_changes = []
        for wd in handlers:
            total_changes.extend(wd.get_used_changes())
            #dont allow ones at the start/end -> handler adds at least a cost of 2
            model.Add(wd.is_used_table[0] == 0)
            model.Add(wd.is_used_table[-1] == 0)
            model.AddSumConstraint([ws.duration for ws in wd.possible_slots],MIN_VAL,MAX_VAL)

        if print_to_file != None:
            with open(print_to_file+"_before.txt","a") as f:
                f.write("######################################## \n")
                f.write(str(model.Proto()))
                f.write("######################################## \n")

        solver = cp_model.CpSolver()
        solver.parameters.max_time_in_seconds = 30
        solver.parameters.num_search_workers = 8

        result = solver.Solve(model)
        assert result == cp_model.FEASIBLE
        model.Proto().solution_hint.Clear()
        for i,field in enumerate(model._CpModel__model.variables):
            model.Proto().solution_hint.vars.append(i)
            model.Proto().solution_hint.values.append(solver.ResponseProto().solution[i])

        if force_solution:
            model.Add(sum(total_changes) == n_d*2)
        model.Minimize(sum(total_changes))
        if print_to_file != None:
            with open(print_to_file+"_after.txt","a") as f:
                f.write("######################################## \n")
                f.write(str(model.Proto()))
                f.write("######################################## \n")
        result = solver.Solve(model)
        print('Problem = ', n_d, '\n', solver.ResponseStats())
        assert result == cp_model.FEASIBLE or result == cp_model.OPTIMAL
        if result == cp_model.OPTIMAL:
            assert solver.Value(sum(total_changes)) == n_d*2
        else:
            assert solver.Value(sum(total_changes)) <= n_d*2*2

    #force right solution
    def test_f_one():
        refine_x_days(1,True)
    def test_f_two():
        refine_x_days(2,True)
    def test_f_three():
        refine_x_days(3,True)
    def test_f_four():
        refine_x_days(4,True)
    def test_f_five():
        refine_x_days(5,True)
    def test_f_six():
        refine_x_days(6,True)
    def test_f_seven():
        refine_x_days(7,True)

    test_f_one()
    test_f_two()
    test_f_three()
    test_f_four()
    test_f_five()
    test_f_six()
    test_f_seven()

[Phibedy]

I don't have CpSolver.ResponseProto(self). Ah you just added it 5 days ago

In [2]: ortools.__version__
Out[2]: '7.0.6546'

The forced testcases work (using my code), but if refine_x_days(NUMBER,False) fail

def test_one():
    refine_x_days(1,False)
def test_two():
    refine_x_days(2,False)
def test_three():
    refine_x_days(3,False)
def test_four():
    refine_x_days(4,False)
def test_five():
    refine_x_days(5,False)
def test_six():
    refine_x_days(6,False)
def test_seven():
    refine_x_days(7,False)
test_one()
test_two()
test_three()
test_four()
test_five()
test_six()
test_seven()

No running the tests with the current master commit still fail on some machines if you run the test several times.

[syxolk]

I just cloned the or-tools repository and checked out the master branch. This is the last commit:

commit ee89733a600fd2d424df1392c85f557c2bc2a045 (HEAD -> master, origin/master)
Author: Laurent Perron <[email protected]>
Date:   Mon Apr 1 13:27:21 2019 +0200

    change CP-SAT internals

Then I compiled it from sources and installed it, resulting in this version:

Finished processing dependencies for ortools==7.0.6571

My OS/Python/or-tools versions:

$ lsb_release -d
Description:    Ubuntu 18.04.2 LTS

$ python3 --version
Python 3.6.7

$  python3 -c "import ortools; print(ortools.__version__)"
7.0.6571

Futher, I am using the modified example of @lperron but readded the test_one/two/... test cases.
See here fore the full code: https://gist.github.com/syxolk/f185b8c1f2ea689d969b6190a5b42f45

I executed it for like 10 times, all looks good (shortened output):

$ pytest minimal_example.py

cpsat/minimal_example.py ..............
============== 14 passed in 1.47 seconds =================

But once in a while (can be after 10 to 20 runs), I get failed test cases like this (shortened output):

$ pytest minimal_example.py

cpsat/minimal_example.py ........F.....
n_d = 2, force_solution = False, print_to_file = None

sometimes other test cases fail:

n_d = 1, force_solution = False, print_to_file = None

In fact, I threw together a short bash script that runs the test cases 100-times:

#!/bin/bash
FAILS=0
for i in {1..100}; do
    if ! pytest minimal_example.py; then
        ((FAILS++))
    fi
done

echo $FAILS

And I got: 11, which means that 11 out of 100 test runs failed.

[lperron]

I just ran it 100 times without any error. I will continue digging in.

[Phibedy]

On my computer (Ubuntu 18.04.2 LTS) it works also 100 times. That is quite interesting, I thought that maybe docker+debian breaks it.

##################Edit

Dockerfile that works on my computer (Ubuntu 18.04.2 LTS) and not on syxolk computer (Ubuntu 18.04.2 LTS).

FROM python:3.6.7-slim
RUN pip install ortools==7.0.6546 pytest
WORKDIR /app
COPY minimal_example.py .
CMD ["pytest", "minimal_example.py"]

• We have the same kernel
• don't know if it matters but syxolk has avx2, my cpu only supports avx. I just tested in on another computer with avx2 and it works.

Code I used for test:
https://github.com/Phibedy/or_tools_roulette

So even if you start a docker image it fails on machines that fail and vice versa.
It seems to be hardware/driver related?

[lperron]

I did manage to replicate it.
I have a problem in protobuf format that returns 2 in sequential mode and 4 in MT mode.
Now, we need to debug it.

Thanks for the report!

[lperron]

More info, the problem happens in single thread when using the core algorithms with linearization disabled.
parameters.optimize_with_core = True
parameters.linearization_level = 0
parameters.num_search_workers = 1

Fix is coming.

[lperron]

This is fixed on master.

Thanks again for the bug report.

[Phibedy]

Thank you, works :+1:

[syxolk]

Thank you for the quick fix! I can confirm that the minimal example didn't have a single failed test case on my machine (with and without multi-threading) in 100 test runs with the current master branch.

btw: What is your release schedule? It's fine for us to compile from source but it's more convenient to just install from pip :)

[lperron]

Most likely a minor release in April.
We would like to fix the cmake build, so timing is unpredictable.
Laurent Perron | Operations Research | [email protected] | (33) 1 42 68 53
00

Le jeu. 4 avr. 2019 à 21:49, Hans [e-mail] a écrit :

Thank you for the quick fix! I can confirm that the minimal example had
not a single failed test case on my machine (with and without
multi-threading) in 100 test runs with the current master branch.

btw: What is your release schedule? It's fine for us to compile from
source but it's more convenient to just install from pip :)

—
You are receiving this because you modified the open/close state.
Reply to this email directly, view it on GitHub
https://github.com/google/or-tools/issues/1162#issuecomment-480038360,
or mute the thread
https://github.com/notifications/unsubscribe-auth/AKj17Y2ap8KEY1IH9q4ZUoO81Ngt6Gmhks5vdldAgaJpZM4cRGpi
.