lemon-project-template-glpk: 33de93886c88 deps/glpk/src/glpapi20.c

lemon-project-template-glpk

view deps/glpk/src/glpapi20.c @ 9:33de93886c88

Import GLPK 4.47

author	Alpar Juttner <alpar@cs.elte.hu>
date	Sun, 06 Nov 2011 20:59:10 +0100
parents
children

line source

1 /* glpapi20.c */

3 /***********************************************************************

4 * This code is part of GLPK (GNU Linear Programming Kit).

5 *

7 * 2009, 2010, 2011 Andrew Makhorin, Department for Applied Informatics,

9 * E-mail: <mao@gnu.org>.

10 *

11 * GLPK is free software: you can redistribute it and/or modify it

12 * under the terms of the GNU General Public License as published by

13 * the Free Software Foundation, either version 3 of the License, or

14 * (at your option) any later version.

15 *

16 * GLPK is distributed in the hope that it will be useful, but WITHOUT

17 * ANY WARRANTY; without even the implied warranty of MERCHANTABILITY

18 * or FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public

19 * License for more details.

20 *

21 * You should have received a copy of the GNU General Public License

22 * along with GLPK. If not, see <http://www.gnu.org/licenses/>.

23 ***********************************************************************/

25 #include "glpnpp.h"

27 int glp_intfeas1(glp_prob *P, int use_bound, int obj_bound)

28 { /* solve integer feasibility problem */

29 NPP *npp = NULL;

30 glp_prob *mip = NULL;

31 int *obj_ind = NULL;

32 double *obj_val = NULL;

33 int obj_row = 0;

34 int i, j, k, obj_len, temp, ret;

35 /* check the problem object */

36 if (P == NULL || P->magic != GLP_PROB_MAGIC)

37 xerror("glp_intfeas1: P = %p; invalid problem object\n",

38 P);

39 if (P->tree != NULL)

40 xerror("glp_intfeas1: operation not allowed\n");

41 /* integer solution is currently undefined */

42 P->mip_stat = GLP_UNDEF;

43 P->mip_obj = 0.0;

44 /* check columns (variables) */

45 for (j = 1; j <= P->n; j++)

46 { GLPCOL *col = P->col[j];

47 #if 0 /* currently binarization is not yet implemented */

48 if (!(col->kind == GLP_IV || col->type == GLP_FX))

49 { xprintf("glp_intfeas1: column %d: non-integer non-fixed var"

50 "iable not allowed\n", j);

51 #else

52 if (!((col->kind == GLP_IV && col->lb == 0.0 && col->ub == 1.0)

53 || col->type == GLP_FX))

54 { xprintf("glp_intfeas1: column %d: non-binary non-fixed vari"

55 "able not allowed\n", j);

56 #endif

57 ret = GLP_EDATA;

58 goto done;

59 }

60 temp = (int)col->lb;

61 if ((double)temp != col->lb)

62 { if (col->type == GLP_FX)

63 xprintf("glp_intfeas1: column %d: fixed value %g is non-"

64 "integer or out of range\n", j, col->lb);

65 else

66 xprintf("glp_intfeas1: column %d: lower bound %g is non-"

67 "integer or out of range\n", j, col->lb);

68 ret = GLP_EDATA;

69 goto done;

70 }

71 temp = (int)col->ub;

72 if ((double)temp != col->ub)

73 { xprintf("glp_intfeas1: column %d: upper bound %g is non-int"

74 "eger or out of range\n", j, col->ub);

75 ret = GLP_EDATA;

76 goto done;

77 }

78 if (col->type == GLP_DB && col->lb > col->ub)

79 { xprintf("glp_intfeas1: column %d: lower bound %g is greater"

80 " than upper bound %g\n", j, col->lb, col->ub);

81 ret = GLP_EBOUND;

82 goto done;

83 }

84 }

85 /* check rows (constraints) */

86 for (i = 1; i <= P->m; i++)

87 { GLPROW *row = P->row[i];

88 GLPAIJ *aij;

89 for (aij = row->ptr; aij != NULL; aij = aij->r_next)

90 { temp = (int)aij->val;

91 if ((double)temp != aij->val)

92 { xprintf("glp_intfeas1: row = %d, column %d: constraint c"

93 "oefficient %g is non-integer or out of range\n",

94 i, aij->col->j, aij->val);

95 ret = GLP_EDATA;

96 goto done;

97 }

98 }

99 temp = (int)row->lb;

100 if ((double)temp != row->lb)

101 { if (row->type == GLP_FX)

102 xprintf("glp_intfeas1: row = %d: fixed value %g is non-i"

103 "nteger or out of range\n", i, row->lb);

104 else

105 xprintf("glp_intfeas1: row = %d: lower bound %g is non-i"

106 "nteger or out of range\n", i, row->lb);

107 ret = GLP_EDATA;

108 goto done;

109 }

110 temp = (int)row->ub;

111 if ((double)temp != row->ub)

112 { xprintf("glp_intfeas1: row = %d: upper bound %g is non-inte"

113 "ger or out of range\n", i, row->ub);

114 ret = GLP_EDATA;

115 goto done;

116 }

117 if (row->type == GLP_DB && row->lb > row->ub)

118 { xprintf("glp_intfeas1: row %d: lower bound %g is greater th"

119 "an upper bound %g\n", i, row->lb, row->ub);

120 ret = GLP_EBOUND;

121 goto done;

122 }

123 }

124 /* check the objective function */

125 temp = (int)P->c0;

126 if ((double)temp != P->c0)

127 { xprintf("glp_intfeas1: objective constant term %g is non-integ"

128 "er or out of range\n", P->c0);

129 ret = GLP_EDATA;

130 goto done;

131 }

132 for (j = 1; j <= P->n; j++)

133 { temp = (int)P->col[j]->coef;

134 if ((double)temp != P->col[j]->coef)

135 { xprintf("glp_intfeas1: column %d: objective coefficient is "

136 "non-integer or out of range\n", j, P->col[j]->coef);

137 ret = GLP_EDATA;

138 goto done;

139 }

140 }

141 /* save the objective function and set it to zero */

142 obj_ind = xcalloc(1+P->n, sizeof(int));

143 obj_val = xcalloc(1+P->n, sizeof(double));

144 obj_len = 0;

145 obj_ind[0] = 0;

146 obj_val[0] = P->c0;

147 P->c0 = 0.0;

148 for (j = 1; j <= P->n; j++)

149 { if (P->col[j]->coef != 0.0)

150 { obj_len++;

151 obj_ind[obj_len] = j;

152 obj_val[obj_len] = P->col[j]->coef;

153 P->col[j]->coef = 0.0;

154 }

155 }

156 /* add inequality to bound the objective function, if required */

157 if (!use_bound)

158 xprintf("Will search for ANY feasible solution\n");

159 else

160 { xprintf("Will search only for solution not worse than %d\n",

161 obj_bound);

162 obj_row = glp_add_rows(P, 1);

163 glp_set_mat_row(P, obj_row, obj_len, obj_ind, obj_val);

164 if (P->dir == GLP_MIN)

165 glp_set_row_bnds(P, obj_row,

166 GLP_UP, 0.0, (double)obj_bound - obj_val[0]);

167 else if (P->dir == GLP_MAX)

168 glp_set_row_bnds(P, obj_row,

169 GLP_LO, (double)obj_bound - obj_val[0], 0.0);

170 else

171 xassert(P != P);

172 }

173 /* create preprocessor workspace */

174 xprintf("Translating to CNF-SAT...\n");

175 xprintf("Original problem has %d row%s, %d column%s, and %d non-z"

176 "ero%s\n", P->m, P->m == 1 ? "" : "s", P->n, P->n == 1 ? "" :

177 "s", P->nnz, P->nnz == 1 ? "" : "s");

178 npp = npp_create_wksp();

179 /* load the original problem into the preprocessor workspace */

180 npp_load_prob(npp, P, GLP_OFF, GLP_MIP, GLP_OFF);

181 /* perform translation to SAT-CNF problem instance */

182 ret = npp_sat_encode_prob(npp);

183 if (ret == 0)

184 ;

185 else if (ret == GLP_ENOPFS)

186 xprintf("PROBLEM HAS NO INTEGER FEASIBLE SOLUTION\n");

187 else if (ret == GLP_ERANGE)

188 xprintf("glp_intfeas1: translation to SAT-CNF failed because o"

189 "f integer overflow\n");

190 else

191 xassert(ret != ret);

192 if (ret != 0)

193 goto done;

194 /* build SAT-CNF problem instance and try to solve it */

195 mip = glp_create_prob();

196 npp_build_prob(npp, mip);

197 ret = glp_minisat1(mip);

198 /* only integer feasible solution can be postprocessed */

199 if (!(mip->mip_stat == GLP_OPT || mip->mip_stat == GLP_FEAS))

200 { P->mip_stat = mip->mip_stat;

201 goto done;

202 }

203 /* postprocess the solution found */

204 npp_postprocess(npp, mip);

205 /* the transformed problem is no longer needed */

206 glp_delete_prob(mip), mip = NULL;

207 /* store solution to the original problem object */

208 npp_unload_sol(npp, P);

209 /* change the solution status to 'integer feasible' */

210 P->mip_stat = GLP_FEAS;

211 /* check integer feasibility */

212 for (i = 1; i <= P->m; i++)

213 { GLPROW *row;

214 GLPAIJ *aij;

215 double sum;

216 row = P->row[i];

217 sum = 0.0;

218 for (aij = row->ptr; aij != NULL; aij = aij->r_next)

219 sum += aij->val * aij->col->mipx;

220 xassert(sum == row->mipx);

221 if (row->type == GLP_LO || row->type == GLP_DB ||

222 row->type == GLP_FX)

223 xassert(sum >= row->lb);

224 if (row->type == GLP_UP || row->type == GLP_DB ||

225 row->type == GLP_FX)

226 xassert(sum <= row->ub);

227 }

228 /* compute value of the original objective function */

229 P->mip_obj = obj_val[0];

230 for (k = 1; k <= obj_len; k++)

231 P->mip_obj += obj_val[k] * P->col[obj_ind[k]]->mipx;

232 xprintf("Objective value = %17.9e\n", P->mip_obj);

233 done: /* delete the transformed problem, if it exists */

234 if (mip != NULL)

235 glp_delete_prob(mip);

236 /* delete the preprocessor workspace, if it exists */

237 if (npp != NULL)

238 npp_delete_wksp(npp);

239 /* remove inequality used to bound the objective function */

240 if (obj_row > 0)

241 { int ind[1+1];

242 ind[1] = obj_row;

243 glp_del_rows(P, 1, ind);

244 }

245 /* restore the original objective function */

246 if (obj_ind != NULL)

247 { P->c0 = obj_val[0];

248 for (k = 1; k <= obj_len; k++)

249 P->col[obj_ind[k]]->coef = obj_val[k];

250 xfree(obj_ind);

251 xfree(obj_val);

252 }

253 return ret;

254 }

255

256 /* eof */