source: trunk/Clp/examples/network.cpp @ 1559

Last change on this file since 1559 was 1559, checked in by stefan, 10 years ago

merge split branch into trunk

  • Property svn:eol-style set to native
  • Property svn:keywords set to Id
File size: 5.7 KB
Line 
1/* $Id: network.cpp 1559 2010-06-05 19:42:36Z stefan $ */
2// Copyright (C) 2002, International Business Machines
3// Corporation and others.  All Rights Reserved.
4
5
6// The point of this example is to show how to create a model without
7// using mps files.
8
9// This reads a network problem created by netgen which can be
10// downloaded from www.netlib.org/lp/generators/netgen
11// This is a generator due to Darwin Klingman
12
13#include "ClpSimplex.hpp"
14#include "ClpFactorization.hpp"
15#include "ClpNetworkMatrix.hpp"
16#include <stdio.h>
17#include <assert.h>
18#include <cmath>
19
20int main(int argc, const char *argv[])
21{
22     int numberColumns;
23     int numberRows;
24
25     FILE * fp;
26     if (argc > 1) {
27          fp = fopen(argv[1], "r");
28          if (!fp) {
29               fprintf(stderr, "Unable to open file %s\n", argv[1]);
30               exit(1);
31          }
32     } else {
33          fp = fopen("input.130", "r");
34          if (!fp) {
35               fprintf(stderr, "Unable to open file input.l30 in Samples directory\n");
36               exit(1);
37          }
38     }
39     int problem;
40     char temp[100];
41     // read and skip
42     fscanf(fp, "%s", temp);
43     assert(!strcmp(temp, "BEGIN"));
44     fscanf(fp, "%*s %*s %d %d %*s %*s %d %*s", &problem, &numberRows,
45            &numberColumns);
46     // scan down to SUPPLY
47     while (fgets(temp, 100, fp)) {
48          if (!strncmp(temp, "SUPPLY", 6))
49               break;
50     }
51     if (strncmp(temp, "SUPPLY", 6)) {
52          fprintf(stderr, "Unable to find SUPPLY\n");
53          exit(2);
54     }
55     // get space for rhs
56     double * lower = new double[numberRows];
57     double * upper = new double[numberRows];
58     int i;
59     for (i = 0; i < numberRows; i++) {
60          lower[i] = 0.0;
61          upper[i] = 0.0;
62     }
63     // ***** Remember to convert to C notation
64     while (fgets(temp, 100, fp)) {
65          int row;
66          int value;
67          if (!strncmp(temp, "ARCS", 4))
68               break;
69          sscanf(temp, "%d %d", &row, &value);
70          upper[row-1] = -value;
71          lower[row-1] = -value;
72     }
73     if (strncmp(temp, "ARCS", 4)) {
74          fprintf(stderr, "Unable to find ARCS\n");
75          exit(2);
76     }
77     // number of columns may be underestimate
78     int * head = new int[2*numberColumns];
79     int * tail = new int[2*numberColumns];
80     int * ub = new int[2*numberColumns];
81     int * cost = new int[2*numberColumns];
82     // ***** Remember to convert to C notation
83     numberColumns = 0;
84     while (fgets(temp, 100, fp)) {
85          int iHead;
86          int iTail;
87          int iUb;
88          int iCost;
89          if (!strncmp(temp, "DEMAND", 6))
90               break;
91          sscanf(temp, "%d %d %d %d", &iHead, &iTail, &iCost, &iUb);
92          iHead--;
93          iTail--;
94          head[numberColumns] = iHead;
95          tail[numberColumns] = iTail;
96          ub[numberColumns] = iUb;
97          cost[numberColumns] = iCost;
98          numberColumns++;
99     }
100     if (strncmp(temp, "DEMAND", 6)) {
101          fprintf(stderr, "Unable to find DEMAND\n");
102          exit(2);
103     }
104     // ***** Remember to convert to C notation
105     while (fgets(temp, 100, fp)) {
106          int row;
107          int value;
108          if (!strncmp(temp, "END", 3))
109               break;
110          sscanf(temp, "%d %d", &row, &value);
111          upper[row-1] = value;
112          lower[row-1] = value;
113     }
114     if (strncmp(temp, "END", 3)) {
115          fprintf(stderr, "Unable to find END\n");
116          exit(2);
117     }
118     printf("Problem %d has %d rows and %d columns\n", problem,
119            numberRows, numberColumns);
120     fclose(fp);
121     ClpSimplex  model;
122     // now build model - we have rhs so build columns - two elements
123     // per column
124
125     double * objective = new double[numberColumns];
126     double * lowerColumn = new double[numberColumns];
127     double * upperColumn = new double[numberColumns];
128
129     double * element = new double [2*numberColumns];
130     int * start = new int[numberColumns+1];
131     int * row = new int[2*numberColumns];
132     start[numberColumns] = 2 * numberColumns;
133     for (i = 0; i < numberColumns; i++) {
134          start[i] = 2 * i;
135          element[2*i] = -1.0;
136          element[2*i+1] = 1.0;
137          row[2*i] = head[i];
138          row[2*i+1] = tail[i];
139          lowerColumn[i] = 0.0;
140          upperColumn[i] = ub[i];
141          objective[i] = cost[i];
142     }
143     // Create Packed Matrix
144     CoinPackedMatrix matrix;
145     int * lengths = NULL;
146     matrix.assignMatrix(true, numberRows, numberColumns,
147                         2 * numberColumns, element, row, start, lengths);
148     ClpNetworkMatrix network(matrix);
149     // load model
150
151     model.loadProblem(network,
152                       lowerColumn, upperColumn, objective,
153                       lower, upper);
154
155     delete [] lower;
156     delete [] upper;
157     delete [] head;
158     delete [] tail;
159     delete [] ub;
160     delete [] cost;
161     delete [] objective;
162     delete [] lowerColumn;
163     delete [] upperColumn;
164     delete [] element;
165     delete [] start;
166     delete [] row;
167
168     /* The confusing flow below is in to exercise both dual and primal
169        when ClpNetworkMatrix storage used.
170
171        For practical use just one call e.g. model.dual(); would be used.
172
173        If network then factorization scheme is changed
174        to be much faster.
175
176        Still not as fast as a real network code, but more flexible
177     */
178     model.factorization()->maximumPivots(200 + model.numberRows() / 100);
179     model.factorization()->maximumPivots(1000);
180     //model.factorization()->maximumPivots(1);
181     if (model.numberRows() < 50)
182          model.messageHandler()->setLogLevel(63);
183     model.dual();
184     model.setOptimizationDirection(-1);
185     //model.messageHandler()->setLogLevel(63);
186     model.primal();
187     model.setOptimizationDirection(1);
188     model.primal();
189     return 0;
190}
Note: See TracBrowser for help on using the repository browser.