5.4.3. C语言的SDP建模与优化¶
在本节中，我们将使用 MindOpt C 语言的 API 来建模以及求解 半定规划问题示例 中的问题。
5.4.3.1. SDP示例1¶
首先，需要引入头文件：
22#include "Mindopt.h"
第一步：创建优化模型：
先建立一空的MindOpt模型。
67 /*------------------------------------------------------------------*/
68 /* Step 1. Create a model and change the parameters. */
69 /*------------------------------------------------------------------*/
70 /* Create an empty model. */
71 MDO_CHECK_CALL(Mdo_createMdl(&model));
第二步：SDP模型输入
接着，我们将目标函数改为maximization，并新增 \(\mathbf{X}\) 对称矩阵变量。
73 /*------------------------------------------------------------------*/
74 /* Step 2. Input model. */
75 /*------------------------------------------------------------------*/
76 /* Change to maximization problem. */
77 MDO_CHECK_CALL(Mdo_setIntAttr(model, MDO_INT_ATTR_MIN_SENSE, MDO_NO));
78
79 /* Add matrix variable. */
80 MDO_CHECK_CALL(Mdo_addSymMat(model, dim_mat, mat_name));
我们利用 Mdo_replaceSymMatObjs() 来输入目标函数中的 \(\mathbf{C}\) 的非零元素。
其中，第一个参数为 模型， 第二个参数为 矩阵的索引值 (此处为0)，第三个元素代表着 矩阵非零的元素的个数 第四个和第五个参数分别代表 矩阵中非零元素其对应的行和列之索引 ，
而最后一个参数则为 非零元的数值的索引 。
82 /* Input objective coefficients. */
83 MDO_CHECK_CALL(Mdo_replaceSymMatObjs(
84 model, 0, C_size, C_row_indices, C_col_indices, C_values));
接着，我们输入第一个约束。我们首先输一空约束；接着再利用 Mdo_replaceSymMatElements() 输入约束中的
\(\mathbf{A}\) 的非零元素。其中，第二和第三个参数分别为： 约束的索引值 (此处为0)以及 矩阵的索引值 (此处为0)，第四个参数代表 矩阵的非零元素的个数，第五个和第六个参数代表 矩阵中非零元素其对应的行和列之索引，而最后一个参数则为 非零元的数值的索引。
86 /* Input first constraint. */
87 MDO_CHECK_CALL(Mdo_addRow(model, 1.0, 1.0, 0, NULL, NULL, "c0"));
88 MDO_CHECK_CALL(Mdo_replaceSymMatElements(model, 0, 0, A_size,
89 A_row_indices_in_cone, A_col_indices_in_cone, A_values_in_cone));
第三步：求解SDP模型
模型输入后，我们接着用 Mdo_solveProb() 来求解问题，并用 Mdo_displayResults() 呈现求解结果。
94 /* Solve the problem. */
95 MDO_CHECK_CALL(Mdo_solveProb(model));
96 Mdo_displayResults(model);
第四步: 取得SDP模型的解
我们利用泛用函数 Mdo_getRealAttr() 来取得最优的的目标函数值
104 MDO_CHECK_CALL(Mdo_getRealAttr(model, MDO_REAL_ATTR_PRIMAL_OBJ_VAL, &val));
最后，我们利用 Mdo_getRealAttrSymMat() 来取得 \(\mathbf{X}\) 的变量值。其中，
第二个参数定义 属性值 (附注: Mdo_getRealAttrSymMat() 为泛用函数，因此必须配合
属性值来定义函数的使用目的; MDO_REAL_ATTR_SYM_MAT_PRIMAL_SOLN定义目的为
获取(对称矩阵)变量的数值)，第三个参数为: 矩阵的索引值 (此处为0)，第四个参数为 矩阵的尺寸 返回值soln则为矩阵中对应元素的索引。
107 MDO_CHECK_CALL(Mdo_getRealAttrSymMat(
108 model, MDO_REAL_ATTR_SYM_MAT_PRIMAL_SOLN, 0, dim_mat * dim_mat, NULL, NULL, soln));
文件链接 MdoSdoEx1.c 提供了完整源代码：
1/**
2 * Description
3 * -----------
4 *
5 * Semidefinite optimization (row-wise input).
6 *
7 * Formulation
8 * -----------
9 *
10 * Maximize
11 * obj: tr(C X)
12 *
13 * Subject To
14 * c0 : tr(A X) = 1
15 * Matrix
16 * C = [ -3 0 1 ] A = [ 3 0 1 ]
17 * [ 0 -2 0 ] [ 0 4 0 ]
18 * [ 1 0 -3 ] [ 1 0 5 ]
19 * End
20 */
21#include <stdio.h>
22#include "Mindopt.h"
23
24/* Macro to check the return code. */
25#define MDO_CHECK_CALL(MDO_CALL) \
26 code = MDO_CALL; \
27 if (code != MDO_OKAY) \
28 { \
29 Mdo_explainResult(model, code, str); \
30 Mdo_freeMdl(&model); \
31 fprintf(stderr, "===================================\n"); \
32 fprintf(stderr, "Error : code <%d>\n", code); \
33 fprintf(stderr, "Reason : %s\n", str); \
34 fprintf(stderr, "===================================\n"); \
35 return (int)code; \
36 }
37
38int main(void)
39{
40 /* Variables. */
41 char str[1024] = { "\0" };
42 MdoMdl * model = NULL;
43 MdoResult code = MDO_OKAY;
44 MdoStatus status = MDO_UNKNOWN;
45
46 /* Input data. */
47 const int num_mats = 1;
48 const int dim_mat = 3 ; /* Dimension of the matrix variables. */
49 const char * mat_name = "X";
50
51 const int C_size = 4;
52 const int C_row_indices[] = { 0, 0, 1, 2 }; /* Row index of a matrix variable. */
53 const int C_col_indices[] = { 0, 2, 1, 2 }; /* Column index of a matrix variable. */
54 const double C_values[] = { -3.0, 1.0, -2.0, -3.0 }; /* Values of a matrix variable. */
55
56 const int A_size = 4;
57 const int A_mat_idx = 1;
58 const int A_row_indices_in_cone[] = { 0, 2, 1, 2 }; /* Row index in a matrix variable. */
59 const int A_col_indices_in_cone[] = { 0, 0, 1, 2 }; /* Column index in a matrix variable. */
60 const double A_values_in_cone[] = { 3.0, 1.0, 4.0, 5.0 }; /* Values of a matrix variable. */
61
62 /* Temp. */
63 double soln[15];
64 double val;
65 int i, j, b;
66
67 /*------------------------------------------------------------------*/
68 /* Step 1. Create a model and change the parameters. */
69 /*------------------------------------------------------------------*/
70 /* Create an empty model. */
71 MDO_CHECK_CALL(Mdo_createMdl(&model));
72
73 /*------------------------------------------------------------------*/
74 /* Step 2. Input model. */
75 /*------------------------------------------------------------------*/
76 /* Change to maximization problem. */
77 MDO_CHECK_CALL(Mdo_setIntAttr(model, MDO_INT_ATTR_MIN_SENSE, MDO_NO));
78
79 /* Add matrix variable. */
80 MDO_CHECK_CALL(Mdo_addSymMat(model, dim_mat, mat_name));
81
82 /* Input objective coefficients. */
83 MDO_CHECK_CALL(Mdo_replaceSymMatObjs(
84 model, 0, C_size, C_row_indices, C_col_indices, C_values));
85
86 /* Input first constraint. */
87 MDO_CHECK_CALL(Mdo_addRow(model, 1.0, 1.0, 0, NULL, NULL, "c0"));
88 MDO_CHECK_CALL(Mdo_replaceSymMatElements(model, 0, 0, A_size,
89 A_row_indices_in_cone, A_col_indices_in_cone, A_values_in_cone));
90
91 /*------------------------------------------------------------------*/
92 /* Step 3. Solve the problem and populate the result. */
93 /*------------------------------------------------------------------*/
94 /* Solve the problem. */
95 MDO_CHECK_CALL(Mdo_solveProb(model));
96 Mdo_displayResults(model);
97
98 switch (Mdo_getStatus(model))
99 {
100 case MDO_UNKNOWN:
101 printf("Optimizer terminated with an UNKNOWN status.\n");
102 break;
103 case MDO_OPTIMAL:
104 MDO_CHECK_CALL(Mdo_getRealAttr(model, MDO_REAL_ATTR_PRIMAL_OBJ_VAL, &val));
105 printf("Optimizer terminated with an OPTIMAL status.\n");
106 printf(" - Primal objective : %e.\n", val);
107 MDO_CHECK_CALL(Mdo_getRealAttrSymMat(
108 model, MDO_REAL_ATTR_SYM_MAT_PRIMAL_SOLN, 0, dim_mat * dim_mat, NULL, NULL, soln));
109 printf("X = \n");
110 for (i = 0; i < dim_mat; ++i)
111 {
112 printf(" (");
113 for (j = 0; j < dim_mat; ++j)
114 {
115 printf(" %+8.6e", soln[i * dim_mat + j]);
116 }
117 printf(" )\n");
118 }
119 break;
120 case MDO_INFEASIBLE:
121 printf("Optimizer terminated with an INFEASIBLE status.\n");
122 break;
123 case MDO_UNBOUNDED:
124 printf("Optimizer terminated with an UNBOUNDED status.\n");
125 break;
126 case MDO_INF_OR_UBD:
127 printf("Optimizer terminated with an INFEASIBLE or UNBOUNDED status.\n");
128 break;
129 }
130
131 /*------------------------------------------------------------------*/
132 /* Step 4. Free the model. */
133 /*------------------------------------------------------------------*/
134 /* Free the model. */
135 Mdo_freeMdl(&model);
136
137 return (int)code;
138}
5.4.3.2. SDP示例2¶
首先，需要引入头文件：
29 Mdo_explainResult(model, code, str); \
第一步：创建MindOpt模型
首先，我们必须先建立一空的MindOpt模型。
95 /*------------------------------------------------------------------*/
96 /* Step 1. Create a model and change the parameters. */
97 /*------------------------------------------------------------------*/
98 /* Create an empty model. */
99 MDO_CHECK_CALL(Mdo_createMdl(&model));
第二步：SDP模型输入
接着，我们将目标函数改为maximization，并新增 \(x_0\) 和 \(x_1\) 两个变量，以及 \(\mathbf{X}_0\) 和 \(\mathbf{X}_1\) 两个对称矩阵变量。
101 /*------------------------------------------------------------------*/
102 /* Step 2. Input model. */
103 /*------------------------------------------------------------------*/
104 /* Change to maximization problem. */
105 MDO_CHECK_CALL(Mdo_setIntAttr(model, MDO_INT_ATTR_MIN_SENSE, MDO_NO));
106
107 /* Add variables. */
108 MDO_CHECK_CALL(Mdo_addCol(model, 0.0, MDO_INFINITY, 0.0, 0, NULL, NULL, "x0", MDO_NO));
109 MDO_CHECK_CALL(Mdo_addCol(model, 0.0, MDO_INFINITY, 0.0, 0, NULL, NULL, "x1", MDO_NO));
110
111 /* Add matrix variables. */
112 MDO_CHECK_CALL(Mdo_addSymMats(model, num_mats, dim_mats, mat_names));
我们利用 Mdo_replaceSymMatObjs() 来输入目标函数中的 \(\mathbf{C}_0\) 的非零元素。
其中，第二个参数为 矩阵的索引值 (此处为 C0_mat_idx)，第三个元素代表着 矩阵非零的元素的个数 第四个和第五个参数分别代表 矩阵中非零元素其对应的行和列之索引 ，
而最后一个参数则为 非零元的数值的索引 。
114 /* Input objective coefficients. */
115 MDO_CHECK_CALL(Mdo_replaceSymMatObjs(
116 model, C0_mat_idx, C0_size, C0_row_indices, C0_col_indices, C0_values));
同理，我们再接着输入目标函数中的 \(\mathbf{C}_1\) 的非零元素。
117 MDO_CHECK_CALL(Mdo_replaceSymMatObjs(
118 model, C1_mat_idx, C1_size, C1_row_indices, C1_col_indices, C1_values));
接着，我们输入第一个约束。我们首先输入空约束，接着再利用 Mdo_replaceSymMatElements()
输入约束中的 \(\mathbf{A}_{00}\) 的非零元素。其中，第二个和第三个参数分别为： 约束的索引值 (此处为0)以及 矩阵的索引值 (此处为A00_mat_idx)，第四个参数为 矩阵的非零元素的个数，第四个和第五个参数 矩阵中非零元素其对应的行和列之索引 ，而最后一个参数则为 非零元的数值的索引 。
120 /* Input first constraint. */
121 MDO_CHECK_CALL(Mdo_addRow(model, 1.0, 1.0, row0_size, row0_idx, row0_val, "c0"));
122 MDO_CHECK_CALL(Mdo_replaceSymMatElements(model, 0, A00_mat_idx, A00_size,
123 A00_row_indices_in_cone, A00_col_indices_in_cone, A00_values_in_cone));
同理，我们再接着输入第二个约束。
125 /* Input second constraint. */
126 MDO_CHECK_CALL(Mdo_addRow(model, 2.0, 2.0, row1_size, row1_idx, row1_val, "c1"));
127 MDO_CHECK_CALL(Mdo_replaceSymMatElements(model, 1, A11_mat_idx, A11_size,
128 A11_row_indices_in_cone, A11_col_indices_in_cone, A11_values_in_cone));
第三步：求解SDP模型
模型输入后，我们接着用 Mdo_solveProb() 来求解问题，并用 Mdo_displayResults() 呈现求解结果。
130 /*------------------------------------------------------------------*/
131 /* Step 3. Solve the problem and populate the result. */
132 /*------------------------------------------------------------------*/
133 /* Solve the problem. */
134 MDO_CHECK_CALL(Mdo_solveProb(model));
135 Mdo_displayResults(model);
第四步: 取得SDP模型的解
我们利用泛用函数 Mdo_getRealAttr() 来取得最优的的目标函数值，并用泛用函数 Mdo_getRealAttrArray()
来取得 \(x_0\) 和 \(x_1\) 两个变量值；此处第二个参数定义了 属性值 (附注: Mdo_getRealAttrArray() 为泛用函数，
因此必须配合属性值来定义函数的使用目的；MDO_REAL_ATTR_PRIMAL_SOLN 定义目的为获取(线性)变量的数值)，
第三和第四个参数分别代表第一个索引和最后一个索引 加上1，
也就是 \(x_j\), 其中 \(0 \leq j < 2\) 。
143 MDO_CHECK_CALL(Mdo_getRealAttr(model, MDO_REAL_ATTR_PRIMAL_OBJ_VAL, &val));
144 printf("Optimizer terminated with an OPTIMAL status.\n");
145 printf(" - Primal objective : %e.\n", val);
146 MDO_CHECK_CALL(Mdo_getRealAttrArray(model, MDO_REAL_ATTR_PRIMAL_SOLN, 0, 2, soln));
最后，我们利用 Mdo_getRealAttrSymMat() 来取得 \(\mathbf{X}_0\) 和 \(\mathbf{X}_1\) 的变量值。其中，
第二个参数定义 属性值 (附注: Mdo_getRealAttrSymMat() 为泛用函数，
因此必须配合属性值来定义函数的使用目的；MDO_REAL_ATTR_SYM_MAT_PRIMAL_SOLN
定义目的为获取(对称矩阵)变量的数值)，第三个参数为： 矩阵的序列数 (此处为b)，
第四个参数为 矩阵的尺寸 ，
返回值soln则为矩阵中对应的元素的索引。
153 MDO_CHECK_CALL(Mdo_getRealAttrSymMat(
154 model, MDO_REAL_ATTR_SYM_MAT_PRIMAL_SOLN, b, dim_mats[b] * dim_mats[b], NULL, NULL, soln));
文件链接 MdoSdoEx2.c 提供了完整源代码：
1/**
2 * Description
3 * -----------
4 *
5 * Semidefinite optimization (row-wise input).
6 *
7 * Formulation
8 * -----------
9 *
10 * Maximize
11 * obj: tr(C0 X0) + tr(C1 X1) + 0 x0 + 0 x1
12 *
13 * Subject To
14 * c0 : tr(A00 X0) + 1 x0 = 1
15 * c1 : tr(A11 X1) + 1 x1 = 2
16 * Bounds
17 * 0 <= x0
18 * 0 <= x1
19 * Matrix
20 * C0 = [ 2 1 ] A00 = [ 3 1 ]
21 * [ 1 2 ] [ 1 3 ]
22 *
23 * C1 = [ 3 0 1 ] A11 = [ 3 0 1 ]
24 * [ 0 2 0 ] [ 0 4 0 ]
25 * [ 1 0 3 ] [ 1 0 5 ]
26 * End
27 */
28#include <stdio.h>
29#include "Mindopt.h"
30
31/* Macro to check the return code. */
32#define MDO_CHECK_CALL(MDO_CALL) \
33 code = MDO_CALL; \
34 if (code != MDO_OKAY) \
35 { \
36 Mdo_explainResult(model, code, str); \
37 Mdo_freeMdl(&model); \
38 fprintf(stderr, "===================================\n"); \
39 fprintf(stderr, "Error : code <%d>\n", code); \
40 fprintf(stderr, "Reason : %s\n", str); \
41 fprintf(stderr, "===================================\n"); \
42 return (int)code; \
43 }
44
45int main(void)
46{
47 /* Variables. */
48 char str[1024] = { "\0" };
49 MdoMdl * model = NULL;
50 MdoResult code = MDO_OKAY;
51 MdoStatus status = MDO_UNKNOWN;
52
53 /* Input data. */
54 const int num_cols = 2;
55 const int num_mats = 2;
56 const int dim_mats[] = { 2, 3 }; /* Dimension of the matrix variables. */
57 const char * mat_names[] = { "X0", "X1" };
58
59 const int row0_size = 1;
60 const int row0_idx[] = { 0 };
61 const double row0_val[] = { 1.0 };
62 const int row1_size = 1;
63 const int row1_idx[] = { 1 };
64 const double row1_val[] = { 1.0 };
65
66 const int C0_size = 3;
67 const int C0_mat_idx = 0;
68 const int C0_row_indices[] = { 0, 1, 1 }; /* Row index of a matrix variable. */
69 const int C0_col_indices[] = { 0, 0, 1 }; /* Column index of a matrix variable. */
70 const double C0_values[] = { 2.0, 1.0, 2.0 }; /* Values of a matrix variable. */
71
72 const int C1_size = 4;
73 const int C1_mat_idx = 1;
74 const int C1_row_indices[] = { 0, 0, 1, 2 }; /* Row index of a matrix variable. */
75 const int C1_col_indices[] = { 0, 2, 1, 2 }; /* Column index of a matrix variable. */
76 const double C1_values[] = { 3.0, 1.0, 2.0, 3.0 }; /* Values of a matrix variable. */
77
78 const int A00_size = 3;
79 const int A00_mat_idx = 0;
80 const int A00_row_indices_in_cone[] = { 0, 1, 1 }; /* Row index in a matrix variable. */
81 const int A00_col_indices_in_cone[] = { 0, 0, 1 }; /* Column index in a matrix variable. */
82 const double A00_values_in_cone[] = { 3.0, 1.0, 3.0 }; /* Values of a matrix variable. */
83
84 const int A11_size = 4;
85 const int A11_mat_idx = 1;
86 const int A11_row_indices_in_cone[] = { 0, 2, 1, 2 }; /* Row index in a matrix variable. */
87 const int A11_col_indices_in_cone[] = { 0, 0, 1, 2 }; /* Column index in a matrix variable. */
88 const double A11_values_in_cone[] = { 3.0, 1.0, 4.0, 5.0 }; /* Values of a matrix variable. */
89
90 /* Temp. */
91 double soln[15];
92 double val;
93 int i, j, b;
94
95 /*------------------------------------------------------------------*/
96 /* Step 1. Create a model and change the parameters. */
97 /*------------------------------------------------------------------*/
98 /* Create an empty model. */
99 MDO_CHECK_CALL(Mdo_createMdl(&model));
100
101 /*------------------------------------------------------------------*/
102 /* Step 2. Input model. */
103 /*------------------------------------------------------------------*/
104 /* Change to maximization problem. */
105 MDO_CHECK_CALL(Mdo_setIntAttr(model, MDO_INT_ATTR_MIN_SENSE, MDO_NO));
106
107 /* Add variables. */
108 MDO_CHECK_CALL(Mdo_addCol(model, 0.0, MDO_INFINITY, 0.0, 0, NULL, NULL, "x0", MDO_NO));
109 MDO_CHECK_CALL(Mdo_addCol(model, 0.0, MDO_INFINITY, 0.0, 0, NULL, NULL, "x1", MDO_NO));
110
111 /* Add matrix variables. */
112 MDO_CHECK_CALL(Mdo_addSymMats(model, num_mats, dim_mats, mat_names));
113
114 /* Input objective coefficients. */
115 MDO_CHECK_CALL(Mdo_replaceSymMatObjs(
116 model, C0_mat_idx, C0_size, C0_row_indices, C0_col_indices, C0_values));
117 MDO_CHECK_CALL(Mdo_replaceSymMatObjs(
118 model, C1_mat_idx, C1_size, C1_row_indices, C1_col_indices, C1_values));
119
120 /* Input first constraint. */
121 MDO_CHECK_CALL(Mdo_addRow(model, 1.0, 1.0, row0_size, row0_idx, row0_val, "c0"));
122 MDO_CHECK_CALL(Mdo_replaceSymMatElements(model, 0, A00_mat_idx, A00_size,
123 A00_row_indices_in_cone, A00_col_indices_in_cone, A00_values_in_cone));
124
125 /* Input second constraint. */
126 MDO_CHECK_CALL(Mdo_addRow(model, 2.0, 2.0, row1_size, row1_idx, row1_val, "c1"));
127 MDO_CHECK_CALL(Mdo_replaceSymMatElements(model, 1, A11_mat_idx, A11_size,
128 A11_row_indices_in_cone, A11_col_indices_in_cone, A11_values_in_cone));
129
130 /*------------------------------------------------------------------*/
131 /* Step 3. Solve the problem and populate the result. */
132 /*------------------------------------------------------------------*/
133 /* Solve the problem. */
134 MDO_CHECK_CALL(Mdo_solveProb(model));
135 Mdo_displayResults(model);
136
137 switch (Mdo_getStatus(model))
138 {
139 case MDO_UNKNOWN:
140 printf("Optimizer terminated with an UNKNOWN status.\n");
141 break;
142 case MDO_OPTIMAL:
143 MDO_CHECK_CALL(Mdo_getRealAttr(model, MDO_REAL_ATTR_PRIMAL_OBJ_VAL, &val));
144 printf("Optimizer terminated with an OPTIMAL status.\n");
145 printf(" - Primal objective : %e.\n", val);
146 MDO_CHECK_CALL(Mdo_getRealAttrArray(model, MDO_REAL_ATTR_PRIMAL_SOLN, 0, 2, soln));
147 for (j = 0; j < 2; ++j)
148 {
149 printf("x[%d] = %+8.6e\n", j, soln[j]);
150 }
151 for (b = 0; b < num_mats; ++b)
152 {
153 MDO_CHECK_CALL(Mdo_getRealAttrSymMat(
154 model, MDO_REAL_ATTR_SYM_MAT_PRIMAL_SOLN, b, dim_mats[b] * dim_mats[b], NULL, NULL, soln));
155 printf("X[%d] = \n", b);
156 for (i = 0; i < dim_mats[b]; ++i)
157 {
158 printf(" (");
159 for (j = 0; j < dim_mats[b]; ++j)
160 {
161 printf(" %+8.6e", soln[i * dim_mats[b] + j]);
162 }
163 printf(" )\n");
164 }
165 }
166 break;
167 case MDO_INFEASIBLE:
168 printf("Optimizer terminated with an INFEASIBLE status.\n");
169 break;
170 case MDO_UNBOUNDED:
171 printf("Optimizer terminated with an UNBOUNDED status.\n");
172 break;
173 case MDO_INF_OR_UBD:
174 printf("Optimizer terminated with an INFEASIBLE or UNBOUNDED status.\n");
175 break;
176 }
177
178 /*------------------------------------------------------------------*/
179 /* Step 4. Free the model. */
180 /*------------------------------------------------------------------*/
181 /* Free the model. */
182 Mdo_freeMdl(&model);
183
184 return (int)code;
185}