Repository: entropicalabs/openqaoa (language: python)

Issue:
Transform a docplex optimization problem into openQAOA Hamiltonian
---
name: [Docplex](https://pypi.org/project/docplex/) converter.
about: Having an optimization problem in Docplex form transforms it into the OpenQAOA Hamiltonian.
title: 'docplex2hamiltonian'
labels: 'feature request'
assignees: @alejomonbar @MaldoAlberto 

---

### Prerequisites
Before raising this feature request, I have made sure:

- [x] It does not lead to duplication.
- [x] It has not already been requested.

### Is your feature request related to a problem? Please describe.

Implementing optimization problems using openQAOA is restricted to certain problems with a predefined structure ( NumberPartition, TSP, Knapsack, etc). For example, if we want to modify the Knapsack problem with a Heuristic equation, the implemented Knapsack problem would not help me to solve my problem and I must create my own Hamiltonian of this problem. Docplex a python based library used to solve optimization problems can be used to codify the optimization problem and then translate it directly into its OpenQAOA Hamiltonian.

It will be something like:

``` python
from docplex.mp.model import Model
from openqaoa.utilities import X_mixer_hamiltonian
from openqaoa.problems.helper_functions import docplex2hamiltonian #The function we want to implement

mdl = Model(name="MyOptimizationProblem") # Docplex model

z = mdl.binary_var_list(10, name="Z") #My binary optimization variables


objective = mdl.sum(z) - 2 * z[0]  + z[3] * z[4]  # Cost function

mdl.minimize(objective) # Adding the objective function to the model 

cost_hamil = docplex2hamiltonian(mdl) # The function we want to implement

mixer_hamil = X_mixer_hamiltonian(num_assets)

circuit_params = QAOACircuitParams(cost_hamil, mixer_hamil, p=1)


variate_params = create_qaoa_variational_params(circuit_params, 'standard', 'ramp')
backend_obj_qiskit = get_qaoa_backend(circuit_params, 'qiskit_shot_simulator', n_shots=10000)
backend_obj_qiskit.expectation_w_uncertainty(variate_params)
backend_obj_qiskit.expectation(variate_params)
optimizer_dict = {'method': 'cobyla', 'maxiter': 10}
variate_params_qiskit = variate_params
optimizer_obj_qiskit = get_optimizer(backend_obj_qiskit, variate_params_qiskit, optimizer_dict)
optimizer_obj_qiskit()
optimizer_obj_qiskit.results_information()
```

### Describe the solution you'd like
Create a class function that translates an optimization problem into its Hamiltonian representation.

### Describe alternatives you've considered
we think Numba can be used instead of Docplex. But, we should start first with Docplex.

### Additional context
<img width="917" alt="Screen Shot 2022-06-17 at 8 07 02 AM" src="https://user-images.githubusercontent.com/11268357/174235813-1bb688aa-facd-4997-a9ff-014bea12e794.png">
<img width="918" alt="Screen Shot 2022-06-17 at 8 06 52 AM" src="https://user-images.githubusercontent.com/11268357/174235821-c31d105f-671c-4a2d-948e-c4e6b48f441c.png">


<img width="686" alt="Screen Shot 2022-06-17 at 7 30 50 AM" src="https://user-images.githubusercontent.com/11268357/174232005-4887f154-f52d-47e9-88c8-76b28cec3942.png">


Dev from docplex2 qubo
This is the new method mentioned in #36.
At this point, we have implemented the conversion of the cost function from a docplex model to the PUBO representation. We will be updating the inequality constraint and the equality constraints soon.


``` python
from docplex.mp.model import Model

from openqaoa.problems.helper_functions import docplex2qubo

mdl = Model("Test") # Docplex model
num_z = 5 # Number of variables
z = mdl.binary_var_list(num_z, name="Z")


objective = mdl.sum(z) - 2 * z[0]  + z[3] * z[4] + 5

mdl.minimize(objective)

sol = mdl.solve()

mdl.print_solution(print_zeros=True)

qubo_problem = docplex2qubo(mdl)

isingOQ = qubo_problem.qubo_to_ising()

from openqaoa.workflows.optimizer import QAOA  
q = QAOA()
q.compile(qubo_problem)
q.optimize() 
print(q.results_information["best probability"])

```
Transform a docplex optimization problem into openQAOA Hamiltonian
---
name: [Docplex](https://pypi.org/project/docplex/) converter.
about: Having an optimization problem in Docplex form transforms it into the OpenQAOA Hamiltonian.
title: 'docplex2hamiltonian'
labels: 'feature request'
assignees: @alejomonbar @MaldoAlberto 

---

### Prerequisites
Before raising this feature request, I have made sure:

- [x] It does not lead to duplication.
- [x] It has not already been requested.

### Is your feature request related to a problem? Please describe.

Implementing optimization problems using openQAOA is restricted to certain problems with a predefined structure ( NumberPartition, TSP, Knapsack, etc). For example, if we want to modify the Knapsack problem with a Heuristic equation, the implemented Knapsack problem would not help me to solve my problem and I must create my own Hamiltonian of this problem. Docplex a python based library used to solve optimization problems can be used to codify the optimization problem and then translate it directly into its OpenQAOA Hamiltonian.

It will be something like:

``` python
from docplex.mp.model import Model
from openqaoa.utilities import X_mixer_hamiltonian
from openqaoa.problems.helper_functions import docplex2hamiltonian #The function we want to implement

mdl = Model(name="MyOptimizationProblem") # Docplex model

z = mdl.binary_var_list(10, name="Z") #My binary optimization variables


objective = mdl.sum(z) - 2 * z[0]  + z[3] * z[4]  # Cost function

mdl.minimize(objective) # Adding the objective function to the model 

cost_hamil = docplex2hamiltonian(mdl) # The function we want to implement

mixer_hamil = X_mixer_hamiltonian(num_assets)

circuit_params = QAOACircuitParams(cost_hamil, mixer_hamil, p=1)


variate_params = create_qaoa_variational_params(circuit_params, 'standard', 'ramp')
backend_obj_qiskit = get_qaoa_backend(circuit_params, 'qiskit_shot_simulator', n_shots=10000)
backend_obj_qiskit.expectation_w_uncertainty(variate_params)
backend_obj_qiskit.expectation(variate_params)
optimizer_dict = {'method': 'cobyla', 'maxiter': 10}
variate_params_qiskit = variate_params
optimizer_obj_qiskit = get_optimizer(backend_obj_qiskit, variate_params_qiskit, optimizer_dict)
optimizer_obj_qiskit()
optimizer_obj_qiskit.results_information()
```

### Describe the solution you'd like
Create a class function that translates an optimization problem into its Hamiltonian representation.

### Describe alternatives you've considered
we think Numba can be used instead of Docplex. But, we should start first with Docplex.

### Additional context
<img width="917" alt="Screen Shot 2022-06-17 at 8 07 02 AM" src="https://user-images.githubusercontent.com/11268357/174235813-1bb688aa-facd-4997-a9ff-014bea12e794.png">
<img width="918" alt="Screen Shot 2022-06-17 at 8 06 52 AM" src="https://user-images.githubusercontent.com/11268357/174235821-c31d105f-671c-4a2d-948e-c4e6b48f441c.png">


<img width="686" alt="Screen Shot 2022-06-17 at 7 30 50 AM" src="https://user-images.githubusercontent.com/11268357/174232005-4887f154-f52d-47e9-88c8-76b28cec3942.png">


Transform a docplex optimization problem into openQAOA Hamiltonian
---
name: [Docplex](https://pypi.org/project/docplex/) converter.
about: Having an optimization problem in Docplex form transforms it into the OpenQAOA Hamiltonian.
title: 'docplex2hamiltonian'
labels: 'feature request'
assignees: @alejomonbar @MaldoAlberto 

---

### Prerequisites
Before raising this feature request, I have made sure:

- [x] It does not lead to duplication.
- [x] It has not already been requested.

### Is your feature request related to a problem? Please describe.

Implementing optimization problems using openQAOA is restricted to certain problems with a predefined structure ( NumberPartition, TSP, Knapsack, etc). For example, if we want to modify the Knapsack problem with a Heuristic equation, the implemented Knapsack problem would not help me to solve my problem and I must create my own Hamiltonian of this problem. Docplex a python based library used to solve optimization problems can be used to codify the optimization problem and then translate it directly into its OpenQAOA Hamiltonian.

It will be something like:

``` python
from docplex.mp.model import Model
from openqaoa.utilities import X_mixer_hamiltonian
from openqaoa.problems.helper_functions import docplex2hamiltonian #The function we want to implement

mdl = Model(name="MyOptimizationProblem") # Docplex model

z = mdl.binary_var_list(10, name="Z") #My binary optimization variables


objective = mdl.sum(z) - 2 * z[0]  + z[3] * z[4]  # Cost function

mdl.minimize(objective) # Adding the objective function to the model 

cost_hamil = docplex2hamiltonian(mdl) # The function we want to implement

mixer_hamil = X_mixer_hamiltonian(num_assets)

circuit_params = QAOACircuitParams(cost_hamil, mixer_hamil, p=1)


variate_params = create_qaoa_variational_params(circuit_params, 'standard', 'ramp')
backend_obj_qiskit = get_qaoa_backend(circuit_params, 'qiskit_shot_simulator', n_shots=10000)
backend_obj_qiskit.expectation_w_uncertainty(variate_params)
backend_obj_qiskit.expectation(variate_params)
optimizer_dict = {'method': 'cobyla', 'maxiter': 10}
variate_params_qiskit = variate_params
optimizer_obj_qiskit = get_optimizer(backend_obj_qiskit, variate_params_qiskit, optimizer_dict)
optimizer_obj_qiskit()
optimizer_obj_qiskit.results_information()
```

### Describe the solution you'd like
Create a class function that translates an optimization problem into its Hamiltonian representation.

### Describe alternatives you've considered
we think Numba can be used instead of Docplex. But, we should start first with Docplex.

### Additional context
<img width="917" alt="Screen Shot 2022-06-17 at 8 07 02 AM" src="https://user-images.githubusercontent.com/11268357/174235813-1bb688aa-facd-4997-a9ff-014bea12e794.png">
<img width="918" alt="Screen Shot 2022-06-17 at 8 06 52 AM" src="https://user-images.githubusercontent.com/11268357/174235821-c31d105f-671c-4a2d-948e-c4e6b48f441c.png">


<img width="686" alt="Screen Shot 2022-06-17 at 7 30 50 AM" src="https://user-images.githubusercontent.com/11268357/174232005-4887f154-f52d-47e9-88c8-76b28cec3942.png">




Discussion / hints:
Thanks, @alejomonbar and @MaldoAlberto !

I like the idea. This way OpenQAOA will also have a simple wrapper to convert qiskit QUBO problems into hamiltonians (cc @vishal-ph ). 

Do you think `docplex2hamiltonian` can be written without needing to install docplex? This way we could make the installation of docplex optional, and not compulsory: if you want to use this feature you'd go ahead and install docplex. But otherwise, the package stays leaner.

Ideally, we `docplex2hamiltonian` would take some python primitives (i.e., dictionaries and strings, and not an object coming from docplex). Do you think this will be possible?

Also, you call `docplex2hamiltonian` a class function but you plan to put it in the utilities. What class are you proposing here?  
@Q-lds and @vishal-ph 

- Do you think docplex2hamiltonian can be written without needing to install docplex? 
   Sure, we will make it optional.

- Ideally, we docplex2hamiltonian would take some python primitives (i.e., dictionaries and strings, and not an object coming from docplex).

I think it needs a docplex document as input. Otherwise, we lose the whole purpose of the method because it requires the flexibility of adjusting the cost function, which is only possible if we work outside the function docplex2hamiltonian with the docplex method. 

- Also, you call docplex2hamiltonian a class function but you plan to put it in the utilities. What class are you proposing here?
Well, I'm not sure at this point what is best. But, let us work on the problem and come up with some options. Because the implementation requires other functions that translate equality and inequality constraints into its cost function representation.
> I think it needs a docplex document as input. Otherwise, we lose the whole purpose of the method because it requires the flexibility of adjusting the cost function, which is only possible if we work outside the function docplex2hamiltonian with the docplex method.

Well, what I was thinking is whether you can easily convert the model into python primitives by hand so that the parameters of `docplex2hamiltonian()` are all python primitives. But ok, let's start with the parameter being the docplex model as you are suggesting.

> Well, I'm not sure at this point what is best. But, let us work on the problem and come up with some options. Because the implementation requires other functions that translate equality and inequality constraints into its cost function representation.

Ok. Intuitively I think this could be an abstract method of the QUBO class. But feel free to propose an implementation and we'll take it from there!

I am really looking forward to seeing how this feature will come along! 


Thanks, @alejomonbar and @MaldoAlberto !

I like the idea. This way OpenQAOA will also have a simple wrapper to convert qiskit QUBO problems into hamiltonians (cc @vishal-ph ). 

Do you think `docplex2hamiltonian` can be written without needing to install docplex? This way we could make the installation of docplex optional, and not compulsory: if you want to use this feature you'd go ahead and install docplex. But otherwise, the package stays leaner.

Ideally, we `docplex2hamiltonian` would take some python primitives (i.e., dictionaries and strings, and not an object coming from docplex). Do you think this will be possible?

Also, you call `docplex2hamiltonian` a class function but you plan to put it in the utilities. What class are you proposing here?  
@Q-lds and @vishal-ph 

- Do you think docplex2hamiltonian can be written without needing to install docplex? 
   Sure, we will make it optional.

- Ideally, we docplex2hamiltonian would take some python primitives (i.e., dictionaries and strings, and not an object coming from docplex).

I think it needs a docplex document as input. Otherwise, we lose the whole purpose of the method because it requires the flexibility of adjusting the cost function, which is only possible if we work outside the function docplex2hamiltonian with the docplex method. 

- Also, you call docplex2hamiltonian a class function but you plan to put it in the utilities. What class are you proposing here?
Well, I'm not sure at this point what is best. But, let us work on the problem and come up with some options. Because the implementation requires other functions that translate equality and inequality constraints into its cost function representation.
> I think it needs a docplex document as input. Otherwise, we lose the whole purpose of the method because it requires the flexibility of adjusting the cost function, which is only possible if we work outside the function docplex2hamiltonian with the docplex method.

Well, what I was thinking is whether you can easily convert the model into python primitives by hand so that the parameters of `docplex2hamiltonian()` are all python primitives. But ok, let's start with the parameter being the docplex model as you are suggesting.

> Well, I'm not sure at this point what is best. But, let us work on the problem and come up with some options. Because the implementation requires other functions that translate equality and inequality constraints into its cost function representation.

Ok. Intuitively I think this could be an abstract method of the QUBO class. But feel free to propose an implementation and we'll take it from there!

I am really looking forward to seeing how this feature will come along! 

Thanks, @alejomonbar and @MaldoAlberto !

I like the idea. This way OpenQAOA will also have a simple wrapper to convert qiskit QUBO problems into hamiltonians (cc @vishal-ph ). 

Do you think `docplex2hamiltonian` can be written without needing to install docplex? This way we could make the installation of docplex optional, and not compulsory: if you want to use this feature you'd go ahead and install docplex. But otherwise, the package stays leaner.

Ideally, we `docplex2hamiltonian` would take some python primitives (i.e., dictionaries and strings, and not an object coming from docplex). Do you think this will be possible?

Also, you call `docplex2hamiltonian` a class function but you plan to put it in the utilities. What class are you proposing here?  
@Q-lds and @vishal-ph 

- Do you think docplex2hamiltonian can be written without needing to install docplex? 
   Sure, we will make it optional.

- Ideally, we docplex2hamiltonian would take some python primitives (i.e., dictionaries and strings, and not an object coming from docplex).

I think it needs a docplex document as input. Otherwise, we lose the whole purpose of the method because it requires the flexibility of adjusting the cost function, which is only possible if we work outside the function docplex2hamiltonian with the docplex method. 

- Also, you call docplex2hamiltonian a class function but you plan to put it in the utilities. What class are you proposing here?
Well, I'm not sure at this point what is best. But, let us work on the problem and come up with some options. Because the implementation requires other functions that translate equality and inequality constraints into its cost function representation.
> I think it needs a docplex document as input. Otherwise, we lose the whole purpose of the method because it requires the flexibility of adjusting the cost function, which is only possible if we work outside the function docplex2hamiltonian with the docplex method.

Well, what I was thinking is whether you can easily convert the model into python primitives by hand so that the parameters of `docplex2hamiltonian()` are all python primitives. But ok, let's start with the parameter being the docplex model as you are suggesting.

> Well, I'm not sure at this point what is best. But, let us work on the problem and come up with some options. Because the implementation requires other functions that translate equality and inequality constraints into its cost function representation.

Ok. Intuitively I think this could be an abstract method of the QUBO class. But feel free to propose an implementation and we'll take it from there!

I am really looking forward to seeing how this feature will come along! 


Fix patch:
diff --git a/openqaoa/problems/converters.py b/openqaoa/problems/converters.py
new file mode 100644
index 000000000..b48bf2c45
--- /dev/null
+++ b/openqaoa/problems/converters.py
@@ -0,0 +1,387 @@
+#   Copyright 2022 Entropica Labs
+#
+#   Licensed under the Apache License, Version 2.0 (the "License");
+#   you may not use this file except in compliance with the License.
+#   You may obtain a copy of the License at
+#
+#       http://www.apache.org/licenses/LICENSE-2.0
+#
+#   Unless required by applicable law or agreed to in writing, software
+#   distributed under the License is distributed on an "AS IS" BASIS,
+#   WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
+#   See the License for the specific language governing permissions and
+#   limitations under the License.
+
+import numpy as np
+from collections import defaultdict
+from .problem import QUBO
+
+
+class FromDocplex2IsingModel:
+    def __init__(self, model, multipliers: float = None):
+
+        """
+        Creates an instance to translate Docplex models to its Ising Model representation
+
+        Parameters
+        ----------
+        model : Docplex model
+            It is an object that has the mathematical expressions (Cost function,
+            Inequality and Equality constraints) of an optimization problem in
+            binary representation.
+
+        """
+        # assign the docplex Model
+        self.model = model.copy()
+
+        # save the index in a dict
+        self.idx_terms = {}
+        for x in self.model.iter_variables():
+            self.idx_terms[x] = x.index
+            if x.vartype.short_name != "binary":
+                TypeError(f"Variable {x.vartype.short_name} is not allowed.")
+        # get doclex qubo and ising model
+        self.qubo_docplex, self.ising_model = self.get_models(multipliers)
+
+    def linear_expr(self, expr):
+        """
+        Adds linear expresions to the objective function stored in self.qubo_dict.
+
+        Parameters
+        ----------
+        expr : model.objective_expr
+            A docplex model attribute.
+
+        """
+        for x, weight in expr.get_linear_part().iter_terms():
+            self.qubo_dict[(self.idx_terms[x],)] += weight
+
+    def quadratic_expr(self, expr):
+        """
+        Adds quadratic terms to the objective function stored in self.qubo_dict.
+
+        Parameters
+        ----------
+        expr : model.objective_expr
+            A docplex model attribute.
+        """
+
+        # save the terms and coeffs of the quadratic part,
+        # considering two index i and j
+        for x, y, weight in expr.iter_quad_triplets():
+            i = self.idx_terms[x]
+            j = self.idx_terms[y]
+            if i == j:  # if this is the term is  Z1**2 for example
+                self.qubo_dict[(i,)] += weight
+            else:
+                self.qubo_dict[(i, j)] += weight
+    
+    @staticmethod
+    def equality_to_penalty(expression, multiplier: float):
+        """
+        Add equality constraints to the cost function using the penality representation.
+        The constraints should be linear.
+        
+        Parameters
+        ----------
+        expression : docplex.mp.linear.LinearExpr
+            docplex equality constraint "expression == 0".
+        multiplier : float
+            Lagrange multiplier of the penalty
+
+        Returns
+        -------
+        penalty : docplex.mp.quad.QuadExpr
+            Penalty that will be added to the cost function.
+
+        """
+        penalty = multiplier * (expression) ** 2
+        return penalty
+    
+    @staticmethod
+    def bounds(linear_expression):
+        """
+        Generates the limits of a linear term
+
+        Parameters
+        ----------
+        linear_exp : generator
+            Iterable term with the quadratic terms of the cost function.
+
+        Returns
+        -------
+        l_bound : float
+            Lower bound of the quadratic term.
+        u_bound : float
+            Upper bound of the quadratic term
+        """
+
+        l_bound = (
+            u_bound
+        ) = linear_expression.constant  # Lower and upper bound of the constraint
+        for term, coeff in linear_expression.iter_terms():
+            l_bound += min(0, coeff)
+            u_bound += max(0, coeff)
+
+        return l_bound, u_bound
+    
+    @staticmethod
+    def quadratic_bounds(iter_exp):
+        """
+        Generates the limits of the quadratic terms
+
+        Parameters
+        ----------
+        iter_exp : generator
+            Iterable term with the quadratic terms of the cost function.
+
+        Returns
+        -------
+        l_bound : float
+            Lower bound of the quadratic term.
+        u_bound : float
+            Upper bound of the quadratic term
+
+        """
+        l_bound = 0
+        u_bound = 0
+        for z1, z2, coeff in iter_exp:
+            l_bound += min(0, coeff)
+            u_bound += max(0, coeff)
+        return l_bound, u_bound
+
+    def inequality_to_equality(self, constraint):
+        """
+        Transform inequality contraints into equality constriants using 
+        slack variables.
+
+        Parameters
+        ----------
+        constraint : docplex.mp.constr.LinearConstraint
+            DESCRIPTION.
+
+        Returns
+        -------
+        new_exp : docplex.mp.linear.LinearExpr
+            The equality constraint representation of the inequality constraint.
+
+        """
+
+        if constraint.sense_string == "LE":  # Less or equal inequality constraint
+            new_exp = constraint.get_right_expr() + -1 * constraint.get_left_expr()
+        elif constraint.sense_string == "GE":  # Great or equal inequality constriant
+            new_exp = constraint.get_left_expr() + -1 * constraint.get_right_expr()
+        else:
+            AttributeError(f"It is not possible to implement constraint {constraint.sense_string}.")
+
+        lower_bound, upper_bound = self.bounds(new_exp)
+        slack_lim = upper_bound  # Slack var limit
+
+        if slack_lim > 0:
+            sign = -1
+        else:
+            sign = 1
+        # Create the slack variables to represent valid terms in the inequality
+        n_slack = int(np.ceil(np.log2(abs(slack_lim + 1))))
+        if n_slack > 0:
+            slack_vars = self.model.binary_var_list(
+                n_slack, name=f"slack_{constraint.name}"
+            )
+            for x in slack_vars:
+                self.idx_terms[x] = x.index
+
+            for nn, var in enumerate(slack_vars[:-1]):
+                new_exp += sign * (2 ** nn) * var
+
+            new_exp += (
+                sign * (slack_lim - 2 ** (n_slack - 1) + 1) * slack_vars[-1]
+            )  # restrict the last term to fit the upper bound
+
+        return new_exp
+
+    def multipliers_generators(self):
+        """
+        Penality term size adapter, this is the Lagrange multiplier of the cost
+        function penalties for every constraint if the multiplier is not indicated
+        by the user. 
+
+        Returns
+        -------
+        float
+            the multiplier resized by the cost function limits.
+
+        """
+        cost_func = self.model.objective_expr
+        l_bound_linear, u_bound_linear = self.bounds(cost_func.get_linear_part())
+        l_bound_quad, u_bound_quad = self.quadratic_bounds(
+            cost_func.iter_quad_triplets()
+        )
+        return 1.0 + (u_bound_linear - l_bound_linear) + (u_bound_quad - l_bound_quad)
+
+    def linear_constraints(self, multipliers=None) -> None:
+        """
+        Adds the constraints of the problem to the objective function. 
+
+        Parameters
+        ----------
+        multiplier : List
+            For each constraint a multiplier, if None it automatically is selected.
+
+        Returns
+        -------
+        None.
+
+        """
+
+        constraints_list = list(self.model.iter_linear_constraints())
+        n_constraints = len(constraints_list)
+
+        if (
+            multipliers == None
+        ):  # Default penalties are choosen from the bounds of the objective func.
+            multipliers = n_constraints * [self.multipliers_generators()]
+
+        elif type(multipliers) in [int, float]:
+            multipliers = n_constraints * [multipliers]
+
+        elif type(multipliers) != list:
+            print(f"{type(multipliers)} is not a accepted format")
+
+        for cn, constraint in enumerate(constraints_list):
+            if (
+                constraint.sense_string == "EQ"
+            ):  # Equality constraint added as a penalty.
+                left_exp = constraint.get_left_expr()
+                right_exp = constraint.get_right_expr()
+                expression = left_exp + -1 * right_exp
+                penalty = self.equality_to_penalty(expression, multipliers[cn])
+            elif constraint.sense_string in [
+                "LE",
+                "GE",
+            ]:  # Inequality constraint added as a penalty with additional slack variables.
+                constraint.name = f"C{cn}"
+                ineq2eq = self.inequality_to_equality(constraint)
+                penalty = self.equality_to_penalty(ineq2eq, multipliers[cn])
+            else:
+                print("Not an accepted constraint!")
+
+            self.linear_expr(penalty)
+            self.quadratic_expr(penalty)
+            self.constant += penalty.constant
+            self.objective_qubo += penalty
+    
+    @staticmethod
+    def qubo_to_ising(n_variables, qubo_terms, qubo_weights):
+        """
+        Converts the terms and weights in QUBO representation ([0,1])
+        to the Ising representation ([-1, 1]). 
+
+        Parameters
+        ----------
+        n_variables : int
+            number of variables.
+        qubo_terms : List
+            List of QUBO variables.
+        qubo_weights : List
+            coefficients of the variables
+            
+
+        Returns
+        -------
+        Ising Model stored on QUBO class
+
+        """
+        ising_terms, ising_weights = [], []
+        linear_terms = np.zeros(n_variables)
+
+        constant_term = 0
+        # Process the given terms and weights
+        for weight, term in zip(qubo_weights, qubo_terms):
+
+            if len(term) == 2:
+                u, v = term
+
+                if u != v:
+                    ising_terms.append([u, v])
+                    ising_weights.append(weight / 4)
+                else:
+                    constant_term += weight / 4
+
+                linear_terms[term[0]] -= weight / 4
+                linear_terms[term[1]] -= weight / 4
+                constant_term += weight / 4
+            elif len(term) == 1:
+                linear_terms[term[0]] -= weight / 2
+                constant_term += weight / 2
+            elif len(term) == 0:
+                constant_term += weight
+            else:
+                TypeError(f"Term {term} is not recognized!")
+
+        for variable, linear_term in enumerate(linear_terms):
+            ising_terms.append([variable])
+            ising_weights.append(linear_term)
+
+        ising_terms.append([])
+        ising_weights.append(constant_term)
+
+        return QUBO(n_variables, ising_terms, ising_weights)
+
+    def get_models(self, multipliers: float = None):
+        """
+        Creates a QUBO docplex model, QUBO dict OQ model, and an Ising Model form
+        a Docplex quadratic program.
+
+        Parameters
+        ----------
+        model : Docplex model
+            It is an object that has the mathematical expressions (Cost function,
+            Inequality and Equality constraints) of an optimization problem.
+
+        Returns
+        -------
+        qubo_docplex, ising_model 
+                 
+
+        """
+        # save a dictionary with the qubo information
+        self.qubo_dict = defaultdict(float)
+
+        # objective sense
+        if self.model.objective_sense.is_minimize():
+            # If it is minimized
+            self.objective_expr = self.model.objective_expr
+        else:
+            # If it is maximized, multiplied by -1 to minimize.
+            self.objective_expr = -1 * self.model.objective_expr
+
+        # Objective QUBO
+        self.objective_qubo = self.model.objective_expr
+        # Obtain the constant from the model
+        self.constant = self.objective_expr.constant
+
+        # Save the terms and coeffs form the linear part
+        self.linear_expr(self.objective_expr)
+
+        # Save the terms and coeffs form the quadratic part
+        self.quadratic_expr(self.objective_expr)
+
+        # Add the linear constraints into the qubo
+        self.linear_constraints(multipliers=multipliers)
+
+        terms = list(self.qubo_dict.keys()) + [
+            []
+        ]  # The right term is for adding the constant part of the QUBO
+
+        weights = list(self.qubo_dict.values()) + [self.constant]
+
+        n_variables = self.model.number_of_variables
+        # QUBO docplex
+        qubo_docplex = self.model.copy()
+        qubo_docplex.remove_objective()
+        qubo_docplex.minimize(self.objective_qubo)
+
+        # convert the docplex terms in an Ising Model
+        ising_model = self.qubo_to_ising(n_variables, terms, weights)
+
+        return qubo_docplex, ising_model
diff --git a/setup.py b/setup.py
index 2ead4e2df..cd2698b69 100644
--- a/setup.py
+++ b/setup.py
@@ -20,7 +20,8 @@
     "nbsphinx==0.8.9",
     "ipython==8.2.0" ,
     "nbconvert==6.5.0",
-    "pytest-cov==3.0.0"
+    "pytest-cov==3.0.0",
+    "docplex==2.23.222"
     ]
 
 setup(


Test patch:
diff --git a/tests/test_converters.py b/tests/test_converters.py
new file mode 100644
index 000000000..b0882249b
--- /dev/null
+++ b/tests/test_converters.py
@@ -0,0 +1,80 @@
+#   Copyright 2022 Entropica Labs
+#
+#   Licensed under the Apache License, Version 2.0 (the "License");
+#   you may not use this file except in compliance with the License.
+#   You may obtain a copy of the License at
+#
+#       http://www.apache.org/licenses/LICENSE-2.0
+#
+#   Unless required by applicable law or agreed to in writing, software
+#   distributed under the License is distributed on an "AS IS" BASIS,
+#   WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
+#   See the License for the specific language governing permissions and
+#   limitations under the License.
+
+import unittest
+
+from openqaoa.problems.problem import QUBO, MaximumCut
+from openqaoa.problems.converters import FromDocplex2IsingModel
+from docplex.mp.model import Model
+import networkx as nx
+
+
+class TestDocplex2IsingClass(unittest.TestCase):
+    
+    """
+    Test the converter from docplex models
+    
+    """
+
+    def test_qubo(self): 
+        """
+        Test the QUBO class is generated from the function FromDocplex2IsingModel
+
+        """
+        # Creating a basic docplex model
+        mdl = Model("Test") # Docplex model
+        num_z = 5 # Number of variables
+        z = mdl.binary_var_list(num_z, name="Z") # docplex variables
+        objective = mdl.sum(z) - 2 * z[0]  + z[3] * z[4] + 5 # objective function
+
+        mdl.minimize(objective) # Optimization
+
+        ising_problem = FromDocplex2IsingModel(mdl).ising_model # Ising model of the Docplex model
+        
+        self.assertIsInstance(ising_problem, QUBO)
+        
+    def test_model_maxcut(self):
+        """
+        Test the Maxcut application of OpenQAOA gives the same result as the 
+        Docplex translation model.
+        """
+        
+        # Graph representing the maxcut problem
+        n_nodes = 5 # number of nodes
+        G = nx.Graph()
+        G.add_nodes_from(range(n_nodes))
+        G.add_edges_from([[0,1], [0,2], [1,2], [4,3], [3,2], [4,0], [2,4]])
+        
+        # Docplex model 
+        mdl = Model(name="Max-cut")
+        x = mdl.binary_var_list(n_nodes, name="x")
+        for w, v in G.edges:
+            G.edges[w, v].setdefault("weight", 1)
+        objective = mdl.sum(G.edges[i, j]["weight"] * x[i] * (1 - x[j]) + G.edges[i, j]["weight"] * x[j] * (1 - x[i]) for i, j in G.edges)
+        mdl.maximize(objective)
+        
+        # Translating the problem to OQ ising Model
+        ModelOQ = FromDocplex2IsingModel(mdl)
+        Ising_model_OQ = ModelOQ.ising_model.asdict()
+        # Using the predefine function of this problem
+        IsingModelDirect = MaximumCut(G).get_qubo_problem().asdict()
+        # Comparing both results in this specific case MaxCut from OpenQAOa gives 
+        # the coefficients omultiplied by two of the DocplexToIsingModel
+        for nn, term in enumerate(IsingModelDirect["terms"]):
+            idx_OQ = Ising_model_OQ["terms"].index(term)
+            self.assertAlmostEqual(2 * Ising_model_OQ["weights"][idx_OQ], IsingModelDirect["weights"][nn])
+
+
+if __name__ == '__main__':
+    unittest.main()