5.4.5. Python 的SDP建模与优化¶
在本节中，我们将使用 MindOpt Python 语言的 API 来建模以及求解 半定规划问题示例 中的问题。
5.4.5.1. SDP示例1¶
首先，需要引入 Python 包：
23from mindoptpy import *
第一步：创建优化模型：
先建立一空的MindOpt模型。
30 # Step 1. Create a model and change the parameters.
31 model = MdoModel()
第二步：SDP模型输入
接着，我们将目标函数改为maximization，并新增 \(\mathbf{X}\) 对称矩阵变量。
34 # Step 2. Input model.
35 # Change to maximization problem.
36 model.set_int_attr(MDO_INT_ATTR.MIN_SENSE, 0)
37
38 # Add matrix variable.
39 model.add_sym_mat(3, "X")
我们利用 mindoptpy.MdoModel.replace_sym_mat_objs() 来输入目标函数中的 \(\mathbf{C}\) 的非零元素。
其中，第一个参数为 矩阵的索引值 (此处为0)，第二和第三个参数分别代表 矩阵中非零元素其对应的行和
列之索引 ，而最后一个参数则为 非零元的数值。
41 # Input objective coefficients.
42 model.replace_sym_mat_objs(0, [ 0, 0, 1, 2 ], [ 0, 2, 1, 2], [ -3.0, 1.0, -2.0, -3.0])
接着，我们输入第一个约束。我们首先输一空约束；接着再利用 mindoptpy.MdoModel.replace_sym_mat_elements() 输入约束中的
\(\mathbf{A}\) 的非零元素。其中，第一和第二个参数分别为： 约束的索引值 (此处为0)以及 矩阵的索引值 (此处为0)，第三第和第四个参数 矩阵中非零元素其对应的行和列之索引，而最后一个参数则为 非零元的数值。
44 # Input first constraint.
45 model.add_cons(1.0, 1.0, name="c0")
46 model.replace_sym_mat_elements(0, 0, [ 0, 2, 1, 2 ], [ 0, 0, 1, 2 ], [ 3.0, 1.0, 4.0, 5.0 ])
第三步：求解SDP模型
模型输入后，我们接着用 mindoptpy.MdoModel.solve_prob() 来求解问题，并用 mindoptpy.MdoModel.display_results() 呈现求解结果。
48 # Step 3. Solve the problem and populate the result.
49 model.solve_prob()
50 model.display_results()
第四步: 取得SDP模型的解
我们利用泛用函数 mindoptpy.MdoModel.get_real_attr() 来取得最优的的目标函数值
54 print("Optimizer terminated with an OPTIMAL status (code {0}).".format(status_code))
55 print(" - Primal objective : {:8.6f}".format(model.get_real_attr(MDO_REAL_ATTR.PRIMAL_OBJ_VAL)))
最后，我们利用 mindoptpy.MdoModel.get_real_attr_sym_mat() 来取得 \(\mathbf{X}\) 的变量值。其中，
第一个参数定义 属性值 (附注: mindoptpy.MdoModel.get_real_attr_sym_mat() 为泛用函数，因此必须配合
属性值来定义函数的使用目的；MDO_REAL_ATTR.SYM_MAT_PRIMAL_SOLN 定义目的为
获取(对称矩阵)变量的数值)，第二个参数为： 矩阵的索引值 (此处为0)，第三第和第四
个参数为 矩阵中需要回传的元素其对应的行和列之索引。返回值soln则为矩阵中对应的元素数值。
57 soln = model.get_real_attr_sym_mat(MDO_REAL_ATTR.SYM_MAT_PRIMAL_SOLN, 0,
58 [i * 3 + j for i in range(3) for j in range(3)],
59 [j * 3 + i for i in range(3) for j in range(3)])
文件链接 mdo_sdo_ex1.py 提供了完整源代码：
1"""
2/**
3 * Description
4 * -----------
5 *
6 * Semidefinite optimization (row-wise input).
7 *
8 * Formulation
9 * -----------
10 *
11 * Maximize
12 * obj: tr(C X)
13 *
14 * Subject To
15 * c0 : tr(A X) = 1
16 * Matrix
17 * C = [ -3 0 1 ] A = [ 3 0 1 ]
18 * [ 0 -2 0 ] [ 0 4 0 ]
19 * [ 1 0 -3 ] [ 1 0 5 ]
20 * End
21 */
22 """
23from mindoptpy import *
24
25
26if __name__ == "__main__":
27
28 MDO_INFINITY = MdoModel.get_infinity()
29
30 # Step 1. Create a model and change the parameters.
31 model = MdoModel()
32
33 try:
34 # Step 2. Input model.
35 # Change to maximization problem.
36 model.set_int_attr(MDO_INT_ATTR.MIN_SENSE, 0)
37
38 # Add matrix variable.
39 model.add_sym_mat(3, "X")
40
41 # Input objective coefficients.
42 model.replace_sym_mat_objs(0, [ 0, 0, 1, 2 ], [ 0, 2, 1, 2], [ -3.0, 1.0, -2.0, -3.0])
43
44 # Input first constraint.
45 model.add_cons(1.0, 1.0, name="c0")
46 model.replace_sym_mat_elements(0, 0, [ 0, 2, 1, 2 ], [ 0, 0, 1, 2 ], [ 3.0, 1.0, 4.0, 5.0 ])
47
48 # Step 3. Solve the problem and populate the result.
49 model.solve_prob()
50 model.display_results()
51
52 status_code, status_msg = model.get_status()
53 if status_msg == "OPTIMAL":
54 print("Optimizer terminated with an OPTIMAL status (code {0}).".format(status_code))
55 print(" - Primal objective : {:8.6f}".format(model.get_real_attr(MDO_REAL_ATTR.PRIMAL_OBJ_VAL)))
56
57 soln = model.get_real_attr_sym_mat(MDO_REAL_ATTR.SYM_MAT_PRIMAL_SOLN, 0,
58 [i * 3 + j for i in range(3) for j in range(3)],
59 [j * 3 + i for i in range(3) for j in range(3)])
60 print("X = ")
61 for i in range(3):
62 print(" (", end="")
63 for j in range(3):
64 print(" {0:8.6f}".format(soln[i * 3 + j]), end=""),
65 print(" )")
66
67 else:
68 print("Optimizer terminated with a(n) {0} status (code {1}).".format(status_msg, status_code))
69
70 except MdoError as e:
71 print("Received Mindopt exception.")
72 print(" - Code : {}".format(e.code))
73 print(" - Reason : {}".format(e.message))
74 except Exception as e:
75 print("Received exception.")
76 print(" - Reason : {}".format(e))
77 finally:
78 # Step 4. Free the model.
79 model.free_mdl()
5.4.5.2. SDP示例2¶
首先，需要引入 Python 包：
30from mindoptpy import *
第一步：创建MindOpt模型
首先，我们必须先建立一空的MindOpt模型。
37 # Step 1. Create a model and change the parameters.
38 model = MdoModel()
第二步：SDP模型输入
接着，我们将目标函数改为maximization，并新增 \(x_0\) 和 \(x_1\) 两个变量，以及 \(\mathbf{X}_0\) 和 \(\mathbf{X}_1\) 两个对称矩阵变量。
41 # Step 2. Input model.
42 # Change to maximization problem.
43 model.set_int_attr(MDO_INT_ATTR.MIN_SENSE, 0)
44
45 # Add variables.
46 x = []
47 x.append(model.add_var(0.0, MDO_INFINITY, 0.0, None, "x0", False))
48 x.append(model.add_var(0.0, MDO_INFINITY, 0.0, None, "x1", False))
49 # Add matrix variables.
50 model.add_sym_mats([ 2, 3 ], [ "X0", "X1" ]);
我们利用 mindoptpy.MdoModel.replace_sym_mat_objs() 来输入目标函数中的 \(\mathbf{C}_0\) 的非零元素。
其中，第一个参数为 矩阵的索引值 (此处为0)，第二和第三个参数分别代表 矩阵中非零元素其对应的行和列之索引 ，
而最后一个参数则为 非零元的数值 。
52 # Input objective coefficients.
53 model.replace_sym_mat_objs(0, [ 0, 1, 1 ], [ 0, 0, 1 ], [ 2.0, 1.0, 2.0 ]);
同理，我们再接着输入目标函数中的 \(\mathbf{C}_1\) 的非零元素。
54 model.replace_sym_mat_objs(1, [ 0, 0, 1, 2 ], [ 0, 2, 1, 2], [ 3.0, 1.0, 2.0, 3.0]);
接着，我们输入第一个约束。我们首先输入 \(x_0=1\) 接着再利用 mindoptpy.MdoModel.replace_sym_mat_elements()
输入约束中的 \(\mathbf{A}_{00}\) 的非零元素。其中，第一和第二个参数分别为： 约束的索引值 (此处为0)以及 矩阵的索引值 (此处为0)，第三第和第四个参数 矩阵中非零元素其对应的行和列之索引 ，而最后一个参数则为 非零元的数值 。
56 # Input first constraint.
57 model.add_cons(1.0 * x[0] == 1.0, "c0");
58 model.replace_sym_mat_elements(0, 0, [ 0, 1, 1 ], [ 0, 0, 1 ], [ 3.0, 1.0, 3.0 ]);
同理，我们再接着输入第二个约束。
60 # Input second constraint.
61 model.add_cons(1.0 * x[1] == 2.0, "c1");
62 model.replace_sym_mat_elements(1, 1, [ 0, 2, 1, 2 ], [ 0, 0, 1, 2 ], [ 3.0, 1.0, 4.0, 5.0 ]);
第三步：求解SDP模型
模型输入后，我们接着用 mindoptpy.MdoModel.solve_prob() 来求解问题，并用 mindoptpy.MdoModel.display_results() 呈现求解结果。
64 # Step 3. Solve the problem and populate the result.
65 model.solve_prob()
66 model.display_results()
第四步: 取得SDP模型的解
我们利用泛用函数 mindoptpy.MdoModel.get_real_attr() 来
取得最优的的目标函数值，并用泛用函数 mindoptpy.MdoModel.get_real_attr_array() 来
取得 \(x_0\) 和 \(x_1\) 两个变量值；此处第一个参数定义了 属性值 (附注: mindoptpy.MdoModel.get_real_attr_array() 为泛用函数，
因此必须配合属性值来定义函数的使用目的；MDO_REAL_ATTR.PRIMAL_SOLN定义目的为获取(线性)变量的数值)，
第二和第三个参数分别代表第一个索引和最后一个索引 加上1，
也就是 \(x_j,\) 其中 \(0 \leq j < 2\)。
71 print(" - Primal objective : {:8.6f}".format(model.get_real_attr(MDO_REAL_ATTR.PRIMAL_OBJ_VAL)))
72 soln = model.get_real_attr_array(MDO_REAL_ATTR.PRIMAL_SOLN, 0, 2)
最后，我们利用 mindoptpy.MdoModel.get_real_attr_sym_mat() 来取得 \(\mathbf{X}_0\) 的变量值。其中，
第一个参数定义 属性值 (附注: mindoptpy.MdoModel.get_real_attr_sym_mat() 为泛用函数，
因此必须配合属性值来定义函数的使用目的；MDO_REAL_ATTR.SYM_MAT_PRIMAL_SOLN
定义目的为获取(对称矩阵)变量的数值)，第二个参数为： 矩阵的索引值 (此处为0)，
第三第和第四个参数为 矩阵中需要回传的元素其对应的行和列之索引 。
返回值soln则为矩阵中对应的元素数值。
76 soln = model.get_real_attr_sym_mat(MDO_REAL_ATTR.SYM_MAT_PRIMAL_SOLN, 0,
77 [i * 2 + j for i in range(2) for j in range(2)],
78 [j * 2 + i for i in range(2) for j in range(2)])
文件链接 mdo_sdo_ex2.py 提供了完整源代码：
1"""
2/**
3 * Description
4 * -----------
5 *
6 * Semidefinite optimization (row-wise input).
7 *
8 * Formulation
9 * -----------
10 *
11 * Maximize
12 * obj: tr(C0 X0) + tr(C1 X1) + 0 x0 + 0 x1
13 *
14 * Subject To
15 * c0 : tr(A00 X0) + 1 x0 = 1
16 * c1 : tr(A11 X1) + 1 x1 = 2
17 * Bounds
18 * 0 <= x0
19 * 0 <= x1
20 * Matrix
21 * C0 = [ 2 1 ] A00 = [ 3 1 ]
22 * [ 1 2 ] [ 1 3 ]
23 *
24 * C1 = [ 3 0 1 ] A11 = [ 3 0 1 ]
25 * [ 0 2 0 ] [ 0 4 0 ]
26 * [ 1 0 3 ] [ 1 0 5 ]
27 * End
28 */
29 """
30from mindoptpy import *
31
32
33if __name__ == "__main__":
34
35 MDO_INFINITY = MdoModel.get_infinity()
36
37 # Step 1. Create a model and change the parameters.
38 model = MdoModel()
39
40 try:
41 # Step 2. Input model.
42 # Change to maximization problem.
43 model.set_int_attr(MDO_INT_ATTR.MIN_SENSE, 0)
44
45 # Add variables.
46 x = []
47 x.append(model.add_var(0.0, MDO_INFINITY, 0.0, None, "x0", False))
48 x.append(model.add_var(0.0, MDO_INFINITY, 0.0, None, "x1", False))
49 # Add matrix variables.
50 model.add_sym_mats([ 2, 3 ], [ "X0", "X1" ]);
51
52 # Input objective coefficients.
53 model.replace_sym_mat_objs(0, [ 0, 1, 1 ], [ 0, 0, 1 ], [ 2.0, 1.0, 2.0 ]);
54 model.replace_sym_mat_objs(1, [ 0, 0, 1, 2 ], [ 0, 2, 1, 2], [ 3.0, 1.0, 2.0, 3.0]);
55
56 # Input first constraint.
57 model.add_cons(1.0 * x[0] == 1.0, "c0");
58 model.replace_sym_mat_elements(0, 0, [ 0, 1, 1 ], [ 0, 0, 1 ], [ 3.0, 1.0, 3.0 ]);
59
60 # Input second constraint.
61 model.add_cons(1.0 * x[1] == 2.0, "c1");
62 model.replace_sym_mat_elements(1, 1, [ 0, 2, 1, 2 ], [ 0, 0, 1, 2 ], [ 3.0, 1.0, 4.0, 5.0 ]);
63
64 # Step 3. Solve the problem and populate the result.
65 model.solve_prob()
66 model.display_results()
67
68 status_code, status_msg = model.get_status()
69 if status_msg == "OPTIMAL":
70 print("Optimizer terminated with an OPTIMAL status (code {0}).".format(status_code))
71 print(" - Primal objective : {:8.6f}".format(model.get_real_attr(MDO_REAL_ATTR.PRIMAL_OBJ_VAL)))
72 soln = model.get_real_attr_array(MDO_REAL_ATTR.PRIMAL_SOLN, 0, 2)
73 for index, value in enumerate(soln):
74 print("x[{0}]={1:8.6f}".format(index, value))
75
76 soln = model.get_real_attr_sym_mat(MDO_REAL_ATTR.SYM_MAT_PRIMAL_SOLN, 0,
77 [i * 2 + j for i in range(2) for j in range(2)],
78 [j * 2 + i for i in range(2) for j in range(2)])
79 print("X[0] = ")
80 for i in range(2):
81 print(" (", end="")
82 for j in range(2):
83 print(" {0:8.6f}".format(soln[i * 2 + j]), end=""),
84 print(" )")
85
86 soln = model.get_real_attr_sym_mat(MDO_REAL_ATTR.SYM_MAT_PRIMAL_SOLN, 1,
87 [i * 3 + j for i in range(3) for j in range(3)],
88 [j * 3 + i for i in range(3) for j in range(3)])
89 print("X[1] = ")
90 for i in range(3):
91 print(" (", end=""),
92 for j in range(3):
93 print(" {0:8.6f}".format(soln[i * 3 + j]), end=""),
94 print(" )")
95
96 else:
97 print("Optimizer terminated with a(n) {0} status (code {1}).".format(status_msg, status_code))
98
99 except MdoError as e:
100 print("Received Mindopt exception.")
101 print(" - Code : {}".format(e.code))
102 print(" - Reason : {}".format(e.message))
103 except Exception as e:
104 print("Received exception.")
105 print(" - Reason : {}".format(e))
106 finally:
107 # Step 4. Free the model.
108 model.free_mdl()