source: trunk/Test/unitTest.cpp @ 115

Last change on this file since 115 was 115, checked in by forrest, 17 years ago

OsiSimplexInterface?

  • Property svn:eol-style set to native
  • Property svn:keywords set to Author Date Id Revision
File size: 35.1 KB
Line 
1// Copyright (C) 2002, International Business Machines
2// Corporation and others.  All Rights Reserved.
3
4#include "CoinPragma.hpp"
5#include <cassert>
6#include <cstdio>
7#include <cmath>
8#include <cfloat>
9#include <string>
10#include <iostream>
11
12#include <time.h>
13
14#include "CoinMpsIO.hpp"
15#include "CoinPackedMatrix.hpp"
16#include "CoinPackedVector.hpp"
17#include "CoinHelperFunctions.hpp"
18
19#include "ClpFactorization.hpp"
20#include "ClpSimplex.hpp"
21#include "ClpDualRowSteepest.hpp"
22#include "ClpDualRowDantzig.hpp"
23#include "ClpPrimalColumnSteepest.hpp"
24#include "ClpPrimalColumnDantzig.hpp"
25#include "ClpParameters.hpp"
26
27#include "Presolve.hpp"
28#ifdef CLP_IDIOT
29#include "Idiot.hpp"
30#endif
31
32
33#include <time.h>
34#ifndef _MSC_VER
35#include <sys/times.h>
36#include <sys/resource.h>
37#include <unistd.h>
38#endif
39
40static double cpuTime()
41{
42  double cpu_temp;
43#if defined(_MSC_VER)
44  unsigned int ticksnow;        /* clock_t is same as int */
45 
46  ticksnow = (unsigned int)clock();
47 
48  cpu_temp = (double)((double)ticksnow/CLOCKS_PER_SEC);
49#else
50  struct rusage usage;
51  getrusage(RUSAGE_SELF,&usage);
52  cpu_temp = usage.ru_utime.tv_sec;
53  cpu_temp += 1.0e-6*((double) usage.ru_utime.tv_usec);
54#endif
55  return cpu_temp;
56}
57
58
59//#############################################################################
60
61#ifdef NDEBUG
62#undef NDEBUG
63#endif
64
65// Function Prototypes. Function definitions is in this file.
66void testingMessage( const char * const msg );
67
68//----------------------------------------------------------------
69// unitTest [-mpsDir=V1] [-netlibDir=V2] [-netlib]
70//
71// where:
72//   -mpsDir: directory containing mps test files
73//       Default value V1="../Mps/Sample"   
74//   -netlibDir: directory containing netlib files
75//       Default value V2="../Mps/Netlib"
76//   -netlib
77//       If specified, then netlib test set run
78//
79// All parameters are optional.
80//----------------------------------------------------------------
81int mainTest (int argc, const char *argv[],bool doDual,
82              ClpSimplex empty, bool doPresolve,int doIdiot)
83{
84  int i;
85
86  // define valid parameter keywords
87  std::set<std::string> definedKeyWords;
88  definedKeyWords.insert("-mpsDir");
89  definedKeyWords.insert("-netlibDir");
90  definedKeyWords.insert("-netlib");
91
92  // Create a map of parameter keys and associated data
93  std::map<std::string,std::string> parms;
94  for ( i=1; i<argc; i++ ) {
95    std::string parm(argv[i]);
96    std::string key,value;
97    unsigned int  eqPos = parm.find('=');
98
99    // Does parm contain and '='
100    if ( eqPos==std::string::npos ) {
101      //Parm does not contain '='
102      key = parm;
103    }
104    else {
105      key=parm.substr(0,eqPos);
106      value=parm.substr(eqPos+1);
107    }
108
109    // Is specifed key valid?
110    if ( definedKeyWords.find(key) == definedKeyWords.end() ) {
111      // invalid key word.
112      // Write help text
113      std::cerr <<"Undefined parameter \"" <<key <<"\".\n";
114      std::cerr <<"Correct usage: \n";
115      std::cerr <<"  unitTest [-mpsDir=V1] [-netlibDir=V2] [-netlib]\n";
116      std::cerr <<"  where:\n";
117      std::cerr <<"    -mpsDir: directory containing mps test files\n";
118      std::cerr <<"        Default value V1=\"../Mps/Sample\"\n";
119      std::cerr <<"    -netlibDir: directory containing netlib files\n";
120      std::cerr <<"        Default value V2=\"../Mps/Netlib\"\n";
121      std::cerr <<"    -netlib\n";
122      std::cerr <<"        If specified, then netlib testset run.\n";
123      return 1;
124    }
125    parms[key]=value;
126  }
127 
128  const char dirsep =  CoinFindDirSeparator();
129  // Set directory containing mps data files.
130  std::string mpsDir;
131  if (parms.find("-mpsDir") != parms.end())
132    mpsDir=parms["-mpsDir"] + dirsep;
133  else 
134    mpsDir = dirsep == '/' ? "../Mps/Sample/" : "..\\Mps\\Sample\\";
135 
136  // Set directory containing netlib data files.
137  std::string netlibDir;
138  if (parms.find("-netlibDir") != parms.end())
139    netlibDir=parms["-netlibDir"] + dirsep;
140  else 
141    netlibDir = dirsep == '/' ? "../Mps/Netlib/" : "..\\Mps\\Netlib\\";
142
143  testingMessage( "Testing ClpSimplex\n" );
144  ClpSimplexUnitTest(mpsDir,netlibDir);
145  if (parms.find("-netlib") != parms.end())
146  {
147    unsigned int m;
148   
149    // Define test problems:
150    //   mps names,
151    //   maximization or minimization,
152    //   Number of rows and columns in problem, and
153    //   objective function value
154    std::vector<std::string> mpsName;
155    std::vector<bool> min;
156    std::vector<int> nRows;
157    std::vector<int> nCols;
158    std::vector<double> objValue;
159    std::vector<double> objValueTol;
160    mpsName.push_back("25fv47");
161    min.push_back(true);
162    nRows.push_back(822);
163    nCols.push_back(1571);
164    objValueTol.push_back(1.E-10);
165    objValue.push_back(5.5018458883E+03);
166
167    mpsName.push_back("80bau3b");min.push_back(true);nRows.push_back(2263);nCols.push_back(9799);objValueTol.push_back(1.e-10);objValue.push_back(9.8722419241E+05);
168    mpsName.push_back("blend");min.push_back(true);nRows.push_back(75);nCols.push_back(83);objValueTol.push_back(1.e-10);objValue.push_back(-3.0812149846e+01);
169    mpsName.push_back("pilotnov");min.push_back(true);nRows.push_back(976);nCols.push_back(2172);objValueTol.push_back(1.e-10);objValue.push_back(-4.4972761882e+03);
170    mpsName.push_back("maros-r7");min.push_back(true);nRows.push_back(3137);nCols.push_back(9408);objValueTol.push_back(1.e-10);objValue.push_back(1.4971851665e+06);
171   
172    mpsName.push_back("pilot");min.push_back(true);nRows.push_back(1442);nCols.push_back(3652);objValueTol.push_back(1.e-5);objValue.push_back(/*-5.5740430007e+02*/-557.48972927292);
173    mpsName.push_back("pilot4");min.push_back(true);nRows.push_back(411);nCols.push_back(1000);objValueTol.push_back(1.e-8);objValue.push_back(-2.5811392641e+03);
174    mpsName.push_back("pilot87");min.push_back(true);nRows.push_back(2031);nCols.push_back(4883);objValueTol.push_back(1.e-4);objValue.push_back(3.0171072827e+02);
175    mpsName.push_back("adlittle");min.push_back(true);nRows.push_back(57);nCols.push_back(97);objValueTol.push_back(1.e-10);objValue.push_back(2.2549496316e+05);
176    mpsName.push_back("afiro");min.push_back(true);nRows.push_back(28);nCols.push_back(32);objValueTol.push_back(1.e-10);objValue.push_back(-4.6475314286e+02);
177    mpsName.push_back("agg");min.push_back(true);nRows.push_back(489);nCols.push_back(163);objValueTol.push_back(1.e-10);objValue.push_back(-3.5991767287e+07);
178    mpsName.push_back("agg2");min.push_back(true);nRows.push_back(517);nCols.push_back(302);objValueTol.push_back(1.e-10);objValue.push_back(-2.0239252356e+07);
179    mpsName.push_back("agg3");min.push_back(true);nRows.push_back(517);nCols.push_back(302);objValueTol.push_back(1.e-10);objValue.push_back(1.0312115935e+07);
180    mpsName.push_back("bandm");min.push_back(true);nRows.push_back(306);nCols.push_back(472);objValueTol.push_back(1.e-10);objValue.push_back(-1.5862801845e+02);
181    mpsName.push_back("beaconfd");min.push_back(true);nRows.push_back(174);nCols.push_back(262);objValueTol.push_back(1.e-10);objValue.push_back(3.3592485807e+04);
182    mpsName.push_back("bnl1");min.push_back(true);nRows.push_back(644);nCols.push_back(1175);objValueTol.push_back(1.e-10);objValue.push_back(1.9776295615E+03);
183    mpsName.push_back("bnl2");min.push_back(true);nRows.push_back(2325);nCols.push_back(3489);objValueTol.push_back(1.e-10);objValue.push_back(1.8112365404e+03);
184    mpsName.push_back("boeing1");min.push_back(true);nRows.push_back(/*351*/352);nCols.push_back(384);objValueTol.push_back(1.e-10);objValue.push_back(-3.3521356751e+02);
185    mpsName.push_back("boeing2");min.push_back(true);nRows.push_back(167);nCols.push_back(143);objValueTol.push_back(1.e-10);objValue.push_back(-3.1501872802e+02);
186    mpsName.push_back("bore3d");min.push_back(true);nRows.push_back(234);nCols.push_back(315);objValueTol.push_back(1.e-10);objValue.push_back(1.3730803942e+03);
187    mpsName.push_back("brandy");min.push_back(true);nRows.push_back(221);nCols.push_back(249);objValueTol.push_back(1.e-10);objValue.push_back(1.5185098965e+03);
188    mpsName.push_back("capri");min.push_back(true);nRows.push_back(272);nCols.push_back(353);objValueTol.push_back(1.e-10);objValue.push_back(2.6900129138e+03);
189    mpsName.push_back("cycle");min.push_back(true);nRows.push_back(1904);nCols.push_back(2857);objValueTol.push_back(1.e-9);objValue.push_back(-5.2263930249e+00);
190    mpsName.push_back("czprob");min.push_back(true);nRows.push_back(930);nCols.push_back(3523);objValueTol.push_back(1.e-10);objValue.push_back(2.1851966989e+06);
191    mpsName.push_back("d2q06c");min.push_back(true);nRows.push_back(2172);nCols.push_back(5167);objValueTol.push_back(1.e-7);objValue.push_back(122784.21557456);
192    mpsName.push_back("d6cube");min.push_back(true);nRows.push_back(416);nCols.push_back(6184);objValueTol.push_back(1.e-7);objValue.push_back(3.1549166667e+02);
193    mpsName.push_back("degen2");min.push_back(true);nRows.push_back(445);nCols.push_back(534);objValueTol.push_back(1.e-10);objValue.push_back(-1.4351780000e+03);
194    mpsName.push_back("degen3");min.push_back(true);nRows.push_back(1504);nCols.push_back(1818);objValueTol.push_back(1.e-10);objValue.push_back(-9.8729400000e+02);
195    mpsName.push_back("dfl001");min.push_back(true);nRows.push_back(6072);nCols.push_back(12230);objValueTol.push_back(1.e-5);objValue.push_back(1.1266396047E+07);
196    mpsName.push_back("e226");min.push_back(true);nRows.push_back(224);nCols.push_back(282);objValueTol.push_back(1.e-10);objValue.push_back(-1.8751929066e+01+7.113); // The correct answer includes -7.113 term. This is a constant in the objective function. See line 1683 of the mps file.
197    mpsName.push_back("etamacro");min.push_back(true);nRows.push_back(401);nCols.push_back(688);objValueTol.push_back(1.e-6);objValue.push_back(-7.5571521774e+02 );
198    mpsName.push_back("fffff800");min.push_back(true);nRows.push_back(525);nCols.push_back(854);objValueTol.push_back(1.e-6);objValue.push_back(5.5567961165e+05);
199    mpsName.push_back("finnis");min.push_back(true);nRows.push_back(498);nCols.push_back(614);objValueTol.push_back(1.e-6);objValue.push_back(1.7279096547e+05);
200    mpsName.push_back("fit1d");min.push_back(true);nRows.push_back(25);nCols.push_back(1026);objValueTol.push_back(1.e-10);objValue.push_back(-9.1463780924e+03);
201    mpsName.push_back("fit1p");min.push_back(true);nRows.push_back(628);nCols.push_back(1677);objValueTol.push_back(1.e-10);objValue.push_back(9.1463780924e+03);
202    mpsName.push_back("fit2d");min.push_back(true);nRows.push_back(26);nCols.push_back(10500);objValueTol.push_back(1.e-10);objValue.push_back(-6.8464293294e+04);
203    mpsName.push_back("fit2p");min.push_back(true);nRows.push_back(3001);nCols.push_back(13525);objValueTol.push_back(1.e-9);objValue.push_back(6.8464293232e+04);
204    mpsName.push_back("forplan");min.push_back(true);nRows.push_back(162);nCols.push_back(421);objValueTol.push_back(1.e-6);objValue.push_back(-6.6421873953e+02);
205    mpsName.push_back("ganges");min.push_back(true);nRows.push_back(1310);nCols.push_back(1681);objValueTol.push_back(1.e-5);objValue.push_back(-1.0958636356e+05);
206    mpsName.push_back("gfrd-pnc");min.push_back(true);nRows.push_back(617);nCols.push_back(1092);objValueTol.push_back(1.e-10);objValue.push_back(6.9022359995e+06);
207    mpsName.push_back("greenbea");min.push_back(true);nRows.push_back(2393);nCols.push_back(5405);objValueTol.push_back(1.e-10);objValue.push_back(/*-7.2462405908e+07*/-72555248.129846);
208    mpsName.push_back("greenbeb");min.push_back(true);nRows.push_back(2393);nCols.push_back(5405);objValueTol.push_back(1.e-10);objValue.push_back(/*-4.3021476065e+06*/-4302260.2612066);
209    mpsName.push_back("grow15");min.push_back(true);nRows.push_back(301);nCols.push_back(645);objValueTol.push_back(1.e-10);objValue.push_back(-1.0687094129e+08);
210    mpsName.push_back("grow22");min.push_back(true);nRows.push_back(441);nCols.push_back(946);objValueTol.push_back(1.e-10);objValue.push_back(-1.6083433648e+08);
211    mpsName.push_back("grow7");min.push_back(true);nRows.push_back(141);nCols.push_back(301);objValueTol.push_back(1.e-10);objValue.push_back(-4.7787811815e+07);
212    mpsName.push_back("israel");min.push_back(true);nRows.push_back(175);nCols.push_back(142);objValueTol.push_back(1.e-10);objValue.push_back(-8.9664482186e+05);
213    mpsName.push_back("kb2");min.push_back(true);nRows.push_back(44);nCols.push_back(41);objValueTol.push_back(1.e-10);objValue.push_back(-1.7499001299e+03);
214    mpsName.push_back("lotfi");min.push_back(true);nRows.push_back(154);nCols.push_back(308);objValueTol.push_back(1.e-10);objValue.push_back(-2.5264706062e+01);
215    mpsName.push_back("maros");min.push_back(true);nRows.push_back(847);nCols.push_back(1443);objValueTol.push_back(1.e-10);objValue.push_back(-5.8063743701e+04);
216    mpsName.push_back("modszk1");min.push_back(true);nRows.push_back(688);nCols.push_back(1620);objValueTol.push_back(1.e-10);objValue.push_back(3.2061972906e+02);
217    mpsName.push_back("nesm");min.push_back(true);nRows.push_back(663);nCols.push_back(2923);objValueTol.push_back(1.e-5);objValue.push_back(1.4076073035e+07);
218    mpsName.push_back("perold");min.push_back(true);nRows.push_back(626);nCols.push_back(1376);objValueTol.push_back(1.e-6);objValue.push_back(-9.3807580773e+03);
219    //mpsName.push_back("qap12");min.push_back(true);nRows.push_back(3193);nCols.push_back(8856);objValueTol.push_back(1.e-6);objValue.push_back(5.2289435056e+02);
220    //mpsName.push_back("qap15");min.push_back(true);nRows.push_back(6331);nCols.push_back(22275);objValueTol.push_back(1.e-10);objValue.push_back(1.0409940410e+03);
221    mpsName.push_back("recipe");min.push_back(true);nRows.push_back(92);nCols.push_back(180);objValueTol.push_back(1.e-10);objValue.push_back(-2.6661600000e+02);
222    mpsName.push_back("sc105");min.push_back(true);nRows.push_back(106);nCols.push_back(103);objValueTol.push_back(1.e-10);objValue.push_back(-5.2202061212e+01);
223    mpsName.push_back("sc205");min.push_back(true);nRows.push_back(206);nCols.push_back(203);objValueTol.push_back(1.e-10);objValue.push_back(-5.2202061212e+01);
224    mpsName.push_back("sc50a");min.push_back(true);nRows.push_back(51);nCols.push_back(48);objValueTol.push_back(1.e-10);objValue.push_back(-6.4575077059e+01);
225    mpsName.push_back("sc50b");min.push_back(true);nRows.push_back(51);nCols.push_back(48);objValueTol.push_back(1.e-10);objValue.push_back(-7.0000000000e+01);
226    mpsName.push_back("scagr25");min.push_back(true);nRows.push_back(472);nCols.push_back(500);objValueTol.push_back(1.e-10);objValue.push_back(-1.4753433061e+07);
227    mpsName.push_back("scagr7");min.push_back(true);nRows.push_back(130);nCols.push_back(140);objValueTol.push_back(1.e-6);objValue.push_back(-2.3313892548e+06);
228    mpsName.push_back("scfxm1");min.push_back(true);nRows.push_back(331);nCols.push_back(457);objValueTol.push_back(1.e-10);objValue.push_back(1.8416759028e+04);
229    mpsName.push_back("scfxm2");min.push_back(true);nRows.push_back(661);nCols.push_back(914);objValueTol.push_back(1.e-10);objValue.push_back(3.6660261565e+04);
230    mpsName.push_back("scfxm3");min.push_back(true);nRows.push_back(991);nCols.push_back(1371);objValueTol.push_back(1.e-10);objValue.push_back(5.4901254550e+04);
231    mpsName.push_back("scorpion");min.push_back(true);nRows.push_back(389);nCols.push_back(358);objValueTol.push_back(1.e-10);objValue.push_back(1.8781248227e+03);
232    mpsName.push_back("scrs8");min.push_back(true);nRows.push_back(491);nCols.push_back(1169);objValueTol.push_back(1.e-5);objValue.push_back(9.0429998619e+02);
233    mpsName.push_back("scsd1");min.push_back(true);nRows.push_back(78);nCols.push_back(760);objValueTol.push_back(1.e-10);objValue.push_back(8.6666666743e+00);
234    mpsName.push_back("scsd6");min.push_back(true);nRows.push_back(148);nCols.push_back(1350);objValueTol.push_back(1.e-10);objValue.push_back(5.0500000078e+01);
235    mpsName.push_back("scsd8");min.push_back(true);nRows.push_back(398);nCols.push_back(2750);objValueTol.push_back(1.e-10);objValue.push_back(9.0499999993e+02);
236    mpsName.push_back("sctap1");min.push_back(true);nRows.push_back(301);nCols.push_back(480);objValueTol.push_back(1.e-10);objValue.push_back(1.4122500000e+03);
237    mpsName.push_back("sctap2");min.push_back(true);nRows.push_back(1091);nCols.push_back(1880);objValueTol.push_back(1.e-10);objValue.push_back(1.7248071429e+03);
238    mpsName.push_back("sctap3");min.push_back(true);nRows.push_back(1481);nCols.push_back(2480);objValueTol.push_back(1.e-10);objValue.push_back(1.4240000000e+03);
239    mpsName.push_back("seba");min.push_back(true);nRows.push_back(516);nCols.push_back(1028);objValueTol.push_back(1.e-10);objValue.push_back(1.5711600000e+04);
240    mpsName.push_back("share1b");min.push_back(true);nRows.push_back(118);nCols.push_back(225);objValueTol.push_back(1.e-10);objValue.push_back(-7.6589318579e+04);
241    mpsName.push_back("share2b");min.push_back(true);nRows.push_back(97);nCols.push_back(79);objValueTol.push_back(1.e-10);objValue.push_back(-4.1573224074e+02);
242    mpsName.push_back("shell");min.push_back(true);nRows.push_back(537);nCols.push_back(1775);objValueTol.push_back(1.e-10);objValue.push_back(1.2088253460e+09);
243    mpsName.push_back("ship04l");min.push_back(true);nRows.push_back(403);nCols.push_back(2118);objValueTol.push_back(1.e-10);objValue.push_back(1.7933245380e+06);
244    mpsName.push_back("ship04s");min.push_back(true);nRows.push_back(403);nCols.push_back(1458);objValueTol.push_back(1.e-10);objValue.push_back(1.7987147004e+06);
245    mpsName.push_back("ship08l");min.push_back(true);nRows.push_back(779);nCols.push_back(4283);objValueTol.push_back(1.e-10);objValue.push_back(1.9090552114e+06);
246    mpsName.push_back("ship08s");min.push_back(true);nRows.push_back(779);nCols.push_back(2387);objValueTol.push_back(1.e-10);objValue.push_back(1.9200982105e+06);
247    mpsName.push_back("ship12l");min.push_back(true);nRows.push_back(1152);nCols.push_back(5427);objValueTol.push_back(1.e-10);objValue.push_back(1.4701879193e+06);
248    mpsName.push_back("ship12s");min.push_back(true);nRows.push_back(1152);nCols.push_back(2763);objValueTol.push_back(1.e-10);objValue.push_back(1.4892361344e+06);
249    mpsName.push_back("sierra");min.push_back(true);nRows.push_back(1228);nCols.push_back(2036);objValueTol.push_back(1.e-10);objValue.push_back(1.5394362184e+07);
250    mpsName.push_back("stair");min.push_back(true);nRows.push_back(357);nCols.push_back(467);objValueTol.push_back(1.e-10);objValue.push_back(-2.5126695119e+02);
251    mpsName.push_back("standata");min.push_back(true);nRows.push_back(360);nCols.push_back(1075);objValueTol.push_back(1.e-10);objValue.push_back(1.2576995000e+03);
252    //mpsName.push_back("standgub");min.push_back(true);nRows.push_back(362);nCols.push_back(1184);objValueTol.push_back(1.e-10);objValue.push_back(1257.6995);
253    mpsName.push_back("standmps");min.push_back(true);nRows.push_back(468);nCols.push_back(1075);objValueTol.push_back(1.e-10);objValue.push_back(1.4060175000E+03); 
254    mpsName.push_back("stocfor1");min.push_back(true);nRows.push_back(118);nCols.push_back(111);objValueTol.push_back(1.e-10);objValue.push_back(-4.1131976219E+04);
255    mpsName.push_back("stocfor2");min.push_back(true);nRows.push_back(2158);nCols.push_back(2031);objValueTol.push_back(1.e-10);objValue.push_back(-3.9024408538e+04);
256    //mpsName.push_back("stocfor3");min.push_back(true);nRows.push_back(16676);nCols.push_back(15695);objValueTol.push_back(1.e-10);objValue.push_back(-3.9976661576e+04);
257    //mpsName.push_back("truss");min.push_back(true);nRows.push_back(1001);nCols.push_back(8806);objValueTol.push_back(1.e-10);objValue.push_back(4.5881584719e+05);
258    mpsName.push_back("tuff");min.push_back(true);nRows.push_back(334);nCols.push_back(587);objValueTol.push_back(1.e-10);objValue.push_back(2.9214776509e-01);
259    mpsName.push_back("vtpbase");min.push_back(true);nRows.push_back(199);nCols.push_back(203);objValueTol.push_back(1.e-10);objValue.push_back(1.2983146246e+05);
260    mpsName.push_back("wood1p");min.push_back(true);nRows.push_back(245);nCols.push_back(2594);objValueTol.push_back(1.e-10);objValue.push_back(1.4429024116e+00);
261    mpsName.push_back("woodw");min.push_back(true);nRows.push_back(1099);nCols.push_back(8405);objValueTol.push_back(1.e-10);objValue.push_back(1.3044763331E+00);
262
263    double timeTaken =0.0;
264  // Loop once for each Mps File
265    for (m=0; m<mpsName.size(); m++ ) {
266      std::cerr <<"  processing mps file: " <<mpsName[m] 
267                <<" (" <<m+1 <<" out of " <<mpsName.size() <<")" <<std::endl;
268   
269      // Read data mps file,
270      std::string fn = netlibDir+mpsName[m];
271      CoinMpsIO mps;
272      mps.readMps(fn.c_str(),"mps");
273      double time1 = cpuTime();
274      ClpSimplex solution=empty;
275      solution.loadProblem(*mps.getMatrixByCol(),mps.getColLower(),
276                           mps.getColUpper(),
277                           mps.getObjCoefficients(),
278                           mps.getRowLower(),mps.getRowUpper());
279
280      solution.setDblParam(ClpObjOffset,mps.objectiveOffset());
281#if 0
282      solution.setOptimizationDirection(-1);
283      {
284        int j;
285        double * obj = solution.objective();
286        int n=solution.numberColumns();
287        for (j=0;j<n;j++)
288          obj[j] *= -1.0;
289      }
290#endif
291      if (doPresolve) {
292#ifdef USE_PRESOLVE
293        Presolve pinfo;
294        ClpSimplex * model2 = pinfo.presolvedModel(solution,1.0e-8);
295        // change from 200 (unless user has changed)
296        if (model2->factorization()->maximumPivots()==200)
297          model2->factorization()->maximumPivots(100+model2->numberRows()/100);
298        if (doDual) {
299          // faster if bounds tightened
300          int numberInfeasibilities = model2->tightenPrimalBounds();
301          if (numberInfeasibilities)
302            std::cout<<"** Analysis indicates model infeasible"
303                     <<std::endl;
304          if (doIdiot<0)
305            model2->crash(1000,2);
306          model2->dual();
307        } else {
308#ifdef CLP_IDIOT
309          if (doIdiot>0) {
310            Idiot info(*model2);
311            info.crash(doIdiot);
312          }
313#endif
314          model2->primal(1);
315        }
316        pinfo.postsolve(true);
317       
318        delete model2;
319#if 1
320        printf("Resolving from postsolved model\n");
321        // later try without (1) and check duals before solve
322        solution.primal(1);
323        if (solution.numberIterations())
324          printf("****** iterated %d\n",solution.numberIterations());
325        /*solution.checkSolution();
326        printf("%g dual %g(%d) Primal %g(%d)\n",
327               solution.objectiveValue(),
328               solution.sumDualInfeasibilities(),
329               solution.numberDualInfeasibilities(),
330               solution.sumPrimalInfeasibilities(),
331               solution.numberPrimalInfeasibilities());*/
332#endif
333        if (0) {
334          Presolve pinfoA;
335          model2 = pinfoA.presolvedModel(solution,1.0e-8);
336
337          printf("Resolving from presolved optimal solution\n");
338          model2->primal(1);
339               
340          delete model2;
341        }
342#else
343        if (doDual) {
344          if (doIdiot<0)
345            solution.crash(1000,2);
346          solution.dual();
347        } else {
348#ifdef CLP_IDIOT
349          if (doIdiot>0) {
350            Idiot info(solution);
351            info.crash(doIdiot);
352          }
353#endif
354          solution.primal(1);
355        }
356#endif
357      } else {
358        if (doDual) {
359          if (doIdiot<0)
360            solution.crash(1000,2);
361          solution.dual();
362        } else {
363#ifdef CLP_IDIOT
364          if (doIdiot>0) {
365            Idiot info(solution);
366            info.crash(doIdiot);
367          }
368#endif
369          solution.primal(1);
370        }
371      }
372      double time2 = cpuTime()-time1;
373      timeTaken += time2;
374      printf("Took %g seconds\n",time2);
375      // Test objective solution value
376      {
377        double soln = solution.objectiveValue();
378        CoinRelFltEq eq(objValueTol[m]);
379        std::cerr <<soln <<",  " <<objValue[m] <<" diff "<<
380          soln-objValue[m]<<std::endl;
381        if(!eq(soln,objValue[m]))
382          printf("** difference fails\n");
383      }
384    }
385    printf("Total time %g seconds\n",timeTaken);
386  }
387  else {
388    testingMessage( "***Skipped Testing on netlib    ***\n" );
389    testingMessage( "***use -netlib to test class***\n" );
390  }
391 
392  testingMessage( "All tests completed successfully\n" );
393  return 0;
394}
395
396 
397// Display message on stdout and stderr
398void testingMessage( const char * const msg )
399{
400  std::cerr <<msg;
401  //cout <<endl <<"*****************************************"
402  //     <<endl <<msg <<endl;
403}
404
405//--------------------------------------------------------------------------
406// test factorization methods and simplex method
407void
408ClpSimplexUnitTest(const std::string & mpsDir,
409                   const std::string & netlibDir)
410{
411 
412  CoinRelFltEq eq(0.000001);
413
414  {
415    ClpSimplex solution;
416 
417    // matrix data
418    //deliberate hiccup of 2 between 0 and 1
419    CoinBigIndex start[5]={0,4,7,8,9};
420    int length[5]={2,3,1,1,1};
421    int rows[11]={0,2,-1,-1,0,1,2,0,1,2};
422    double elements[11]={7.0,2.0,1.0e10,1.0e10,-2.0,1.0,-2.0,1,1,1};
423    CoinPackedMatrix matrix(true,3,5,8,elements,rows,start,length);
424   
425    // rim data
426    double objective[7]={-4.0,1.0,0.0,0.0,0.0,0.0,0.0};
427    double rowLower[5]={14.0,3.0,3.0,1.0e10,1.0e10};
428    double rowUpper[5]={14.0,3.0,3.0,-1.0e10,-1.0e10};
429    double colLower[7]={0.0,0.0,0.0,0.0,0.0,0.0,0.0};
430    double colUpper[7]={100.0,100.0,100.0,100.0,100.0,100.0,100.0};
431   
432    // basis 1
433    int rowBasis1[3]={-1,-1,-1};
434    int colBasis1[5]={1,1,-1,-1,1};
435    solution.loadProblem(matrix,colLower,colUpper,objective,
436                         rowLower,rowUpper);
437    int i;
438    solution.createStatus();
439    for (i=0;i<3;i++) {
440      if (rowBasis1[i]<0) {
441        solution.setRowStatus(i,ClpSimplex::atLowerBound);
442      } else {
443        solution.setRowStatus(i,ClpSimplex::basic);
444      }
445    }
446    for (i=0;i<5;i++) {
447      if (colBasis1[i]<0) {
448        solution.setColumnStatus(i,ClpSimplex::atLowerBound);
449      } else {
450        solution.setColumnStatus(i,ClpSimplex::basic);
451      }
452    }
453    solution.setLogLevel(3+4+8+16+32);
454    solution.primal();
455    for (i=0;i<3;i++) {
456      if (rowBasis1[i]<0) {
457        solution.setRowStatus(i,ClpSimplex::atLowerBound);
458      } else {
459        solution.setRowStatus(i,ClpSimplex::basic);
460      }
461    }
462    for (i=0;i<5;i++) {
463      if (colBasis1[i]<0) {
464        solution.setColumnStatus(i,ClpSimplex::atLowerBound);
465      } else {
466        solution.setColumnStatus(i,ClpSimplex::basic);
467      }
468    }
469    // intricate stuff does not work with scaling
470    solution.scaling(0);
471    assert(!solution.factorize ( ));
472    const double * colsol = solution.primalColumnSolution();
473    const double * rowsol = solution.primalRowSolution();
474    solution.getSolution(rowsol,colsol);
475    double colsol1[5]={20.0/7.0,3.0,0.0,0.0,23.0/7.0};
476    for (i=0;i<5;i++) {
477      assert(eq(colsol[i],colsol1[i]));
478    }
479    // now feed in again without actually doing factorization
480    ClpFactorization factorization2 = *solution.factorization();
481    ClpSimplex solution2 = solution;
482    solution2.setFactorization(factorization2);
483    solution2.createStatus();
484    for (i=0;i<3;i++) {
485      if (rowBasis1[i]<0) {
486        solution2.setRowStatus(i,ClpSimplex::atLowerBound);
487      } else {
488        solution2.setRowStatus(i,ClpSimplex::basic);
489      }
490    }
491    for (i=0;i<5;i++) {
492      if (colBasis1[i]<0) {
493        solution2.setColumnStatus(i,ClpSimplex::atLowerBound);
494      } else {
495        solution2.setColumnStatus(i,ClpSimplex::basic);
496      }
497    }
498    // intricate stuff does not work with scaling
499    solution2.scaling(0);
500    solution2.getSolution(rowsol,colsol);
501    colsol = solution2.primalColumnSolution();
502    rowsol = solution2.primalRowSolution();
503    for (i=0;i<5;i++) {
504      assert(eq(colsol[i],colsol1[i]));
505    }
506    solution2.setDualBound(0.1);
507    solution2.dual();
508    objective[2]=-1.0;
509    objective[3]=-0.5;
510    objective[4]=10.0;
511    solution.dual();
512    for (i=0;i<3;i++) {
513      rowLower[i]=-1.0e50;
514      colUpper[i+2]=0.0;
515    }
516    solution.setLogLevel(3);
517    solution.dual();
518    double rowObjective[]={1.0,0.5,-10.0};
519    solution.loadProblem(matrix,colLower,colUpper,objective,
520                         rowLower,rowUpper,rowObjective);
521    solution.dual();
522  }
523  {   
524    CoinMpsIO m;
525    std::string fn = mpsDir+"exmip1";
526    m.readMps(fn.c_str(),"mps");
527    ClpSimplex solution;
528    solution.loadProblem(*m.getMatrixByCol(),m.getColLower(),m.getColUpper(),
529                         m.getObjCoefficients(),
530                         m.getRowLower(),m.getRowUpper());
531    solution.dual();
532  }
533  // test steepest edge
534  {   
535    CoinMpsIO m;
536    std::string fn = netlibDir+"finnis";
537    m.readMps(fn.c_str(),"mps");
538    ClpModel model;
539    model.loadProblem(*m.getMatrixByCol(),m.getColLower(),
540                    m.getColUpper(),
541                    m.getObjCoefficients(),
542                    m.getRowLower(),m.getRowUpper());
543    ClpSimplex solution(model);
544
545    solution.scaling(1); 
546    solution.setDualBound(1.0e8);
547    //solution.factorization()->maximumPivots(1);
548    //solution.setLogLevel(3);
549    solution.setDualTolerance(1.0e-7);
550    // set objective sense,
551    ClpDualRowSteepest steep;
552    solution.setDualRowPivotAlgorithm(steep);
553    solution.setDblParam(ClpObjOffset,m.objectiveOffset());
554    solution.dual();
555  }
556  // test normal solution
557  {   
558    CoinMpsIO m;
559    std::string fn = netlibDir+"afiro";
560    m.readMps(fn.c_str(),"mps");
561    ClpSimplex solution;
562    ClpModel model;
563    // do twice - without and with scaling
564    int iPass;
565    for (iPass=0;iPass<2;iPass++) {
566      // explicit row objective for testing
567      int nr = m.getNumRows();
568      double * rowObj = new double[nr];
569      CoinFillN(rowObj,nr,0.0);
570      model.loadProblem(*m.getMatrixByCol(),m.getColLower(),m.getColUpper(),
571                      m.getObjCoefficients(),
572                      m.getRowLower(),m.getRowUpper(),rowObj);
573      delete [] rowObj;
574      solution = ClpSimplex(model);
575      if (iPass) {
576        solution.scaling();
577      }
578      solution.dual();
579      solution.dual();
580      // test optimal
581      assert (solution.status()==0);
582      int numberColumns = solution.numberColumns();
583      int numberRows = solution.numberRows();
584      CoinPackedVector colsol(numberColumns,solution.primalColumnSolution());
585      double * objective = solution.objective();
586      double objValue = colsol.dotProduct(objective);
587      CoinRelFltEq eq(1.0e-8);
588      assert(eq(objValue,-4.6475314286e+02));
589      double * lower = solution.columnLower();
590      double * upper = solution.columnUpper();
591      double * sol = solution.primalColumnSolution();
592      double * result = new double[numberColumns];
593      CoinFillN ( result, numberColumns,0.0);
594      solution.matrix()->transposeTimes(solution.dualRowSolution(), result);
595      int iRow , iColumn;
596      // see if feasible and dual feasible
597      for (iColumn=0;iColumn<numberColumns;iColumn++) {
598        double value = sol[iColumn];
599        assert(value<upper[iColumn]+1.0e-8);
600        assert(value>lower[iColumn]-1.0e-8);
601        value = objective[iColumn]-result[iColumn];
602        assert (value>-1.0e-5);
603        if (sol[iColumn]>1.0e-5)
604          assert (value<1.0e-5);
605      }
606      delete [] result;
607      result = new double[numberRows];
608      CoinFillN ( result, numberRows,0.0);
609      solution.matrix()->times(colsol, result);
610      lower = solution.rowLower();
611      upper = solution.rowUpper();
612      sol = solution.primalRowSolution();
613      for (iRow=0;iRow<numberRows;iRow++) {
614        double value = result[iRow];
615        assert(eq(value,sol[iRow]));
616        assert(value<upper[iRow]+1.0e-8);
617        assert(value>lower[iRow]-1.0e-8);
618      }
619      delete [] result;
620      // test row objective
621      double * rowObjective = solution.rowObjective();
622      CoinDisjointCopyN(solution.dualRowSolution(),numberRows,rowObjective);
623      CoinDisjointCopyN(solution.dualColumnSolution(),numberColumns,objective);
624      // this sets up all slack basis
625      solution.createStatus();
626      solution.dual();
627      CoinFillN(rowObjective,numberRows,0.0);
628      CoinDisjointCopyN(m.getObjCoefficients(),numberColumns,objective);
629      solution.dual();
630    }
631  }
632  // test unbounded
633  {   
634    CoinMpsIO m;
635    std::string fn = netlibDir+"brandy";
636    m.readMps(fn.c_str(),"mps");
637    ClpSimplex solution;
638    // do twice - without and with scaling
639    int iPass;
640    for (iPass=0;iPass<2;iPass++) {
641      solution.loadProblem(*m.getMatrixByCol(),m.getColLower(),m.getColUpper(),
642                      m.getObjCoefficients(),
643                      m.getRowLower(),m.getRowUpper());
644      if (iPass)
645        solution.scaling();
646      solution.setOptimizationDirection(-1);
647      // test unbounded and ray
648#ifdef DUAL
649      solution.setDualBound(100.0);
650      solution.dual();
651#else
652      solution.primal();
653#endif
654      assert (solution.status()==2);
655      int numberColumns = solution.numberColumns();
656      int numberRows = solution.numberRows();
657      double * lower = solution.columnLower();
658      double * upper = solution.columnUpper();
659      double * sol = solution.primalColumnSolution();
660      double * ray = solution.unboundedRay();
661      double * objective = solution.objective();
662      double objChange=0.0;
663      int iRow , iColumn;
664      // make sure feasible and columns form ray
665      for (iColumn=0;iColumn<numberColumns;iColumn++) {
666        double value = sol[iColumn];
667        assert(value<upper[iColumn]+1.0e-8);
668        assert(value>lower[iColumn]-1.0e-8);
669        value = ray[iColumn];
670        if (value>0.0)
671          assert(upper[iColumn]>1.0e30);
672        else if (value<0.0)
673          assert(lower[iColumn]<-1.0e30);
674        objChange += value*objective[iColumn];
675      }
676      // make sure increasing objective
677      assert(objChange>0.0);
678      double * result = new double[numberRows];
679      CoinFillN ( result, numberRows,0.0);
680      solution.matrix()->times(sol, result);
681      lower = solution.rowLower();
682      upper = solution.rowUpper();
683      sol = solution.primalRowSolution();
684      for (iRow=0;iRow<numberRows;iRow++) {
685        double value = result[iRow];
686        assert(eq(value,sol[iRow]));
687        assert(value<upper[iRow]+1.0e-8);
688        assert(value>lower[iRow]-1.0e-8);
689      }
690      CoinFillN ( result, numberRows,0.0);
691      solution.matrix()->times(ray, result);
692      // there may be small differences (especially if scaled)
693      for (iRow=0;iRow<numberRows;iRow++) {
694        double value = result[iRow];
695        if (value>1.0e-8)
696          assert(upper[iRow]>1.0e30);
697        else if (value<-1.0e-8)
698          assert(lower[iRow]<-1.0e30);
699      }
700      delete [] result;
701      delete [] ray;
702    }
703  }
704  // test infeasible
705  {   
706    CoinMpsIO m;
707    std::string fn = netlibDir+"brandy";
708    m.readMps(fn.c_str(),"mps");
709    ClpSimplex solution;
710    // do twice - without and with scaling
711    int iPass;
712    for (iPass=0;iPass<2;iPass++) {
713      solution.loadProblem(*m.getMatrixByCol(),m.getColLower(),m.getColUpper(),
714                      m.getObjCoefficients(),
715                      m.getRowLower(),m.getRowUpper());
716      if (iPass)
717        solution.scaling();
718      // test infeasible and ray
719      solution.columnUpper()[0]=0.0;
720#ifdef DUAL
721      solution.setDualBound(100.0);
722      solution.dual();
723#else
724      solution.primal();
725#endif
726      assert (solution.status()==1);
727      int numberColumns = solution.numberColumns();
728      int numberRows = solution.numberRows();
729      double * lower = solution.rowLower();
730      double * upper = solution.rowUpper();
731      double * ray = solution.infeasibilityRay();
732      assert(ray);
733      // construct proof of infeasibility
734      int iRow , iColumn;
735      double lo=0.0,up=0.0;
736      int nl=0,nu=0;
737      for (iRow=0;iRow<numberRows;iRow++) {
738        if (lower[iRow]>-1.0e20) {
739          lo += ray[iRow]*lower[iRow];
740        } else {
741          if (ray[iRow]>1.0e-8) 
742            nl++;
743        }
744        if (upper[iRow]<1.0e20) {
745          up += ray[iRow]*upper[iRow];
746        } else {
747          if (ray[iRow]>1.0e-8) 
748            nu++;
749        }
750      }
751      if (nl)
752        lo=-1.0e100;
753      if (nu)
754        up=1.0e100;
755      double * result = new double[numberColumns];
756      double lo2=0.0,up2=0.0;
757      CoinFillN ( result, numberColumns,0.0);
758      solution.matrix()->transposeTimes(ray, result);
759      lower = solution.columnLower();
760      upper = solution.columnUpper();
761      nl=nu=0;
762      for (iColumn=0;iColumn<numberColumns;iColumn++) {
763        if (result[iColumn]>1.0e-8) {
764          if (lower[iColumn]>-1.0e20)
765            lo2 += result[iColumn]*lower[iColumn];
766          else
767            nl++;
768          if (upper[iColumn]<1.0e20)
769            up2 += result[iColumn]*upper[iColumn];
770          else
771            nu++;
772        } else if (result[iColumn]<-1.0e-8) {
773          if (lower[iColumn]>-1.0e20)
774            up2 += result[iColumn]*lower[iColumn];
775          else
776            nu++;
777          if (upper[iColumn]<1.0e20)
778            lo2 += result[iColumn]*upper[iColumn];
779          else
780            nl++;
781        }
782      }
783      if (nl)
784        lo2=-1.0e100;
785      if (nu)
786        up2=1.0e100;
787      // make sure inconsistency
788      assert(lo2>up||up2<lo);
789      delete [] result;
790      delete [] ray;
791    }
792  }
793 
794}
Note: See TracBrowser for help on using the repository browser.