Add feasibility checks to optimizers/converters (#1199)

* Draft for issue#1134

* resolved issue in quadratic_problem.py

* incorporated changes as suggested by @adekusar-drl

* fix linting issues

* added is_feasible method in all optimizers, added changes suggested by @ woodsp-ibm

* removed cast import

* made changes as per points 1 to 3

* made suggested changes by @adekusar-drl

* some cleaning

* added reno

* fixed reno

* Update releasenotes/notes/feasibility-check-b99605f771e745b7.yaml

Co-authored-by: Steve Wood <[email protected]>

* Update releasenotes/notes/feasibility-check-b99605f771e745b7.yaml

Co-authored-by: Steve Wood <[email protected]>

* code review

Co-authored-by: Manoel Marques <[email protected]>
Co-authored-by: Takashi Imamichi <[email protected]>
Co-authored-by: Anton Dekusar <[email protected]>
Co-authored-by: Anton Dekusar <[email protected]>
Co-authored-by: Steve Wood <[email protected]>

qiskit/optimization/algorithms/admm_optimizer.py

@@ -21,7 +21,8 @@
 from qiskit.aqua.algorithms import NumPyMinimumEigensolver
 
 from .minimum_eigen_optimizer import MinimumEigenOptimizer
-from .optimization_algorithm import OptimizationAlgorithm, OptimizationResult
+from .optimization_algorithm import (OptimizationResultStatus, OptimizationAlgorithm,
+                                     OptimizationResult)
 from .slsqp_optimizer import SlsqpOptimizer
 from ..problems.constraint import Constraint
 from ..problems.linear_constraint import LinearConstraint
@@ -180,15 +181,16 @@ class ADMMOptimizationResult(OptimizationResult):
     """ ADMMOptimization Result."""
 
     def __init__(self, x: np.ndarray, fval: float, variables: List[Variable],
-                 state: ADMMState) -> None:
+                 state: ADMMState, status: OptimizationResultStatus) -> None:
         """
         Args:
             x: the optimal value found by ADMM.
             fval: the optimal function value.
             variables: the list of variables of the optimization problem.
             state: the internal computation state of ADMM.
+            status: Termination status of an optimization algorithm
         """
-        super().__init__(x=x, fval=fval, variables=variables, raw_results=state)
+        super().__init__(x=x, fval=fval, variables=variables, status=status, raw_results=state)
 
     @property
     def state(self) -> ADMMState:
@@ -375,13 +377,15 @@ def solve(self, problem: QuadraticProgram) -> ADMMOptimizationResult:
         objective_value = objective_value * sense
 
         # convert back integer to binary
-        base_result = OptimizationResult(solution, objective_value, original_variables)
+        base_result = OptimizationResult(solution, objective_value, original_variables,
+                                         OptimizationResultStatus.SUCCESS)
         base_result = int2bin.interpret(base_result)
 
         # third parameter is our internal state of computations.
         result = ADMMOptimizationResult(x=base_result.x, fval=base_result.fval,
                                         variables=base_result.variables,
-                                        state=self._state)
+                                        state=self._state,
+                                        status=self._get_feasibility_status(problem, base_result.x))
 
         # debug
         self._log.debug("solution=%s, objective=%s at iteration=%s",

qiskit/optimization/algorithms/cplex_optimizer.py

@@ -136,9 +136,9 @@ def solve(self, problem: QuadraticProgram) -> OptimizationResult:
         sol = cplex.solution
 
         # create results
-        result = OptimizationResult(x=sol.get_values(),
-                                    fval=sol.get_objective_value(),
+        result = OptimizationResult(x=sol.get_values(), fval=sol.get_objective_value(),
                                     variables=problem.variables,
+                                    status=self._get_feasibility_status(problem, sol.get_values()),
                                     raw_results=sol)
 
         # return solution

qiskit/optimization/algorithms/grover_optimizer.py

@@ -24,9 +24,10 @@
 from qiskit.aqua.algorithms.amplitude_amplifiers.grover import Grover
 from qiskit.aqua.components.initial_states import Custom
 from qiskit.aqua.components.oracles import CustomCircuitOracle
-from qiskit.circuit.library import QuadraticForm
 from qiskit.providers import BaseBackend
-from .optimization_algorithm import OptimizationAlgorithm, OptimizationResult
+from qiskit.circuit.library import QuadraticForm
+from .optimization_algorithm import (OptimizationResultStatus, OptimizationAlgorithm,
+                                     OptimizationResult)
 from ..converters.quadratic_program_to_qubo import QuadraticProgramToQubo
 from ..problems import Variable
 from ..problems.quadratic_program import QuadraticProgram
@@ -261,15 +262,17 @@ def solve(self, problem: QuadraticProgram) -> OptimizationResult:
 
         # Compute function value
         fval = problem_init.objective.evaluate(opt_x)
-        result = OptimizationResult(x=opt_x, fval=fval, variables=problem_.variables)
+        result = OptimizationResult(x=opt_x, fval=fval, variables=problem_.variables,
+                                    status=OptimizationResultStatus.SUCCESS)
 
         # cast binaries back to integers
         result = self._qubo_converter.interpret(result)
 
         return GroverOptimizationResult(x=result.x, fval=result.fval, variables=result.variables,
                                         operation_counts=operation_count, n_input_qubits=n_key,
                                         n_output_qubits=n_value, intermediate_fval=fval,
-                                        threshold=threshold)
+                                        threshold=threshold,
+                                        status=self._get_feasibility_status(problem, result.x))
 
     def _measure(self, circuit: QuantumCircuit) -> str:
         """Get probabilities from the given backend, and picks a random outcome."""
@@ -333,7 +336,8 @@ class GroverOptimizationResult(OptimizationResult):
 
     def __init__(self, x: Union[List[float], np.ndarray], fval: float, variables: List[Variable],
                  operation_counts: Dict[int, Dict[str, int]], n_input_qubits: int,
-                 n_output_qubits: int, intermediate_fval: float, threshold: float) -> None:
+                 n_output_qubits: int, intermediate_fval: float, threshold: float,
+                 status: OptimizationResultStatus) -> None:
         """
         Constructs a result object with the specific Grover properties.
 
@@ -347,8 +351,9 @@ def __init__(self, x: Union[List[float], np.ndarray], fval: float, variables: Li
             intermediate_fval: The intermediate value of the objective function of the solution,
                 that is expected to be identical with ``fval``.
             threshold: The threshold of Grover algorithm.
+            status: the termination status of the optimization algorithm.
         """
-        super().__init__(x, fval, variables, None)
+        super().__init__(x, fval, variables, status, None)
         self._operation_counts = operation_counts
         self._n_input_qubits = n_input_qubits
         self._n_output_qubits = n_output_qubits

qiskit/optimization/algorithms/minimum_eigen_optimizer.py

@@ -17,27 +17,28 @@
 from qiskit.aqua.algorithms import MinimumEigensolver, MinimumEigensolverResult
 from qiskit.aqua.operators import StateFn, DictStateFn
 
-from .optimization_algorithm import OptimizationAlgorithm, OptimizationResult
+from .optimization_algorithm import (OptimizationResultStatus, OptimizationAlgorithm,
+                                     OptimizationResult)
 from ..converters.quadratic_program_to_qubo import QuadraticProgramToQubo
 from ..problems.quadratic_program import QuadraticProgram, Variable
 
 
 class MinimumEigenOptimizationResult(OptimizationResult):
     """ Minimum Eigen Optimizer Result."""
 
-    def __init__(self, x: Union[List[float], np.ndarray], fval: float,
-                 variables: List[Variable],
-                 samples: List[Tuple[str, float, float]],
+    def __init__(self, x: Union[List[float], np.ndarray], fval: float, variables: List[Variable],
+                 status: OptimizationResultStatus, samples: List[Tuple[str, float, float]],
                  min_eigen_solver_result: Optional[MinimumEigensolverResult] = None) -> None:
         """
         Args:
             x: the optimal value found by ``MinimumEigensolver``.
             fval: the optimal function value.
             variables: the list of variables of the optimization problem.
+            status: the termination status of the optimization algorithm.
             samples: the basis state as bitstring, the QUBO value, and the probability of sampling.
             min_eigen_solver_result: the result obtained from the underlying algorithm.
         """
-        super().__init__(x, fval, variables, None)
+        super().__init__(x, fval, variables, status, None)
         self._samples = samples
         self._min_eigen_solver_result = min_eigen_solver_result
 
@@ -178,12 +179,15 @@ def solve(self, problem: QuadraticProgram) -> MinimumEigenOptimizationResult:
             samples = [(x_str, offset, 1.0)]
 
         # translate result back to integers
-        result = OptimizationResult(x=x, fval=fval, variables=problem_.variables)
+        result = OptimizationResult(x=x, fval=fval, variables=problem_.variables,
+                                    status=OptimizationResultStatus.SUCCESS)
         result = self._qubo_converter.interpret(result)
+
         return MinimumEigenOptimizationResult(x=result.x, fval=result.fval,
                                               variables=result.variables,
-                                              samples=samples,
-                                              min_eigen_solver_result=eigen_result)
+                                              status=self._get_feasibility_status(problem,
+                                                                                  result.x),
+                                              samples=samples, min_eigen_solver_result=eigen_result)
 
 
 def _eigenvector_to_solutions(eigenvector: Union[dict, np.ndarray, StateFn],

qiskit/optimization/algorithms/multistart_optimizer.py

@@ -23,8 +23,7 @@
 from scipy.stats import uniform
 
 from qiskit.optimization import QuadraticProgram, INFINITY
-from qiskit.optimization.algorithms.optimization_algorithm import (OptimizationAlgorithm,
-                                                                   OptimizationResult)
+from .optimization_algorithm import OptimizationAlgorithm, OptimizationResult
 
 logger = logging.getLogger(__name__)
 
@@ -93,6 +92,7 @@ def multi_start_solve(self, minimize: Callable[[np.array], Tuple[np.array, Any]]
                 rest_sol = rest
 
         return OptimizationResult(x=x_sol, fval=fval_sol, variables=problem.variables,
+                                  status=self._get_feasibility_status(problem, x_sol),
                                   raw_results=rest_sol)
 
     @property

qiskit/optimization/algorithms/optimization_algorithm.py

@@ -22,6 +22,19 @@
 from ..problems.quadratic_program import QuadraticProgram, Variable
 
 
+class OptimizationResultStatus(Enum):
+    """Termination status of an optimization algorithm."""
+
+    SUCCESS = 0
+    """the optimization algorithm succeeded to find an optimal solution."""
+
+    FAILURE = 1
+    """the optimization algorithm ended in a failure."""
+
+    INFEASIBLE = 2
+    """the optimization algorithm obtained an infeasible solution."""
+
+
 class OptimizationAlgorithm(ABC):
     """An abstract class for optimization algorithms in Qiskit's optimization module."""
 
@@ -84,18 +97,21 @@ def _verify_compatibility(self, problem: QuadraticProgram) -> None:
         if msg:
             raise QiskitOptimizationError('Incompatible problem: {}'.format(msg))
 
+    def _get_feasibility_status(self, problem: QuadraticProgram,
+                                x: Union[List[float], np.ndarray]) -> OptimizationResultStatus:
+        """Returns whether the input result is feasible or not for the given problem.
 
-class OptimizationResultStatus(Enum):
-    """Termination status of an optimization algorithm."""
-
-    SUCCESS = 0
-    """the optimization algorithm succeeded to find an optimal solution."""
+        Args:
+            problem: Problem to verify.
+            x: the input result list.
 
-    FAILURE = 1
-    """the optimization algorithm ended in a failure."""
+        Returns:
+            The status of the result.
+        """
+        is_feasible = problem.is_feasible(x)
 
-    INFEASIBLE = 2
-    """the optimization algorithm obtained an infeasible solution."""
+        return OptimizationResultStatus.SUCCESS if is_feasible \
+            else OptimizationResultStatus.INFEASIBLE
 
 
 class OptimizationResult:
@@ -140,8 +156,8 @@ class OptimizationResult:
 
     def __init__(self, x: Union[List[float], np.ndarray], fval: float,
                  variables: List[Variable],
-                 raw_results: Optional[Any] = None,
-                 status: OptimizationResultStatus = OptimizationResultStatus.SUCCESS) -> None:
+                 status: OptimizationResultStatus,
+                 raw_results: Optional[Any] = None) -> None:
         """
         Args:
             x: the optimal value found in the optimization.

qiskit/optimization/algorithms/recursive_minimum_eigen_optimizer.py

@@ -21,7 +21,8 @@
 from qiskit.aqua.algorithms import NumPyMinimumEigensolver
 from qiskit.aqua.utils.validation import validate_min
 
-from .optimization_algorithm import OptimizationAlgorithm, OptimizationResult
+from .optimization_algorithm import (OptimizationResultStatus, OptimizationAlgorithm,
+                                     OptimizationResult)
 from .minimum_eigen_optimizer import MinimumEigenOptimizer, MinimumEigenOptimizationResult
 from ..converters.quadratic_program_to_qubo import QuadraticProgramToQubo
 from ..exceptions import QiskitOptimizationError
@@ -50,9 +51,8 @@ class IntermediateResult(Enum):
 
 class RecursiveMinimumEigenOptimizationResult(OptimizationResult):
     """Recursive Eigen Optimizer Result."""
-    def __init__(self, x: Union[List[float], np.ndarray], fval: float,
-                 variables: List[Variable],
-                 replacements: Dict[str, Tuple[str, int]],
+    def __init__(self, x: Union[List[float], np.ndarray], fval: float, variables: List[Variable],
+                 status: OptimizationResultStatus, replacements: Dict[str, Tuple[str, int]],
                  history: Tuple[List[MinimumEigenOptimizationResult], OptimizationResult]) -> None:
         """
         Constructs an instance of the result class.
@@ -61,6 +61,7 @@ def __init__(self, x: Union[List[float], np.ndarray], fval: float,
             x: the optimal value found in the optimization.
             fval: the optimal function value.
             variables: the list of variables of the optimization problem.
+            status: the termination status of the optimization algorithm.
             replacements: a dictionary of substituted variables. Key is a variable being
                 substituted, value is a tuple of substituting variable and a weight, either 1 or -1.
             history: a tuple containing intermediate results. The first element is a list of
@@ -70,7 +71,7 @@ def __init__(self, x: Union[List[float], np.ndarray], fval: float,
                 :class:`~qiskit.optimization.algorithm.OptimizationResult` obtained at the last step
                 via `min_num_vars_optimizer`.
         """
-        super().__init__(x, fval, variables, None)
+        super().__init__(x, fval, variables, status, None)
         self._replacements = replacements
         self._history = history
 
@@ -281,13 +282,16 @@ def find_value(x, replacements, var_values):
         # construct result
         x_v = [var_values[x_aux.name] for x_aux in problem_ref.variables]
         fval = result.fval
-        result = OptimizationResult(x=x_v, fval=fval, variables=problem_ref.variables)
+        result = OptimizationResult(x=x_v, fval=fval, variables=problem_ref.variables,
+                                    status=OptimizationResultStatus.SUCCESS)
         result = self._qubo_converter.interpret(result)
 
         return RecursiveMinimumEigenOptimizationResult(x=result.x, fval=result.fval,
                                                        variables=result.variables,
                                                        replacements=replacements,
-                                                       history=history)
+                                                       history=history,
+                                                       status=(self._get_feasibility_status
+                                                               (problem, result.x)))
 
     def _find_strongest_correlation(self, correlations):
 

qiskit/optimization/algorithms/slsqp_optimizer.py

@@ -18,7 +18,7 @@
 from scipy.optimize import fmin_slsqp
 
 from .multistart_optimizer import MultiStartOptimizer
-from .optimization_algorithm import OptimizationResult
+from .optimization_algorithm import OptimizationResultStatus, OptimizationResult
 from ..exceptions import QiskitOptimizationError
 from ..problems import Variable
 from ..problems.constraint import Constraint
@@ -32,8 +32,9 @@ class SlsqpOptimizationResult(OptimizationResult):
     SLSQP optimization result, defines additional properties that may be returned by the optimizer.
     """
     def __init__(self, x: Union[List[float], np.ndarray], fval: float, variables: List[Variable],
-                 fx: Optional[np.ndarray] = None, its: Optional[int] = None,
-                 imode: Optional[int] = None, smode: Optional[str] = None) -> None:
+                 status: OptimizationResultStatus, fx: Optional[np.ndarray] = None,
+                 its: Optional[int] = None, imode: Optional[int] = None,
+                 smode: Optional[str] = None) -> None:
         """
         Constructs a result object with properties specific to SLSQP.
 
@@ -46,8 +47,9 @@ def __init__(self, x: Union[List[float], np.ndarray], fval: float, variables: Li
             imode: The exit mode from the optimizer
                 (see the documentation of ``scipy.optimize.fmin_slsqp``).
             smode: Message describing the exit mode from the optimizer.
+            status: the termination status of the optimization algorithm.
         """
-        super().__init__(x, fval, variables, None)
+        super().__init__(x, fval, variables, status, None)
         self._fx = fx
         self._its = its
         self._imode = imode
@@ -217,8 +219,10 @@ def _minimize(x_0: np.array) -> Tuple[np.array, Any]:
         result = self.multi_start_solve(_minimize, problem)
 
         if self._full_output:
-            return SlsqpOptimizationResult(result.x, result.fval, result.variables,
+            return SlsqpOptimizationResult(x=result.x, fval=result.fval, variables=result.variables,
+                                           status=self._get_feasibility_status(problem, result.x),
                                            fx=result.raw_results[0], its=result.raw_results[1],
                                            imode=result.raw_results[2], smode=result.raw_results[3])
         else:
-            return SlsqpOptimizationResult(result.x, result.fval, result.variables)
+            return SlsqpOptimizationResult(x=result.x, fval=result.fval, variables=result.variables,
+                                           status=self._get_feasibility_status(problem, result.x))

qiskit/optimization/converters/inequality_to_equality.py

@@ -352,7 +352,7 @@ def interpret(self, result: OptimizationResult) -> OptimizationResult:
         names = [x.name for x in self._dst.variables]
         new_x = self._interpret_var(names, result.x)
         return OptimizationResult(x=new_x, fval=result.fval, variables=self._src.variables,
-                                  raw_results=result.raw_results, status=result.status)
+                                  status=result.status, raw_results=result.raw_results)
 
     def _interpret_var(self, names, vals) -> List[int]:
         # interpret slack variables

qiskit/optimization/converters/integer_to_binary.py

@@ -219,7 +219,7 @@ def interpret(self, result: OptimizationResult) -> OptimizationResult:
         """
         new_x = self._interpret_var(result.x)
         return OptimizationResult(x=new_x, fval=result.fval, variables=self._src.variables,
-                                  raw_results=result.raw_results)
+                                  status=result.status, raw_results=result.raw_results)
 
     def _interpret_var(self, vals: Union[List[float], np.ndarray]) -> List[float]:
         # interpret integer values

qiskit/optimization/converters/linear_equality_to_penalty.py

@@ -190,8 +190,8 @@ def interpret(self, result: OptimizationResult) -> OptimizationResult:
             new_status = OptimizationResultStatus.INFEASIBLE
 
         return OptimizationResult(x=result.x, fval=substituted_qp.objective.constant,
-                                  variables=self._src.variables, raw_results=result.raw_results,
-                                  status=new_status)
+                                  variables=self._src.variables, status=new_status,
+                                  raw_results=result.raw_results)
 
     @property
     def penalty(self) -> Optional[float]: