5.6.2. C 的MIQCP建模和优化¶
在本节中,我们将使用 MindOpt C API,以按行输入的形式来建模以及求解 MIQCP题示例 中的问题。
首先,引入头文件:
27#include "Mindopt.h"
并创建优化模型:
78 CHECK_RESULT(MDOemptyenv(&env));
79 CHECK_RESULT(MDOstartenv(env));
80 CHECK_RESULT(MDOnewmodel(env, &model, MODEL_NAME, 0, NULL, NULL, NULL, NULL, NULL));
接下来,我们通过 MDOsetIntAttr()
将目标函数设置为 最小化,并调用 MDOaddvar()
来添加四个优化变量,定义其下界、上界、名称和类型,以及其在目标函数中线性项的系数(关于 MDOsetIntAttr()
和 MDOaddvar()
的详细使用方式,请参考 属性):
86 /* Change to minimization problem. */
87 CHECK_RESULT(MDOsetintattr(model, MODEL_SENSE, MDO_MINIMIZE));
88
89 /* Add variables. */
90 CHECK_RESULT(MDOaddvar(model, 0, NULL, NULL, 1.0, 0, 10.0, MDO_INTEGER, "x0"));
91 CHECK_RESULT(MDOaddvar(model, 0, NULL, NULL, 1.0, 0, MDO_INFINITY, MDO_CONTINUOUS, "x1"));
92 CHECK_RESULT(MDOaddvar(model, 0, NULL, NULL, 1.0, 0, MDO_INFINITY, MDO_CONTINUOUS, "x2"));
93 CHECK_RESULT(MDOaddvar(model, 0, NULL, NULL, 1.0, 0, MDO_INFINITY, MDO_CONTINUOUS, "x3"));
Note
在函数 MDOaddvar()
中一个参数位是 vtype
,设置为 MDO_INTEGER
代表这个变量是整形变量。
矩阵的非零元随后将按 列 输入;因此,MDOaddvar()
中,与约束矩阵相关联的参数 size
、 indices
、 value
分别用 0
、 NULL
、 NULL
代替(换句话说,此时问题无约束)。
接下来我们将添加二次规划中的目标函数的二次项系数。我们使用以下三列数组来定义:其中 qo_col1
和 qo_col2
分别记录二次项中所有非零项的两个变量索引,而 qo_values
是与之相对应的非零系数值。
49 /* Quadratic part in objective: 1/2 [ x0^2 + x1^2 + x2^2 + x3^2 + x0 x1] */
50 int qo_nnz = 5;
51 int qo_col1[] = { 0, 1, 2, 3, 0 };
52 int qo_col2[] = { 0, 1, 2, 3, 1 };
53 double qo_values[] = { 0.5, 0.5, 0.5, 0.5, 0.5 };
我们调用 MDOaddqpterms()
设置目标的二次项:
95 /* Add quadratic objective term. */
96 CHECK_RESULT(MDOaddqpterms(model, qo_nnz, qo_col1, qo_col2, qo_values));
以下我们将开始添加模型中的两个二次约束。其中的线性项与添加线性约束的方式类似,我们使用以下数组来定义线性项;其中, row1_idx
和 row2_idx
分别表示第一和第二个约束中非零元素的位置(索引),而 row1_val
和 row2_val
则是与之相对应的非零数值。
55 /* Linear part in the first constraint: 1 x0 + 1 x1 + 2 x2 + 3 x3 */
56 int row1_nnz = 4;
57 int row1_idx[] = { 0, 1, 2, 3 };
58 double row1_val[] = { 1.0, 1.0, 2.0, 3.0 };
65 /* Linear part in the second constraint: 1 x0 - 1 x2 + 6 x3 */
66 int row2_nnz = 3;
67 int row2_idx[] = { 0, 2, 3 };
68 double row2_val[] = { 1.0, -1.0, 6.0 };
二次约束中的二次项与目标函数中的二次项类似:第一个二次约束中我们用 qc1_col1
与 qc1_col2
分别记录二次项中所有非零项的两个变量索引,而 qc1_values
是与之相对应的非零系数值。第二个二次约束类似。
59 /* Quadratic part in the first constraint: - 1/2 [ x0^2 + x1^2 + x2^2 + x3^2 + x0 x1] */
60 int qc1_nnz = 5;
61 int qc1_col1[] = { 0, 1, 2, 3, 0 };
62 int qc1_col2[] = { 0, 1, 2, 3, 1 };
63 double qc1_values[] = { -0.5, -0.5, -0.5, -0.5, -0.5 };
69 /* Quadratic part in the second constraint: 1/2 [x1^2] */
70 int qc2_nnz = 1;
71 int qc2_col1[] = { 1 };
72 int qc2_col2[] = { 1 };
73 double qc2_values[] = { 0.5 };
我们调用 MDOaddqconstr()
来添加二次约束:
98 /* Add quadratic constraints. */
99 CHECK_RESULT(MDOaddqconstr(model, row1_nnz, row1_idx, row1_val, qc1_nnz, qc1_col1, qc1_col2, qc1_values, MDO_GREATER_EQUAL, 1.0, "c0"));
100 CHECK_RESULT(MDOaddqconstr(model, row2_nnz, row2_idx, row2_val, qc2_nnz, qc2_col1, qc2_col2, qc2_values, MDO_LESS_EQUAL, 1.0, "c1"));
问题输入完成后,再调用 MDOoptimize()
求解优化问题:
105 /* Solve the problem. */
106 CHECK_RESULT(MDOoptimize(model));
然后,我们可以通过获取属性值的方式来获取对应的优化目标值objective和变量的取值:
108 CHECK_RESULT(MDOgetintattr(model, STATUS, &status));
109 if (status == MDO_OPTIMAL)
110 {
111 CHECK_RESULT(MDOgetdblattr(model, OBJ_VAL, &obj));
112 printf("The optimal objective value is: %f\n", obj);
113 for (int i = 0; i < 4; ++i)
114 {
115 CHECK_RESULT(MDOgetdblattrelement(model, X, i, &x));
116 printf("x[%d] = %f\n", i, x);
117 }
118 }
119 else
120 {
121 printf("No feasible solution.\n");
122 }
最后,调用 MDOfreemodel()
和 MDOfreeenv()
来释放模型:
30#define RELEASE_MEMORY \
31 MDOfreemodel(model); \
32 MDOfreeenv(env);
127 RELEASE_MEMORY;
示例 MdoMIQCPEx1.c 提供了完整源代码:
1/**
2 * Description
3 * -----------
4 * Mixed Integer Quadratically constrained quadratic optimization (row-wise input).
5 * Formulation
6 * -----------
7 *
8 * Minimize
9 * obj: 1 x0 + 1 x1 + 1 x2 + 1 x3
10 * + 1/2 [ x0^2 + x1^2 + x2^2 + x3^2 + x0 x1]
11 *
12 * Subject To
13 * c0 : 1 x0 + 1 x1 + 2 x2 + 3 x3 - 1/2 [ x0^2 + x1^2 + x2^2 + x3^2 + x0 x1] >= 1
14 * c1 : 1 x0 - 1 x2 + 6 x3 + 1/2 [x1^2] <= 1
15 * Bounds
16 * 0 <= x0 <= 10
17 * 0 <= x1
18 * 0 <= x2
19 * 0 <= x3
20 * Integer
21 * x0
22 * End
23 */
24
25#include <stdio.h>
26#include <stdlib.h>
27#include "Mindopt.h"
28
29/* Macro to check the return code */
30#define RELEASE_MEMORY \
31 MDOfreemodel(model); \
32 MDOfreeenv(env);
33#define CHECK_RESULT(code) { int res = code; if (res != 0) { fprintf(stderr, "Bad code: %d\n", res); RELEASE_MEMORY; return (res); } }
34#define MODEL_NAME "QCP_01"
35#define MODEL_SENSE "ModelSense"
36#define STATUS "Status"
37#define OBJ_VAL "ObjVal"
38#define X "X"
39
40int main(void)
41{
42 /* Variables. */
43 MDOenv *env;
44 MDOmodel *model;
45 double obj, x;
46 int status, i;
47
48 /* Prepare model data. */
49 /* Quadratic part in objective: 1/2 [ x0^2 + x1^2 + x2^2 + x3^2 + x0 x1] */
50 int qo_nnz = 5;
51 int qo_col1[] = { 0, 1, 2, 3, 0 };
52 int qo_col2[] = { 0, 1, 2, 3, 1 };
53 double qo_values[] = { 0.5, 0.5, 0.5, 0.5, 0.5 };
54
55 /* Linear part in the first constraint: 1 x0 + 1 x1 + 2 x2 + 3 x3 */
56 int row1_nnz = 4;
57 int row1_idx[] = { 0, 1, 2, 3 };
58 double row1_val[] = { 1.0, 1.0, 2.0, 3.0 };
59 /* Quadratic part in the first constraint: - 1/2 [ x0^2 + x1^2 + x2^2 + x3^2 + x0 x1] */
60 int qc1_nnz = 5;
61 int qc1_col1[] = { 0, 1, 2, 3, 0 };
62 int qc1_col2[] = { 0, 1, 2, 3, 1 };
63 double qc1_values[] = { -0.5, -0.5, -0.5, -0.5, -0.5 };
64
65 /* Linear part in the second constraint: 1 x0 - 1 x2 + 6 x3 */
66 int row2_nnz = 3;
67 int row2_idx[] = { 0, 2, 3 };
68 double row2_val[] = { 1.0, -1.0, 6.0 };
69 /* Quadratic part in the second constraint: 1/2 [x1^2] */
70 int qc2_nnz = 1;
71 int qc2_col1[] = { 1 };
72 int qc2_col2[] = { 1 };
73 double qc2_values[] = { 0.5 };
74
75 /*------------------------------------------------------------------*/
76 /* Step 1. Create environment and model. */
77 /*------------------------------------------------------------------*/
78 CHECK_RESULT(MDOemptyenv(&env));
79 CHECK_RESULT(MDOstartenv(env));
80 CHECK_RESULT(MDOnewmodel(env, &model, MODEL_NAME, 0, NULL, NULL, NULL, NULL, NULL));
81
82
83 /*------------------------------------------------------------------*/
84 /* Step 2. Input model. */
85 /*------------------------------------------------------------------*/
86 /* Change to minimization problem. */
87 CHECK_RESULT(MDOsetintattr(model, MODEL_SENSE, MDO_MINIMIZE));
88
89 /* Add variables. */
90 CHECK_RESULT(MDOaddvar(model, 0, NULL, NULL, 1.0, 0, 10.0, MDO_INTEGER, "x0"));
91 CHECK_RESULT(MDOaddvar(model, 0, NULL, NULL, 1.0, 0, MDO_INFINITY, MDO_CONTINUOUS, "x1"));
92 CHECK_RESULT(MDOaddvar(model, 0, NULL, NULL, 1.0, 0, MDO_INFINITY, MDO_CONTINUOUS, "x2"));
93 CHECK_RESULT(MDOaddvar(model, 0, NULL, NULL, 1.0, 0, MDO_INFINITY, MDO_CONTINUOUS, "x3"));
94
95 /* Add quadratic objective term. */
96 CHECK_RESULT(MDOaddqpterms(model, qo_nnz, qo_col1, qo_col2, qo_values));
97
98 /* Add quadratic constraints. */
99 CHECK_RESULT(MDOaddqconstr(model, row1_nnz, row1_idx, row1_val, qc1_nnz, qc1_col1, qc1_col2, qc1_values, MDO_GREATER_EQUAL, 1.0, "c0"));
100 CHECK_RESULT(MDOaddqconstr(model, row2_nnz, row2_idx, row2_val, qc2_nnz, qc2_col1, qc2_col2, qc2_values, MDO_LESS_EQUAL, 1.0, "c1"));
101
102 /*------------------------------------------------------------------*/
103 /* Step 3. Solve the problem and populate optimization result. */
104 /*------------------------------------------------------------------*/
105 /* Solve the problem. */
106 CHECK_RESULT(MDOoptimize(model));
107
108 CHECK_RESULT(MDOgetintattr(model, STATUS, &status));
109 if (status == MDO_OPTIMAL)
110 {
111 CHECK_RESULT(MDOgetdblattr(model, OBJ_VAL, &obj));
112 printf("The optimal objective value is: %f\n", obj);
113 for (int i = 0; i < 4; ++i)
114 {
115 CHECK_RESULT(MDOgetdblattrelement(model, X, i, &x));
116 printf("x[%d] = %f\n", i, x);
117 }
118 }
119 else
120 {
121 printf("No feasible solution.\n");
122 }
123
124 /*------------------------------------------------------------------*/
125 /* Step 4. Free the model. */
126 /*------------------------------------------------------------------*/
127 RELEASE_MEMORY;
128
129 return 0;
130}