glpnet03.c

著名的大规模线性规划求解器源码GLPK.C语言版本,可以修剪.内有详细帮助文档.

            j = iran(csa, 1, nsksr);
            itail[nsort] = k;
            li = lsinks[i];
            ihead[nsort] = li;
            ipred[li] += ksp;
            li = lsinks[j];
            ipred[li] += ks - ksp;
            n = iran(csa, 1, nsrchn);
            k = lsorc;
            for (ii = 1; ii <= n; ii++)
               k = ipred[k];
         }
         li = lsinks[1];
         ipred[li] += isup[lsorc] - ks * nsksr;
         nskel += nsort;
         /* Sort the arcs in the chain from source lsorc using itail as
            sort key. */
         sort(csa);
         /* Print this part of skeleton and create the arcs for these
            nodes. */
         i = 1;
         itail[nsort+1] = 0;
L300:    for (j = nftsor; j <= nodes; j++)
            iflag[j] = 0;
         ktl = nonsor - 1;
         it = itail[i];
         iflag[it] = 1;
L320:    ih = ihead[i];
         iflag[ih] = 1;
         narcs++;
         ktl--;
         /* Determine if this skeleton arc should be capacitated. */
         icap = itsup;
         jcap = iran(csa, 1, 100);
         if (jcap <= ipcap)
         {  icap = isup[lsorc];
            if (mincap > icap) icap = mincap;
         }
         /* Determine if this skeleton arc should have the maximum
            cost. */
         icost = maxcst;
         jcost = iran(csa, 1, 100);
         if (jcost > iphic)
            icost = iran(csa, mincst, maxcst);
         if (G == NULL)
            xprintf("%6s%6d%6d%2s%10d%10d\n", "", it, ih, "", icost,
               icap);
         else
         {  glp_arc *a = glp_add_arc(G, it, ih);
            if (a_cap >= 0)
            {  double temp = (double)icap;
               memcpy((char *)a->data + a_cap, &temp, sizeof(double));
            }
            if (a_cost >= 0)
            {  double temp = (double)icost;
               memcpy((char *)a->data + a_cost, &temp, sizeof(double));
            }
         }
         i++;
         if (itail[i] == it) goto L320;
         pickj(csa, it);
         if (i <= nsort) goto L300;
      }
      /* Create arcs from the transshipment sinks. */
      if (ntsink != 0)
      {  for (i = nfsink; i <= ltsink; i++)
         {  for (j = nftsor; j <= nodes; j++)
               iflag[j] = 0;
            ktl = nonsor - 1;
            iflag[i] = 1;
            pickj(csa, i);
         }
      }
L390: /* Print the demand records and end record. */
      if (G == NULL)
      {  xprintf("DEMAND\n");
         for (i = nfsink; i <= nodes; i++)
            xprintf("%6s%6d%18s%10d\n", "", i, "", ipred[i]);
         xprintf("END\n");
      }
      else
      {  if (v_rhs >= 0)
         {  for (i = nfsink; i <= nodes; i++)
            {  double temp = - (double)ipred[i];
               glp_vertex *v = G->v[i];
               memcpy((char *)v->data + v_rhs, &temp, sizeof(double));
            }
         }
      }
      /* Free working arrays. */
      xfree(ipred);
      xfree(ihead);
      xfree(itail);
      xfree(iflag);
      xfree(isup);
      xfree(lsinks);
      /* The instance has been successfully generated. */
      ret = 0;
done: return ret;
}

/***********************************************************************
*  The routine cresup randomly distributes the total supply among the
*  source nodes. */

static void cresup(struct csa *csa)
{     int i, j, ks, ksp;
      xassert(itsup > nsorc);
      ks = itsup / nsorc;
      for (i = 1; i <= nsorc; i++)
         isup[i] = 0;
      for (i = 1; i <= nsorc; i++)
      {  ksp = iran(csa, 1, ks);
         j = iran(csa, 1, nsorc);
         isup[i] += ksp;
         isup[j] += ks - ksp;
      }
      j = iran(csa, 1, nsorc);
      isup[j] += itsup - ks * nsorc;
      return;
}

/***********************************************************************
*  The routine chain adds node lpick to the end of the chain with source
*  node lsorc. */

static void chain(struct csa *csa, int lpick, int lsorc)
{     int i, j, k, l, m;
      k = 0;
      m = ist;
      for (i = 1; i <= lpick; i++)
      {  l = k;
         k = m;
         m = ipred[k];
      }
      ipred[l] = m;
      j = ipred[lsorc];
      ipred[k] = j;
      ipred[lsorc] = k;
      return;
}

/***********************************************************************
*  The routine chnarc puts the arcs in the chain from source lsorc into
*  the ihead and itail arrays for sorting. */

static void chnarc(struct csa *csa, int lsorc)
{     int ito, ifrom;
      nsort = 0;
      ito = ipred[lsorc];
L10:  if (ito == lsorc) return;
      nsort++;
      ifrom = ipred[ito];
      ihead[nsort] = ito;
      itail[nsort] = ifrom;
      ito = ifrom;
      goto L10;
}

/***********************************************************************
*  The routine sort sorts the nsort arcs in the ihead and itail arrays.
*  ihead is used as the sort key (i.e. forward star sort order). */

static void sort(struct csa *csa)
{     int i, j, k, l, m, n, it;
      n = nsort;
      m = n;
L10:  m /= 2;
      if (m == 0) return;
      k = n - m;
      j = 1;
L20:  i = j;
L30:  l = i + m;
      if (itail[i] <= itail[l]) goto L40;
      it = itail[i];
      itail[i] = itail[l];
      itail[l] = it;
      it = ihead[i];
      ihead[i] = ihead[l];
      ihead[l] = it;
      i -= m;
      if (i >= 1) goto L30;
L40:  j++;
      if (j <= k) goto L20;
      goto L10;
}

/***********************************************************************
*  The routine pickj creates a random number of arcs out of node 'it'.
*  Various parameters are dynamically adjusted in an attempt to ensure
*  that the generated network has the correct number of arcs. */

static void pickj(struct csa *csa, int it)
{     int j, k, l, nn, nupbnd, icap, jcap, icost;
      if ((nodlft - 1) * 2 > iarcs - narcs - 1)
      {  nodlft--;
         return;
      }
      if ((iarcs - narcs + nonsor - ktl - 1) / nodlft - nonsor + 1 >= 0)
         k = nonsor;
      else
      {  nupbnd = (iarcs - narcs - nodlft) / nodlft * 2;
L40:     k = iran(csa, 1, nupbnd);
         if (nodlft == 1) k = iarcs - narcs;
         if ((nodlft - 1) * (nonsor - 1) < iarcs - narcs - k) goto L40;
      }
      nodlft--;
      for (j = 1; j <= k; j++)
      {  nn = iran(csa, 1, ktl);
         ktl--;
         for (l = nftsor; l <= nodes; l++)
         {  if (iflag[l] != 1)
            {  nn--;
               if (nn == 0) goto L70;
            }
         }
         return;
L70:     iflag[l] = 1;
         icap = itsup;
         jcap = iran(csa, 1, 100);
         if (jcap <= ipcap)
            icap = iran(csa, mincap, maxcap);
         icost = iran(csa, mincst, maxcst);
         if (G == NULL)
            xprintf("%6s%6d%6d%2s%10d%10d\n", "", it, l, "", icost,
               icap);
         else
         {  glp_arc *a = glp_add_arc(G, it, l);
            if (a_cap >= 0)
            {  double temp = (double)icap;
               memcpy((char *)a->data + a_cap, &temp, sizeof(double));
            }
            if (a_cost >= 0)
            {  double temp = (double)icost;
               memcpy((char *)a->data + a_cost, &temp, sizeof(double));
            }
         }
         narcs++;
      }
      return;
}

/***********************************************************************
*  The routine assign generate assignment problems. It defines the unit
*  supplies, builds a skeleton, then calls pickj to create the arcs. */

static void assign(struct csa *csa)
{     int i, it, nn, l, ll, icost;
      if (G == NULL)
         xprintf("SUPPLY\n");
      for (i = 1; i <= nsorc; i++)
      {  isup[i] = 1;
         iflag[i] = 0;
         if (G == NULL)
            xprintf("%6s%6d%18s%10d\n", "", i, "", isup[i]);
         else
         {  if (v_rhs >= 0)
            {  double temp = (double)isup[i];
               glp_vertex *v = G->v[i];
               memcpy((char *)v->data + v_rhs, &temp, sizeof(double));
            }
         }
      }
      if (G == NULL)
         xprintf("ARCS\n");
      for (i = nfsink; i <= nodes; i++)
         ipred[i] = 1;
      for (it = 1; it <= nsorc; it++)
      {  for (i = nfsink; i <= nodes; i++)
            iflag[i] = 0;
         ktl = nsink - 1;
         nn = iran(csa, 1, nsink - it + 1);
         for (l = 1; l <= nsorc; l++)
         {  if (iflag[l] != 1)
            {  nn--;
               if (nn == 0) break;
            }
         }
         narcs++;
         ll = nsorc + l;
         icost = iran(csa, mincst, maxcst);
         if (G == NULL)
            xprintf("%6s%6d%6d%2s%10d%10d\n", "", it, ll, "", icost,
               isup[1]);
         else
         {  glp_arc *a = glp_add_arc(G, it, ll);
            if (a_cap >= 0)
            {  double temp = (double)isup[1];
               memcpy((char *)a->data + a_cap, &temp, sizeof(double));
            }
            if (a_cost >= 0)
            {  double temp = (double)icost;
               memcpy((char *)a->data + a_cost, &temp, sizeof(double));
            }
         }
         iflag[l] = 1;
         iflag[ll] = 1;
         pickj(csa, it);
      }
      return;
}

/***********************************************************************
*  Portable congruential (uniform) random number generator:
*
*     next_value = ((7**5) * previous_value) modulo ((2**31)-1)
*
*  This generator consists of three routines:
*
*  (1) setran - initializes constants and seed
*  (2) iran   - generates an integer random number
*  (3) rran   - generates a real random number
*
*  The generator requires a machine with at least 32 bits of precision.
*  The seed (iseed) must be in the range [1,(2**31)-1]. */

static void setran(struct csa *csa, int iseed)
{     xassert(iseed >= 1);
      mult = 16807;
      modul = 2147483647;
      i15 = 1 << 15;
      i16 = 1 << 16;
      jran = iseed;
      return;
}

/***********************************************************************
*  The routine iran generates an integer random number between ilow and
*  ihigh. If ilow > ihigh then iran returns ihigh. */

static int iran(struct csa *csa, int ilow, int ihigh)
{     int ixhi, ixlo, ixalo, leftlo, ixahi, ifulhi, irtlo, iover,
         irthi, j;
      ixhi = jran / i16;
      ixlo = jran - ixhi * i16;
      ixalo = ixlo * mult;
      leftlo = ixalo / i16;
      ixahi = ixhi * mult;
      ifulhi = ixahi + leftlo;
      irtlo = ixalo - leftlo * i16;
      iover = ifulhi / i15;
      irthi = ifulhi - iover * i15;
      jran = ((irtlo - modul) + irthi * i16) + iover;
      if (jran < 0) jran += modul;
      j = ihigh - ilow + 1;
      if (j > 0)
         return jran % j + ilow;
      else
         return ihigh;
}

/**********************************************************************/

#if 0
static int scan(char card[80+1], int pos, int len)
{     char buf[10+1];
      memcpy(buf, &card[pos-1], len);
      buf[len] = '\0';
      return atoi(buf);
}

int main(void)
{     int parm[1+15];
      char card[80+1];
      xassert(fgets(card, sizeof(card), stdin) == card);
      parm[1] = scan(card, 1, 8);
      parm[2] = scan(card, 9, 8);
      xassert(fgets(card, sizeof(card), stdin) == card);
      parm[3] = scan(card, 1, 5);
      parm[4] = scan(card, 6, 5);
      parm[5] = scan(card, 11, 5);
      parm[6] = scan(card, 16, 5);
      parm[7] = scan(card, 21, 5);
      parm[8] = scan(card, 26, 5);
      parm[9] = scan(card, 31, 10);
      parm[10] = scan(card, 41, 5);
      parm[11] = scan(card, 46, 5);
      parm[12] = scan(card, 51, 5);
      parm[13] = scan(card, 56, 5);
      parm[14] = scan(card, 61, 10);
      parm[15] = scan(card, 71, 10);
      glp_netgen(NULL, 0, 0, 0, parm);
      return 0;
}
#endif

/* eof */