lemon-project-template-glpk

annotate deps/glpk/src/glpnet03.c @ 11:4fc6ad2fb8a6

Test GLPK in src/main.cc
author Alpar Juttner <alpar@cs.elte.hu>
date Sun, 06 Nov 2011 21:43:29 +0100
parents
children
rev   line source
alpar@9 1 /* glpnet03.c (Klingman's network problem generator) */
alpar@9 2
alpar@9 3 /***********************************************************************
alpar@9 4 * This code is part of GLPK (GNU Linear Programming Kit).
alpar@9 5 *
alpar@9 6 * This code is the result of translation of the Fortran program NETGEN
alpar@9 7 * developed by Dr. Darwin Klingman, which is publically available from
alpar@9 8 * NETLIB at <http://www.netlib.org/lp/generators>.
alpar@9 9 *
alpar@9 10 * The translation was made by Andrew Makhorin <mao@gnu.org>.
alpar@9 11 *
alpar@9 12 * GLPK is free software: you can redistribute it and/or modify it
alpar@9 13 * under the terms of the GNU General Public License as published by
alpar@9 14 * the Free Software Foundation, either version 3 of the License, or
alpar@9 15 * (at your option) any later version.
alpar@9 16 *
alpar@9 17 * GLPK is distributed in the hope that it will be useful, but WITHOUT
alpar@9 18 * ANY WARRANTY; without even the implied warranty of MERCHANTABILITY
alpar@9 19 * or FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public
alpar@9 20 * License for more details.
alpar@9 21 *
alpar@9 22 * You should have received a copy of the GNU General Public License
alpar@9 23 * along with GLPK. If not, see <http://www.gnu.org/licenses/>.
alpar@9 24 ***********************************************************************/
alpar@9 25
alpar@9 26 #include "glpapi.h"
alpar@9 27
alpar@9 28 /***********************************************************************
alpar@9 29 * NAME
alpar@9 30 *
alpar@9 31 * glp_netgen - Klingman's network problem generator
alpar@9 32 *
alpar@9 33 * SYNOPSIS
alpar@9 34 *
alpar@9 35 * int glp_netgen(glp_graph *G, int v_rhs, int a_cap, int a_cost,
alpar@9 36 * const int parm[1+15]);
alpar@9 37 *
alpar@9 38 * DESCRIPTION
alpar@9 39 *
alpar@9 40 * The routine glp_netgen is a network problem generator developed by
alpar@9 41 * Dr. Darwin Klingman. It can create capacitated and uncapacitated
alpar@9 42 * minimum cost flow (or transshipment), transportation, and assignment
alpar@9 43 * problems.
alpar@9 44 *
alpar@9 45 * The parameter G specifies the graph object, to which the generated
alpar@9 46 * problem data have to be stored. Note that on entry the graph object
alpar@9 47 * is erased with the routine glp_erase_graph.
alpar@9 48 *
alpar@9 49 * The parameter v_rhs specifies an offset of the field of type double
alpar@9 50 * in the vertex data block, to which the routine stores the supply or
alpar@9 51 * demand value. If v_rhs < 0, the value is not stored.
alpar@9 52 *
alpar@9 53 * The parameter a_cap specifies an offset of the field of type double
alpar@9 54 * in the arc data block, to which the routine stores the arc capacity.
alpar@9 55 * If a_cap < 0, the capacity is not stored.
alpar@9 56 *
alpar@9 57 * The parameter a_cost specifies an offset of the field of type double
alpar@9 58 * in the arc data block, to which the routine stores the per-unit cost
alpar@9 59 * if the arc flow. If a_cost < 0, the cost is not stored.
alpar@9 60 *
alpar@9 61 * The array parm contains description of the network to be generated:
alpar@9 62 *
alpar@9 63 * parm[0] not used
alpar@9 64 * parm[1] (iseed) 8-digit positive random number seed
alpar@9 65 * parm[2] (nprob) 8-digit problem id number
alpar@9 66 * parm[3] (nodes) total number of nodes
alpar@9 67 * parm[4] (nsorc) total number of source nodes (including
alpar@9 68 * transshipment nodes)
alpar@9 69 * parm[5] (nsink) total number of sink nodes (including
alpar@9 70 * transshipment nodes)
alpar@9 71 * parm[6] (iarcs) number of arcs
alpar@9 72 * parm[7] (mincst) minimum cost for arcs
alpar@9 73 * parm[8] (maxcst) maximum cost for arcs
alpar@9 74 * parm[9] (itsup) total supply
alpar@9 75 * parm[10] (ntsorc) number of transshipment source nodes
alpar@9 76 * parm[11] (ntsink) number of transshipment sink nodes
alpar@9 77 * parm[12] (iphic) percentage of skeleton arcs to be given
alpar@9 78 * the maximum cost
alpar@9 79 * parm[13] (ipcap) percentage of arcs to be capacitated
alpar@9 80 * parm[14] (mincap) minimum upper bound for capacitated arcs
alpar@9 81 * parm[15] (maxcap) maximum upper bound for capacitated arcs
alpar@9 82 *
alpar@9 83 * The routine generates a transportation problem if:
alpar@9 84 *
alpar@9 85 * nsorc + nsink = nodes, ntsorc = 0, and ntsink = 0.
alpar@9 86 *
alpar@9 87 * The routine generates an assignment problem if the requirements for
alpar@9 88 * a transportation problem are met and:
alpar@9 89 *
alpar@9 90 * nsorc = nsink and itsup = nsorc.
alpar@9 91 *
alpar@9 92 * RETURNS
alpar@9 93 *
alpar@9 94 * If the instance was successfully generated, the routine glp_netgen
alpar@9 95 * returns zero; otherwise, if specified parameters are inconsistent,
alpar@9 96 * the routine returns a non-zero error code.
alpar@9 97 *
alpar@9 98 * REFERENCES
alpar@9 99 *
alpar@9 100 * D.Klingman, A.Napier, and J.Stutz. NETGEN: A program for generating
alpar@9 101 * large scale capacitated assignment, transportation, and minimum cost
alpar@9 102 * flow networks. Management Science 20 (1974), 814-20. */
alpar@9 103
alpar@9 104 struct csa
alpar@9 105 { /* common storage area */
alpar@9 106 glp_graph *G;
alpar@9 107 int v_rhs, a_cap, a_cost;
alpar@9 108 int nodes, iarcs, mincst, maxcst, itsup, nsorc, nsink, nonsor,
alpar@9 109 nfsink, narcs, nsort, nftsor, ipcap, mincap, maxcap, ktl,
alpar@9 110 nodlft, *ipred, *ihead, *itail, *iflag, *isup, *lsinks, mult,
alpar@9 111 modul, i15, i16, jran;
alpar@9 112 };
alpar@9 113
alpar@9 114 #define G (csa->G)
alpar@9 115 #define v_rhs (csa->v_rhs)
alpar@9 116 #define a_cap (csa->a_cap)
alpar@9 117 #define a_cost (csa->a_cost)
alpar@9 118 #define nodes (csa->nodes)
alpar@9 119 #define iarcs (csa->iarcs)
alpar@9 120 #define mincst (csa->mincst)
alpar@9 121 #define maxcst (csa->maxcst)
alpar@9 122 #define itsup (csa->itsup)
alpar@9 123 #define nsorc (csa->nsorc)
alpar@9 124 #define nsink (csa->nsink)
alpar@9 125 #define nonsor (csa->nonsor)
alpar@9 126 #define nfsink (csa->nfsink)
alpar@9 127 #define narcs (csa->narcs)
alpar@9 128 #define nsort (csa->nsort)
alpar@9 129 #define nftsor (csa->nftsor)
alpar@9 130 #define ipcap (csa->ipcap)
alpar@9 131 #define mincap (csa->mincap)
alpar@9 132 #define maxcap (csa->maxcap)
alpar@9 133 #define ktl (csa->ktl)
alpar@9 134 #define nodlft (csa->nodlft)
alpar@9 135 #if 0
alpar@9 136 /* spent a day to find out this bug */
alpar@9 137 #define ist (csa->ist)
alpar@9 138 #else
alpar@9 139 #define ist (ipred[0])
alpar@9 140 #endif
alpar@9 141 #define ipred (csa->ipred)
alpar@9 142 #define ihead (csa->ihead)
alpar@9 143 #define itail (csa->itail)
alpar@9 144 #define iflag (csa->iflag)
alpar@9 145 #define isup (csa->isup)
alpar@9 146 #define lsinks (csa->lsinks)
alpar@9 147 #define mult (csa->mult)
alpar@9 148 #define modul (csa->modul)
alpar@9 149 #define i15 (csa->i15)
alpar@9 150 #define i16 (csa->i16)
alpar@9 151 #define jran (csa->jran)
alpar@9 152
alpar@9 153 static void cresup(struct csa *csa);
alpar@9 154 static void chain(struct csa *csa, int lpick, int lsorc);
alpar@9 155 static void chnarc(struct csa *csa, int lsorc);
alpar@9 156 static void sort(struct csa *csa);
alpar@9 157 static void pickj(struct csa *csa, int it);
alpar@9 158 static void assign(struct csa *csa);
alpar@9 159 static void setran(struct csa *csa, int iseed);
alpar@9 160 static int iran(struct csa *csa, int ilow, int ihigh);
alpar@9 161
alpar@9 162 int glp_netgen(glp_graph *G_, int _v_rhs, int _a_cap, int _a_cost,
alpar@9 163 const int parm[1+15])
alpar@9 164 { struct csa _csa, *csa = &_csa;
alpar@9 165 int iseed, nprob, ntsorc, ntsink, iphic, i, nskel, nltr, ltsink,
alpar@9 166 ntrans, npsink, nftr, npsorc, ntravl, ntrrem, lsorc, lpick,
alpar@9 167 nsksr, nsrchn, j, item, l, ks, k, ksp, li, n, ii, it, ih, icap,
alpar@9 168 jcap, icost, jcost, ret;
alpar@9 169 G = G_;
alpar@9 170 v_rhs = _v_rhs;
alpar@9 171 a_cap = _a_cap;
alpar@9 172 a_cost = _a_cost;
alpar@9 173 if (G != NULL)
alpar@9 174 { if (v_rhs >= 0 && v_rhs > G->v_size - (int)sizeof(double))
alpar@9 175 xerror("glp_netgen: v_rhs = %d; invalid offset\n", v_rhs);
alpar@9 176 if (a_cap >= 0 && a_cap > G->a_size - (int)sizeof(double))
alpar@9 177 xerror("glp_netgen: a_cap = %d; invalid offset\n", a_cap);
alpar@9 178 if (a_cost >= 0 && a_cost > G->a_size - (int)sizeof(double))
alpar@9 179 xerror("glp_netgen: a_cost = %d; invalid offset\n", a_cost);
alpar@9 180 }
alpar@9 181 /* Input the user's random number seed and fix it if
alpar@9 182 non-positive. */
alpar@9 183 iseed = parm[1];
alpar@9 184 nprob = parm[2];
alpar@9 185 if (iseed <= 0) iseed = 13502460;
alpar@9 186 setran(csa, iseed);
alpar@9 187 /* Input the user's problem characteristics. */
alpar@9 188 nodes = parm[3];
alpar@9 189 nsorc = parm[4];
alpar@9 190 nsink = parm[5];
alpar@9 191 iarcs = parm[6];
alpar@9 192 mincst = parm[7];
alpar@9 193 maxcst = parm[8];
alpar@9 194 itsup = parm[9];
alpar@9 195 ntsorc = parm[10];
alpar@9 196 ntsink = parm[11];
alpar@9 197 iphic = parm[12];
alpar@9 198 ipcap = parm[13];
alpar@9 199 mincap = parm[14];
alpar@9 200 maxcap = parm[15];
alpar@9 201 /* Check the size of the problem. */
alpar@9 202 if (!(10 <= nodes && nodes <= 100000))
alpar@9 203 { ret = 1;
alpar@9 204 goto done;
alpar@9 205 }
alpar@9 206 /* Check user supplied parameters for consistency. */
alpar@9 207 if (!(nsorc >= 0 && nsink >= 0 && nsorc + nsink <= nodes))
alpar@9 208 { ret = 2;
alpar@9 209 goto done;
alpar@9 210 }
alpar@9 211 if (iarcs < 0)
alpar@9 212 { ret = 3;
alpar@9 213 goto done;
alpar@9 214 }
alpar@9 215 if (mincst > maxcst)
alpar@9 216 { ret = 4;
alpar@9 217 goto done;
alpar@9 218 }
alpar@9 219 if (itsup < 0)
alpar@9 220 { ret = 5;
alpar@9 221 goto done;
alpar@9 222 }
alpar@9 223 if (!(0 <= ntsorc && ntsorc <= nsorc))
alpar@9 224 { ret = 6;
alpar@9 225 goto done;
alpar@9 226 }
alpar@9 227 if (!(0 <= ntsink && ntsink <= nsink))
alpar@9 228 { ret = 7;
alpar@9 229 goto done;
alpar@9 230 }
alpar@9 231 if (!(0 <= iphic && iphic <= 100))
alpar@9 232 { ret = 8;
alpar@9 233 goto done;
alpar@9 234 }
alpar@9 235 if (!(0 <= ipcap && ipcap <= 100))
alpar@9 236 { ret = 9;
alpar@9 237 goto done;
alpar@9 238 }
alpar@9 239 if (mincap > maxcap)
alpar@9 240 { ret = 10;
alpar@9 241 goto done;
alpar@9 242 }
alpar@9 243 /* Initailize the graph object. */
alpar@9 244 if (G != NULL)
alpar@9 245 { glp_erase_graph(G, G->v_size, G->a_size);
alpar@9 246 glp_add_vertices(G, nodes);
alpar@9 247 if (v_rhs >= 0)
alpar@9 248 { double zero = 0.0;
alpar@9 249 for (i = 1; i <= nodes; i++)
alpar@9 250 { glp_vertex *v = G->v[i];
alpar@9 251 memcpy((char *)v->data + v_rhs, &zero, sizeof(double));
alpar@9 252 }
alpar@9 253 }
alpar@9 254 }
alpar@9 255 /* Allocate working arrays. */
alpar@9 256 ipred = xcalloc(1+nodes, sizeof(int));
alpar@9 257 ihead = xcalloc(1+nodes, sizeof(int));
alpar@9 258 itail = xcalloc(1+nodes, sizeof(int));
alpar@9 259 iflag = xcalloc(1+nodes, sizeof(int));
alpar@9 260 isup = xcalloc(1+nodes, sizeof(int));
alpar@9 261 lsinks = xcalloc(1+nodes, sizeof(int));
alpar@9 262 /* Print the problem documentation records. */
alpar@9 263 if (G == NULL)
alpar@9 264 { xprintf("BEGIN\n");
alpar@9 265 xprintf("NETGEN PROBLEM%8d%10s%10d NODES AND%10d ARCS\n",
alpar@9 266 nprob, "", nodes, iarcs);
alpar@9 267 xprintf("USER:%11d%11d%11d%11d%11d%11d\nDATA:%11d%11d%11d%11d%"
alpar@9 268 "11d%11d\n", iseed, nsorc, nsink, mincst,
alpar@9 269 maxcst, itsup, ntsorc, ntsink, iphic, ipcap,
alpar@9 270 mincap, maxcap);
alpar@9 271 }
alpar@9 272 else
alpar@9 273 glp_set_graph_name(G, "NETGEN");
alpar@9 274 /* Set various constants used in the program. */
alpar@9 275 narcs = 0;
alpar@9 276 nskel = 0;
alpar@9 277 nltr = nodes - nsink;
alpar@9 278 ltsink = nltr + ntsink;
alpar@9 279 ntrans = nltr - nsorc;
alpar@9 280 nfsink = nltr + 1;
alpar@9 281 nonsor = nodes - nsorc + ntsorc;
alpar@9 282 npsink = nsink - ntsink;
alpar@9 283 nodlft = nodes - nsink + ntsink;
alpar@9 284 nftr = nsorc + 1;
alpar@9 285 nftsor = nsorc - ntsorc + 1;
alpar@9 286 npsorc = nsorc - ntsorc;
alpar@9 287 /* Randomly distribute the supply among the source nodes. */
alpar@9 288 if (npsorc + npsink == nodes && npsorc == npsink &&
alpar@9 289 itsup == nsorc)
alpar@9 290 { assign(csa);
alpar@9 291 nskel = nsorc;
alpar@9 292 goto L390;
alpar@9 293 }
alpar@9 294 cresup(csa);
alpar@9 295 /* Print the supply records. */
alpar@9 296 if (G == NULL)
alpar@9 297 { xprintf("SUPPLY\n");
alpar@9 298 for (i = 1; i <= nsorc; i++)
alpar@9 299 xprintf("%6s%6d%18s%10d\n", "", i, "", isup[i]);
alpar@9 300 xprintf("ARCS\n");
alpar@9 301 }
alpar@9 302 else
alpar@9 303 { if (v_rhs >= 0)
alpar@9 304 { for (i = 1; i <= nsorc; i++)
alpar@9 305 { double temp = (double)isup[i];
alpar@9 306 glp_vertex *v = G->v[i];
alpar@9 307 memcpy((char *)v->data + v_rhs, &temp, sizeof(double));
alpar@9 308 }
alpar@9 309 }
alpar@9 310 }
alpar@9 311 /* Make the sources point to themselves in ipred array. */
alpar@9 312 for (i = 1; i <= nsorc; i++)
alpar@9 313 ipred[i] = i;
alpar@9 314 if (ntrans == 0) goto L170;
alpar@9 315 /* Chain the transshipment nodes together in the ipred array. */
alpar@9 316 ist = nftr;
alpar@9 317 ipred[nltr] = 0;
alpar@9 318 for (i = nftr; i < nltr; i++)
alpar@9 319 ipred[i] = i+1;
alpar@9 320 /* Form even length chains for 60 percent of the transshipments.*/
alpar@9 321 ntravl = 6 * ntrans / 10;
alpar@9 322 ntrrem = ntrans - ntravl;
alpar@9 323 L140: lsorc = 1;
alpar@9 324 while (ntravl != 0)
alpar@9 325 { lpick = iran(csa, 1, ntravl + ntrrem);
alpar@9 326 ntravl--;
alpar@9 327 chain(csa, lpick, lsorc);
alpar@9 328 if (lsorc == nsorc) goto L140;
alpar@9 329 lsorc++;
alpar@9 330 }
alpar@9 331 /* Add the remaining transshipments to the chains. */
alpar@9 332 while (ntrrem != 0)
alpar@9 333 {
alpar@9 334 lpick = iran(csa, 1, ntrrem);
alpar@9 335 ntrrem--;
alpar@9 336 lsorc = iran(csa, 1, nsorc);
alpar@9 337 chain(csa, lpick, lsorc);
alpar@9 338 }
alpar@9 339 L170: /* Set all demands equal to zero. */
alpar@9 340 for (i = nfsink; i <= nodes; i++)
alpar@9 341 ipred[i] = 0;
alpar@9 342 /* The following loop takes one chain at a time (through the use
alpar@9 343 of logic contained in the loop and calls to other routines) and
alpar@9 344 creates the remaining network arcs. */
alpar@9 345 for (lsorc = 1; lsorc <= nsorc; lsorc++)
alpar@9 346 { chnarc(csa, lsorc);
alpar@9 347 for (i = nfsink; i <= nodes; i++)
alpar@9 348 iflag[i] = 0;
alpar@9 349 /* Choose the number of sinks to be hooked up to the current
alpar@9 350 chain. */
alpar@9 351 if (ntrans != 0)
alpar@9 352 nsksr = (nsort * 2 * nsink) / ntrans;
alpar@9 353 else
alpar@9 354 nsksr = nsink / nsorc + 1;
alpar@9 355 if (nsksr < 2) nsksr = 2;
alpar@9 356 if (nsksr > nsink) nsksr = nsink;
alpar@9 357 nsrchn = nsort;
alpar@9 358 /* Randomly pick nsksr sinks and put their names in lsinks. */
alpar@9 359 ktl = nsink;
alpar@9 360 for (j = 1; j <= nsksr; j++)
alpar@9 361 { item = iran(csa, 1, ktl);
alpar@9 362 ktl--;
alpar@9 363 for (l = nfsink; l <= nodes; l++)
alpar@9 364 { if (iflag[l] != 1)
alpar@9 365 { item--;
alpar@9 366 if (item == 0) goto L230;
alpar@9 367 }
alpar@9 368 }
alpar@9 369 break;
alpar@9 370 L230: lsinks[j] = l;
alpar@9 371 iflag[l] = 1;
alpar@9 372 }
alpar@9 373 /* If last source chain, add all sinks with zero demand to
alpar@9 374 lsinks list. */
alpar@9 375 if (lsorc == nsorc)
alpar@9 376 { for (j = nfsink; j <= nodes; j++)
alpar@9 377 { if (ipred[j] == 0 && iflag[j] != 1)
alpar@9 378 { nsksr++;
alpar@9 379 lsinks[nsksr] = j;
alpar@9 380 iflag[j] = 1;
alpar@9 381 }
alpar@9 382 }
alpar@9 383 }
alpar@9 384 /* Create demands for group of sinks in lsinks. */
alpar@9 385 ks = isup[lsorc] / nsksr;
alpar@9 386 k = ipred[lsorc];
alpar@9 387 for (i = 1; i <= nsksr; i++)
alpar@9 388 { nsort++;
alpar@9 389 ksp = iran(csa, 1, ks);
alpar@9 390 j = iran(csa, 1, nsksr);
alpar@9 391 itail[nsort] = k;
alpar@9 392 li = lsinks[i];
alpar@9 393 ihead[nsort] = li;
alpar@9 394 ipred[li] += ksp;
alpar@9 395 li = lsinks[j];
alpar@9 396 ipred[li] += ks - ksp;
alpar@9 397 n = iran(csa, 1, nsrchn);
alpar@9 398 k = lsorc;
alpar@9 399 for (ii = 1; ii <= n; ii++)
alpar@9 400 k = ipred[k];
alpar@9 401 }
alpar@9 402 li = lsinks[1];
alpar@9 403 ipred[li] += isup[lsorc] - ks * nsksr;
alpar@9 404 nskel += nsort;
alpar@9 405 /* Sort the arcs in the chain from source lsorc using itail as
alpar@9 406 sort key. */
alpar@9 407 sort(csa);
alpar@9 408 /* Print this part of skeleton and create the arcs for these
alpar@9 409 nodes. */
alpar@9 410 i = 1;
alpar@9 411 itail[nsort+1] = 0;
alpar@9 412 L300: for (j = nftsor; j <= nodes; j++)
alpar@9 413 iflag[j] = 0;
alpar@9 414 ktl = nonsor - 1;
alpar@9 415 it = itail[i];
alpar@9 416 iflag[it] = 1;
alpar@9 417 L320: ih = ihead[i];
alpar@9 418 iflag[ih] = 1;
alpar@9 419 narcs++;
alpar@9 420 ktl--;
alpar@9 421 /* Determine if this skeleton arc should be capacitated. */
alpar@9 422 icap = itsup;
alpar@9 423 jcap = iran(csa, 1, 100);
alpar@9 424 if (jcap <= ipcap)
alpar@9 425 { icap = isup[lsorc];
alpar@9 426 if (mincap > icap) icap = mincap;
alpar@9 427 }
alpar@9 428 /* Determine if this skeleton arc should have the maximum
alpar@9 429 cost. */
alpar@9 430 icost = maxcst;
alpar@9 431 jcost = iran(csa, 1, 100);
alpar@9 432 if (jcost > iphic)
alpar@9 433 icost = iran(csa, mincst, maxcst);
alpar@9 434 if (G == NULL)
alpar@9 435 xprintf("%6s%6d%6d%2s%10d%10d\n", "", it, ih, "", icost,
alpar@9 436 icap);
alpar@9 437 else
alpar@9 438 { glp_arc *a = glp_add_arc(G, it, ih);
alpar@9 439 if (a_cap >= 0)
alpar@9 440 { double temp = (double)icap;
alpar@9 441 memcpy((char *)a->data + a_cap, &temp, sizeof(double));
alpar@9 442 }
alpar@9 443 if (a_cost >= 0)
alpar@9 444 { double temp = (double)icost;
alpar@9 445 memcpy((char *)a->data + a_cost, &temp, sizeof(double));
alpar@9 446 }
alpar@9 447 }
alpar@9 448 i++;
alpar@9 449 if (itail[i] == it) goto L320;
alpar@9 450 pickj(csa, it);
alpar@9 451 if (i <= nsort) goto L300;
alpar@9 452 }
alpar@9 453 /* Create arcs from the transshipment sinks. */
alpar@9 454 if (ntsink != 0)
alpar@9 455 { for (i = nfsink; i <= ltsink; i++)
alpar@9 456 { for (j = nftsor; j <= nodes; j++)
alpar@9 457 iflag[j] = 0;
alpar@9 458 ktl = nonsor - 1;
alpar@9 459 iflag[i] = 1;
alpar@9 460 pickj(csa, i);
alpar@9 461 }
alpar@9 462 }
alpar@9 463 L390: /* Print the demand records and end record. */
alpar@9 464 if (G == NULL)
alpar@9 465 { xprintf("DEMAND\n");
alpar@9 466 for (i = nfsink; i <= nodes; i++)
alpar@9 467 xprintf("%6s%6d%18s%10d\n", "", i, "", ipred[i]);
alpar@9 468 xprintf("END\n");
alpar@9 469 }
alpar@9 470 else
alpar@9 471 { if (v_rhs >= 0)
alpar@9 472 { for (i = nfsink; i <= nodes; i++)
alpar@9 473 { double temp = - (double)ipred[i];
alpar@9 474 glp_vertex *v = G->v[i];
alpar@9 475 memcpy((char *)v->data + v_rhs, &temp, sizeof(double));
alpar@9 476 }
alpar@9 477 }
alpar@9 478 }
alpar@9 479 /* Free working arrays. */
alpar@9 480 xfree(ipred);
alpar@9 481 xfree(ihead);
alpar@9 482 xfree(itail);
alpar@9 483 xfree(iflag);
alpar@9 484 xfree(isup);
alpar@9 485 xfree(lsinks);
alpar@9 486 /* The instance has been successfully generated. */
alpar@9 487 ret = 0;
alpar@9 488 done: return ret;
alpar@9 489 }
alpar@9 490
alpar@9 491 /***********************************************************************
alpar@9 492 * The routine cresup randomly distributes the total supply among the
alpar@9 493 * source nodes. */
alpar@9 494
alpar@9 495 static void cresup(struct csa *csa)
alpar@9 496 { int i, j, ks, ksp;
alpar@9 497 xassert(itsup > nsorc);
alpar@9 498 ks = itsup / nsorc;
alpar@9 499 for (i = 1; i <= nsorc; i++)
alpar@9 500 isup[i] = 0;
alpar@9 501 for (i = 1; i <= nsorc; i++)
alpar@9 502 { ksp = iran(csa, 1, ks);
alpar@9 503 j = iran(csa, 1, nsorc);
alpar@9 504 isup[i] += ksp;
alpar@9 505 isup[j] += ks - ksp;
alpar@9 506 }
alpar@9 507 j = iran(csa, 1, nsorc);
alpar@9 508 isup[j] += itsup - ks * nsorc;
alpar@9 509 return;
alpar@9 510 }
alpar@9 511
alpar@9 512 /***********************************************************************
alpar@9 513 * The routine chain adds node lpick to the end of the chain with source
alpar@9 514 * node lsorc. */
alpar@9 515
alpar@9 516 static void chain(struct csa *csa, int lpick, int lsorc)
alpar@9 517 { int i, j, k, l, m;
alpar@9 518 k = 0;
alpar@9 519 m = ist;
alpar@9 520 for (i = 1; i <= lpick; i++)
alpar@9 521 { l = k;
alpar@9 522 k = m;
alpar@9 523 m = ipred[k];
alpar@9 524 }
alpar@9 525 ipred[l] = m;
alpar@9 526 j = ipred[lsorc];
alpar@9 527 ipred[k] = j;
alpar@9 528 ipred[lsorc] = k;
alpar@9 529 return;
alpar@9 530 }
alpar@9 531
alpar@9 532 /***********************************************************************
alpar@9 533 * The routine chnarc puts the arcs in the chain from source lsorc into
alpar@9 534 * the ihead and itail arrays for sorting. */
alpar@9 535
alpar@9 536 static void chnarc(struct csa *csa, int lsorc)
alpar@9 537 { int ito, ifrom;
alpar@9 538 nsort = 0;
alpar@9 539 ito = ipred[lsorc];
alpar@9 540 L10: if (ito == lsorc) return;
alpar@9 541 nsort++;
alpar@9 542 ifrom = ipred[ito];
alpar@9 543 ihead[nsort] = ito;
alpar@9 544 itail[nsort] = ifrom;
alpar@9 545 ito = ifrom;
alpar@9 546 goto L10;
alpar@9 547 }
alpar@9 548
alpar@9 549 /***********************************************************************
alpar@9 550 * The routine sort sorts the nsort arcs in the ihead and itail arrays.
alpar@9 551 * ihead is used as the sort key (i.e. forward star sort order). */
alpar@9 552
alpar@9 553 static void sort(struct csa *csa)
alpar@9 554 { int i, j, k, l, m, n, it;
alpar@9 555 n = nsort;
alpar@9 556 m = n;
alpar@9 557 L10: m /= 2;
alpar@9 558 if (m == 0) return;
alpar@9 559 k = n - m;
alpar@9 560 j = 1;
alpar@9 561 L20: i = j;
alpar@9 562 L30: l = i + m;
alpar@9 563 if (itail[i] <= itail[l]) goto L40;
alpar@9 564 it = itail[i];
alpar@9 565 itail[i] = itail[l];
alpar@9 566 itail[l] = it;
alpar@9 567 it = ihead[i];
alpar@9 568 ihead[i] = ihead[l];
alpar@9 569 ihead[l] = it;
alpar@9 570 i -= m;
alpar@9 571 if (i >= 1) goto L30;
alpar@9 572 L40: j++;
alpar@9 573 if (j <= k) goto L20;
alpar@9 574 goto L10;
alpar@9 575 }
alpar@9 576
alpar@9 577 /***********************************************************************
alpar@9 578 * The routine pickj creates a random number of arcs out of node 'it'.
alpar@9 579 * Various parameters are dynamically adjusted in an attempt to ensure
alpar@9 580 * that the generated network has the correct number of arcs. */
alpar@9 581
alpar@9 582 static void pickj(struct csa *csa, int it)
alpar@9 583 { int j, k, l, nn, nupbnd, icap, jcap, icost;
alpar@9 584 if ((nodlft - 1) * 2 > iarcs - narcs - 1)
alpar@9 585 { nodlft--;
alpar@9 586 return;
alpar@9 587 }
alpar@9 588 if ((iarcs - narcs + nonsor - ktl - 1) / nodlft - nonsor + 1 >= 0)
alpar@9 589 k = nonsor;
alpar@9 590 else
alpar@9 591 { nupbnd = (iarcs - narcs - nodlft) / nodlft * 2;
alpar@9 592 L40: k = iran(csa, 1, nupbnd);
alpar@9 593 if (nodlft == 1) k = iarcs - narcs;
alpar@9 594 if ((nodlft - 1) * (nonsor - 1) < iarcs - narcs - k) goto L40;
alpar@9 595 }
alpar@9 596 nodlft--;
alpar@9 597 for (j = 1; j <= k; j++)
alpar@9 598 { nn = iran(csa, 1, ktl);
alpar@9 599 ktl--;
alpar@9 600 for (l = nftsor; l <= nodes; l++)
alpar@9 601 { if (iflag[l] != 1)
alpar@9 602 { nn--;
alpar@9 603 if (nn == 0) goto L70;
alpar@9 604 }
alpar@9 605 }
alpar@9 606 return;
alpar@9 607 L70: iflag[l] = 1;
alpar@9 608 icap = itsup;
alpar@9 609 jcap = iran(csa, 1, 100);
alpar@9 610 if (jcap <= ipcap)
alpar@9 611 icap = iran(csa, mincap, maxcap);
alpar@9 612 icost = iran(csa, mincst, maxcst);
alpar@9 613 if (G == NULL)
alpar@9 614 xprintf("%6s%6d%6d%2s%10d%10d\n", "", it, l, "", icost,
alpar@9 615 icap);
alpar@9 616 else
alpar@9 617 { glp_arc *a = glp_add_arc(G, it, l);
alpar@9 618 if (a_cap >= 0)
alpar@9 619 { double temp = (double)icap;
alpar@9 620 memcpy((char *)a->data + a_cap, &temp, sizeof(double));
alpar@9 621 }
alpar@9 622 if (a_cost >= 0)
alpar@9 623 { double temp = (double)icost;
alpar@9 624 memcpy((char *)a->data + a_cost, &temp, sizeof(double));
alpar@9 625 }
alpar@9 626 }
alpar@9 627 narcs++;
alpar@9 628 }
alpar@9 629 return;
alpar@9 630 }
alpar@9 631
alpar@9 632 /***********************************************************************
alpar@9 633 * The routine assign generate assignment problems. It defines the unit
alpar@9 634 * supplies, builds a skeleton, then calls pickj to create the arcs. */
alpar@9 635
alpar@9 636 static void assign(struct csa *csa)
alpar@9 637 { int i, it, nn, l, ll, icost;
alpar@9 638 if (G == NULL)
alpar@9 639 xprintf("SUPPLY\n");
alpar@9 640 for (i = 1; i <= nsorc; i++)
alpar@9 641 { isup[i] = 1;
alpar@9 642 iflag[i] = 0;
alpar@9 643 if (G == NULL)
alpar@9 644 xprintf("%6s%6d%18s%10d\n", "", i, "", isup[i]);
alpar@9 645 else
alpar@9 646 { if (v_rhs >= 0)
alpar@9 647 { double temp = (double)isup[i];
alpar@9 648 glp_vertex *v = G->v[i];
alpar@9 649 memcpy((char *)v->data + v_rhs, &temp, sizeof(double));
alpar@9 650 }
alpar@9 651 }
alpar@9 652 }
alpar@9 653 if (G == NULL)
alpar@9 654 xprintf("ARCS\n");
alpar@9 655 for (i = nfsink; i <= nodes; i++)
alpar@9 656 ipred[i] = 1;
alpar@9 657 for (it = 1; it <= nsorc; it++)
alpar@9 658 { for (i = nfsink; i <= nodes; i++)
alpar@9 659 iflag[i] = 0;
alpar@9 660 ktl = nsink - 1;
alpar@9 661 nn = iran(csa, 1, nsink - it + 1);
alpar@9 662 for (l = 1; l <= nsorc; l++)
alpar@9 663 { if (iflag[l] != 1)
alpar@9 664 { nn--;
alpar@9 665 if (nn == 0) break;
alpar@9 666 }
alpar@9 667 }
alpar@9 668 narcs++;
alpar@9 669 ll = nsorc + l;
alpar@9 670 icost = iran(csa, mincst, maxcst);
alpar@9 671 if (G == NULL)
alpar@9 672 xprintf("%6s%6d%6d%2s%10d%10d\n", "", it, ll, "", icost,
alpar@9 673 isup[1]);
alpar@9 674 else
alpar@9 675 { glp_arc *a = glp_add_arc(G, it, ll);
alpar@9 676 if (a_cap >= 0)
alpar@9 677 { double temp = (double)isup[1];
alpar@9 678 memcpy((char *)a->data + a_cap, &temp, sizeof(double));
alpar@9 679 }
alpar@9 680 if (a_cost >= 0)
alpar@9 681 { double temp = (double)icost;
alpar@9 682 memcpy((char *)a->data + a_cost, &temp, sizeof(double));
alpar@9 683 }
alpar@9 684 }
alpar@9 685 iflag[l] = 1;
alpar@9 686 iflag[ll] = 1;
alpar@9 687 pickj(csa, it);
alpar@9 688 }
alpar@9 689 return;
alpar@9 690 }
alpar@9 691
alpar@9 692 /***********************************************************************
alpar@9 693 * Portable congruential (uniform) random number generator:
alpar@9 694 *
alpar@9 695 * next_value = ((7**5) * previous_value) modulo ((2**31)-1)
alpar@9 696 *
alpar@9 697 * This generator consists of three routines:
alpar@9 698 *
alpar@9 699 * (1) setran - initializes constants and seed
alpar@9 700 * (2) iran - generates an integer random number
alpar@9 701 * (3) rran - generates a real random number
alpar@9 702 *
alpar@9 703 * The generator requires a machine with at least 32 bits of precision.
alpar@9 704 * The seed (iseed) must be in the range [1,(2**31)-1]. */
alpar@9 705
alpar@9 706 static void setran(struct csa *csa, int iseed)
alpar@9 707 { xassert(iseed >= 1);
alpar@9 708 mult = 16807;
alpar@9 709 modul = 2147483647;
alpar@9 710 i15 = 1 << 15;
alpar@9 711 i16 = 1 << 16;
alpar@9 712 jran = iseed;
alpar@9 713 return;
alpar@9 714 }
alpar@9 715
alpar@9 716 /***********************************************************************
alpar@9 717 * The routine iran generates an integer random number between ilow and
alpar@9 718 * ihigh. If ilow > ihigh then iran returns ihigh. */
alpar@9 719
alpar@9 720 static int iran(struct csa *csa, int ilow, int ihigh)
alpar@9 721 { int ixhi, ixlo, ixalo, leftlo, ixahi, ifulhi, irtlo, iover,
alpar@9 722 irthi, j;
alpar@9 723 ixhi = jran / i16;
alpar@9 724 ixlo = jran - ixhi * i16;
alpar@9 725 ixalo = ixlo * mult;
alpar@9 726 leftlo = ixalo / i16;
alpar@9 727 ixahi = ixhi * mult;
alpar@9 728 ifulhi = ixahi + leftlo;
alpar@9 729 irtlo = ixalo - leftlo * i16;
alpar@9 730 iover = ifulhi / i15;
alpar@9 731 irthi = ifulhi - iover * i15;
alpar@9 732 jran = ((irtlo - modul) + irthi * i16) + iover;
alpar@9 733 if (jran < 0) jran += modul;
alpar@9 734 j = ihigh - ilow + 1;
alpar@9 735 if (j > 0)
alpar@9 736 return jran % j + ilow;
alpar@9 737 else
alpar@9 738 return ihigh;
alpar@9 739 }
alpar@9 740
alpar@9 741 /**********************************************************************/
alpar@9 742
alpar@9 743 #if 0
alpar@9 744 static int scan(char card[80+1], int pos, int len)
alpar@9 745 { char buf[10+1];
alpar@9 746 memcpy(buf, &card[pos-1], len);
alpar@9 747 buf[len] = '\0';
alpar@9 748 return atoi(buf);
alpar@9 749 }
alpar@9 750
alpar@9 751 int main(void)
alpar@9 752 { int parm[1+15];
alpar@9 753 char card[80+1];
alpar@9 754 xassert(fgets(card, sizeof(card), stdin) == card);
alpar@9 755 parm[1] = scan(card, 1, 8);
alpar@9 756 parm[2] = scan(card, 9, 8);
alpar@9 757 xassert(fgets(card, sizeof(card), stdin) == card);
alpar@9 758 parm[3] = scan(card, 1, 5);
alpar@9 759 parm[4] = scan(card, 6, 5);
alpar@9 760 parm[5] = scan(card, 11, 5);
alpar@9 761 parm[6] = scan(card, 16, 5);
alpar@9 762 parm[7] = scan(card, 21, 5);
alpar@9 763 parm[8] = scan(card, 26, 5);
alpar@9 764 parm[9] = scan(card, 31, 10);
alpar@9 765 parm[10] = scan(card, 41, 5);
alpar@9 766 parm[11] = scan(card, 46, 5);
alpar@9 767 parm[12] = scan(card, 51, 5);
alpar@9 768 parm[13] = scan(card, 56, 5);
alpar@9 769 parm[14] = scan(card, 61, 10);
alpar@9 770 parm[15] = scan(card, 71, 10);
alpar@9 771 glp_netgen(NULL, 0, 0, 0, parm);
alpar@9 772 return 0;
alpar@9 773 }
alpar@9 774 #endif
alpar@9 775
alpar@9 776 /* eof */