search

IBMDecisionOptimization / Decision-Optimization-with-CPLEX-samples

20 stars 9 forks source link

CPLEX MIP Limit Solutions returning bad solver status = error #3

Closed bhimabi closed 3 years ago

bhimabi commented 3 years ago

Code: Python CPLEX Version: CPLEX Studio 12.10

I wanted my MIP solver to exit when it finds its first non-zero solution (i.e. 2nd solution).

So i set option 'mip_limit_solutions' = 2 to run my code. It does run and find a solution but when it tries to returns the output it gives the error. I tried to increase the limit to 999(i.e. solver exits due to finding optimal instead of this limit) and it finds the global optimal and returns the solution. Hence, I suspect the error related to the solution limit. I have included the terminal output and a code snippet below. I have not included the entire model as it solves perfectly if that option is not used.

In case the formatting below is bad, here is a screenshot of the terminal output.Terminal Error Output

Detecting symmetries...

  0     2     2466.7967    37     1475.9071     2465.8441      280   67.07%

Elapsed time = 0.50 sec. (218.90 ticks, tree = 0.01 MB, solutions = 1)

  • 228 155 integral 0 1476.1593 2315.5068 3073 56.86%

Implied bound cuts applied: 1

Flow cuts applied: 4

Mixed integer rounding cuts applied: 19

Gomory fractional cuts applied: 1

Root node processing (before b&c):

Real time = 0.50 sec. (218.93 ticks)

Sequential b&c:

Real time = 0.16 sec. (71.10 ticks) ------------ Total (root+branch&cut) = 0.66 sec. (290.03 ticks)

    0.67 seconds required for solver

Traceback (most recent call last):

File "Castro MILP.py", line 324, in 

X.solve(report_timing=True, tee=True, warmstart=True)

File "\venv\lib\site-packages\pyomo\solvers\plugins\solvers\persistent_solver.py", line 527, in solve default_variable_value=self._default_variable_value)

File "\venv\lib\site-packages\pyomo\core\base\PyomoModel.py", line 242, in load_from % str(results.solver.status))

ValueError: Cannot load a SolverResults object with bad status: error

X.options['threads'] = 1
X.options['mip_limits_solutions'] = 2
X.options['timelimit'] = 1
X.options['emphasis_mip'] = 1
X.solve(report_timing=True, tee=True, warmstart=True)
vlkong commented 3 years ago

Hello,

In the traces you show, it looks like cplex is returning no solution , and pyomo does not handle this correctly.

I'd be happy to help you if you can provide an example using only cplex or docplex, but I have no expertise in pyomo code.

bhimabi commented 3 years ago

I shall try converting this into a cplex/docplex version and try it again. Is there any difference between the 2(cplex vs docplex) in terms of computational time?

vlkong commented 3 years ago

CPLEX is the matrix api that maps directly to CPLEX native.

docplex is an object oriented modeling layer that does the mapping to the CPLEX native API for you. docplex uses CPLEX to solve, so for the exact same model, the solve times should be the same.

bhimabi commented 3 years ago

Sorry. Ive been busy recently and just started on the conversion. I am almost done with the conversion. I wanted to know is there any way to fix a variable's value without fixing its lower and upper bound?

Context: I have a MILP model where I fix all the binary variables making it into a LP. I keep resolving it while fixing the binary variables to different values each time. Since the model is identical between solves, I wrote the model as an MILP before fixing the binaries.

bhimabi commented 3 years ago

I tinkered around with docplex and had a few other questions.

I also managed to get the solver times through "Model.solve_details.time". This however, provides the total time including the model building time. Is there a way to just access the cplex solve durations?

How does resolve() work as compared to solve()? Does it just look for updates and solve the model after that (thus resulting in faster solve times)?

Does docplex/cplex warm start LPs be default? If not would it be possible to used the previous iteration's solution as the initial basis? Would set_lp_start_basis() be the right approach?

Thanks

PhilippeCouronne commented 3 years ago

Hello,

I will try to answer your( many) questions:

  1. Model.solve_details.time returns CPLEX solve time. Model build time is not measured by DOcplex. You'll have to write Python code for this (use time.time() for this). If you have performance problems for the build phase, we have a notebook explaining how to avoid traps and pitfalls here: https://github.com/IBMDecisionOptimization/docplex-examples/blob/master/examples/mp/jupyter/efficient.ipynb

  2. repeated solves: By default, DOcplex keeps the state of the previous solve, so that adding a constraint, or modifying a constraint in place will use the last solution as warm start. If this is not the desired behavior add the clean_before_solve=True keyword argument to solve to start a new solve from fresh. Same for LPs.

m.solve(clean_before_solve=True)

to start a new solve from scratch (default is False)

  1. Fixing variable values: you can either set bounds (lb+ub) or add a constraint:

    m.add(x ==3)

or a ranged constraint:

m.add_range_constraint( 3 -001 <= x <= 3 +001)

The third solution has the advantage of letting you adjust the range precision, while the second uses the standard cplex rhs precision.

Let me know if this helps or if you have other questions.

bhimabi commented 3 years ago

That solved a number of issues. But that got me thinking of the purpose of resolve().

After adding in the new 'm.add(x ==3)' constraints, running resolve() did nothing while the solve() function actually worked. Thus, what is the difference between the 2 and when would one use the resolve function?

Secondly, I have created my variables as var_dict

T_n = m.continuous_var_dict(nn, name='T_n')

If I wanted to operate on each variable(change lb/ub/etc) or extract all the variable solution values as an array/dictionary, is there a built in function for that? Currently I am using a for loop to iterate through each variable. I am guessing using built in functionality(if any) would be faster than normal for loops.

Lastly, I was trying to access the solve status to know if a model is optimal or infeasible. I was looking at the Cplex manual status code website linked below.

According to this, optimal should have a value of 1 and infeasible a value of 3. However, it am getting a value of 2 for my optimal solution. Does docplex use a different status code system? In[2]: m.solve_status Out[2]: <JobSolveStatus.OPTIMAL_SOLUTION: 2>

PhilippeCouronne commented 3 years ago

Hello,   The resolve() function is marked as internal, and should not be used for "solving again", actually it resolves (as in resolution) internal computed expressions. It will be renamed into "_resolve". Only solve() is to be used for solving. Sorry for the misunderstanding.

As for variable dicts, There is a method get_value_dict on solution objects, which return a dictionary of variable -> values , is this what you were looking for?

JobSolveStatus is an enumerated value which was used to indicate the result of cloud solve. What you are looking for is the CPLEX solve status which you can get from solve_details, e.g.

mdl.solve_details.status

The exact CPLEX code, as you mentioned, can also be obtained form solve details by

mdl.solve_details.status_code

Actually, Model.get_solve_status() will be deprecated in next versions, as cloud solve is also deprecated. You should use Model.solve_details.status (or status_code) instead.

bhimabi commented 3 years ago

That get_value_dict was perfect. However, is there a similar function to set the lb/ub of all the variables(or a set of them) to a dictionary value? Currently I have set all variables to have a lb and ub of 0. I then generate a list of variables which need to be equal to 1 and loop through all of them to set their lb and ub to 1. After this, I reset the lb and ub back to 0. I repeat this 1000s of time in my code and this step is eating up a significant portion of my run time. I wanted to know if there is a more efficient built in method/function?

On a related note, I was cloning the model periodically to have a copy of the model and its state stored. When I needed to get the solution, I was not able to use get_value_dict since the dict of vars I previously had would have been for the master model from which the clone was created. 1 alternative was to manually loop through the various dicts and clone() each variable. Once again, I wanted to know if there is a more efficient built in method/function?

I was also testing out the warm start/advance basis. Solution 2 is just marginally different from solution 1. Thus, the 2 solutions are expected to be close/similar.

<Fix binaries to solution 1>
m.solve(clean_before_solve=True)
<Fix binaries to solution 2>
m.solve(clean_before_solve=False)
<Fix binaries to solution 2>
m.solve(clean_before_solve=True)

It would be expected that warm starting would speed up performance(i.e. the middle run would be faster than the last run). But when i compared the Model.solve_details.time from each of the runs, the solve times were rather similar and for a sizable portion of the test runs, the clean start seemed to perform faster. Just wanted to check if this could be an implementation error on my part. If not, it could be something intrinsically related to my model which I shall proceed to investigate. On a related note, "Model.solution.solve_details.nb_iterations" does this refer to the number of pivots/simplex iterations performed or does this refer to a different iteration count?

PhilippeCouronne commented 3 years ago

I will try to answer your (many) questions, feel free to reply if I miss something.

  1. solve_details is actually stored in Model. th enum,ber of iterations refers either to number of MIP iterations (if problem is a mip), or number of LP iterations (if lp). I'm using LP as it may refer to different types of algorithms: number of barrier iterations if barrier was used, for example, or number of simplex iterations.

  2. MIP start/warm start: as you probably know, the impact of mip start can greatly vary , depending on the problem. To make sure the MIP start was actually used, check the CPLEX log for messages of the type:

MIP start 'm1' defined initial solution with objective 210500.0000.

  1. cloning: yes, variables belong to one model, so get_value_dict() might not be the best option here. You might compute a dictionary from indices (or names) to values, which is independent from models. To build a solution from such a dict, use Model.get_var_by_name or Model.get_var_by_index to convert indices/names to variables.

  2. performance of bounds change: there is currently no batch operation for bounds change, but this makes sense, and I will consider this for a future release. In addition, performance can most probably be improved, as time is spent checking new bound values, and this can be made optional. We will definitely improve this, but there is no simple workaround for now. This said, there are alternate ways for flipping bounds of a collection of binary variables. For example, to set a collection of binary variables to zero, add a constraint stating their sum is zero:

mdl.add(mdl.sum_vars(xs) == 0)

Or to filp them to 1, state their sum is equal to the number of variables:

mdl.add(mdl.sum_vars(xs) == len(xs))

You may then remove the constraint later to lift the constraint. Another variant is to trigger one of the above constraint with an extra binary "switch" variable, using an equivalence:

mdl.add_equivalence(switch, (mdl.sum_vars(xs) == 0))

if switch equals 1 then the constraint is satisfied, if =0, then it is not. Equivalence constraints implement two-way logical connections, if you need only one-way implication, sue Model.add_indicator instead:

mdl.add_indicator(switch, (mdl.sum_vars(xs) == 0))

if switch equals 1, the constraint is satisfied, if =0, the constraint is not guaranteed to be satisfied (but it may, by chance)

Hope this helps.

  Philippe.
bhimabi commented 3 years ago

If I am running an LP, what do you look out for in the CPLEX log to check for warm start? The log message you said is for MIPs.

I also realised an interesting phenomenon. Since I have binary values which I will be fixing, I could define them as both binary or continuous variables. However, the performance(runtime) of DOCPLEX varies due to that. My run currently consists of fixing said 'binary' variables. solving model and unfixing the variables. I have group the time taken for fixing and unfixing together. When I define the variables that get fixed as binary CPLEX takes a total of 170 seconds and fixing/unfixing takes 110 seconds. When I define the variables that get fixed as continuous CPLEX takes a total of 90 seconds and fixing/unfixing takes 190 seconds. Do you happen to know why there is a such a large difference in these processes just due to how the variables are defined? Intuitively I would have expected no differnece since it gets solved as an LP and hence there should be no big difference.

NOTE: Those timings are averages from multiple runs and datasets. Thus, I suspect there is a deeper meaning for this

PhilippeCouronne commented 3 years ago

The difference in fixing/unfixing probably comes from Docplex spending too much time in checking new bounds: for binary variables it checks that new bounds are inside [0..1], for continuous variables, checks are simpler. As I said in a previous reply, this issue is now registered and will be addressed in a forthcoming version (by which I mean performance will be improved).

If I understand well, your binary variables come out as 0 or 1, even if you redefine them as continuous? Maybe your matrix is unimodular, that happens for classical problems, e.g. network problems (https://en.wikipedia.org/wiki/Unimodular_matrix) Solving a LP problem is very different from solving a MIP, in particular, presolve is different.

To warm start a LP, you need different information from MIP starts. The right technique is to export/read a basis file (or set basis statuses by code). First, export basis statuses as a BAS file, and then read it back. See documentation of Model.export_basis and Model.read_basis_file. I am attaching a tiny example to demonstrate this.

lp_tutorial.zip

bhimabi commented 3 years ago

Response to your 2nd paragraph of your reply: I currently have an MIP model(every binary value is feasible). Instead of running it through a branch and bound algo. I am using my own algorithm to determine the value of my W variable which is binary. I fix W values (using ub and lb) according to my algorithm and then solve the model. Effectively, the model has become an LP. Since the binaries are fixed, THEORETICALLY there should be no difference in performance regardless of how the W variables are defined(binary or continuous). However, like you mentioned in your first paragraph, there are certain bounds checks by CPLEX which are affecting the performance. My test do indicate that defining them as continuous are indeed much faster. [Above section does not need a reply unless there is a different suggestion you have based on this information]

Response to final paragraph: Does this mean, the clean start you had mentioned only applies to MIP? I shall try using the export and read basis and test its performance. I shall edit this comment once I have done the tests.

PhilippeCouronne commented 3 years ago

Once you define a model with binary (or integer) variables, CPLEX qualifies it as a MIP model and uses MIP presolve and branch& bound, regardless of the bounds you set on those binaries.

MIP presolve is quite different from LP presolve, and this alone may explain (part of) performance differences. You can see problem type in Model.print_information()

As for clean_before_solve: no it applies both to MIP and LP. In particular it cleans presolve results, so that presolve is redone from scratch.

You might also be interested in class docplex.mp.relax_linear.LinearRelaxer: it computes a linear relaxation of a MIP model and returns a new model. This might help you get a linear version of a given MIP model. To exchange solutions between models, export solutions as index-based dictionaries from the relaxed model and use Model.get_var_by_index() to retrieve the corresponding MIP model variable (as both have same variables, with different types).

bhimabi commented 3 years ago

I was trying to implement the basis tests and I was getting a few errors. It did not appear everytime I ran the same code. Occassionally the error would appear in the first iteration while sometimes it appears in the middle(a few 1000 iterations later) and in other cases it would never appear. The error message I recieved varied between the 5 errors i mentioned below. I will just list the error codes below and will post the entire error at the end of this comment to keep things neat. All the errors have the occur while invoking the basis commands and have identical tracebacks.

I have a snippet of my code where this error originates from for your reference. In all the errors, the error originated from the 2nd time the basis file was input into the model after solving (marked with <== ).

def lp(chromosome, args, model, plot=False):
    # Getting values to which the 'binaries' will be fixed to
    w = convert2winn_(chromosome, args)

    # Disabling my slack variables in my model to see if it is feasible
    model.SslackPLUS.set_ub(0)
    model.SslackMINUS.set_ub(0)
    model.UslackPLUS.set_ub(0)
    model.Hslack.set_ub(0)

    # Fixing 'binaries'vars
    for slot in w:
        model.W_inn_[slot].set_lb(1)
        model.W_inn_[slot].set_ub(1)

    model.m.solve()

    # If feasible
    if model.m.solve_status.value == 2:
        model.m.read_basis_file(model.m.export_basis())

        # Recalculate same model with SLIGHTLY different objective function
        model.m.maximize(model.z)
        model.m.solve()

        # update model with basis for future runs
        model.m.read_basis_file(model.m.export_basis())          <==

    # If infeasible
    else:
        # Enable slack to set up un-slacked LP
        model.SslackPLUS.set_ub(1e+20)
        model.SslackMINUS.set_ub(1e+20)
        model.UslackPLUS.set_ub(1e+20)
        model.Hslack.set_ub(1e+20)

        model.m.solve()
        model.m.read_basis_file(model.m.export_basis())          <==
        .
        .
        .
        .

Traceback (most recent call last):
  File "SAv12_cluster.py", line 992, in <module>
    Sol_new = move(Sol, Data, Model)
  File "SAv12_cluster.py", line 307, in move
    new_c.fitness, new_c.mag_penalty, new_c.raw_fitness, new_c.feasible = lp(new_c, args, model_object)
  File "SAv12_cluster.py", line 179, in lp
    model.m.read_basis_file(model.m.export_basis())
  File "/home/asubbara/env/lib/python3.6/site-packages/docplex/mp/model.py", line 634, in read_basis_file
    cpx_read_fn=lambda cpx_, path_: cpx_.start.read_basis(path_))
  File "/home/asubbara/env/lib/python3.6/site-packages/docplex/mp/model.py", line 616, in _read_cplex_file
    cpx_read_fn(cpx, path)
  File "/home/asubbara/env/lib/python3.6/site-packages/docplex/mp/model.py", line 634, in <lambda>
    cpx_read_fn=lambda cpx_, path_: cpx_.start.read_basis(path_))
  File "/home/asubbara/env/lib/python3.6/site-packages/cplex/_internal/_subinterfaces.py", line 9107, in read_basis
    CPX_PROC.readcopybase(self._env._e, self._cplex._lp, filename)
  File "/home/asubbara/env/lib/python3.6/site-packages/cplex/_internal/_procedural.py", line 3371, in readcopybase
    check_status(env, status)
  File "/home/asubbara/env/lib/python3.6/site-packages/cplex/_internal/_procedural.py", line 236, in __call__
    raise CplexSolverError(error_string, env, status)
cplex.exceptions.errors.CplexSolverError: CPLEX Error  1434: Line 355: Couldn't convert 'XU' to a number.
PhilippeCouronne commented 3 years ago

Assuming your model is indeed a LP (check with Model.print_information()), this might be a bug. At least, it looks like the basis file is incorrect. To investigate further, we'd need a model file (preferably SAV , or LP) and the basis file (.bas extension) that creates the error. If we can reproduce with few lines, e.g.

m = ModelReader.read(sav_file)
m.solve()
bp = m.export_basis()
... change bounds
m.read_basis_file(bp)

Then we can investigate.

The documentation for basis files is here: :\https://www.ibm.com/support/knowledgecenter/SSSA5P_20.1.0/ilog.odms.cplex.help/CPLEX/FileFormats/topics/BAS.html

bhimabi commented 3 years ago

I attached the bas and sav files below. When I was running the code in debug mode, when I call the export_basis function, my model is reporting an optimal solution.

Note: This error happens when I run the model with no slacks and I receive an infeasible solution. The model then enables slacks and finds a feasible solution. Only after finding a feasible solution, do i export and read the basis file.

I also noticed that you separated your export and read lines. I combined them into one. I doubt that should be the source of the error.

LP-Castro.zip

PhilippeCouronne commented 3 years ago

I have a problem with the zip file: my corporate antivirus prevents me from opening it, and reports it as infected (!). You might retry with two separate zips, you never know. Nevertheless, I have an idea to try: there is a predicate Model.has_basis() which tells whether the solution indeed has a basis available. the export_basis does not check for this and this might be the cause of the problem. If the model has no basis, there is no way export_basis can work. Normally, export_basis should check for this, so this is very strange. Unfortunately, I cannot go further as my computer cannot open the archive.

The following code prints basis status for each variable.

    def print_basis(mdl):
        assert mdl.has_basis()

        all_vars = list(mdl.iter_variables())
        all_bs = mdl.var_basis_statuses(all_vars)

        for dv, bas in zip(all_vars, all_bs):
            print(f"- variable {dv.lp_name} has status {bas.name}")
bhimabi commented 3 years ago

While on debugigng mode, I did check if there was a basis and there was one. The model also ran and returned a feasible solution.

Here are the files agains. 2nd time is the charm lp-castro (2).zip LP-Castro.zip

PhilippeCouronne commented 3 years ago

Thanks. This time I was able to open and read both files. My best guess is ."bas" extension is maybe mis-interpreted as a BASIC file. Anyway, the good news is both files are correct and run fine both with CPLEX 12.10 and the latest version 20.1. Here's the code I executed:

from docplex.mp.model_reader import read_model
from docplex.mp.environment import Environment

def run_basis_bug(filename, basis_filename, **kwargs):
    Environment().print_information()
    mdl = read_model(filename, **kwargs)
    mdl.print_information()
    print(mdl.cplex_matrix_stats)
    mdl.read_basis_file(basis_filename)
    sol = mdl.solve(log_output=True)
    sol.display()
    return sol

And the output is:

* system is: Windows 64bit
* Python version 3.7.7, located at: C:\python\anaconda2020.02\envs\docplex_dev37\python.exe
* docplex is present, version is 2.19.0
* CPLEX library is present, version is 12.10.0.0, located at: C:\OPTIM\cplex_distrib\cplex1210R0\python\3.7\x64_win64
* pandas is present, version is 1.0.3
Model: LP-Castro
 - number of variables: 309
   - binary=0, integer=0, continuous=309
 - number of constraints: 370
   - linear=370
 - parameters: defaults
 - objective: maximize
 - problem type is: LP
Problem name         : 
Objective sense      : Maximize
Variables            :     309  [Nneg: 175,  Fix: 99,  Free: 35]
Objective nonzeros   :       1
Linear constraints   :     370  [Less: 132,  Greater: 162,  Equal: 76]
  Nonzeros           :    1089
  RHS nonzeros       :      36

Variables            : Min LB: 0.000000         Max UB: 1.000000       
Objective nonzeros   : Min   : 1.000000         Max   : 1.000000       
Linear constraints   :
  Nonzeros           : Min   : 0.0009000000     Max   : 7500.000       
  RHS nonzeros       : Min   : 11.00000         Max   : 500.0000       

Version identifier: 12.10.0.0 | 2019-09-19 | b4d7cc16e3
CPXPARAM_Read_DataCheck                          1
CPXPARAM_RandomSeed                              201909284
Using devex.
solution for: LP-Castro
objective: -1552.167
x2 = 1.000
x29 = 1.000
x33 = 1.000
x38 = 1.000
x93 = 1.000
x101 = 5.675
x102 = 11.000
x103 = 494.000
x104 = 494.000
x105 = 494.000
x106 = 492.000
x107 = 485.750
x108 = 485.750
x109 = -5.000
x110 = 1.000
x111 = 1.000
x113 = 8.000
x114 = 14.250
x116 = 5.000
x117 = 5.000
x118 = -5.000
x119 = -5.000
x120 = -5.000
x135 = 5.000
x137 = 6.000
x164 = 8.000
x168 = 6.250
x173 = 5.000
x228 = 5.000
x274 = 5.000
x283 = 5.000
x284 = 5.000
x285 = 5.000
x303 = 3.611
x305 = 250.000
x306 = 253.388
x307 = 1805.556
x308 = -1555.556
x309 = -1552.167

Second, basis file does not make any difference as , without it solving the LP takes only one iteration with dual. I am probably missing something but I am unable to reproduce the stack you sent. So far, everything looks OK. This said, I will be on vacation for two weeks; I might check support from time to time, but full support will start again beggining of January.

Hope this helps

PS: I am running the latest version of DOcplex (2.19). I suggest you also upgrade to this version.

bhimabi commented 3 years ago

Hi

I am running the latest cplex and docplex 2.19. The odd issue is that I am running this iteration 1000s of times. It only messes up after several runs. This means it works for the first n-number of runs and then it starts messing up. I also noticed that when i run this on my windows laptop this error never appears. This has been only on my linux.

As for your suggestions:

My lp run is defined in a function which takes in the model. After that iteration, I import the new basis so that the next time the lp function is called, it has a basis imported for it to start the simplex algorithm.

bhimabi commented 3 years ago

Hi,

I have been experiencing inconsistant outputs with the LP model on an unrelated issue. I created another issue regarding this to keep issues separate [link]. Perhaps this and my problem of reading basis files failing could be related.

Regarding your other comment:

Second, basis file does not make any difference as , without it solving the LP takes only one iteration with dual. I am probably missing something but I am unable to reproduce the stack you sent. So far, everything looks OK.

The datasets I am running are very tiny instances and I would be expanding these to very large datasets which would increase the size of the LP by over 1000 times. I was trying to see if using the basis file would result in any appreciable speed up in those cases. If you know it from experience/fact that there is no appreciable speed up, I shan't continuing trying to implement this.