5.5.2. C 的MIQP建模和优化

在本节中,我们将使用 MindOpt C API,以按行输入的形式来建模以及求解 MIQP题示例 中的问题。

首先,引入头文件:

29#include "Mindopt.h"

创建优化模型:

93    /*------------------------------------------------------------------*/
94    CHECK_RESULT(MDOemptyenv(&env));
95    CHECK_RESULT(MDOstartenv(env));

接下来,我们通过 MDOsetintattr() 将目标函数设置为 最小化,并调用 MDOaddvar() 来添加四个优化变量。(更多API和使用方式,请参考 C API):

105    /* Add variables. */
106    CHECK_RESULT(MDOaddvar(model, 0, NULL, NULL, 1.0, 0,         10.0, MDO_INTEGER, "x0"));
107    CHECK_RESULT(MDOaddvar(model, 0, NULL, NULL, 2.0, 0, MDO_INFINITY, MDO_CONTINUOUS, "x1"));
108    CHECK_RESULT(MDOaddvar(model, 0, NULL, NULL, 1.0, 0, MDO_INFINITY, MDO_CONTINUOUS, "x2"));

备注

在函数 MDOaddvar() 中一个参数位是 vtype,设置为 MDO_INTEGER 代表这个变量是整形变量。

接下来我们将添加二次规划中的目标函数的二次项系数。我们使用以下三列数组来定义:其中 qo_col1qo_col2 分别记录二次项中所有非零项的两个变量索引,而 qo_values 是与之相对应的非零系数值。

68     */
69    int qo_col1[] = 
70    {
71        0, 
72        1,   1,
73                  2,
74                       3  
75    };
76    int qo_col2[] =
77    {
78        0,
79        0,   1,
80                  2,
81                       3
82    };
83    double qo_values[] =
84    {
85        1.0,
86        0.5, 1.0,
87                  1.0, 
88                       1.0

我们调用 MDOaddqpterms() 设置目标的二次项:

110    /* Add constraints.

调用 MDOaddconstr() 来输入约束:

116    /* Add quadratic objective term. */

问题输入完成后,再调用 MDOoptimize() 求解优化问题。

124    CHECK_RESULT(MDOoptimize(model));

然后,我们可以通过获取属性值的方式来获取求解后的最优值 (optimal value) 和最优解 (optimal solution).

128    /* For MIP(MILP,MIQP, MIQCP) problems, if the solving process       */
129    /* terminates early due to reasons such as timeout or interruption, */
130    /* the model status will indicate termination by timeout (or        */
131    /* interruption, etc.). However, suboptimal solutions may still     */
132    /* exist, making it necessary to check the SolCount property.       */
133    /*------------------------------------------------------------------*/
134    CHECK_RESULT(MDOgetintattr(m, STATUS, &status));
135    CHECK_RESULT(MDOgetintattr(m, SOL_COUNT, &solcount));
136    if (status == MDO_OPTIMAL || status == MDO_SUB_OPTIMAL || solcount != 0)
137    {

最后,调用 MDOfreemodel()MDOfreeenv() 来释放模型:

30/* Macro to check the return code */
31#define RELEASE_MEMORY  \
146    else 

示例 MdoMIQPEx1.c 提供了完整源代码:

  1/**
  2 *  Description
  3 *  -----------
  4 *
  5 *  Mixed Integer Quadratic optimization (row-wise input).
  6 *
  7 *  Formulation
  8
  9 *  -----------
 10 *
 11 *  Minimize
 12 *    obj: 1 x0 + 1 x1 + 1 x2 + 1 x3 
 13 *         + 1/2 [ x0^2 + x1^2 + x2^2 + x3^2 + x0 x1]
 14 *  Subject To
 15 *   c0 : 1 x0 + 1 x1 + 2 x2 + 3 x3 >= 1
 16 *   c1 : 1 x0 - 1 x2 + 6 x3 = 1
 17 *  Bounds
 18 *    0 <= x0 <= 10
 19 *    0 <= x1
 20 *    0 <= x2
 21 *    0 <= x3
 22 *  Integers
 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  "MIQCP_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    /* Model data. */
 52    int    row1_idx[] = { 0,   1,   2,   3   };
 53    double row1_val[] = { 1.0, 1.0, 2.0, 3.0 };
 54    int    row2_idx[] = { 0,    2,   3   };
 55    double row2_val[] = { 1.0, -1.0, 6.0 };
 56
 57    /* Quadratic objective matrix Q.
 58     * 
 59     *  Note.
 60     *  1. The objective function is defined as c^Tx + 1/2 x^TQx, where Q is stored with coordinate format.
 61     *  2. Q will be scaled by 1/2 internally.
 62     *  3. To ensure the symmetricity of Q, user needs to input only the lower triangular part.
 63     * 
 64     * Q = [ 1.0  0.5  0    0   ]
 65     *     [ 0.5  1.0  0    0   ]
 66     *     [ 0.0  0.0  1.0  0   ]
 67     *     [ 0    0    0    1.0 ]
 68     */
 69    int qo_col1[] = 
 70    {
 71        0, 
 72        1,   1,
 73                  2,
 74                       3  
 75    };
 76    int qo_col2[] =
 77    {
 78        0,
 79        0,   1,
 80                  2,
 81                       3
 82    };
 83    double qo_values[] =
 84    {
 85        1.0,
 86        0.5, 1.0,
 87                  1.0, 
 88                       1.0
 89    };
 90
 91     /*------------------------------------------------------------------*/
 92    /* Step 1. Create environment and model.                            */
 93    /*------------------------------------------------------------------*/
 94    CHECK_RESULT(MDOemptyenv(&env));
 95    CHECK_RESULT(MDOstartenv(env));
 96    CHECK_RESULT(MDOnewmodel(env, &model, MODEL_NAME, 0, NULL, NULL, NULL, NULL, NULL));
 97
 98
 99    /*------------------------------------------------------------------*/
100    /* Step 2. Input model.                                             */
101    /*------------------------------------------------------------------*/
102    /* Change to minimization problem. */
103    CHECK_RESULT(MDOsetintattr(model, MODEL_SENSE, MDO_MINIMIZE));
104
105    /* Add variables. */
106    CHECK_RESULT(MDOaddvar(model, 0, NULL, NULL, 1.0, 0,         10.0, MDO_INTEGER, "x0"));
107    CHECK_RESULT(MDOaddvar(model, 0, NULL, NULL, 2.0, 0, MDO_INFINITY, MDO_CONTINUOUS, "x1"));
108    CHECK_RESULT(MDOaddvar(model, 0, NULL, NULL, 1.0, 0, MDO_INFINITY, MDO_CONTINUOUS, "x2"));
109    CHECK_RESULT(MDOaddvar(model, 0, NULL, NULL, 1.0, 0, MDO_INFINITY, MDO_CONTINUOUS, "x3"));
110
111    /* Add constraints.
112     * Note that the nonzero elements are inputted in a row-wise order here.
113     */
114    CHECK_RESULT(MDOaddconstr(model, 4, row1_idx, row1_val, MDO_GREATER_EQUAL, 1.0, "c0"));
115    CHECK_RESULT(MDOaddconstr(model, 3, row2_idx, row2_val, MDO_EQUAL,         1.0, "c1"));
116
117    /* Add quadratic objective term. */
118    CHECK_RESULT(MDOaddqpterms(model, 5, qo_col1, qo_col2, qo_values));
119    
120    /*------------------------------------------------------------------*/
121    /* Step 3. Solve the problem and populate optimization result.                */
122    /*------------------------------------------------------------------*/
123    /* Solve the problem. */
124    CHECK_RESULT(MDOoptimize(model));
125
126    /*------------------------------------------------------------------*/
127    /* Step 4. Retrive model status and objective.                      */
128    /* For MIP(MILP,MIQP, MIQCP) problems, if the solving process       */
129    /* terminates early due to reasons such as timeout or interruption, */
130    /* the model status will indicate termination by timeout (or        */
131    /* interruption, etc.). However, suboptimal solutions may still     */
132    /* exist, making it necessary to check the SolCount property.       */
133    /*------------------------------------------------------------------*/
134    CHECK_RESULT(MDOgetintattr(m, STATUS, &status));
135    CHECK_RESULT(MDOgetintattr(m, SOL_COUNT, &solcount));
136    if (status == MDO_OPTIMAL || status == MDO_SUB_OPTIMAL || solcount != 0)
137    {
138        CHECK_RESULT(MDOgetdblattr(model, OBJ_VAL, &obj));
139        printf("The optimal objective value is: %f\n", obj);
140        for (int i = 0; i < 4; ++i) 
141        {
142            CHECK_RESULT(MDOgetdblattrelement(model, X, i, &x));
143            printf("x[%d] = %f\n", i, x);
144        }
145    } 
146    else 
147    {
148        printf("No feasible solution.\n");
149    }
150 
151    /*------------------------------------------------------------------*/
152    /* Step 4. Free the model.                                          */
153    /*------------------------------------------------------------------*/
154    RELEASE_MEMORY;
155       
156    return 0;
157}