add reduced cost in python interface

Author: ZedongPeng

python/README.md

@@ -74,6 +74,7 @@ print("Status:", m.Status)
 print("Objective:", m.ObjVal)
 print("Primal solution:", m.X)
 print("Dual solution:", m.Pi)
+print("Reduced cost:", m.RC)
 ```
 
 ## Modeling
@@ -190,6 +191,7 @@ After calling `m.optimize()`, the solver stores results in a set of read-only at
 | `RelGap` | float | Relative primal-dual gap. |
 | `X` | numpy.ndarray | Primal solution vector \(x\). May be `None` if no feasible solution was found. |
 | `Pi` | numpy.ndarray | Dual solution vector (Lagrange multipliers). |
+| `RC` | numpy.ndarray | Reduced Cost vector. |
 | `IterCount` | int | Number of iterations performed. |
 | `Runtime` | float | Total wall-clock runtime in seconds. |
 | `RescalingTime` | float | Time spent on preprocessing and rescaling (seconds). |
@@ -214,6 +216,7 @@ print("Iterations:", m.IterCount, " Runtime (s):", m.Runtime)
 # Access solutions
 print("Primal solution:", m.X)
 print("Dual solution:", m.Pi)
+print("Reduced cost:", m.RC)
 
 # Check residuals
 print("Primal residual:", m.RelPrimalResidual)

python/cupdlpx/model.py

@@ -128,6 +128,7 @@ def __init__(
         # initialize solution attributes
         self._x: Optional[np.ndarray] = None # primal solution
         self._y: Optional[np.ndarray] = None # dual solution
+        self._rc: Optional[np.ndarray] = None # reduced costs
         self._objval: Optional[float] = None # objective value
         self._dualobj: Optional[float] = None # dual objective value
         self._gap: Optional[float] = None # primal-dual gap
@@ -376,6 +377,7 @@ def optimize(self):
         # solutions
         self._x = np.asarray(info.get("X")) if info.get("X") is not None else None
         self._y = np.asarray(info.get("Pi")) if info.get("Pi") is not None else None
+        self._rc = np.asarray(info.get("RC")) if info.get("RC") is not None else None
         # objectives & gaps
         primal_obj_eff = info.get("PrimalObj")
         dual_obj_eff   = info.get("DualObj")
@@ -404,7 +406,7 @@ def _clear_solution_cache(self) -> None:
         """
         Clear cached solution attributes.
         """
-        self._x = self._y = None
+        self._x = self._y = self._rc = None
         self._objval = self._dualobj = None
         self._gap = self._rel_gap = None
         self._status = None
@@ -423,6 +425,10 @@ def X(self) -> Optional[np.ndarray]:
     @property
     def Pi(self) -> Optional[np.ndarray]:
         return self._y
+    
+    @property
+    def RC(self) -> Optional[np.ndarray]:
+        return self._rc
 
     @property
     def ObjVal(self) -> Optional[float]:

python_bindings/_core_bindings.cpp

@@ -518,16 +518,19 @@ static py::dict solve_once(
     const int m_out = res->num_constraints;
     py::array_t<double> x({n_out});
     py::array_t<double> y({m_out});
+    py::array_t<double> rc({n_out});
     {
-        auto xb = x.request(), yb = y.request();
+        auto xb = x.request(), yb = y.request(), rcb = rc.request();
         std::memcpy(xb.ptr, res->primal_solution, sizeof(double) * n_out);
         std::memcpy(yb.ptr, res->dual_solution, sizeof(double) * m_out);
+        std::memcpy(rcb.ptr, res->reduced_costs, sizeof(double) * n_out);
     }
     // build info dict
     py::dict info;
     // solution
     info["X"] = x;
     info["Pi"] = y;
+    info["RC"] = rc;
     // objectives and gaps
     info["PrimalObj"] = res->primal_objective_value;
     info["DualObj"] = res->dual_objective_value;