src/glpnet04.c
changeset 2 4c8956a7bdf4
equal deleted inserted replaced
-1:000000000000 0:71ceecb470ba
       
     1 /* glpnet04.c (grid-like network problem generator) */
       
     2 
       
     3 /***********************************************************************
       
     4 *  This code is part of GLPK (GNU Linear Programming Kit).
       
     5 *
       
     6 *  This code is a modified version of the program GRIDGEN, a grid-like
       
     7 *  network problem generator developed by Yusin Lee and Jim Orlin.
       
     8 *  The original code is publically available on the DIMACS ftp site at:
       
     9 *  <ftp://dimacs.rutgers.edu/pub/netflow/generators/network/gridgen>.
       
    10 *
       
    11 *  All changes concern only the program interface, so this modified
       
    12 *  version produces exactly the same instances as the original version.
       
    13 *
       
    14 *  Changes were made by Andrew Makhorin <mao@gnu.org>.
       
    15 *
       
    16 *  GLPK is free software: you can redistribute it and/or modify it
       
    17 *  under the terms of the GNU General Public License as published by
       
    18 *  the Free Software Foundation, either version 3 of the License, or
       
    19 *  (at your option) any later version.
       
    20 *
       
    21 *  GLPK is distributed in the hope that it will be useful, but WITHOUT
       
    22 *  ANY WARRANTY; without even the implied warranty of MERCHANTABILITY
       
    23 *  or FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public
       
    24 *  License for more details.
       
    25 *
       
    26 *  You should have received a copy of the GNU General Public License
       
    27 *  along with GLPK. If not, see <http://www.gnu.org/licenses/>.
       
    28 ***********************************************************************/
       
    29 
       
    30 #include "glpapi.h"
       
    31 
       
    32 /***********************************************************************
       
    33 *  NAME
       
    34 *
       
    35 *  glp_gridgen - grid-like network problem generator
       
    36 *
       
    37 *  SYNOPSIS
       
    38 *
       
    39 *  int glp_gridgen(glp_graph *G, int v_rhs, int a_cap, int a_cost,
       
    40 *     const int parm[1+14]);
       
    41 *
       
    42 *  DESCRIPTION
       
    43 *
       
    44 *  The routine glp_gridgen is a grid-like network problem generator
       
    45 *  developed by Yusin Lee and Jim Orlin.
       
    46 *
       
    47 *  The parameter G specifies the graph object, to which the generated
       
    48 *  problem data have to be stored. Note that on entry the graph object
       
    49 *  is erased with the routine glp_erase_graph.
       
    50 *
       
    51 *  The parameter v_rhs specifies an offset of the field of type double
       
    52 *  in the vertex data block, to which the routine stores the supply or
       
    53 *  demand value. If v_rhs < 0, the value is not stored.
       
    54 *
       
    55 *  The parameter a_cap specifies an offset of the field of type double
       
    56 *  in the arc data block, to which the routine stores the arc capacity.
       
    57 *  If a_cap < 0, the capacity is not stored.
       
    58 *
       
    59 *  The parameter a_cost specifies an offset of the field of type double
       
    60 *  in the arc data block, to which the routine stores the per-unit cost
       
    61 *  if the arc flow. If a_cost < 0, the cost is not stored.
       
    62 *
       
    63 *  The array parm contains description of the network to be generated:
       
    64 *
       
    65 *  parm[0]  not used
       
    66 *  parm[1]  two-ways arcs indicator:
       
    67 *           1 - if links in both direction should be generated
       
    68 *           0 - otherwise
       
    69 *  parm[2]  random number seed (a positive integer)
       
    70 *  parm[3]  number of nodes (the number of nodes generated might be
       
    71 *           slightly different to make the network a grid)
       
    72 *  parm[4]  grid width
       
    73 *  parm[5]  number of sources
       
    74 *  parm[6]  number of sinks
       
    75 *  parm[7]  average degree
       
    76 *  parm[8]  total flow
       
    77 *  parm[9]  distribution of arc costs:
       
    78 *           1 - uniform
       
    79 *           2 - exponential
       
    80 *  parm[10] lower bound for arc cost (uniform)
       
    81 *           100 * lambda (exponential)
       
    82 *  parm[11] upper bound for arc cost (uniform)
       
    83 *           not used (exponential)
       
    84 *  parm[12] distribution of arc capacities:
       
    85 *           1 - uniform
       
    86 *           2 - exponential
       
    87 *  parm[13] lower bound for arc capacity (uniform)
       
    88 *           100 * lambda (exponential)
       
    89 *  parm[14] upper bound for arc capacity (uniform)
       
    90 *           not used (exponential)
       
    91 *
       
    92 *  RETURNS
       
    93 *
       
    94 *  If the instance was successfully generated, the routine glp_gridgen
       
    95 *  returns zero; otherwise, if specified parameters are inconsistent,
       
    96 *  the routine returns a non-zero error code.
       
    97 *
       
    98 *  COMMENTS
       
    99 *
       
   100 *  This network generator generates a grid-like network plus a super
       
   101 *  node. In additional to the arcs connecting the nodes in the grid,
       
   102 *  there is an arc from each supply node to the super node and from the
       
   103 *  super node to each demand node to guarantee feasiblity. These arcs
       
   104 *  have very high costs and very big capacities.
       
   105 *
       
   106 *  The idea of this network generator is as follows: First, a grid of
       
   107 *  n1 * n2 is generated. For example, 5 * 3. The nodes are numbered as
       
   108 *  1 to 15, and the supernode is numbered as n1*n2+1. Then arcs between
       
   109 *  adjacent nodes are generated. For these arcs, the user is allowed to
       
   110 *  specify either to generate two-way arcs or one-way arcs. If two-way
       
   111 *  arcs are to be generated, two arcs, one in each direction, will be
       
   112 *  generated between each adjacent node pairs. Otherwise, only one arc
       
   113 *  will be generated. If this is the case, the arcs will be generated
       
   114 *  in alterntive directions as shown below.
       
   115 *
       
   116 *      1 ---> 2 ---> 3 ---> 4 ---> 5
       
   117 *      |      ^      |      ^      |
       
   118 *      |      |      |      |      |
       
   119 *      V      |      V      |      V
       
   120 *      6 <--- 7 <--- 8 <--- 9 <--- 10
       
   121 *      |      ^      |      ^      |
       
   122 *      |      |      |      |      |
       
   123 *      V      |      V      |      V
       
   124 *     11 --->12 --->13 --->14 ---> 15
       
   125 *
       
   126 *  Then the arcs between the super node and the source/sink nodes are
       
   127 *  added as mentioned before. If the number of arcs still doesn't reach
       
   128 *  the requirement, additional arcs will be added by uniformly picking
       
   129 *  random node pairs. There is no checking to prevent multiple arcs
       
   130 *  between any pair of nodes. However, there will be no self-arcs (arcs
       
   131 *  that poins back to its tail node) in the network.
       
   132 *
       
   133 *  The source and sink nodes are selected uniformly in the network, and
       
   134 *  the imbalances of each source/sink node are also assigned by uniform
       
   135 *  distribution. */
       
   136 
       
   137 struct stat_para
       
   138 {     /* structure for statistical distributions */
       
   139       int distribution;
       
   140       /* the distribution: */
       
   141 #define UNIFORM      1  /* uniform distribution */
       
   142 #define EXPONENTIAL  2  /* exponential distribution */
       
   143       double parameter[5];
       
   144       /* the parameters of the distribution */
       
   145 };
       
   146 
       
   147 struct arcs
       
   148 {     int from;
       
   149       /* the FROM node of that arc */
       
   150       int to;
       
   151       /* the TO node of that arc */
       
   152       int cost;
       
   153       /* original cost of that arc */
       
   154       int u;
       
   155       /* capacity of the arc */
       
   156 };
       
   157 
       
   158 struct imbalance
       
   159 {     int node;
       
   160       /* Node ID */
       
   161       int supply;
       
   162       /* Supply of that node */
       
   163 };
       
   164 
       
   165 struct csa
       
   166 {     /* common storage area */
       
   167       glp_graph *G;
       
   168       int v_rhs, a_cap, a_cost;
       
   169       int seed;
       
   170       /* random number seed */
       
   171       int seed_original;
       
   172       /* the original seed from input */
       
   173       int two_way;
       
   174       /* 0: generate arcs in both direction for the basic grid, except
       
   175          for the arcs to/from the super node.  1: o/w */
       
   176       int n_node;
       
   177       /* total number of nodes in the network, numbered 1 to n_node,
       
   178          including the super node, which is the last one */
       
   179       int n_arc;
       
   180       /* total number of arcs in the network, counting EVERY arc. */
       
   181       int n_grid_arc;
       
   182       /* number of arcs in the basic grid, including the arcs to/from
       
   183          the super node */
       
   184       int n_source, n_sink;
       
   185       /* number of source and sink nodes */
       
   186       int avg_degree;
       
   187       /* average degree, arcs to and from the super node are counted */
       
   188       int t_supply;
       
   189       /* total supply in the network */
       
   190       int n1, n2;
       
   191       /* the two edges of the network grid.  n1 >= n2 */
       
   192       struct imbalance *source_list, *sink_list;
       
   193       /* head of the array of source/sink nodes */
       
   194       struct stat_para arc_costs;
       
   195       /* the distribution of arc costs */
       
   196       struct stat_para capacities;
       
   197       /* distribution of the capacities of the arcs */
       
   198       struct arcs *arc_list;
       
   199       /* head of the arc list array.  Arcs in this array are in the
       
   200          order of grid_arcs, arcs to/from super node, and other arcs */
       
   201 };
       
   202 
       
   203 #define G (csa->G)
       
   204 #define v_rhs (csa->v_rhs)
       
   205 #define a_cap (csa->a_cap)
       
   206 #define a_cost (csa->a_cost)
       
   207 #define seed (csa->seed)
       
   208 #define seed_original (csa->seed_original)
       
   209 #define two_way (csa->two_way)
       
   210 #define n_node (csa->n_node)
       
   211 #define n_arc (csa->n_arc)
       
   212 #define n_grid_arc (csa->n_grid_arc)
       
   213 #define n_source (csa->n_source)
       
   214 #define n_sink (csa->n_sink)
       
   215 #define avg_degree (csa->avg_degree)
       
   216 #define t_supply (csa->t_supply)
       
   217 #define n1 (csa->n1)
       
   218 #define n2 (csa->n2)
       
   219 #define source_list (csa->source_list)
       
   220 #define sink_list (csa->sink_list)
       
   221 #define arc_costs (csa->arc_costs)
       
   222 #define capacities (csa->capacities)
       
   223 #define arc_list (csa->arc_list)
       
   224 
       
   225 static void assign_capacities(struct csa *csa);
       
   226 static void assign_costs(struct csa *csa);
       
   227 static void assign_imbalance(struct csa *csa);
       
   228 static int exponential(struct csa *csa, double lambda[1]);
       
   229 static struct arcs *gen_additional_arcs(struct csa *csa, struct arcs
       
   230       *arc_ptr);
       
   231 static struct arcs *gen_basic_grid(struct csa *csa, struct arcs
       
   232       *arc_ptr);
       
   233 static void gen_more_arcs(struct csa *csa, struct arcs *arc_ptr);
       
   234 static void generate(struct csa *csa);
       
   235 static void output(struct csa *csa);
       
   236 static double randy(struct csa *csa);
       
   237 static void select_source_sinks(struct csa *csa);
       
   238 static int uniform(struct csa *csa, double a[2]);
       
   239 
       
   240 int glp_gridgen(glp_graph *G_, int _v_rhs, int _a_cap, int _a_cost,
       
   241       const int parm[1+14])
       
   242 {     struct csa _csa, *csa = &_csa;
       
   243       int n, ret;
       
   244       G = G_;
       
   245       v_rhs = _v_rhs;
       
   246       a_cap = _a_cap;
       
   247       a_cost = _a_cost;
       
   248       if (G != NULL)
       
   249       {  if (v_rhs >= 0 && v_rhs > G->v_size - (int)sizeof(double))
       
   250             xerror("glp_gridgen: v_rhs = %d; invalid offset\n", v_rhs);
       
   251          if (a_cap >= 0 && a_cap > G->a_size - (int)sizeof(double))
       
   252             xerror("glp_gridgen: a_cap = %d; invalid offset\n", a_cap);
       
   253          if (a_cost >= 0 && a_cost > G->a_size - (int)sizeof(double))
       
   254             xerror("glp_gridgen: a_cost = %d; invalid offset\n", a_cost)
       
   255                ;
       
   256       }
       
   257       /* Check the parameters for consistency. */
       
   258       if (!(parm[1] == 0 || parm[1] == 1))
       
   259       {  ret = 1;
       
   260          goto done;
       
   261       }
       
   262       if (parm[2] < 1)
       
   263       {  ret = 2;
       
   264          goto done;
       
   265       }
       
   266       if (!(10 <= parm[3] && parm[3] <= 40000))
       
   267       {  ret = 3;
       
   268          goto done;
       
   269       }
       
   270       if (!(1 <= parm[4] && parm[4] <= 40000))
       
   271       {  ret = 4;
       
   272          goto done;
       
   273       }
       
   274       if (!(parm[5] >= 0 && parm[6] >= 0 && parm[5] + parm[6] <=
       
   275          parm[3]))
       
   276       {  ret = 5;
       
   277          goto done;
       
   278       }
       
   279       if (!(1 <= parm[7] && parm[7] <= parm[3]))
       
   280       {  ret = 6;
       
   281          goto done;
       
   282       }
       
   283       if (parm[8] < 0)
       
   284       {  ret = 7;
       
   285          goto done;
       
   286       }
       
   287       if (!(parm[9] == 1 || parm[9] == 2))
       
   288       {  ret = 8;
       
   289          goto done;
       
   290       }
       
   291       if (parm[9] == 1 && parm[10] > parm[11] ||
       
   292           parm[9] == 2 && parm[10] < 1)
       
   293       {  ret = 9;
       
   294          goto done;
       
   295       }
       
   296       if (!(parm[12] == 1 || parm[12] == 2))
       
   297       {  ret = 10;
       
   298          goto done;
       
   299       }
       
   300       if (parm[12] == 1 && !(0 <= parm[13] && parm[13] <= parm[14]) ||
       
   301           parm[12] == 2 && parm[13] < 1)
       
   302       {  ret = 11;
       
   303          goto done;
       
   304       }
       
   305       /* Initialize the graph object. */
       
   306       if (G != NULL)
       
   307       {  glp_erase_graph(G, G->v_size, G->a_size);
       
   308          glp_set_graph_name(G, "GRIDGEN");
       
   309       }
       
   310       /* Copy the generator parameters. */
       
   311       two_way = parm[1];
       
   312       seed_original = seed = parm[2];
       
   313       n_node = parm[3];
       
   314       n = parm[4];
       
   315       n_source = parm[5];
       
   316       n_sink = parm[6];
       
   317       avg_degree = parm[7];
       
   318       t_supply = parm[8];
       
   319       arc_costs.distribution = parm[9];
       
   320       if (parm[9] == 1)
       
   321       {  arc_costs.parameter[0] = parm[10];
       
   322          arc_costs.parameter[1] = parm[11];
       
   323       }
       
   324       else
       
   325       {  arc_costs.parameter[0] = (double)parm[10] / 100.0;
       
   326          arc_costs.parameter[1] = 0.0;
       
   327       }
       
   328       capacities.distribution = parm[12];
       
   329       if (parm[12] == 1)
       
   330       {  capacities.parameter[0] = parm[13];
       
   331          capacities.parameter[1] = parm[14];
       
   332       }
       
   333       else
       
   334       {  capacities.parameter[0] = (double)parm[13] / 100.0;
       
   335          capacities.parameter[1] = 0.0;
       
   336       }
       
   337       /* Calculate the edge lengths of the grid according to the
       
   338          input. */
       
   339       if (n * n >= n_node)
       
   340       {  n1 = n;
       
   341          n2 = (int)((double)n_node / (double)n + 0.5);
       
   342       }
       
   343       else
       
   344       {  n2 = n;
       
   345          n1 = (int)((double)n_node / (double)n + 0.5);
       
   346       }
       
   347       /* Recalculate the total number of nodes and plus 1 for the super
       
   348          node. */
       
   349       n_node = n1 * n2 + 1;
       
   350       n_arc = n_node * avg_degree;
       
   351       n_grid_arc = (two_way + 1) * ((n1 - 1) * n2 + (n2 - 1) * n1) +
       
   352          n_source + n_sink;
       
   353       if (n_grid_arc > n_arc) n_arc = n_grid_arc;
       
   354       arc_list = xcalloc(n_arc, sizeof(struct arcs));
       
   355       source_list = xcalloc(n_source, sizeof(struct imbalance));
       
   356       sink_list = xcalloc(n_sink, sizeof(struct imbalance));
       
   357       /* Generate a random network. */
       
   358       generate(csa);
       
   359       /* Output the network. */
       
   360       output(csa);
       
   361       /* Free all allocated memory. */
       
   362       xfree(arc_list);
       
   363       xfree(source_list);
       
   364       xfree(sink_list);
       
   365       /* The instance has been successfully generated. */
       
   366       ret = 0;
       
   367 done: return ret;
       
   368 }
       
   369 
       
   370 #undef random
       
   371 
       
   372 static void assign_capacities(struct csa *csa)
       
   373 {     /* Assign a capacity to each arc. */
       
   374       struct arcs *arc_ptr = arc_list;
       
   375       int (*random)(struct csa *csa, double *);
       
   376       int i;
       
   377       /* Determine the random number generator to use. */
       
   378       switch (arc_costs.distribution)
       
   379       {  case UNIFORM:
       
   380             random = uniform;
       
   381             break;
       
   382          case EXPONENTIAL:
       
   383             random = exponential;
       
   384             break;
       
   385          default:
       
   386             xassert(csa != csa);
       
   387       }
       
   388       /* Assign capacities to grid arcs. */
       
   389       for (i = n_source + n_sink; i < n_grid_arc; i++, arc_ptr++)
       
   390          arc_ptr->u = random(csa, capacities.parameter);
       
   391       i = i - n_source - n_sink;
       
   392       /* Assign capacities to arcs to/from supernode. */
       
   393       for (; i < n_grid_arc; i++, arc_ptr++)
       
   394          arc_ptr->u = t_supply;
       
   395       /* Assign capacities to all other arcs. */
       
   396       for (; i < n_arc; i++, arc_ptr++)
       
   397          arc_ptr->u = random(csa, capacities.parameter);
       
   398       return;
       
   399 }
       
   400 
       
   401 static void assign_costs(struct csa *csa)
       
   402 {     /* Assign a cost to each arc. */
       
   403       struct arcs *arc_ptr = arc_list;
       
   404       int (*random)(struct csa *csa, double *);
       
   405       int i;
       
   406       /* A high cost assigned to arcs to/from the supernode. */
       
   407       int high_cost;
       
   408       /* The maximum cost assigned to arcs in the base grid. */
       
   409       int max_cost = 0;
       
   410       /* Determine the random number generator to use. */
       
   411       switch (arc_costs.distribution)
       
   412       {  case UNIFORM:
       
   413             random = uniform;
       
   414             break;
       
   415          case EXPONENTIAL:
       
   416             random = exponential;
       
   417             break;
       
   418          default:
       
   419             xassert(csa != csa);
       
   420       }
       
   421       /* Assign costs to arcs in the base grid. */
       
   422       for (i = n_source + n_sink; i < n_grid_arc; i++, arc_ptr++)
       
   423       {  arc_ptr->cost = random(csa, arc_costs.parameter);
       
   424          if (max_cost < arc_ptr->cost) max_cost = arc_ptr->cost;
       
   425       }
       
   426       i = i - n_source - n_sink;
       
   427       /* Assign costs to arcs to/from the super node. */
       
   428       high_cost = max_cost * 2;
       
   429       for (; i < n_grid_arc; i++, arc_ptr++)
       
   430          arc_ptr->cost = high_cost;
       
   431       /* Assign costs to all other arcs. */
       
   432       for (; i < n_arc; i++, arc_ptr++)
       
   433          arc_ptr->cost = random(csa, arc_costs.parameter);
       
   434       return;
       
   435 }
       
   436 
       
   437 static void assign_imbalance(struct csa *csa)
       
   438 {     /* Assign an imbalance to each node. */
       
   439       int total, i;
       
   440       double avg;
       
   441       struct imbalance *ptr;
       
   442       /* assign the supply nodes */
       
   443       avg = 2.0 * t_supply / n_source;
       
   444       do
       
   445       {  for (i = 1, total = t_supply, ptr = source_list + 1;
       
   446             i < n_source; i++, ptr++)
       
   447          {  ptr->supply = (int)(randy(csa) * avg + 0.5);
       
   448             total -= ptr->supply;
       
   449          }
       
   450          source_list->supply = total;
       
   451       }
       
   452       /* redo all if the assignment "overshooted" */
       
   453       while (total <= 0);
       
   454       /* assign the demand nodes */
       
   455       avg = -2.0 * t_supply / n_sink;
       
   456       do
       
   457       {  for (i = 1, total = t_supply, ptr = sink_list + 1;
       
   458             i < n_sink; i++, ptr++)
       
   459          {  ptr->supply = (int)(randy(csa) * avg - 0.5);
       
   460             total += ptr->supply;
       
   461          }
       
   462          sink_list->supply = - total;
       
   463       }
       
   464       while (total <= 0);
       
   465       return;
       
   466 }
       
   467 
       
   468 static int exponential(struct csa *csa, double lambda[1])
       
   469 {     /* Returns an "exponentially distributed" integer with parameter
       
   470          lambda. */
       
   471       return ((int)(- lambda[0] * log((double)randy(csa)) + 0.5));
       
   472 }
       
   473 
       
   474 static struct arcs *gen_additional_arcs(struct csa *csa, struct arcs
       
   475       *arc_ptr)
       
   476 {     /* Generate an arc from each source to the supernode and from
       
   477          supernode to each sink. */
       
   478       int i;
       
   479       for (i = 0; i < n_source; i++, arc_ptr++)
       
   480       {  arc_ptr->from = source_list[i].node;
       
   481          arc_ptr->to = n_node;
       
   482       }
       
   483       for (i = 0; i < n_sink; i++, arc_ptr++)
       
   484       {  arc_ptr->to = sink_list[i].node;
       
   485          arc_ptr->from = n_node;
       
   486       }
       
   487       return arc_ptr;
       
   488 }
       
   489 
       
   490 static struct arcs *gen_basic_grid(struct csa *csa, struct arcs
       
   491       *arc_ptr)
       
   492 {     /* Generate the basic grid. */
       
   493       int direction = 1, i, j, k;
       
   494       if (two_way)
       
   495       {  /* Generate an arc in each direction. */
       
   496          for (i = 1; i < n_node; i += n1)
       
   497          {  for (j = i, k = j + n1 - 1; j < k; j++)
       
   498             {  arc_ptr->from = j;
       
   499                arc_ptr->to = j + 1;
       
   500                arc_ptr++;
       
   501                arc_ptr->from = j + 1;
       
   502                arc_ptr->to = j;
       
   503                arc_ptr++;
       
   504             }
       
   505          }
       
   506          for (i = 1; i <= n1; i++)
       
   507          {  for (j = i + n1; j < n_node; j += n1)
       
   508             {  arc_ptr->from = j;
       
   509                arc_ptr->to = j - n1;
       
   510                arc_ptr++;
       
   511                arc_ptr->from = j - n1;
       
   512                arc_ptr->to = j;
       
   513                arc_ptr++;
       
   514             }
       
   515          }
       
   516       }
       
   517       else
       
   518       {  /* Generate one arc in each direction. */
       
   519          for (i = 1; i < n_node; i += n1)
       
   520          {  if (direction == 1)
       
   521                j = i;
       
   522             else
       
   523                j = i + 1;
       
   524             for (k = j + n1 - 1; j < k; j++)
       
   525             {  arc_ptr->from = j;
       
   526                arc_ptr->to = j + direction;
       
   527                arc_ptr++;
       
   528             }
       
   529             direction = - direction;
       
   530          }
       
   531          for (i = 1; i <= n1; i++)
       
   532          {  j = i + n1;
       
   533             if (direction == 1)
       
   534             {  for (; j < n_node; j += n1)
       
   535                {  arc_ptr->from = j - n1;
       
   536                   arc_ptr->to = j;
       
   537                   arc_ptr++;
       
   538                }
       
   539             }
       
   540             else
       
   541             {  for (; j < n_node; j += n1)
       
   542                {  arc_ptr->from = j - n1;
       
   543                   arc_ptr->to = j;
       
   544                   arc_ptr++;
       
   545                }
       
   546             }
       
   547             direction = - direction;
       
   548          }
       
   549       }
       
   550       return arc_ptr;
       
   551 }
       
   552 
       
   553 static void gen_more_arcs(struct csa *csa, struct arcs *arc_ptr)
       
   554 {     /* Generate random arcs to meet the specified density. */
       
   555       int i;
       
   556       double ab[2];
       
   557       ab[0] = 0.9;
       
   558       ab[1] = n_node - 0.99;  /* upper limit is n_node-1 because the
       
   559                                  supernode cannot be selected */
       
   560       for (i = n_grid_arc; i < n_arc; i++, arc_ptr++)
       
   561       {  arc_ptr->from = uniform(csa, ab);
       
   562          arc_ptr->to = uniform(csa, ab);
       
   563          if (arc_ptr->from == arc_ptr->to)
       
   564          {  arc_ptr--;
       
   565             i--;
       
   566          }
       
   567       }
       
   568       return;
       
   569 }
       
   570 
       
   571 static void generate(struct csa *csa)
       
   572 {     /* Generate a random network. */
       
   573       struct arcs *arc_ptr = arc_list;
       
   574       arc_ptr = gen_basic_grid(csa, arc_ptr);
       
   575       select_source_sinks(csa);
       
   576       arc_ptr = gen_additional_arcs(csa, arc_ptr);
       
   577       gen_more_arcs(csa, arc_ptr);
       
   578       assign_costs(csa);
       
   579       assign_capacities(csa);
       
   580       assign_imbalance(csa);
       
   581       return;
       
   582 }
       
   583 
       
   584 static void output(struct csa *csa)
       
   585 {     /* Output the network in DIMACS format. */
       
   586       struct arcs *arc_ptr;
       
   587       struct imbalance *imb_ptr;
       
   588       int i;
       
   589       if (G != NULL) goto skip;
       
   590       /* Output "c", "p" records. */
       
   591       xprintf("c generated by GRIDGEN\n");
       
   592       xprintf("c seed %d\n", seed_original);
       
   593       xprintf("c nodes %d\n", n_node);
       
   594       xprintf("c grid size %d X %d\n", n1, n2);
       
   595       xprintf("c sources %d sinks %d\n", n_source, n_sink);
       
   596       xprintf("c avg. degree %d\n", avg_degree);
       
   597       xprintf("c supply %d\n", t_supply);
       
   598       switch (arc_costs.distribution)
       
   599       {  case UNIFORM:
       
   600             xprintf("c arc costs: UNIFORM distr. min %d max %d\n",
       
   601                (int)arc_costs.parameter[0],
       
   602                (int)arc_costs.parameter[1]);
       
   603             break;
       
   604          case EXPONENTIAL:
       
   605             xprintf("c arc costs: EXPONENTIAL distr. lambda %d\n",
       
   606                (int)arc_costs.parameter[0]);
       
   607             break;
       
   608          default:
       
   609             xassert(csa != csa);
       
   610       }
       
   611       switch (capacities.distribution)
       
   612       {  case UNIFORM:
       
   613             xprintf("c arc caps :  UNIFORM distr. min %d max %d\n",
       
   614                (int)capacities.parameter[0],
       
   615                (int)capacities.parameter[1]);
       
   616             break;
       
   617          case EXPONENTIAL:
       
   618             xprintf("c arc caps :  EXPONENTIAL distr. %d lambda %d\n",
       
   619                (int)capacities.parameter[0]);
       
   620             break;
       
   621          default:
       
   622             xassert(csa != csa);
       
   623       }
       
   624 skip: if (G == NULL)
       
   625          xprintf("p min %d %d\n", n_node, n_arc);
       
   626       else
       
   627       {  glp_add_vertices(G, n_node);
       
   628          if (v_rhs >= 0)
       
   629          {  double zero = 0.0;
       
   630             for (i = 1; i <= n_node; i++)
       
   631             {  glp_vertex *v = G->v[i];
       
   632                memcpy((char *)v->data + v_rhs, &zero, sizeof(double));
       
   633             }
       
   634          }
       
   635       }
       
   636       /* Output "n node supply". */
       
   637       for (i = 0, imb_ptr = source_list; i < n_source; i++, imb_ptr++)
       
   638       {  if (G == NULL)
       
   639             xprintf("n %d %d\n", imb_ptr->node, imb_ptr->supply);
       
   640          else
       
   641          {  if (v_rhs >= 0)
       
   642             {  double temp = (double)imb_ptr->supply;
       
   643                glp_vertex *v = G->v[imb_ptr->node];
       
   644                memcpy((char *)v->data + v_rhs, &temp, sizeof(double));
       
   645             }
       
   646          }
       
   647       }
       
   648       for (i = 0, imb_ptr = sink_list; i < n_sink; i++, imb_ptr++)
       
   649       {  if (G == NULL)
       
   650             xprintf("n %d %d\n", imb_ptr->node, imb_ptr->supply);
       
   651          else
       
   652          {  if (v_rhs >= 0)
       
   653             {  double temp = (double)imb_ptr->supply;
       
   654                glp_vertex *v = G->v[imb_ptr->node];
       
   655                memcpy((char *)v->data + v_rhs, &temp, sizeof(double));
       
   656             }
       
   657          }
       
   658       }
       
   659       /* Output "a from to lowcap=0 hicap cost". */
       
   660       for (i = 0, arc_ptr = arc_list; i < n_arc; i++, arc_ptr++)
       
   661       {  if (G == NULL)
       
   662             xprintf("a %d %d 0 %d %d\n", arc_ptr->from, arc_ptr->to,
       
   663                arc_ptr->u, arc_ptr->cost);
       
   664          else
       
   665          {  glp_arc *a = glp_add_arc(G, arc_ptr->from, arc_ptr->to);
       
   666             if (a_cap >= 0)
       
   667             {  double temp = (double)arc_ptr->u;
       
   668                memcpy((char *)a->data + a_cap, &temp, sizeof(double));
       
   669             }
       
   670             if (a_cost >= 0)
       
   671             {  double temp = (double)arc_ptr->cost;
       
   672                memcpy((char *)a->data + a_cost, &temp, sizeof(double));
       
   673             }
       
   674          }
       
   675       }
       
   676       return;
       
   677 }
       
   678 
       
   679 static double randy(struct csa *csa)
       
   680 {     /* Returns a random number between 0.0 and 1.0.
       
   681          See Ward Cheney & David Kincaid, "Numerical Mathematics and
       
   682          Computing," 2Ed, pp. 335. */
       
   683       seed = 16807 * seed % 2147483647;
       
   684       if (seed < 0) seed = - seed;
       
   685       return seed * 4.6566128752459e-10;
       
   686 }
       
   687 
       
   688 static void select_source_sinks(struct csa *csa)
       
   689 {     /* Randomly select the source nodes and sink nodes. */
       
   690       int i, *int_ptr;
       
   691       int *temp_list;   /* a temporary list of nodes */
       
   692       struct imbalance *ptr;
       
   693       double ab[2];     /* parameter for random number generator */
       
   694       ab[0] = 0.9;
       
   695       ab[1] = n_node - 0.99;  /* upper limit is n_node-1 because the
       
   696                                  supernode cannot be selected */
       
   697       temp_list = xcalloc(n_node, sizeof(int));
       
   698       for (i = 0, int_ptr = temp_list; i < n_node; i++, int_ptr++)
       
   699          *int_ptr = 0;
       
   700       /* Select the source nodes. */
       
   701       for (i = 0, ptr = source_list; i < n_source; i++, ptr++)
       
   702       {  ptr->node = uniform(csa, ab);
       
   703          if (temp_list[ptr->node] == 1) /* check for duplicates */
       
   704          {  ptr--;
       
   705             i--;
       
   706          }
       
   707          else
       
   708             temp_list[ptr->node] = 1;
       
   709       }
       
   710       /* Select the sink nodes. */
       
   711       for (i = 0, ptr = sink_list; i < n_sink; i++, ptr++)
       
   712       {  ptr->node = uniform(csa, ab);
       
   713          if (temp_list[ptr->node] == 1)
       
   714          {  ptr--;
       
   715             i--;
       
   716          }
       
   717          else
       
   718             temp_list[ptr->node] = 1;
       
   719       }
       
   720       xfree(temp_list);
       
   721       return;
       
   722 }
       
   723 
       
   724 int uniform(struct csa *csa, double a[2])
       
   725 {     /* Generates an integer uniformly selected from [a[0],a[1]]. */
       
   726       return (int)((a[1] - a[0]) * randy(csa) + a[0] + 0.5);
       
   727 }
       
   728 
       
   729 /**********************************************************************/
       
   730 
       
   731 #if 0
       
   732 int main(void)
       
   733 {     int parm[1+14];
       
   734       double temp;
       
   735       scanf("%d", &parm[1]);
       
   736       scanf("%d", &parm[2]);
       
   737       scanf("%d", &parm[3]);
       
   738       scanf("%d", &parm[4]);
       
   739       scanf("%d", &parm[5]);
       
   740       scanf("%d", &parm[6]);
       
   741       scanf("%d", &parm[7]);
       
   742       scanf("%d", &parm[8]);
       
   743       scanf("%d", &parm[9]);
       
   744       if (parm[9] == 1)
       
   745       {  scanf("%d", &parm[10]);
       
   746          scanf("%d", &parm[11]);
       
   747       }
       
   748       else
       
   749       {  scanf("%le", &temp);
       
   750          parm[10] = (int)(100.0 * temp + .5);
       
   751          parm[11] = 0;
       
   752       }
       
   753       scanf("%d", &parm[12]);
       
   754       if (parm[12] == 1)
       
   755       {  scanf("%d", &parm[13]);
       
   756          scanf("%d", &parm[14]);
       
   757       }
       
   758       else
       
   759       {  scanf("%le", &temp);
       
   760          parm[13] = (int)(100.0 * temp + .5);
       
   761          parm[14] = 0;
       
   762       }
       
   763       glp_gridgen(NULL, 0, 0, 0, parm);
       
   764       return 0;
       
   765 }
       
   766 #endif
       
   767 
       
   768 /* eof */