Skip to content

Commit 9c5e2fa

Browse files
Release Managervbraun
Release Manager
authored and
vbraun
committed
Trac #20506: dual variables handling in SDP solver(s)
At present only primal variables are made available in `SemidefiniteProgram()`, although the dual variables are typically also provided by backends. This ticket added access to dual (and slack) variables. URL: http://trac.sagemath.org/20506 Reported by: dimpase Ticket author(s): Dima Pasechnik Reviewer(s): Matthias Koeppe
2 parents 92dc886 + 572b630 commit 9c5e2fa

File tree

6 files changed

+397
-45
lines changed

6 files changed

+397
-45
lines changed

‎src/sage/numerical/backends/cvxopt_backend.pyx

+1-1
Original file line numberDiff line numberDiff line change
@@ -90,7 +90,7 @@ cdef class CVXOPTBackend(GenericBackend):
9090
"reltol":1e-6,
9191
"feastol":1e-7,
9292
"refinement":0 }
93-
93+
self.answer = {}
9494
if maximization:
9595
self.set_sense(+1)
9696
else:

‎src/sage/numerical/backends/cvxopt_sdp_backend.pyx

+129-5
Original file line numberDiff line numberDiff line change
@@ -61,7 +61,7 @@ cdef class CVXOPTSDPBackend(GenericSDPBackend):
6161
"reltol":1e-6,
6262
"feastol":1e-7,
6363
"refinement":0 }
64-
64+
self.answer = {}
6565
if maximization:
6666
self.set_sense(+1)
6767
else:
@@ -255,8 +255,8 @@ cdef class CVXOPTSDPBackend(GenericSDPBackend):
255255
INPUT:
256256
257257
- ``coefficients`` an iterable with ``(c,v)`` pairs where ``c``
258-
is a variable index (integer) and ``v`` is a value (real
259-
value). The pairs come sorted by indices. If c is -1 it
258+
is a variable index (integer) and ``v`` is a value (matrix).
259+
The pairs come sorted by indices. If c is -1 it
260260
represents the constant coefficient.
261261
262262
- ``name`` - an optional name for this row (default: ``None``)
@@ -435,8 +435,6 @@ cdef class CVXOPTSDPBackend(GenericSDPBackend):
435435
436436
EXAMPLE::
437437
438-
sage: from sage.numerical.backends.generic_sdp_backend import get_solver
439-
sage: p = get_solver(solver = "cvxopt")
440438
sage: p = SemidefiniteProgram(solver = "cvxopt", maximization=False)
441439
sage: x = p.new_variable()
442440
sage: p.set_objective(x[0] - x[1] + x[2])
@@ -462,6 +460,32 @@ cdef class CVXOPTSDPBackend(GenericSDPBackend):
462460
i+=1
463461
return sum
464462

463+
cpdef _get_answer(self):
464+
"""
465+
return the complete output dict of the solver
466+
467+
Mainly for testing purposes
468+
469+
TESTS::
470+
471+
sage: p = SemidefiniteProgram(maximization = False, solver='cvxopt')
472+
sage: x = p.new_variable()
473+
sage: p.set_objective(x[0] - x[1])
474+
sage: a1 = matrix([[1, 2.], [2., 3.]])
475+
sage: a2 = matrix([[3, 4.], [4., 5.]])
476+
sage: a3 = matrix([[5, 6.], [6., 7.]])
477+
sage: b1 = matrix([[1, 1.], [1., 1.]])
478+
sage: b2 = matrix([[2, 2.], [2., 2.]])
479+
sage: b3 = matrix([[3, 3.], [3., 3.]])
480+
sage: p.add_constraint(a1*x[0] + a2*x[1] <= a3)
481+
sage: p.add_constraint(b1*x[0] + b2*x[1] <= b3)
482+
sage: p.solve(); # tol 1e-08
483+
-3.0
484+
sage: p.get_backend()._get_answer()
485+
{...}
486+
"""
487+
return self.answer
488+
465489
cpdef get_variable_value(self, int variable):
466490
"""
467491
Return the value of a variable given by the solver.
@@ -617,8 +641,108 @@ cdef class CVXOPTSDPBackend(GenericSDPBackend):
617641
matrices.append(m)
618642
return (indices, matrices)
619643

644+
cpdef dual_variable(self, int i, sparse=False):
645+
"""
646+
The `i`-th dual variable
647+
648+
Available after self.solve() is called, otherwise the result is undefined
649+
650+
- ``index`` (integer) -- the constraint's id.
651+
652+
OUTPUT:
653+
654+
The matrix of the `i`-th dual variable
655+
656+
EXAMPLE::
657+
658+
sage: p = SemidefiniteProgram(maximization = False, solver='cvxopt')
659+
sage: x = p.new_variable()
660+
sage: p.set_objective(x[0] - x[1])
661+
sage: a1 = matrix([[1, 2.], [2., 3.]])
662+
sage: a2 = matrix([[3, 4.], [4., 5.]])
663+
sage: a3 = matrix([[5, 6.], [6., 7.]])
664+
sage: b1 = matrix([[1, 1.], [1., 1.]])
665+
sage: b2 = matrix([[2, 2.], [2., 2.]])
666+
sage: b3 = matrix([[3, 3.], [3., 3.]])
667+
sage: p.add_constraint(a1*x[0] + a2*x[1] <= a3)
668+
sage: p.add_constraint(b1*x[0] + b2*x[1] <= b3)
669+
sage: p.solve() # tol 1e-08
670+
-3.0
671+
sage: B=p.get_backend()
672+
sage: x=p.get_values(x).values()
673+
sage: -(a3*B.dual_variable(0)).trace()-(b3*B.dual_variable(1)).trace() # tol 1e-07
674+
-3.0
675+
sage: g = sum((B.slack(j)*B.dual_variable(j)).trace() for j in range(2)); g # tol 1.5e-08
676+
0.0
677+
678+
679+
TESTS::
680+
681+
sage: B.dual_variable(7)
682+
...
683+
Traceback (most recent call last):
684+
...
685+
IndexError: list index out of range
686+
sage: abs(g - B._get_answer()['gap']) # tol 1e-22
687+
0.0
688+
689+
"""
690+
cdef int n
691+
n = self.answer['zs'][i].size[0]
692+
assert(n == self.answer['zs'][i].size[1]) # must be square matrix
693+
return Matrix(n, n, list(self.answer['zs'][i]), sparse=sparse)
694+
695+
cpdef slack(self, int i, sparse=False):
696+
"""
697+
Slack of the `i`-th constraint
698+
699+
Available after self.solve() is called, otherwise the result is undefined
700+
701+
- ``index`` (integer) -- the constraint's id.
702+
703+
OUTPUT:
704+
705+
The matrix of the slack of the `i`-th constraint
706+
707+
EXAMPLE::
708+
709+
sage: p = SemidefiniteProgram(maximization = False, solver='cvxopt')
710+
sage: x = p.new_variable()
711+
sage: p.set_objective(x[0] - x[1])
712+
sage: a1 = matrix([[1, 2.], [2., 3.]])
713+
sage: a2 = matrix([[3, 4.], [4., 5.]])
714+
sage: a3 = matrix([[5, 6.], [6., 7.]])
715+
sage: b1 = matrix([[1, 1.], [1., 1.]])
716+
sage: b2 = matrix([[2, 2.], [2., 2.]])
717+
sage: b3 = matrix([[3, 3.], [3., 3.]])
718+
sage: p.add_constraint(a1*x[0] + a2*x[1] <= a3)
719+
sage: p.add_constraint(b1*x[0] + b2*x[1] <= b3)
720+
sage: p.solve() # tol 1e-08
721+
-3.0
722+
sage: B = p.get_backend()
723+
sage: B1 = B.slack(1); B1 # tol 1e-08
724+
[0.0 0.0]
725+
[0.0 0.0]
726+
sage: B1.is_positive_definite()
727+
True
728+
sage: x = p.get_values(x).values()
729+
sage: x[0]*b1 + x[1]*b2 - b3 + B1 # tol 1e-09
730+
[0.0 0.0]
731+
[0.0 0.0]
732+
733+
TESTS::
620734
735+
sage: B.slack(7)
736+
...
737+
Traceback (most recent call last):
738+
...
739+
IndexError: list index out of range
621740
741+
"""
742+
cdef int n
743+
n = self.answer['ss'][i].size[0]
744+
assert(n == self.answer['ss'][i].size[1]) # must be square matrix
745+
return Matrix(n, n, list(self.answer['ss'][i]), sparse=sparse)
622746

623747
cpdef row_name(self, int index):
624748
"""

‎src/sage/numerical/backends/generic_sdp_backend.pxd

+2
Original file line numberDiff line numberDiff line change
@@ -18,6 +18,8 @@ cdef class GenericSDPBackend:
1818
cpdef int solve(self) except -1
1919
cpdef get_objective_value(self)
2020
cpdef get_variable_value(self, int variable)
21+
cpdef dual_variable(self, int variable, sparse=*)
22+
cpdef slack(self, int variable, sparse=*)
2123
cpdef bint is_maximization(self)
2224
cpdef row(self, int i)
2325
cpdef int ncols(self)

‎src/sage/numerical/backends/generic_sdp_backend.pyx

+93
Original file line numberDiff line numberDiff line change
@@ -482,6 +482,99 @@ cdef class GenericSDPBackend:
482482
"""
483483
raise NotImplementedError()
484484

485+
cpdef dual_variable(self, int i, sparse=False):
486+
"""
487+
The `i`-th dual variable
488+
489+
Available after self.solve() is called, otherwise the result is undefined
490+
491+
- ``index`` (integer) -- the constraint's id.
492+
493+
OUTPUT:
494+
495+
The matrix of the `i`-th dual variable
496+
497+
EXAMPLE::
498+
499+
sage: p = SemidefiniteProgram(maximization = False,solver = "Nonexistent_LP_solver") # optional - Nonexistent_LP_solver
500+
sage: x = p.new_variable() # optional - Nonexistent_LP_solver
501+
sage: p.set_objective(x[0] - x[1]) # optional - Nonexistent_LP_solver
502+
sage: a1 = matrix([[1, 2.], [2., 3.]]) # optional - Nonexistent_LP_solver
503+
sage: a2 = matrix([[3, 4.], [4., 5.]]) # optional - Nonexistent_LP_solver
504+
sage: a3 = matrix([[5, 6.], [6., 7.]]) # optional - Nonexistent_LP_solver
505+
sage: b1 = matrix([[1, 1.], [1., 1.]]) # optional - Nonexistent_LP_solver
506+
sage: b2 = matrix([[2, 2.], [2., 2.]]) # optional - Nonexistent_LP_solver
507+
sage: b3 = matrix([[3, 3.], [3., 3.]]) # optional - Nonexistent_LP_solver
508+
sage: p.add_constraint(a1*x[0] + a2*x[1] <= a3) # optional - Nonexistent_LP_solver
509+
sage: p.add_constraint(b1*x[0] + b2*x[1] <= b3) # optional - Nonexistent_LP_solver
510+
sage: p.solve() # optional - Nonexistent_LP_solver # tol ???
511+
-3.0
512+
sage: B=p.get_backend() # optional - Nonexistent_LP_solver
513+
sage: x=p.get_values(x).values() # optional - Nonexistent_LP_solver
514+
sage: -(a3*B.dual_variable(0)).trace()-(b3*B.dual_variable(1)).trace() # optional - Nonexistent_LP_solver # tol ???
515+
-3.0
516+
sage: g = sum((B.slack(j)*B.dual_variable(j)).trace() for j in range(2)); g # optional - Nonexistent_LP_solver # tol ???
517+
0.0
518+
519+
TESTS::
520+
521+
sage: B.dual_variable(7) # optional - Nonexistent_LP_solver
522+
...
523+
Traceback (most recent call last):
524+
...
525+
IndexError: list index out of range
526+
sage: abs(g - B._get_answer()['gap']) # optional - Nonexistent_LP_solver # tol 1e-22
527+
0.0
528+
"""
529+
raise NotImplementedError()
530+
531+
cpdef slack(self, int i, sparse=False):
532+
"""
533+
Slack of the `i`-th constraint
534+
535+
Available after self.solve() is called, otherwise the result is undefined
536+
537+
- ``index`` (integer) -- the constraint's id.
538+
539+
OUTPUT:
540+
541+
The matrix of the slack of the `i`-th constraint
542+
543+
EXAMPLE::
544+
545+
sage: p = SemidefiniteProgram(maximization = False,solver = "Nonexistent_LP_solver") # optional - Nonexistent_LP_solver
546+
sage: x = p.new_variable() # optional - Nonexistent_LP_solver
547+
sage: p.set_objective(x[0] - x[1]) # optional - Nonexistent_LP_solver
548+
sage: a1 = matrix([[1, 2.], [2., 3.]]) # optional - Nonexistent_LP_solver
549+
sage: a2 = matrix([[3, 4.], [4., 5.]]) # optional - Nonexistent_LP_solver
550+
sage: a3 = matrix([[5, 6.], [6., 7.]]) # optional - Nonexistent_LP_solver
551+
sage: b1 = matrix([[1, 1.], [1., 1.]]) # optional - Nonexistent_LP_solver
552+
sage: b2 = matrix([[2, 2.], [2., 2.]]) # optional - Nonexistent_LP_solver
553+
sage: b3 = matrix([[3, 3.], [3., 3.]]) # optional - Nonexistent_LP_solver
554+
sage: p.add_constraint(a1*x[0] + a2*x[1] <= a3) # optional - Nonexistent_LP_solver
555+
sage: p.add_constraint(b1*x[0] + b2*x[1] <= b3) # optional - Nonexistent_LP_solver
556+
sage: p.solve() # optional - Nonexistent_LP_solver # tol ???
557+
-3.0
558+
sage: B=p.get_backend() # optional - Nonexistent_LP_solver
559+
sage: B1 = B.slack(1); B1 # optional - Nonexistent_LP_solver # tol ???
560+
[0.0 0.0]
561+
[0.0 0.0]
562+
sage: B1.is_positive_definite() # optional - Nonexistent_LP_solver
563+
True
564+
sage: x = p.get_values(x).values() # optional - Nonexistent_LP_solver
565+
sage: x[0]*b1 + x[1]*b2 - b3 + B1 # optional - Nonexistent_LP_solver # tol ???
566+
[0.0 0.0]
567+
[0.0 0.0]
568+
569+
TESTS::
570+
571+
sage: B.slack(7) # optional - Nonexistent_LP_solver
572+
...
573+
Traceback (most recent call last):
574+
...
575+
IndexError: list index out of range
576+
"""
577+
raise NotImplementedError()
485578

486579
cpdef solver_parameter(self, name, value = None):
487580
"""

‎src/sage/numerical/sdp.pxd

+2
Original file line numberDiff line numberDiff line change
@@ -17,6 +17,8 @@ cdef class SemidefiniteProgram(SageObject):
1717
cpdef int number_of_variables(self)
1818
cdef list _constraints
1919
cpdef sum(self, L)
20+
cpdef dual_variable(self, int i, sparse=*)
21+
cpdef slack(self, int i, sparse=*)
2022

2123

2224
cdef class SDPVariable(Element):

0 commit comments

Comments
 (0)