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