-
Notifications
You must be signed in to change notification settings - Fork 1.2k
/
sos_examples.h
154 lines (129 loc) · 4.33 KB
/
sos_examples.h
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
#pragma once
#include "drake/solvers/mathematical_program.h"
#include "drake/solvers/mathematical_program_result.h"
namespace drake {
namespace solvers {
/**
* This is example 3.35 from Semidefinite Optimization and Convex Algebraic
* Geometry by G. Blekherman, P. Parrilo and R. Thomas. Solve a semidefinite
* programming problem to verify that the univariate quartic polynomial p(x)
* = x⁴+4x³+6x²+4x+5 is sum-of-squares (sos).
*/
class UnivariateQuarticSos {
public:
DRAKE_NO_COPY_NO_MOVE_NO_ASSIGN(UnivariateQuarticSos)
UnivariateQuarticSos();
const MathematicalProgram& prog() const { return prog_; }
void CheckResult(const MathematicalProgramResult& result, double tol) const;
private:
MathematicalProgram prog_;
symbolic::Polynomial p_;
MatrixXDecisionVariable gram_;
VectorX<symbolic::Monomial> monomial_basis_;
};
/**
* This is example 3.38 from Semidefinite Optimization and Convex Algebraic
* Geometry by G. Blekherman, P. Parrilo and R. Thomas. Solve a semidefinite
* programming problem to verify that the bivariate quartic polynomial p(x, y) =
* 2x⁴+5y⁴−2x²y²+2x³y+2x+2 is sum-of-squares (sos).
*/
class BivariateQuarticSos {
public:
DRAKE_NO_COPY_NO_MOVE_NO_ASSIGN(BivariateQuarticSos)
BivariateQuarticSos();
const MathematicalProgram& prog() const { return prog_; }
void CheckResult(const MathematicalProgramResult& result, double tol) const;
private:
MathematicalProgram prog_;
symbolic::Polynomial p_;
MatrixXDecisionVariable gram_;
VectorX<symbolic::Monomial> monomial_basis_;
};
/**
* This is example 3.50 from Semidefinite Optimization and Convex Algebraic
* Geometry by G. Blekherman, P. Parrilo and R. Thomas. Solve a simple sos
* program
* max a + b
* s.t x⁴ + ax + 2+b is sos
* (a-b+1)x² + bx + 1 is sos
*/
class SimpleSos1 {
public:
DRAKE_NO_COPY_NO_MOVE_NO_ASSIGN(SimpleSos1)
SimpleSos1();
const MathematicalProgram& prog() const { return prog_; }
void CheckResult(const MathematicalProgramResult& result, double tol) const;
private:
MathematicalProgram prog_;
symbolic::Variable a_;
symbolic::Variable b_;
symbolic::Variable x_;
symbolic::Polynomial p1_;
symbolic::Polynomial p2_;
MatrixXDecisionVariable gram1_;
MatrixXDecisionVariable gram2_;
VectorX<symbolic::Monomial> monomial_basis1_;
VectorX<symbolic::Monomial> monomial_basis2_;
};
/**
* Prove that the Motzkin polynomial m(x, y) = x⁴y² + x²y⁴ + 1 − 3x²y² is
* always non-negative.
* One ceritificate for the proof is the existence of a polynomial r(x, y)
* satisfying r(x, y) being sos and r(x, y) > 0 for all x, y ≠ 0, such that
* r(x, y) * m(x, y) is sos.
* So we solve the following problem
* find r(x, y)
* s.t r(x, y) - x² − y² is sos
* r(x, y) * m(x, y) is sos
*/
class MotzkinPolynomial {
public:
DRAKE_NO_COPY_NO_MOVE_NO_ASSIGN(MotzkinPolynomial)
MotzkinPolynomial();
const MathematicalProgram& prog() const { return prog_; }
void CheckResult(const MathematicalProgramResult& result, double tol) const;
private:
MathematicalProgram prog_;
symbolic::Variable x_;
symbolic::Variable y_;
symbolic::Polynomial m_;
symbolic::Polynomial r_;
MatrixXDecisionVariable gram1_;
MatrixXDecisionVariable gram2_;
VectorX<symbolic::Monomial> monomial_basis1_;
VectorX<symbolic::Monomial> monomial_basis2_;
};
/**
* Solve the following optimization problem for a univariate polynomial
* max a + b + c
* s.t p(x) = x⁴+ax³+bx²+c+1>=0 for all x >= 0
* p(1) = 1
* According to theorem 3.71 in Semidefinite Optimization and Convex Algebraic
* Geometry by G. Blekherman, P. Parrilo and R. Thomas, this is equivalent to
* the following SOS problem
* max a + b + c
* s.t p(x) = s(x) + x * t(x)
* p(1) = 1
* s(x) is sos
* t(x) is sos
*/
class UnivariateNonnegative1 {
public:
DRAKE_NO_COPY_NO_MOVE_NO_ASSIGN(UnivariateNonnegative1)
UnivariateNonnegative1();
const MathematicalProgram& prog() const { return prog_; }
void CheckResult(const MathematicalProgramResult& result, double tol) const;
private:
MathematicalProgram prog_;
symbolic::Variable a_;
symbolic::Variable b_;
symbolic::Variable c_;
symbolic::Variable x_;
symbolic::Polynomial p_;
symbolic::Polynomial s_;
symbolic::Polynomial t_;
MatrixXDecisionVariable gram_s_;
MatrixXDecisionVariable gram_t_;
};
} // namespace solvers
} // namespace drake