How to avoid violated constraints with pyqubo?
Dear all,
I have written up a problem using the pyqubo package v1.0.7 , and the docs at https://pyqubo.readthedocs.io/en/latest/getting_started.html#solve-qubo-by-dimod-sampler
Starting from the MIP model and then for each constraint f(X) = C, with X some model vars and C some scalar constant and f some function linear each of the X variables, I have added P * (f(X)-C) ** 2 to the QUBO objective function, for example via the call:
Ham += P * Constraint(( x[j] - x[i] + W * (1-p[i,j] + q[i,j]) - s3x[(i,j)] - (w[i]+w[j]) / 2.0 ) ** 2,
label='eq3_{:d}_{:d}'.format(i,j))
However, I get 51 violated constraints via:
decoded_sample.constraints(only_broken=True)
How can I fix the constraints? The penalties P, I have all put at values ranging from P=0.75 to 1.5 to 500 times the max MIP objective value, but that does not help me to get rid of the violated constraints.
This happens both on the sampler = neal.SimulatedAnnealingSampler() as well as if I use a real dwave machine via sampler = DWaveCliqueSampler(). By the way, with DWaveCliqueSampler() I can only get 119 qubits. Is that normal? I thought > 5000 qubits were available?...
Do I do something wrong? Does this simply mean that my problem is infeasible? Maybe a good practice is to
check the problem first with a MIP solver like CBC/CPLEX/XPRESS/Gurobi...
thanks for the help and best regards,
Peter
Details: for example I get:
{'eq10L_6': (False, 0.25),
'eq10L_7': (False, 1.0),
'eq10R_0': (False, 4.0),
'eq10R_2': (False, 1.0),
'eq10R_3': (False, 1.0),
'eq10R_4': (False, 16.0),
...
- I increased the number of samples but apart from taking longer and somewhat less violated constraints sometimes, I never reached zero violated constraints, which should be possible, right?
- Also without the Constraint() around the expression, I get invalid solutions.
Comments
Hi Peter,
Welcome to the Leap community, and thank you for taking the time to post your question here! I have a few questions that might help untangle what's going wrong.
1. What size/density of QUBO are you hoping to solve?
I noticed that you are using the DWaveCliqueSampler to solve your QUBO. This sampler is designed for very dense problems, i.e. where each node is connected to all the others. The D-Wave Advantage QPU has about 5000 qubits as you noted, but they are not fully-connected. We call the process of mapping your problem topology to the QPU topology "embedding" - more on this here: https://docs.dwavesys.com/docs/latest/c_gs_7.html.
The graph size of 119 you found is the largest fully-connected network when mapped to the Advantage QPU. For sparser networks, i.e. not all variables interact with each other, I would recommend a different sampler, for example the EmbeddingComposite. If your problem is larger than the QPU, the LeapHybridSampler would be a good option - it combines quantum and classical resources to solve very large and complex problems. Here is some documentation on available solvers in Ocean: https://docs.ocean.dwavesys.com/projects/system/en/stable/reference/samplers.html
2. How many samples (num_reads) are you using?
One thing to note is that the QPU is a probabilistic device, so it will quickly find good solutions, but statistically it is unlikely to find a "perfect" solution. Increasing the number of samples increases your chance of finding a solution that meets your constraints. More on this here: https://support.dwavesys.com/hc/en-us/articles/360003681974
As you mentioned, adjusting energy scales (penalty weights) can help prioritize meeting constraints over other objectives. However, if there are too many overlapping or conflicting constraints, this might become infeasible. The process of embedding mentioned above can muddy the waters; some of the energy scale fidelity is used to produce coherent results (avoiding "chain breaks", which are explained in the embedding article above).
3. Would you be able to post some code snippets of the full constraints you are using, the resulting QUBO, or maybe a GitHub link?
This might help myself or other Leap users to get a better sense of the size and complexity of your problem. As mentioned above, this could inform what the best approach to solving your QUBO might be. If possible, it is always helpful to try solving a scaled-down version of your problem to better understand all the nuances.
All the best,
Luca
Dear Luca,
and Dear All / Community,
Luca, thanks very much for your elaborate answer. :)
As a response on your questions Luca:
1.I have much bigger problems than I have set (in the cfg table at the top). So I would have loved to go up to using 5000 qubits. (by increasing 'max_number' values in that cfg python dictionary).
2. I have run 10 to 100 samples.
3. I will share my full running code here, and so maybe somebody can get the returned samples to satisfy all (MIP-original) constraints, which is my main objective, of course. :)
It is an (actually running) notebook, but I see no way to attach file here, so I have just put it online at:
https://github.com/PeterSels/scratch/blob/master/CovidProofPositioningUseYourOwnDevKey.ipynb
You can see it directly as a notebook there (includeing the output I got) , and I think you can download it from there as well and then run it yourself.
It should be directly runnable, but you have to just substitute the question marks in
with your own ocean leap DEV-id from dwave.
So I get some invalid constraints (as you can see in the end in out[16]). 15 of them in my run (but this varies from run to run). I am very curious if someone knows how to adapt it so that it returns valid solutions, not violating any constraints. :)
thanks and best regards,
Peter
Note that
has to be changed to
if you want to execute on a real online dwave quantum annealing machine i.o. on a simulator on your own personal computer.
best,
Peter
Hi Peter,
Thanks for sharing your work on Github. I'm not familiar with the paper that was referenced, but it's an interesting project. It's quite helpful to see the comparison with other approaches/algorithms presented in one piece.
I think the size of problem you are working with is a good starting point. Hopefully with some adjustments you will be able to get good results and scale up from there.
Following the questions you answered above (and having browsed your jupyter notebook) I would suggest starting by (1) swapping out DWaveCliqueSampler for a more general-purpose solver, such as EmbeddingComposite and (2) increasing your number of samples, for example up to 1,000 reads just to be sure you are maximizing your chances of finding a good solution.
Cheers,
Luca
Hi Luca,
Thanks for your answer. It was helpful and I could solve some coding things.
I updated the notebook on the github link at
https://github.com/PeterSels/scratch/blob/master/CovidProofPositioningUseYourOwnDevKey.ipynb
just now. It nows also runs a larger example.
But I have some problems left.
(1) The problem (already in the 'mip_qubo' mode) of balancing the penalties so that all objective terms P * (f(x)-c)^2 should become zero. This is a maths thing, or a mip thing already and not yet to do with quantum computing.
(2) In mode m= 'ann_hw', EmbeddingComposite seems to hang
(3) LeapHybridSampler returns but don't know how to specify num_reads (and it returns overlapping rectangles, which is not surprising, given that (a) I don't take the lowest energy from more than 1 read sample."
Could you advise on any of the above?
thanks and best regards,
Peter
Hi Peter,
Our apologies for the delayed response.
I am glad to hear that you are now able to submit larger problems using LeapHybridSampler.
Here are my responses to the problems you still have:
(1) I might suggest that you try to scale down your problem and simplify the constraints to verify if the model works.
(2) For a very large problem, EmbeddingComposite takes a long time to find an embedding. If the problem is too big to fit on the QPU, it will generate an error: "no Embedding found". In either case, it provides an output. Did you notice a different behavior? It will be helpful if you can separate out the call and isolate to observe it.
(3) Unfortunately, LeapHybridSampler does not have a tunable num_reads parameter. It only returns the lowest energy solution it can find within a time frame. The only modifiable parameter this solver has is the time_limit.
Please let us know if you have any further questions.
Please sign in to leave a comment.