source: branches/devel-1/Test/unitTest.cpp @ 30

Last change on this file since 30 was 30, checked in by forrest, 19 years ago

Conditionally do presolve on ifdef

  • Property svn:eol-style set to native
  • Property svn:keywords set to Author Date Id Revision
File size: 33.4 KB
Line 
1// Copyright (C) 2002, International Business Machines
2// Corporation and others.  All Rights Reserved.
3#if defined(_MSC_VER)
4// Turn off compiler warning about long names
5#  pragma warning(disable:4786)
6#endif
7
8#include <cassert>
9#include <cstdio>
10#include <cmath>
11#include <cfloat>
12#include <string>
13#include <iostream>
14
15#include <time.h>
16#include <sys/times.h>
17#include <sys/resource.h>
18#include <unistd.h>
19
20#include "CoinMpsIO.hpp"
21#include "CoinPackedMatrix.hpp"
22#include "CoinPackedVector.hpp"
23#include "CoinHelperFunctions.hpp"
24
25#include "ClpFactorization.hpp"
26#include "ClpSimplex.hpp"
27#include "ClpDualRowSteepest.hpp"
28#include "ClpDualRowDantzig.hpp"
29#include "ClpPrimalColumnSteepest.hpp"
30#include "ClpPrimalColumnDantzig.hpp"
31#include "ClpParameters.hpp"
32
33#include "Presolve.hpp"
34
35//#############################################################################
36
37#ifdef NDEBUG
38#undef NDEBUG
39#endif
40
41// Function Prototypes. Function definitions is in this file.
42void testingMessage( const char * const msg );
43
44//----------------------------------------------------------------
45// unitTest [-mpsDir=V1] [-netlibDir=V2] [-netlib]
46//
47// where:
48//   -mpsDir: directory containing mps test files
49//       Default value V1="../Mps/Sample"   
50//   -netlibDir: directory containing netlib files
51//       Default value V2="../Mps/Netlib"
52//   -netlib
53//       If specified, then netlib test set run
54//
55// All parameters are optional.
56//----------------------------------------------------------------
57
58int mainTest (int argc, const char *argv[],bool doDual,
59              ClpSimplex empty, bool doPresolve)
60{
61  int i;
62
63  // define valid parameter keywords
64  std::set<std::string> definedKeyWords;
65  definedKeyWords.insert("-mpsDir");
66  definedKeyWords.insert("-netlibDir");
67  definedKeyWords.insert("-netlib");
68
69  // Create a map of parameter keys and associated data
70  std::map<std::string,std::string> parms;
71  for ( i=1; i<argc; i++ ) {
72    std::string parm(argv[i]);
73    std::string key,value;
74    unsigned int  eqPos = parm.find('=');
75
76    // Does parm contain and '='
77    if ( eqPos==std::string::npos ) {
78      //Parm does not contain '='
79      key = parm;
80    }
81    else {
82      key=parm.substr(0,eqPos);
83      value=parm.substr(eqPos+1);
84    }
85
86    // Is specifed key valid?
87    if ( definedKeyWords.find(key) == definedKeyWords.end() ) {
88      // invalid key word.
89      // Write help text
90      std::cerr <<"Undefined parameter \"" <<key <<"\".\n";
91      std::cerr <<"Correct usage: \n";
92      std::cerr <<"  unitTest [-mpsDir=V1] [-netlibDir=V2] [-netlib]\n";
93      std::cerr <<"  where:\n";
94      std::cerr <<"    -mpsDir: directory containing mps test files\n";
95      std::cerr <<"        Default value V1=\"../Mps/Sample\"\n";
96      std::cerr <<"    -netlibDir: directory containing netlib files\n";
97      std::cerr <<"        Default value V2=\"../Mps/Netlib\"\n";
98      std::cerr <<"    -netlib\n";
99      std::cerr <<"        If specified, then netlib testset run.\n";
100      return 1;
101    }
102    parms[key]=value;
103  }
104 
105  const char dirsep =  CoinFindDirSeparator();
106  // Set directory containing mps data files.
107  std::string mpsDir;
108  if (parms.find("-mpsDir") != parms.end())
109    mpsDir=parms["-mpsDir"] + dirsep;
110  else 
111    mpsDir = dirsep == '/' ? "../Mps/Sample/" : "..\\Mps\\Sample\\";
112 
113  // Set directory containing netlib data files.
114  std::string netlibDir;
115  if (parms.find("-netlibDir") != parms.end())
116    netlibDir=parms["-netlibDir"] + dirsep;
117  else 
118    netlibDir = dirsep == '/' ? "../Mps/Netlib/" : "..\\Mps\\Netlib\\";
119
120  testingMessage( "Testing ClpSimplex\n" );
121  ClpSimplexUnitTest(mpsDir,netlibDir);
122  if (parms.find("-netlib") != parms.end())
123  {
124    unsigned int m;
125   
126    // Define test problems:
127    //   mps names,
128    //   maximization or minimization,
129    //   Number of rows and columns in problem, and
130    //   objective function value
131    std::vector<std::string> mpsName;
132    std::vector<bool> min;
133    std::vector<int> nRows;
134    std::vector<int> nCols;
135    std::vector<double> objValue;
136    std::vector<double> objValueTol;
137    mpsName.push_back("25fv47");
138    min.push_back(true);
139    nRows.push_back(822);
140    nCols.push_back(1571);
141    objValueTol.push_back(1.E-10);
142    objValue.push_back(5.5018458883E+03);
143
144#if 1
145    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);
146    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);
147    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);
148   
149    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);
150    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);
151    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);
152    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);
153    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);
154    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);
155    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);
156    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);
157    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);
158    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);
159    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);
160    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);
161    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);
162    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);
163    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);
164    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);
165    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);
166    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);
167    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);
168    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);
169    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);
170    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);
171    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);
172    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);
173    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);
174    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.
175    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 );
176    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);
177    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);
178    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);
179    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);
180    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);
181    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);
182    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);
183    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);
184    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);
185#endif
186    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);
187    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);
188    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);
189    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);
190    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);
191    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);
192    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);
193    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);
194    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);
195    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);
196    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);
197    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);
198    //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);
199    //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);
200    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);
201    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);
202    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);
203    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);
204    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);
205    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);
206    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);
207    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);
208    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);
209    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);
210    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);
211    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);
212    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);
213    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);
214    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);
215    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);
216    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);
217    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);
218    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);
219    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);
220    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);
221    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);
222    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);
223    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);
224    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);
225    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);
226    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);
227    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);
228    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);
229    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);
230    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);
231    //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);
232    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); 
233    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);
234    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);
235    //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);
236    //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);
237    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);
238    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);
239    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);
240    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);
241  // Loop once for each Mps File
242    for (m=0; m<mpsName.size(); m++ ) {
243      std::cerr <<"  processing mps file: " <<mpsName[m] 
244                <<" (" <<m+1 <<" out of " <<mpsName.size() <<")" <<std::endl;
245   
246      // Read data mps file,
247      std::string fn = netlibDir+mpsName[m];
248      CoinMpsIO mps;
249      mps.readMps(fn.c_str(),"mps");
250      ClpSimplex solution=empty;
251      solution.loadProblem(*mps.getMatrixByCol(),mps.getColLower(),
252                           mps.getColUpper(),
253                           mps.getObjCoefficients(),
254                           mps.getRowLower(),mps.getRowUpper());
255
256      solution.setDblParam(ClpObjOffset,mps.objectiveOffset());
257      if (doPresolve) {
258#ifdef USE_PRESOLVE
259        Presolve pinfo;
260        ClpSimplex * model2 = pinfo.presolvedModel(solution,1.0e-8);
261        if (doDual)
262          model2->dual();
263        else
264          model2->primal();
265        pinfo.postsolve(true);
266       
267        delete model2;
268        printf("Resolving from postsolved model\n");
269        // later try without (1) and check duals before solve
270        solution.primal(1);
271        if (solution.numberIterations())
272          printf("****** iterated %d\n",solution.numberIterations());
273        solution.checkSolution();
274        printf("%g dual %g(%d) Primal %g(%d)\n",
275               solution.objectiveValue(),
276               solution.sumDualInfeasibilities(),
277               solution.numberDualInfeasibilities(),
278               solution.sumPrimalInfeasibilities(),
279               solution.numberPrimalInfeasibilities());
280        {
281          Presolve pinfoA;
282          model2 = pinfoA.presolvedModel(solution,1.0e-8);
283
284          printf("Resolving from presolved optimal solution\n");
285          model2->primal(1);
286               
287          delete model2;
288        }
289#else
290        if (doDual) {
291          solution.dual();
292        } else {
293          solution.primal();
294        }
295#endif
296      } else {
297        if (doDual) {
298          solution.dual();
299        } else {
300          solution.primal();
301        }
302      }
303      // Test objective solution value
304      {
305        double soln = solution.objectiveValue();
306        CoinRelFltEq eq(objValueTol[m]);
307        std::cerr <<soln <<",  " <<objValue[m] <<" diff "<<
308          soln-objValue[m]<<std::endl;
309        if(!eq(soln,objValue[m]))
310          printf("** difference fails\n");
311      }
312    }
313   
314  }
315  else {
316    testingMessage( "***Skipped Testing on netlib    ***\n" );
317    testingMessage( "***use -netlib to test class***\n" );
318  }
319
320  testingMessage( "All tests completed successfully\n" );
321  return 0;
322}
323
324 
325// Display message on stdout and stderr
326void testingMessage( const char * const msg )
327{
328  std::cerr <<msg;
329  //cout <<endl <<"*****************************************"
330  //     <<endl <<msg <<endl;
331}
332
333//--------------------------------------------------------------------------
334// test factorization methods and simplex method
335void
336ClpSimplexUnitTest(const std::string & mpsDir,
337                   const std::string & netlibDir)
338{
339 
340  CoinRelFltEq eq(0.000001);
341
342  {
343    ClpSimplex solution;
344 
345    // matrix data
346    //deliberate hiccup of 2 between 0 and 1
347    int start[5]={0,4,7,8,9};
348    int length[5]={2,3,1,1,1};
349    int rows[11]={0,2,-1,-1,0,1,2,0,1,2};
350    double elements[11]={7.0,2.0,1.0e10,1.0e10,-2.0,1.0,-2.0,1,1,1};
351    CoinPackedMatrix matrix(true,3,5,8,elements,rows,start,length);
352   
353    // rim data
354    double objective[7]={-4.0,1.0,0.0,0.0,0.0,0.0,0.0};
355    double rowLower[5]={14.0,3.0,3.0,1.0e10,1.0e10};
356    double rowUpper[5]={14.0,3.0,3.0,-1.0e10,-1.0e10};
357    double colLower[7]={0.0,0.0,0.0,0.0,0.0,0.0,0.0};
358    double colUpper[7]={100.0,100.0,100.0,100.0,100.0,100.0,100.0};
359   
360    // basis 1
361    int rowBasis1[3]={-1,-1,-1};
362    int colBasis1[5]={1,1,-1,-1,1};
363    solution.loadProblem(matrix,colLower,colUpper,objective,
364                         rowLower,rowUpper);
365    int i;
366    solution.createStatus();
367    for (i=0;i<3;i++) {
368      if (rowBasis1[i]<0) {
369        solution.setRowStatus(i,ClpSimplex::atLowerBound);
370      } else {
371        solution.setRowStatus(i,ClpSimplex::basic);
372      }
373    }
374    for (i=0;i<5;i++) {
375      if (colBasis1[i]<0) {
376        solution.setColumnStatus(i,ClpSimplex::atLowerBound);
377      } else {
378        solution.setColumnStatus(i,ClpSimplex::basic);
379      }
380    }
381    solution.setLogLevel(3+4+8+16+32);
382    solution.primal();
383    for (i=0;i<3;i++) {
384      if (rowBasis1[i]<0) {
385        solution.setRowStatus(i,ClpSimplex::atLowerBound);
386      } else {
387        solution.setRowStatus(i,ClpSimplex::basic);
388      }
389    }
390    for (i=0;i<5;i++) {
391      if (colBasis1[i]<0) {
392        solution.setColumnStatus(i,ClpSimplex::atLowerBound);
393      } else {
394        solution.setColumnStatus(i,ClpSimplex::basic);
395      }
396    }
397    // intricate stuff does not work with scaling
398    solution.scaling(0);
399    assert(!solution.factorize ( ));
400    const double * colsol = solution.primalColumnSolution();
401    const double * rowsol = solution.primalRowSolution();
402    solution.getSolution(rowsol,colsol);
403    double colsol1[5]={20.0/7.0,3.0,0.0,0.0,23.0/7.0};
404    for (i=0;i<5;i++) {
405      assert(eq(colsol[i],colsol1[i]));
406    }
407    // now feed in again without actually doing factorization
408    ClpFactorization factorization2 = *solution.factorization();
409    ClpSimplex solution2 = solution;
410    solution2.setFactorization(factorization2);
411    solution2.createStatus();
412    for (i=0;i<3;i++) {
413      if (rowBasis1[i]<0) {
414        solution2.setRowStatus(i,ClpSimplex::atLowerBound);
415      } else {
416        solution2.setRowStatus(i,ClpSimplex::basic);
417      }
418    }
419    for (i=0;i<5;i++) {
420      if (colBasis1[i]<0) {
421        solution2.setColumnStatus(i,ClpSimplex::atLowerBound);
422      } else {
423        solution2.setColumnStatus(i,ClpSimplex::basic);
424      }
425    }
426    // intricate stuff does not work with scaling
427    solution2.scaling(0);
428    solution2.getSolution(rowsol,colsol);
429    colsol = solution2.primalColumnSolution();
430    rowsol = solution2.primalRowSolution();
431    for (i=0;i<5;i++) {
432      assert(eq(colsol[i],colsol1[i]));
433    }
434    solution2.setDualBound(0.1);
435    solution2.dual();
436    objective[2]=-1.0;
437    objective[3]=-0.5;
438    objective[4]=10.0;
439    solution.dual();
440    for (i=0;i<3;i++) {
441      rowLower[i]=-1.0e50;
442      colUpper[i+2]=0.0;
443    }
444    solution.setLogLevel(3);
445    solution.dual();
446    double rowObjective[]={1.0,0.5,-10.0};
447    solution.loadProblem(matrix,colLower,colUpper,objective,
448                         rowLower,rowUpper,rowObjective);
449    solution.dual();
450  }
451  {   
452    CoinMpsIO m;
453    std::string fn = mpsDir+"exmip1";
454    m.readMps(fn.c_str(),"mps");
455    ClpSimplex solution;
456    solution.loadProblem(*m.getMatrixByCol(),m.getColLower(),m.getColUpper(),
457                         m.getObjCoefficients(),
458                         m.getRowLower(),m.getRowUpper());
459    solution.dual();
460  }
461  // test steepest edge
462  {   
463    CoinMpsIO m;
464    std::string fn = netlibDir+"finnis";
465    m.readMps(fn.c_str(),"mps");
466    ClpModel model;
467    model.loadProblem(*m.getMatrixByCol(),m.getColLower(),
468                    m.getColUpper(),
469                    m.getObjCoefficients(),
470                    m.getRowLower(),m.getRowUpper());
471    ClpSimplex solution(model);
472
473    solution.scaling(1); 
474    solution.setDualBound(1.0e8);
475    //solution.factorization()->maximumPivots(1);
476    //solution.setLogLevel(3);
477    solution.setDualTolerance(1.0e-7);
478    // set objective sense,
479    ClpDualRowSteepest steep;
480    solution.setDualRowPivotAlgorithm(steep);
481    solution.setDblParam(ClpObjOffset,m.objectiveOffset());
482    solution.dual();
483  }
484  // test normal solution
485  {   
486    CoinMpsIO m;
487    std::string fn = netlibDir+"afiro";
488    m.readMps(fn.c_str(),"mps");
489    ClpSimplex solution;
490    ClpModel model;
491    // do twice - without and with scaling
492    int iPass;
493    for (iPass=0;iPass<2;iPass++) {
494      // explicit row objective for testing
495      int nr = m.getNumRows();
496      double * rowObj = new double[nr];
497      CoinFillN(rowObj,nr,0.0);
498      model.loadProblem(*m.getMatrixByCol(),m.getColLower(),m.getColUpper(),
499                      m.getObjCoefficients(),
500                      m.getRowLower(),m.getRowUpper(),rowObj);
501      delete [] rowObj;
502      solution = ClpSimplex(model);
503      if (iPass) {
504        solution.scaling();
505      }
506      solution.dual();
507      solution.dual();
508      // test optimal
509      assert (solution.status()==0);
510      int numberColumns = solution.numberColumns();
511      int numberRows = solution.numberRows();
512      CoinPackedVector colsol(numberColumns,solution.primalColumnSolution());
513      double * objective = solution.objective();
514      double objValue = colsol.dotProduct(objective);
515      CoinRelFltEq eq(1.0e-8);
516      assert(eq(objValue,-4.6475314286e+02));
517      double * lower = solution.columnLower();
518      double * upper = solution.columnUpper();
519      double * sol = solution.primalColumnSolution();
520      double * result = new double[numberColumns];
521      CoinFillN ( result, numberColumns,0.0);
522      solution.matrix()->transposeTimes(solution.dualRowSolution(), result);
523      int iRow , iColumn;
524      // see if feasible and dual feasible
525      for (iColumn=0;iColumn<numberColumns;iColumn++) {
526        double value = sol[iColumn];
527        assert(value<upper[iColumn]+1.0e-8);
528        assert(value>lower[iColumn]-1.0e-8);
529        value = objective[iColumn]-result[iColumn];
530        assert (value>-1.0e-5);
531        if (sol[iColumn]>1.0e-5)
532          assert (value<1.0e-5);
533      }
534      delete [] result;
535      result = new double[numberRows];
536      CoinFillN ( result, numberRows,0.0);
537      solution.matrix()->times(colsol, result);
538      lower = solution.rowLower();
539      upper = solution.rowUpper();
540      sol = solution.primalRowSolution();
541      for (iRow=0;iRow<numberRows;iRow++) {
542        double value = result[iRow];
543        assert(eq(value,sol[iRow]));
544        assert(value<upper[iRow]+1.0e-8);
545        assert(value>lower[iRow]-1.0e-8);
546      }
547      delete [] result;
548      // test row objective
549      double * rowObjective = solution.rowObjective();
550      CoinDisjointCopyN(solution.dualRowSolution(),numberRows,rowObjective);
551      CoinDisjointCopyN(solution.dualColumnSolution(),numberColumns,objective);
552      // this sets up all slack basis
553      solution.createStatus();
554      solution.dual();
555      CoinFillN(rowObjective,numberRows,0.0);
556      CoinDisjointCopyN(m.getObjCoefficients(),numberColumns,objective);
557      solution.dual();
558    }
559  }
560  // test unbounded
561  {   
562    CoinMpsIO m;
563    std::string fn = netlibDir+"brandy";
564    m.readMps(fn.c_str(),"mps");
565    ClpSimplex solution;
566    // do twice - without and with scaling
567    int iPass;
568    for (iPass=0;iPass<2;iPass++) {
569      solution.loadProblem(*m.getMatrixByCol(),m.getColLower(),m.getColUpper(),
570                      m.getObjCoefficients(),
571                      m.getRowLower(),m.getRowUpper());
572      if (iPass)
573        solution.scaling();
574      solution.setOptimizationDirection(-1);
575      // test unbounded and ray
576#ifdef DUAL
577      solution.setDualBound(100.0);
578      solution.dual();
579#else
580      solution.primal();
581#endif
582      assert (solution.status()==2);
583      int numberColumns = solution.numberColumns();
584      int numberRows = solution.numberRows();
585      double * lower = solution.columnLower();
586      double * upper = solution.columnUpper();
587      double * sol = solution.primalColumnSolution();
588      double * ray = solution.unboundedRay();
589      double * objective = solution.objective();
590      double objChange=0.0;
591      int iRow , iColumn;
592      // make sure feasible and columns form ray
593      for (iColumn=0;iColumn<numberColumns;iColumn++) {
594        double value = sol[iColumn];
595        assert(value<upper[iColumn]+1.0e-8);
596        assert(value>lower[iColumn]-1.0e-8);
597        value = ray[iColumn];
598        if (value>0.0)
599          assert(upper[iColumn]>1.0e30);
600        else if (value<0.0)
601          assert(lower[iColumn]<-1.0e30);
602        objChange += value*objective[iColumn];
603      }
604      // make sure increasing objective
605      assert(objChange>0.0);
606      double * result = new double[numberRows];
607      CoinFillN ( result, numberRows,0.0);
608      solution.matrix()->times(sol, result);
609      lower = solution.rowLower();
610      upper = solution.rowUpper();
611      sol = solution.primalRowSolution();
612      for (iRow=0;iRow<numberRows;iRow++) {
613        double value = result[iRow];
614        assert(eq(value,sol[iRow]));
615        assert(value<upper[iRow]+1.0e-8);
616        assert(value>lower[iRow]-1.0e-8);
617      }
618      CoinFillN ( result, numberRows,0.0);
619      solution.matrix()->times(ray, result);
620      // there may be small differences (especially if scaled)
621      for (iRow=0;iRow<numberRows;iRow++) {
622        double value = result[iRow];
623        if (value>1.0e-8)
624          assert(upper[iRow]>1.0e30);
625        else if (value<-1.0e-8)
626          assert(lower[iRow]<-1.0e30);
627      }
628      delete [] result;
629      delete [] ray;
630    }
631  }
632  // test infeasible
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      // test infeasible and ray
647      solution.columnUpper()[0]=0.0;
648#ifdef DUAL
649      solution.setDualBound(100.0);
650      solution.dual();
651#else
652      solution.primal();
653#endif
654      assert (solution.status()==1);
655      int numberColumns = solution.numberColumns();
656      int numberRows = solution.numberRows();
657      double * lower = solution.rowLower();
658      double * upper = solution.rowUpper();
659      double * ray = solution.infeasibilityRay();
660      assert(ray);
661      // construct proof of infeasibility
662      int iRow , iColumn;
663      double lo=0.0,up=0.0;
664      int nl=0,nu=0;
665      for (iRow=0;iRow<numberRows;iRow++) {
666        if (lower[iRow]>-1.0e20) {
667          lo += ray[iRow]*lower[iRow];
668        } else {
669          if (ray[iRow]>1.0e-8) 
670            nl++;
671        }
672        if (upper[iRow]<1.0e20) {
673          up += ray[iRow]*upper[iRow];
674        } else {
675          if (ray[iRow]>1.0e-8) 
676            nu++;
677        }
678      }
679      if (nl)
680        lo=-1.0e100;
681      if (nu)
682        up=1.0e100;
683      double * result = new double[numberColumns];
684      double lo2=0.0,up2=0.0;
685      CoinFillN ( result, numberColumns,0.0);
686      solution.matrix()->transposeTimes(ray, result);
687      lower = solution.columnLower();
688      upper = solution.columnUpper();
689      nl=nu=0;
690      for (iColumn=0;iColumn<numberColumns;iColumn++) {
691        if (result[iColumn]>1.0e-8) {
692          if (lower[iColumn]>-1.0e20)
693            lo2 += result[iColumn]*lower[iColumn];
694          else
695            nl++;
696          if (upper[iColumn]<1.0e20)
697            up2 += result[iColumn]*upper[iColumn];
698          else
699            nu++;
700        } else if (result[iColumn]<-1.0e-8) {
701          if (lower[iColumn]>-1.0e20)
702            up2 += result[iColumn]*lower[iColumn];
703          else
704            nu++;
705          if (upper[iColumn]<1.0e20)
706            lo2 += result[iColumn]*upper[iColumn];
707          else
708            nl++;
709        }
710      }
711      if (nl)
712        lo2=-1.0e100;
713      if (nu)
714        up2=1.0e100;
715      // make sure inconsistency
716      assert(lo2>up||up2<lo);
717      delete [] result;
718      delete [] ray;
719    }
720  }
721 
722}
Note: See TracBrowser for help on using the repository browser.