5.4.2. C的SDP建模与优化¶
在本节中,我们将使用 MindOpt 的 C API 来建模以及求解 半定规划问题示例 中的问题。
5.4.2.1. SDP示例1¶
首先,需要引入头文件:
26#include "Mindopt.h"
第一步:创建优化模型:
先建立一空的MindOpt模型。
68 /* Create an empty model. */
69 CHECK_RESULT(MDOemptyenv(&env));
70 CHECK_RESULT(MDOstartenv(env));
71 CHECK_RESULT(MDOnewmodel(env, &m, MODEL_NAME, 0, NULL, NULL, NULL, NULL, NULL));
第二步:SDP模型输入
接着,我们将目标函数改为 最大化。
77 /* Change to maximization problem. */
78 CHECK_RESULT(MDOsetintattr(m, MODEL_SENSE, MDO_MAXIMIZE));
- 我们利用
MDOaddpsdvar()
来创建半正定矩阵变量,并可以同时设定目标函数中其对应的矩阵系数: 第一个参数为 模型,
第二个参数为 矩阵变量的维数 (此处为3),
第三个元素代表着 矩阵系数的非零元素个数,
第四个和第五个参数分别代表 矩阵系数中非零元素对应的向量化索引与取值。
79 /* Add matrix variable and corresponding objective coefficients. */
80 CHECK_RESULT(MDOaddpsdvar(m, dim_mat, C_size, C_nnz_indices, C_nnz_values, NULL));
- 接着,我们输入约束。我们利用
MDOaddpsdconstr()
创建带半正定矩阵变量的约束: 第一个参数为 模型,
第二个参数为该约束中 线性项的个数 (此处为0),
第三个参数为该约束中 矩阵项的个数,
第四个参数代表 所有矩阵系数的总非零元素的个数,
第五个和第六个参数代表该约束中 线性项对应的索引与取值,
第七个参数代表该约束中 矩阵变量的索引,
第八个参数代表该约束中 各矩阵系数的初始索引,
第九个和第十个参数代表 所有矩阵系数的非零元素的向量化索引与取值,
第十一个参数代表该 约束的意义,MDO_EQUAL (‘=’)代表等式约束,
第十二个参数代表 约束的右侧值,
最后一个参数代表该约束中 名称。
82 /* Add Constraints. */
83 CHECK_RESULT(MDOaddpsdconstr(m, 0, num_psd_var, A_size, NULL, NULL, psd_var_indices, psd_coeff_start_indices, A_nnz_indices, A_nnz_values, MDO_EQUAL, 1.0, NULL));
第三步:求解SDP模型
模型输入后,我们接着用 MDOoptimize()
来求解问题。
88 /* Solve the problem. */
89 CHECK_RESULT(MDOoptimize(m));
第四步: 取得SDP模型的解
我们利用 MDOgetintattr()
得到求解状态,即 Status 属性;利用泛用函数 MDOgetdblattr()
来取得最优的的目标函数值,即 PrimalObjVal 属性。
90 CHECK_RESULT(MDOgetintattr(m, STATUS, &status));
91 if (status == MDO_OPTIMAL)
92 {
93 CHECK_RESULT(MDOgetdblattr(m, OBJ_VAL, &obj));
94 printf("The objective value is %f\n", obj);
95 }
96 else
97 {
98 printf("No feasible solution exists\n");
99 }
第五步: 释放模型
我们利用 MDOfreemodel()
和 MDOfreeenv()
来释放模型。
30#define RELEASE_MEMORY \ 31 MDOfreemodel(m); \ 32 MDOfreeenv(env);105 RELEASE_MEMORY;
示例 mdo_sdo_ex1.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 * Subject To
13 * c0 : tr(A X) = 1
14 * Bounds
15 * X is p.s.d.
16 *
17 * Matrix
18 * C = [ -3 0 1 ] A = [ 3 0 1 ]
19 * [ 0 -2 0 ] [ 0 4 0 ]
20 * [ 1 0 -3 ] [ 1 0 5 ]
21 * End
22 */
23
24#include <stdio.h>
25#include <stdlib.h>
26#include "Mindopt.h"
27
28/* Macro to check the return code. */
29#define CHECK_RESULT(code) { int res = code; if (res != 0) { fprintf(stderr, "Bad code: %d\n", res); exit(res); } }
30#define RELEASE_MEMORY \
31 MDOfreemodel(m); \
32 MDOfreeenv(env);
33#define MODEL_NAME "SDP_01"
34#define MODEL_SENSE "ModelSense"
35#define STATUS "Status"
36#define OBJ_VAL "PrimalObjVal"
37
38int main(void)
39{
40 /* Variables. */
41 MDOenv *env;
42 MDOmodel *m;
43 double obj;
44 int status;
45
46 /* Model data. */
47 // Variable dims.
48 int num_mats = 1;
49 int dim_mat = 3 ;
50
51 // Objective coefficients. (Only the right upper part in row-major sparse format)
52 int C_size = 4;
53 int C_nnz_indices[] = { 0, 2, 4, 8 };
54 double C_nnz_values[] = { -3.0, 1.0, -2.0, -3.0 };
55
56 // Constraint coefficients. (Only the right upper part in row-major sparse format)
57 int num_psd_var = 1;
58 int psd_var_indices[] = { 0 };
59 int psd_coeff_start_indices[] = { 0 };
60
61 int A_size = 4;
62 int A_nnz_indices[] = { 0, 2, 4, 8 };
63 double A_nnz_values[] = { 3.0, 1.0, 4.0, 5.0 };
64
65 /*------------------------------------------------------------------*/
66 /* Step 1. Create a model and change the parameters. */
67 /*------------------------------------------------------------------*/
68 /* Create an empty model. */
69 CHECK_RESULT(MDOemptyenv(&env));
70 CHECK_RESULT(MDOstartenv(env));
71 CHECK_RESULT(MDOnewmodel(env, &m, MODEL_NAME, 0, NULL, NULL, NULL, NULL, NULL));
72
73 /*------------------------------------------------------------------*/
74 /* Step 2. Input model. */
75 /*------------------------------------------------------------------*/
76 /* Change to maximization problem. */
77 CHECK_RESULT(MDOsetintattr(m, MODEL_SENSE, MDO_MAXIMIZE));
78
79 /* Add matrix variable and corresponding objective coefficients. */
80 CHECK_RESULT(MDOaddpsdvar(m, dim_mat, C_size, C_nnz_indices, C_nnz_values, NULL));
81
82 /* Add Constraints. */
83 CHECK_RESULT(MDOaddpsdconstr(m, 0, num_psd_var, A_size, NULL, NULL, psd_var_indices, psd_coeff_start_indices, A_nnz_indices, A_nnz_values, MDO_EQUAL, 1.0, NULL));
84
85 /*------------------------------------------------------------------*/
86 /* Step 3. Solve the problem and populate the result. */
87 /*------------------------------------------------------------------*/
88 /* Solve the problem. */
89 CHECK_RESULT(MDOoptimize(m));
90 CHECK_RESULT(MDOgetintattr(m, STATUS, &status));
91 if (status == MDO_OPTIMAL)
92 {
93 CHECK_RESULT(MDOgetdblattr(m, OBJ_VAL, &obj));
94 printf("The objective value is %f\n", obj);
95 }
96 else
97 {
98 printf("No feasible solution exists\n");
99 }
100
101 /*------------------------------------------------------------------*/
102 /* Step 4. Free the model. */
103 /*------------------------------------------------------------------*/
104 /* Free the model. */
105 RELEASE_MEMORY;
106
107 return 0;
108}
5.4.2.2. SDP示例2¶
首先,需要引入头文件:
32#include "Mindopt.h"
第一步:创建MindOpt模型
我们先建立一空的MindOpt模型。
94 /* Create an empty model. */
95 CHECK_RESULT(MDOemptyenv(&env));
96 CHECK_RESULT(MDOstartenv(env));
97 CHECK_RESULT(MDOnewmodel(env, &m, MODEL_NAME, 0, NULL, NULL, NULL, NULL, NULL));
第二步:SDP模型输入
接着,我们将目标函数改为 最大化。
102 /* Change to maximization problem. */
103 CHECK_RESULT(MDOsetintattr(m, MODEL_SENSE, MDO_MAXIMIZE));
我们利用 MDOaddvar()
创建两个线性变量,并设定他们的上下界与目标函数中的系数。
105 /* Add linear variable and corresponding objective coefficients. */
106 CHECK_RESULT(MDOaddvar(m, 0, NULL, NULL, 0.0, 0.0, MDO_INFINITY, MDO_CONTINUOUS, NULL));
107 CHECK_RESULT(MDOaddvar(m, 0, NULL, NULL, 0.0, 0.0, MDO_INFINITY, MDO_CONTINUOUS, NULL));
- 我们再利用
MDOaddpsdvar()
来创建两个半正定矩阵变量,并同时设定目标函数中其对应的矩阵系数: 第一个参数为 模型,
第二个参数为 矩阵变量的维数 (此处为3),
第三个元素代表着 矩阵系数的非零元素个数,
第四个和第五个参数分别代表 矩阵系数中非零元素对应的向量化索引与取值。
109 /* Add psd variable and corresponding objective coefficients. */
110 CHECK_RESULT(MDOaddpsdvar(m, dim_mats[0], C0_size, C0_nnz_indices, C0_nnz_values, NULL));
111 CHECK_RESULT(MDOaddpsdvar(m, dim_mats[1], C1_size, C1_nnz_indices, C1_nnz_values, NULL));
- 接着,我们输入约束。我们利用
MDOaddpsdconstr()
创建带半正定矩阵变量的约束: 第一个参数为 模型,
第二个参数为该约束中 线性项的个数 (此处为0),
第三个参数为该约束中 矩阵项的个数,
第四个参数代表 所有矩阵系数的总非零元素的个数,
第五个和第六个参数代表该约束中 线性项对应的索引与取值,
第七个参数代表该约束中 矩阵变量的索引,
第八个参数代表该约束中 各矩阵系数的初始索引,
第九个和第十个参数代表 所有矩阵系数的非零元素的向量化索引与取值,
第十一个参数代表该 约束的意义,MDO_EQUAL(‘=’)表示等式约束,
第十二个参数代表 约束的右侧值,
最后一个参数代表该约束中 名称。
113 /* Add constraints coefficients. */
114 CHECK_RESULT(MDOaddpsdconstr(m, row0_size, num_psd_var[0], A0_size, row0_idx, row0_val, A0_psd_var_indices, A0_psd_coeff_start_indices, A0_nnz_indices, A0_nnz_values, MDO_EQUAL, 1.0, NULL));
115 CHECK_RESULT(MDOaddpsdconstr(m, row1_size, num_psd_var[1], A1_size, row1_idx, row1_val, A1_psd_var_indices, A1_psd_coeff_start_indices, A1_nnz_indices, A1_nnz_values, MDO_EQUAL, 2.0, NULL));
第三步:求解SDP模型
模型输入后,我们接着用 MDOoptimize()
来求解问题。
120 /* Solve the problem. */
121 CHECK_RESULT(MDOoptimize(m));
第四步: 取得SDP模型的解
我们利用 MDOgetintattr()
得到求解状态,即 Status 属性;利用泛用函数 MDOgetdblattr()
来取得最优的的目标函数值,即 PrimalObjVal 属性。
122 CHECK_RESULT(MDOgetintattr(m, STATUS, &status));
123 if (status == MDO_OPTIMAL)
124 {
125 CHECK_RESULT(MDOgetdblattr(m, OBJ_VAL, &obj));
126 printf("The objective value is %f\n", obj);
127 }
128 else
129 {
130 printf("No feasible solution exists\n");
131 }
第五步: 释放模型
我们利用 MDOfreemodel()
和 MDOfreeenv()
来释放模型。
36#define RELEASE_MEMORY \ 37 MDOfreemodel(m); \ 38 MDOfreeenv(env);137 RELEASE_MEMORY;
示例 mdo_sdo_ex2.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 * Subject To
13 * c0 : tr(A00 X0) + 1 x0 = 1
14 * c1 : tr(A11 X1) + 1 x1 = 2
15 * Bounds
16 * 0 <= x0
17 * 0 <= x1
18 * X0, X1 are p.s.d.
19 *
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
30#include <stdio.h>
31#include <stdlib.h>
32#include "Mindopt.h"
33
34/* Macro to check the return code. */
35#define CHECK_RESULT(code) { int res = code; if (res != 0) { fprintf(stderr, "Bad code: %d\n", res); exit(res); } }
36#define RELEASE_MEMORY \
37 MDOfreemodel(m); \
38 MDOfreeenv(env);
39#define MODEL_NAME "SDP_02"
40#define MODEL_SENSE "ModelSense"
41#define STATUS "Status"
42#define OBJ_VAL "PrimalObjVal"
43
44int main(void)
45{
46 /* Variables. */
47 MDOenv *env;
48 MDOmodel *m;
49 double obj;
50 int status;
51
52 /* Model data. */
53 // linear variables.
54 int num_cols = 2;
55
56 // psd variables.
57 int num_mats = 2;
58 int dim_mats[] = { 2, 3 };
59
60 // Objective coefficients. (Only the right upper part in row-major sparse format)
61 int C0_size = 3;
62 int C0_nnz_indices[] = { 0, 1, 3 };
63 double C0_nnz_values[] = { 2.0, 1.0, 2.0 };
64
65 int C1_size = 4;
66 int C1_nnz_indices[] = { 0, 2, 4, 8 };
67 double C1_nnz_values[] = { 3.0, 1.0, 2.0, 3.0 };
68
69 // Constraints (linear).
70 int row0_size = 1;
71 int row0_idx[] = { 0 };
72 double row0_val[] = { 1.0 };
73 int row1_size = 1;
74 int row1_idx[] = { 1 };
75 double row1_val[] = { 1.0 };
76
77 // Constraints (psd). (Only the right upper part in row-major sparse format)
78 int num_psd_var[] = { 1, 1 };
79 int A0_psd_var_indices[] = { 0 };
80 int A0_psd_coeff_start_indices[] = { 0 };
81 int A0_size = 3;
82 int A0_nnz_indices[] = { 0, 1, 3 };
83 double A0_nnz_values[] = { 3.0, 1.0, 3.0 };
84
85 int A1_psd_var_indices[] = { 1 };
86 int A1_psd_coeff_start_indices[] = { 0 };
87 int A1_size = 4;
88 int A1_nnz_indices[] = { 0, 2, 4, 8 };
89 double A1_nnz_values[] = { 3.0, 1.0, 4.0, 5.0 };
90
91 /*------------------------------------------------------------------*/
92 /* Step 1. Create a model and change the parameters. */
93 /*------------------------------------------------------------------*/
94 /* Create an empty model. */
95 CHECK_RESULT(MDOemptyenv(&env));
96 CHECK_RESULT(MDOstartenv(env));
97 CHECK_RESULT(MDOnewmodel(env, &m, MODEL_NAME, 0, NULL, NULL, NULL, NULL, NULL));
98
99 /*------------------------------------------------------------------*/
100 /* Step 2. Input model. */
101 /*------------------------------------------------------------------*/
102 /* Change to maximization problem. */
103 CHECK_RESULT(MDOsetintattr(m, MODEL_SENSE, MDO_MAXIMIZE));
104
105 /* Add linear variable and corresponding objective coefficients. */
106 CHECK_RESULT(MDOaddvar(m, 0, NULL, NULL, 0.0, 0.0, MDO_INFINITY, MDO_CONTINUOUS, NULL));
107 CHECK_RESULT(MDOaddvar(m, 0, NULL, NULL, 0.0, 0.0, MDO_INFINITY, MDO_CONTINUOUS, NULL));
108
109 /* Add psd variable and corresponding objective coefficients. */
110 CHECK_RESULT(MDOaddpsdvar(m, dim_mats[0], C0_size, C0_nnz_indices, C0_nnz_values, NULL));
111 CHECK_RESULT(MDOaddpsdvar(m, dim_mats[1], C1_size, C1_nnz_indices, C1_nnz_values, NULL));
112
113 /* Add constraints coefficients. */
114 CHECK_RESULT(MDOaddpsdconstr(m, row0_size, num_psd_var[0], A0_size, row0_idx, row0_val, A0_psd_var_indices, A0_psd_coeff_start_indices, A0_nnz_indices, A0_nnz_values, MDO_EQUAL, 1.0, NULL));
115 CHECK_RESULT(MDOaddpsdconstr(m, row1_size, num_psd_var[1], A1_size, row1_idx, row1_val, A1_psd_var_indices, A1_psd_coeff_start_indices, A1_nnz_indices, A1_nnz_values, MDO_EQUAL, 2.0, NULL));
116
117 /*------------------------------------------------------------------*/
118 /* Step 3. Solve the problem and populate the result. */
119 /*------------------------------------------------------------------*/
120 /* Solve the problem. */
121 CHECK_RESULT(MDOoptimize(m));
122 CHECK_RESULT(MDOgetintattr(m, STATUS, &status));
123 if (status == MDO_OPTIMAL)
124 {
125 CHECK_RESULT(MDOgetdblattr(m, OBJ_VAL, &obj));
126 printf("The objective value is %f\n", obj);
127 }
128 else
129 {
130 printf("No feasible solution exists\n");
131 }
132
133 /*------------------------------------------------------------------*/
134 /* Step 4. Free the model. */
135 /*------------------------------------------------------------------*/
136 /* Free the model. */
137 RELEASE_MEMORY;
138
139 return 0;
140}