source: trunk/Clp/test/OsiClpSolverInterfaceTest.cpp @ 2022

Last change on this file since 2022 was 1848, checked in by stefan, 9 years ago

merge chgset 1847 from stable/1.14
File size: 60.5 KB
Line 
1// $Id$
2// Copyright (C) 2002, International Business Machines
3// Corporation and others.  All Rights Reserved.
4// This code is licensed under the terms of the Eclipse Public License (EPL).
5
6#include "CoinPragma.hpp"
7
8//#include <cassert>
9//#include <cstdlib>
10//#include <cstdio>
11//#include <iostream>
12
13#include "OsiClpSolverInterface.hpp"
14#include "OsiUnitTests.hpp"
15#include "OsiCuts.hpp"
16#include "OsiRowCut.hpp"
17#include "OsiColCut.hpp"
18#include "CoinMessage.hpp"
19#include "ClpMessage.hpp"
20#include "ClpFactorization.hpp"
21#include "CoinModel.hpp"
22#include "CoinIndexedVector.hpp"
23#include "ClpPlusMinusOneMatrix.hpp"
24
25//#############################################################################
26
27class OsiClpMessageTest :
28   public CoinMessageHandler {
29
30public:
31  virtual int print() ;
32  OsiClpMessageTest();
33};
34
35OsiClpMessageTest::OsiClpMessageTest() : CoinMessageHandler()
36{
37}
38int
39OsiClpMessageTest::print()
40{
41  if (currentMessage().externalNumber()==0&&currentSource()=="Clp") 
42    std::cout<<"This is not actually an advertisement by Dash Associates - just my feeble attempt to test message handling and language - JJHF"<<std::endl;
43  else if (currentMessage().externalNumber()==5&&currentSource()=="Osi") 
44    std::cout<<"End of search trapped"<<std::endl;
45  return CoinMessageHandler::print();
46}
47
48//--------------------------------------------------------------------------
49// test EKKsolution methods.
50void
51OsiClpSolverInterfaceUnitTest(const std::string & mpsDir, const std::string & netlibDir)
52{
53  // Test default constructor
54  {
55    OsiClpSolverInterface m;
56    OSIUNITTEST_ASSERT_ERROR(m.rowsense_           == NULL, {}, "clp", "default constructor");
57    OSIUNITTEST_ASSERT_ERROR(m.rhs_                == NULL, {}, "clp", "default constructor");
58    OSIUNITTEST_ASSERT_ERROR(m.rowrange_           == NULL, {}, "clp", "default constructor");
59    OSIUNITTEST_ASSERT_ERROR(m.matrixByRow_        == NULL, {}, "clp", "default constructor");
60    OSIUNITTEST_ASSERT_ERROR(m.ws_                 == NULL, {}, "clp", "default constructor");
61    OSIUNITTEST_ASSERT_ERROR(m.itlimOrig_       == 9999999, {}, "clp", "default constructor");
62    OSIUNITTEST_ASSERT_ERROR(m.lastAlgorithm_      == 0,    {}, "clp", "default constructor");
63    OSIUNITTEST_ASSERT_ERROR(m.integerInformation_ == NULL, {}, "clp", "default constructor");
64  }
65
66  {   
67    CoinRelFltEq eq;
68    OsiClpSolverInterface m;
69    std::string fn = mpsDir+"exmip1";
70    m.readMps(fn.c_str(),"mps");
71   
72    {
73      OsiClpSolverInterface im;   
74     
75      OSIUNITTEST_ASSERT_ERROR(im.getNumCols()  == 0,    {}, "clp", "default constructor");
76      OSIUNITTEST_ASSERT_ERROR(im.getModelPtr() != NULL, {}, "clp", "default constructor");
77
78      // Test reset
79      im.reset();
80      OSIUNITTEST_ASSERT_ERROR(im.rowsense_           == NULL, {}, "clp", "reset");
81      OSIUNITTEST_ASSERT_ERROR(im.rhs_                == NULL, {}, "clp", "reset");
82      OSIUNITTEST_ASSERT_ERROR(im.rowrange_           == NULL, {}, "clp", "reset");
83      OSIUNITTEST_ASSERT_ERROR(im.matrixByRow_        == NULL, {}, "clp", "reset");
84      OSIUNITTEST_ASSERT_ERROR(im.ws_                 == NULL, {}, "clp", "reset");
85      OSIUNITTEST_ASSERT_ERROR(im.itlimOrig_       == 9999999, {}, "clp", "reset");
86      OSIUNITTEST_ASSERT_ERROR(im.lastAlgorithm_      ==    0, {}, "clp", "reset");
87      OSIUNITTEST_ASSERT_ERROR(im.integerInformation_ == NULL, {}, "clp", "reset");
88    }
89   
90    // Test copy constructor and assignment operator
91    {
92      OsiClpSolverInterface lhs;
93      {     
94        OsiClpSolverInterface im(m);       
95       
96        OsiClpSolverInterface imC1(im);
97        OSIUNITTEST_ASSERT_ERROR(imC1.getModelPtr() != im.getModelPtr(), {}, "clp", "copy constructor");
98        OSIUNITTEST_ASSERT_ERROR(imC1.getNumCols()  == im.getNumCols(),  {}, "clp", "copy constructor");
99        OSIUNITTEST_ASSERT_ERROR(imC1.getNumRows()  == im.getNumRows(),  {}, "clp", "copy constructor");
100       
101        OsiClpSolverInterface imC2(im);
102        OSIUNITTEST_ASSERT_ERROR(imC2.getModelPtr() != im.getModelPtr(), {}, "clp", "copy constructor");
103        OSIUNITTEST_ASSERT_ERROR(imC2.getNumCols()  == im.getNumCols(),  {}, "clp", "copy constructor");
104        OSIUNITTEST_ASSERT_ERROR(imC2.getNumRows()  == im.getNumRows(),  {}, "clp", "copy constructor");
105
106        OSIUNITTEST_ASSERT_ERROR(imC1.getModelPtr() != imC2.getModelPtr(), {}, "clp", "copy constructor");
107       
108        lhs = imC2;
109      }
110
111      // Test that lhs has correct values even though rhs has gone out of scope
112      OSIUNITTEST_ASSERT_ERROR(lhs.getModelPtr() != m.getModelPtr(), {}, "clp", "assignment operator");
113      OSIUNITTEST_ASSERT_ERROR(lhs.getNumCols()  == m.getNumCols(),  {}, "clp", "copy constructor");
114      OSIUNITTEST_ASSERT_ERROR(lhs.getNumRows()  == m.getNumRows(),  {}, "clp", "copy constructor");
115    }
116
117    // Test clone
118    {
119      OsiClpSolverInterface clpSi(m);
120      OsiSolverInterface * siPtr = &clpSi;
121      OsiSolverInterface * siClone = siPtr->clone();
122      OsiClpSolverInterface * clpClone = dynamic_cast<OsiClpSolverInterface*>(siClone);
123     
124      OSIUNITTEST_ASSERT_ERROR(clpClone != NULL, {}, "clp", "clone");
125      OSIUNITTEST_ASSERT_ERROR(clpClone->getModelPtr() != clpSi.getModelPtr(), {}, "clp", "clone");
126      OSIUNITTEST_ASSERT_ERROR(clpClone->getNumRows() == clpSi.getNumRows(),   {}, "clp", "clone");
127      OSIUNITTEST_ASSERT_ERROR(clpClone->getNumCols() == m.getNumCols(),       {}, "clp", "clone");
128
129      delete siClone;
130    }
131 
132    // Test infinity
133    {
134      OsiClpSolverInterface si;
135      OSIUNITTEST_ASSERT_ERROR(si.getInfinity() == OsiClpInfinity, {}, "clp", "infinity");
136    }
137   
138    // Test some catches
139    if (!OsiClpHasNDEBUG())
140    {
141      OsiClpSolverInterface solver;
142      try {
143        solver.setObjCoeff(0,0.0);
144        OSIUNITTEST_ADD_OUTCOME("clp", "setObjCoeff on empty model", "should throw exception", OsiUnitTest::TestOutcome::ERROR, false);
145      }
146      catch (CoinError e) {
147        if (OsiUnitTest::verbosity >= 1)
148          std::cout<<"Correct throw from setObjCoeff on empty model"<<std::endl;
149      }
150
151      std::string fn = mpsDir+"exmip1";
152      solver.readMps(fn.c_str(),"mps");
153      OSIUNITTEST_CATCH_ERROR(solver.setObjCoeff(0,0.0), {}, "clp", "setObjCoeff on nonempty model");
154
155      try {
156        int index[]={0,20};
157        double value[]={0.0,0.0,0.0,0.0};
158        solver.setColSetBounds(index,index+2,value);
159        OSIUNITTEST_ADD_OUTCOME("clp", "setColSetBounds on cols not in model", "should throw exception", OsiUnitTest::TestOutcome::ERROR, false);
160      }
161      catch (CoinError e) {
162        if (OsiUnitTest::verbosity >= 1)
163          std::cout<<"Correct throw from setObjCoeff on empty model"<<std::endl;
164      }
165    }
166   
167    {
168      OsiClpSolverInterface clpSi(m);
169      int nc = clpSi.getNumCols();
170      int nr = clpSi.getNumRows();
171      const double * cl = clpSi.getColLower();
172      const double * cu = clpSi.getColUpper();
173      const double * rl = clpSi.getRowLower();
174      const double * ru = clpSi.getRowUpper();
175      OSIUNITTEST_ASSERT_ERROR(nc == 8, return, "clp", "read and copy exmip1");
176      OSIUNITTEST_ASSERT_ERROR(nr == 5, return, "clp", "read and copy exmip1");
177      OSIUNITTEST_ASSERT_ERROR(eq(cl[0],2.5), {}, "clp", "read and copy exmip1");
178      OSIUNITTEST_ASSERT_ERROR(eq(cl[1],0.0), {}, "clp", "read and copy exmip1");
179      OSIUNITTEST_ASSERT_ERROR(eq(cu[1],4.1), {}, "clp", "read and copy exmip1");
180      OSIUNITTEST_ASSERT_ERROR(eq(cu[2],1.0), {}, "clp", "read and copy exmip1");
181      OSIUNITTEST_ASSERT_ERROR(eq(rl[0],2.5), {}, "clp", "read and copy exmip1");
182      OSIUNITTEST_ASSERT_ERROR(eq(rl[4],3.0), {}, "clp", "read and copy exmip1");
183      OSIUNITTEST_ASSERT_ERROR(eq(ru[1],2.1), {}, "clp", "read and copy exmip1");
184      OSIUNITTEST_ASSERT_ERROR(eq(ru[4],15.), {}, "clp", "read and copy exmip1");
185     
186      const double * cs = clpSi.getColSolution();
187      OSIUNITTEST_ASSERT_ERROR(eq(cs[0],2.5), {}, "clp", "read and copy exmip1");
188      OSIUNITTEST_ASSERT_ERROR(eq(cs[7],0.0), {}, "clp", "read and copy exmip1");
189     
190      OSIUNITTEST_ASSERT_ERROR(!eq(cl[3],1.2345), {}, "clp", "set col lower");
191      clpSi.setColLower( 3, 1.2345 );
192      OSIUNITTEST_ASSERT_ERROR( eq(cl[3],1.2345), {}, "clp", "set col lower");
193     
194      OSIUNITTEST_ASSERT_ERROR(!eq(clpSi.getColUpper()[4],10.2345), {}, "clp", "set col upper");
195      clpSi.setColUpper( 4, 10.2345 );
196      OSIUNITTEST_ASSERT_ERROR( eq(clpSi.getColUpper()[4],10.2345), {}, "clp", "set col upper");
197
198      double objValue = clpSi.getObjValue();
199      OSIUNITTEST_ASSERT_ERROR(eq(objValue,3.5), {}, "clp", "getObjValue() before solve");
200
201      OSIUNITTEST_ASSERT_ERROR(eq(clpSi.getObjCoefficients()[0], 1.0), {}, "clp", "read and copy exmip1");
202      OSIUNITTEST_ASSERT_ERROR(eq(clpSi.getObjCoefficients()[1], 0.0), {}, "clp", "read and copy exmip1");
203      OSIUNITTEST_ASSERT_ERROR(eq(clpSi.getObjCoefficients()[2], 0.0), {}, "clp", "read and copy exmip1");
204      OSIUNITTEST_ASSERT_ERROR(eq(clpSi.getObjCoefficients()[3], 0.0), {}, "clp", "read and copy exmip1");
205      OSIUNITTEST_ASSERT_ERROR(eq(clpSi.getObjCoefficients()[4], 2.0), {}, "clp", "read and copy exmip1");
206      OSIUNITTEST_ASSERT_ERROR(eq(clpSi.getObjCoefficients()[5], 0.0), {}, "clp", "read and copy exmip1");
207      OSIUNITTEST_ASSERT_ERROR(eq(clpSi.getObjCoefficients()[6], 0.0), {}, "clp", "read and copy exmip1");
208      OSIUNITTEST_ASSERT_ERROR(eq(clpSi.getObjCoefficients()[7],-1.0), {}, "clp", "read and copy exmip1");
209    }
210   
211    // Test matrixByRow method
212    {
213      const OsiClpSolverInterface si(m);
214      const CoinPackedMatrix * smP = si.getMatrixByRow();
215
216      OSIUNITTEST_ASSERT_ERROR(smP->getMajorDim()    ==  5, return, "clp", "getMatrixByRow: major dim");
217      OSIUNITTEST_ASSERT_ERROR(smP->getNumElements() == 14, return, "clp", "getMatrixByRow: num elements");
218
219      CoinRelFltEq eq;
220      const double * ev = smP->getElements();
221      OSIUNITTEST_ASSERT_ERROR(eq(ev[0],   3.0), {}, "clp", "getMatrixByRow: elements");
222      OSIUNITTEST_ASSERT_ERROR(eq(ev[1],   1.0), {}, "clp", "getMatrixByRow: elements");
223      OSIUNITTEST_ASSERT_ERROR(eq(ev[2],  -2.0), {}, "clp", "getMatrixByRow: elements");
224      OSIUNITTEST_ASSERT_ERROR(eq(ev[3],  -1.0), {}, "clp", "getMatrixByRow: elements");
225      OSIUNITTEST_ASSERT_ERROR(eq(ev[4],  -1.0), {}, "clp", "getMatrixByRow: elements");
226      OSIUNITTEST_ASSERT_ERROR(eq(ev[5],   2.0), {}, "clp", "getMatrixByRow: elements");
227      OSIUNITTEST_ASSERT_ERROR(eq(ev[6],   1.1), {}, "clp", "getMatrixByRow: elements");
228      OSIUNITTEST_ASSERT_ERROR(eq(ev[7],   1.0), {}, "clp", "getMatrixByRow: elements");
229      OSIUNITTEST_ASSERT_ERROR(eq(ev[8],   1.0), {}, "clp", "getMatrixByRow: elements");
230      OSIUNITTEST_ASSERT_ERROR(eq(ev[9],   2.8), {}, "clp", "getMatrixByRow: elements");
231      OSIUNITTEST_ASSERT_ERROR(eq(ev[10], -1.2), {}, "clp", "getMatrixByRow: elements");
232      OSIUNITTEST_ASSERT_ERROR(eq(ev[11],  5.6), {}, "clp", "getMatrixByRow: elements");
233      OSIUNITTEST_ASSERT_ERROR(eq(ev[12],  1.0), {}, "clp", "getMatrixByRow: elements");
234      OSIUNITTEST_ASSERT_ERROR(eq(ev[13],  1.9), {}, "clp", "getMatrixByRow: elements");
235     
236      const int * mi = smP->getVectorStarts();
237      OSIUNITTEST_ASSERT_ERROR(mi[0] ==  0, {}, "clp", "getMatrixByRow: vector starts");
238      OSIUNITTEST_ASSERT_ERROR(mi[1] ==  5, {}, "clp", "getMatrixByRow: vector starts");
239      OSIUNITTEST_ASSERT_ERROR(mi[2] ==  7, {}, "clp", "getMatrixByRow: vector starts");
240      OSIUNITTEST_ASSERT_ERROR(mi[3] ==  9, {}, "clp", "getMatrixByRow: vector starts");
241      OSIUNITTEST_ASSERT_ERROR(mi[4] == 11, {}, "clp", "getMatrixByRow: vector starts");
242      OSIUNITTEST_ASSERT_ERROR(mi[5] == 14, {}, "clp", "getMatrixByRow: vector starts");
243     
244      const int * ei = smP->getIndices();
245      OSIUNITTEST_ASSERT_ERROR(ei[ 0] == 0, {}, "clp", "getMatrixByRow: indices");
246      OSIUNITTEST_ASSERT_ERROR(ei[ 1] == 1, {}, "clp", "getMatrixByRow: indices");
247      OSIUNITTEST_ASSERT_ERROR(ei[ 2] == 3, {}, "clp", "getMatrixByRow: indices");
248      OSIUNITTEST_ASSERT_ERROR(ei[ 3] == 4, {}, "clp", "getMatrixByRow: indices");
249      OSIUNITTEST_ASSERT_ERROR(ei[ 4] == 7, {}, "clp", "getMatrixByRow: indices");
250      OSIUNITTEST_ASSERT_ERROR(ei[ 5] == 1, {}, "clp", "getMatrixByRow: indices");
251      OSIUNITTEST_ASSERT_ERROR(ei[ 6] == 2, {}, "clp", "getMatrixByRow: indices");
252      OSIUNITTEST_ASSERT_ERROR(ei[ 7] == 2, {}, "clp", "getMatrixByRow: indices");
253      OSIUNITTEST_ASSERT_ERROR(ei[ 8] == 5, {}, "clp", "getMatrixByRow: indices");
254      OSIUNITTEST_ASSERT_ERROR(ei[ 9] == 3, {}, "clp", "getMatrixByRow: indices");
255      OSIUNITTEST_ASSERT_ERROR(ei[10] == 6, {}, "clp", "getMatrixByRow: indices");
256      OSIUNITTEST_ASSERT_ERROR(ei[11] == 0, {}, "clp", "getMatrixByRow: indices");
257      OSIUNITTEST_ASSERT_ERROR(ei[12] == 4, {}, "clp", "getMatrixByRow: indices");
258      OSIUNITTEST_ASSERT_ERROR(ei[13] == 7, {}, "clp", "getMatrixByRow: indices");
259    }
260
261    // Test adding several cuts
262    {
263      OsiClpSolverInterface fim;
264      std::string fn = mpsDir+"exmip1";
265      fim.readMps(fn.c_str(),"mps");
266      // exmip1.mps has 2 integer variables with index 2 & 3
267      fim.initialSolve();
268      OsiRowCut cuts[3];
269
270      // Generate one ineffective cut plus two trivial cuts
271      int c;
272      int nc = fim.getNumCols();
273      int *inx = new int[nc];
274      for (c=0;c<nc;c++) inx[c]=c;
275      double *el = new double[nc];
276      for (c=0;c<nc;c++) el[c]=1.0e-50+((double)c)*((double)c);
277     
278      cuts[0].setRow(nc,inx,el);
279      cuts[0].setLb(-100.);
280      cuts[0].setUb(500.);
281      cuts[0].setEffectiveness(22);
282      el[4]=0.0; // to get inf later
283     
284      for (c=2;c<4;c++) {
285        el[0]=1.0;
286        inx[0]=c;
287        cuts[c-1].setRow(1,inx,el);
288        cuts[c-1].setLb(1.);
289        cuts[c-1].setUb(100.);
290        cuts[c-1].setEffectiveness(c);
291      }
292      fim.writeMps("x1.mps");
293      fim.applyRowCuts(3,cuts);
294      fim.writeMps("x2.mps");
295      // resolve - should get message about zero elements
296      fim.resolve();
297      fim.writeMps("x3.mps");
298      // check integer solution
299      const double * cs = fim.getColSolution();
300      CoinRelFltEq eq;
301      OSIUNITTEST_ASSERT_ERROR(eq(cs[2], 1.0), {}, "clp", "add cuts");
302      OSIUNITTEST_ASSERT_ERROR(eq(cs[3], 1.0), {}, "clp", "add cuts");
303      // check will find invalid matrix
304      el[0]=1.0/el[4];
305      inx[0]=0;
306      cuts[0].setRow(nc,inx,el);
307      cuts[0].setLb(-100.);
308      cuts[0].setUb(500.);
309      cuts[0].setEffectiveness(22);
310      fim.applyRowCut(cuts[0]);
311      // resolve - should get message about zero elements
312      fim.resolve();
313      OSIUNITTEST_ASSERT_WARNING(fim.isAbandoned(), {}, "clp", "add cuts");
314      delete[]el;
315      delete[]inx;
316    }
317
318    // Test using bad basis
319    {
320      OsiClpSolverInterface fim;
321      std::string fn = mpsDir+"exmip1";
322      fim.readMps(fn.c_str(),"mps");
323      fim.initialSolve();
324      int numberRows = fim.getModelPtr()->numberRows();
325      int * rowStatus = new int[numberRows];
326      int numberColumns = fim.getModelPtr()->numberColumns();
327      int * columnStatus = new int[numberColumns];
328      fim.getBasisStatus(columnStatus,rowStatus);
329      for (int i=0; i<numberRows; i++)
330        rowStatus[i]=2;
331      fim.setBasisStatus(columnStatus,rowStatus);
332      fim.resolve();
333      delete [] rowStatus;
334      delete [] columnStatus;
335    }
336
337    // Test matrixByCol method
338    {
339      const OsiClpSolverInterface si(m);
340      const CoinPackedMatrix * smP = si.getMatrixByCol();
341     
342      OSIUNITTEST_ASSERT_ERROR(smP->getMajorDim()    ==  8, return, "clp", "getMatrixByCol: major dim");
343      OSIUNITTEST_ASSERT_ERROR(smP->getMinorDim()    ==  5, return, "clp", "getMatrixByCol: minor dim");
344      OSIUNITTEST_ASSERT_ERROR(smP->getNumElements() == 14, return, "clp", "getMatrixByCol: number of elements");
345      OSIUNITTEST_ASSERT_ERROR(smP->getSizeVectorStarts() == 9, return, "clp", "getMatrixByCol: vector starts size");
346
347      CoinRelFltEq eq;
348      const double * ev = smP->getElements();
349      OSIUNITTEST_ASSERT_ERROR(eq(ev[ 0], 3.0), {}, "clp", "getMatrixByCol: elements");
350      OSIUNITTEST_ASSERT_ERROR(eq(ev[ 1], 5.6), {}, "clp", "getMatrixByCol: elements");
351      OSIUNITTEST_ASSERT_ERROR(eq(ev[ 2], 1.0), {}, "clp", "getMatrixByCol: elements");
352      OSIUNITTEST_ASSERT_ERROR(eq(ev[ 3], 2.0), {}, "clp", "getMatrixByCol: elements");
353      OSIUNITTEST_ASSERT_ERROR(eq(ev[ 4], 1.1), {}, "clp", "getMatrixByCol: elements");
354      OSIUNITTEST_ASSERT_ERROR(eq(ev[ 5], 1.0), {}, "clp", "getMatrixByCol: elements");
355      OSIUNITTEST_ASSERT_ERROR(eq(ev[ 6],-2.0), {}, "clp", "getMatrixByCol: elements");
356      OSIUNITTEST_ASSERT_ERROR(eq(ev[ 7], 2.8), {}, "clp", "getMatrixByCol: elements");
357      OSIUNITTEST_ASSERT_ERROR(eq(ev[ 8],-1.0), {}, "clp", "getMatrixByCol: elements");
358      OSIUNITTEST_ASSERT_ERROR(eq(ev[ 9], 1.0), {}, "clp", "getMatrixByCol: elements");
359      OSIUNITTEST_ASSERT_ERROR(eq(ev[10], 1.0), {}, "clp", "getMatrixByCol: elements");
360      OSIUNITTEST_ASSERT_ERROR(eq(ev[11],-1.2), {}, "clp", "getMatrixByCol: elements");
361      OSIUNITTEST_ASSERT_ERROR(eq(ev[12],-1.0), {}, "clp", "getMatrixByCol: elements");
362      OSIUNITTEST_ASSERT_ERROR(eq(ev[13], 1.9), {}, "clp", "getMatrixByCol: elements");
363     
364      const CoinBigIndex * mi = smP->getVectorStarts();
365      OSIUNITTEST_ASSERT_ERROR(mi[0] ==  0, {}, "clp", "getMatrixByCol: vector starts");
366      OSIUNITTEST_ASSERT_ERROR(mi[1] ==  2, {}, "clp", "getMatrixByCol: vector starts");
367      OSIUNITTEST_ASSERT_ERROR(mi[2] ==  4, {}, "clp", "getMatrixByCol: vector starts");
368      OSIUNITTEST_ASSERT_ERROR(mi[3] ==  6, {}, "clp", "getMatrixByCol: vector starts");
369      OSIUNITTEST_ASSERT_ERROR(mi[4] ==  8, {}, "clp", "getMatrixByCol: vector starts");
370      OSIUNITTEST_ASSERT_ERROR(mi[5] == 10, {}, "clp", "getMatrixByCol: vector starts");
371      OSIUNITTEST_ASSERT_ERROR(mi[6] == 11, {}, "clp", "getMatrixByCol: vector starts");
372      OSIUNITTEST_ASSERT_ERROR(mi[7] == 12, {}, "clp", "getMatrixByCol: vector starts");
373      OSIUNITTEST_ASSERT_ERROR(mi[8] == 14, {}, "clp", "getMatrixByCol: vector starts");
374     
375      const int * ei = smP->getIndices();
376      OSIUNITTEST_ASSERT_ERROR(ei[ 0] == 0, {}, "clp", "getMatrixByCol: indices");
377      OSIUNITTEST_ASSERT_ERROR(ei[ 1] == 4, {}, "clp", "getMatrixByCol: indices");
378      OSIUNITTEST_ASSERT_ERROR(ei[ 2] == 0, {}, "clp", "getMatrixByCol: indices");
379      OSIUNITTEST_ASSERT_ERROR(ei[ 3] == 1, {}, "clp", "getMatrixByCol: indices");
380      OSIUNITTEST_ASSERT_ERROR(ei[ 4] == 1, {}, "clp", "getMatrixByCol: indices");
381      OSIUNITTEST_ASSERT_ERROR(ei[ 5] == 2, {}, "clp", "getMatrixByCol: indices");
382      OSIUNITTEST_ASSERT_ERROR(ei[ 6] == 0, {}, "clp", "getMatrixByCol: indices");
383      OSIUNITTEST_ASSERT_ERROR(ei[ 7] == 3, {}, "clp", "getMatrixByCol: indices");
384      OSIUNITTEST_ASSERT_ERROR(ei[ 8] == 0, {}, "clp", "getMatrixByCol: indices");
385      OSIUNITTEST_ASSERT_ERROR(ei[ 9] == 4, {}, "clp", "getMatrixByCol: indices");
386      OSIUNITTEST_ASSERT_ERROR(ei[10] == 2, {}, "clp", "getMatrixByCol: indices");
387      OSIUNITTEST_ASSERT_ERROR(ei[11] == 3, {}, "clp", "getMatrixByCol: indices");
388      OSIUNITTEST_ASSERT_ERROR(ei[12] == 0, {}, "clp", "getMatrixByCol: indices");
389      OSIUNITTEST_ASSERT_ERROR(ei[13] == 4, {}, "clp", "getMatrixByCol: indices");
390    }
391
392    //--------------
393    // Test rowsense, rhs, rowrange, matrixByRow
394    {
395      OsiClpSolverInterface lhs;
396      {     
397        OsiClpSolverInterface siC1(m);
398        OSIUNITTEST_ASSERT_WARNING(siC1.rowrange_ == NULL, {}, "clp", "row range");
399        OSIUNITTEST_ASSERT_WARNING(siC1.rowsense_ == NULL, {}, "clp", "row sense");
400        OSIUNITTEST_ASSERT_WARNING(siC1.rhs_ == NULL, {}, "clp", "right hand side");
401        OSIUNITTEST_ASSERT_WARNING(siC1.matrixByRow_ == NULL, {}, "clp", "matrix by row");
402
403        const char   * siC1rs  = siC1.getRowSense();
404        OSIUNITTEST_ASSERT_ERROR(siC1rs[0] == 'G', {}, "clp", "row sense");
405        OSIUNITTEST_ASSERT_ERROR(siC1rs[1] == 'L', {}, "clp", "row sense");
406        OSIUNITTEST_ASSERT_ERROR(siC1rs[2] == 'E', {}, "clp", "row sense");
407        OSIUNITTEST_ASSERT_ERROR(siC1rs[3] == 'R', {}, "clp", "row sense");
408        OSIUNITTEST_ASSERT_ERROR(siC1rs[4] == 'R', {}, "clp", "row sense");
409
410        const double * siC1rhs = siC1.getRightHandSide();
411        OSIUNITTEST_ASSERT_ERROR(eq(siC1rhs[0],2.5), {}, "clp", "right hand side");
412        OSIUNITTEST_ASSERT_ERROR(eq(siC1rhs[1],2.1), {}, "clp", "right hand side");
413        OSIUNITTEST_ASSERT_ERROR(eq(siC1rhs[2],4.0), {}, "clp", "right hand side");
414        OSIUNITTEST_ASSERT_ERROR(eq(siC1rhs[3],5.0), {}, "clp", "right hand side");
415        OSIUNITTEST_ASSERT_ERROR(eq(siC1rhs[4],15.), {}, "clp", "right hand side");
416
417        const double * siC1rr  = siC1.getRowRange();
418        OSIUNITTEST_ASSERT_ERROR(eq(siC1rr[0],0.0), {}, "clp", "row range");
419        OSIUNITTEST_ASSERT_ERROR(eq(siC1rr[1],0.0), {}, "clp", "row range");
420        OSIUNITTEST_ASSERT_ERROR(eq(siC1rr[2],0.0), {}, "clp", "row range");
421        OSIUNITTEST_ASSERT_ERROR(eq(siC1rr[3],5.0-1.8), {}, "clp", "row range");
422        OSIUNITTEST_ASSERT_ERROR(eq(siC1rr[4],15.0-3.0), {}, "clp", "row range");
423
424        const CoinPackedMatrix * siC1mbr = siC1.getMatrixByRow();
425        OSIUNITTEST_ASSERT_ERROR(siC1mbr != NULL, {}, "clp", "matrix by row");
426        OSIUNITTEST_ASSERT_ERROR(siC1mbr->getMajorDim()    ==  5, return, "clp", "matrix by row: major dim");
427        OSIUNITTEST_ASSERT_ERROR(siC1mbr->getNumElements() == 14, return, "clp", "matrix by row: num elements");
428
429        const double * ev = siC1mbr->getElements();
430        OSIUNITTEST_ASSERT_ERROR(eq(ev[ 0], 3.0), {}, "clp", "matrix by row: elements");
431        OSIUNITTEST_ASSERT_ERROR(eq(ev[ 1], 1.0), {}, "clp", "matrix by row: elements");
432        OSIUNITTEST_ASSERT_ERROR(eq(ev[ 2],-2.0), {}, "clp", "matrix by row: elements");
433        OSIUNITTEST_ASSERT_ERROR(eq(ev[ 3],-1.0), {}, "clp", "matrix by row: elements");
434        OSIUNITTEST_ASSERT_ERROR(eq(ev[ 4],-1.0), {}, "clp", "matrix by row: elements");
435        OSIUNITTEST_ASSERT_ERROR(eq(ev[ 5], 2.0), {}, "clp", "matrix by row: elements");
436        OSIUNITTEST_ASSERT_ERROR(eq(ev[ 6], 1.1), {}, "clp", "matrix by row: elements");
437        OSIUNITTEST_ASSERT_ERROR(eq(ev[ 7], 1.0), {}, "clp", "matrix by row: elements");
438        OSIUNITTEST_ASSERT_ERROR(eq(ev[ 8], 1.0), {}, "clp", "matrix by row: elements");
439        OSIUNITTEST_ASSERT_ERROR(eq(ev[ 9], 2.8), {}, "clp", "matrix by row: elements");
440        OSIUNITTEST_ASSERT_ERROR(eq(ev[10],-1.2), {}, "clp", "matrix by row: elements");
441        OSIUNITTEST_ASSERT_ERROR(eq(ev[11], 5.6), {}, "clp", "matrix by row: elements");
442        OSIUNITTEST_ASSERT_ERROR(eq(ev[12], 1.0), {}, "clp", "matrix by row: elements");
443        OSIUNITTEST_ASSERT_ERROR(eq(ev[13], 1.9), {}, "clp", "matrix by row: elements");
444
445        const CoinBigIndex * mi = siC1mbr->getVectorStarts();
446        OSIUNITTEST_ASSERT_ERROR(mi[0] ==  0, {}, "clp", "matrix by row: vector starts");
447        OSIUNITTEST_ASSERT_ERROR(mi[1] ==  5, {}, "clp", "matrix by row: vector starts");
448        OSIUNITTEST_ASSERT_ERROR(mi[2] ==  7, {}, "clp", "matrix by row: vector starts");
449        OSIUNITTEST_ASSERT_ERROR(mi[3] ==  9, {}, "clp", "matrix by row: vector starts");
450        OSIUNITTEST_ASSERT_ERROR(mi[4] == 11, {}, "clp", "matrix by row: vector starts");
451        OSIUNITTEST_ASSERT_ERROR(mi[5] == 14, {}, "clp", "matrix by row: vector starts");
452
453        const int * ei = siC1mbr->getIndices();
454        OSIUNITTEST_ASSERT_ERROR(ei[ 0] == 0, {}, "clp", "matrix by row: indices");
455        OSIUNITTEST_ASSERT_ERROR(ei[ 1] == 1, {}, "clp", "matrix by row: indices");
456        OSIUNITTEST_ASSERT_ERROR(ei[ 2] == 3, {}, "clp", "matrix by row: indices");
457        OSIUNITTEST_ASSERT_ERROR(ei[ 3] == 4, {}, "clp", "matrix by row: indices");
458        OSIUNITTEST_ASSERT_ERROR(ei[ 4] == 7, {}, "clp", "matrix by row: indices");
459        OSIUNITTEST_ASSERT_ERROR(ei[ 5] == 1, {}, "clp", "matrix by row: indices");
460        OSIUNITTEST_ASSERT_ERROR(ei[ 6] == 2, {}, "clp", "matrix by row: indices");
461        OSIUNITTEST_ASSERT_ERROR(ei[ 7] == 2, {}, "clp", "matrix by row: indices");
462        OSIUNITTEST_ASSERT_ERROR(ei[ 8] == 5, {}, "clp", "matrix by row: indices");
463        OSIUNITTEST_ASSERT_ERROR(ei[ 9] == 3, {}, "clp", "matrix by row: indices");
464        OSIUNITTEST_ASSERT_ERROR(ei[10] == 6, {}, "clp", "matrix by row: indices");
465        OSIUNITTEST_ASSERT_ERROR(ei[11] == 0, {}, "clp", "matrix by row: indices");
466        OSIUNITTEST_ASSERT_ERROR(ei[12] == 4, {}, "clp", "matrix by row: indices");
467        OSIUNITTEST_ASSERT_ERROR(ei[13] == 7, {}, "clp", "matrix by row: indices");
468
469        OSIUNITTEST_ASSERT_WARNING(siC1rs  == siC1.getRowSense(), {}, "clp", "row sense");
470        OSIUNITTEST_ASSERT_WARNING(siC1rhs == siC1.getRightHandSide(), {}, "clp", "right hand side");
471        OSIUNITTEST_ASSERT_WARNING(siC1rr  == siC1.getRowRange(), {}, "clp", "row range");
472
473        // Change CLP Model by adding free row
474        OsiRowCut rc;
475        rc.setLb(-COIN_DBL_MAX);
476        rc.setUb( COIN_DBL_MAX);
477        OsiCuts cuts;
478        cuts.insert(rc);
479        siC1.applyCuts(cuts);
480
481        // Since model was changed, test that cached data is now freed.
482        OSIUNITTEST_ASSERT_ERROR(siC1.rowrange_ == NULL, {}, "clp", "free cached data after adding row");
483        OSIUNITTEST_ASSERT_ERROR(siC1.rowsense_ == NULL, {}, "clp", "free cached data after adding row");
484        OSIUNITTEST_ASSERT_ERROR(siC1.rhs_ == NULL, {}, "clp", "free cached data after adding row");
485        // siC1.matrixByRow_ is updated, so it does not have to be NULL
486
487        siC1rs  = siC1.getRowSense();
488        OSIUNITTEST_ASSERT_ERROR(siC1rs[0] == 'G', {}, "clp", "row sense after adding row");
489        OSIUNITTEST_ASSERT_ERROR(siC1rs[1] == 'L', {}, "clp", "row sense after adding row");
490        OSIUNITTEST_ASSERT_ERROR(siC1rs[2] == 'E', {}, "clp", "row sense after adding row");
491        OSIUNITTEST_ASSERT_ERROR(siC1rs[3] == 'R', {}, "clp", "row sense after adding row");
492        OSIUNITTEST_ASSERT_ERROR(siC1rs[4] == 'R', {}, "clp", "row sense after adding row");
493        OSIUNITTEST_ASSERT_ERROR(siC1rs[5] == 'N', {}, "clp", "row sense after adding row");
494
495        siC1rhs = siC1.getRightHandSide();
496        OSIUNITTEST_ASSERT_ERROR(eq(siC1rhs[0],2.5), {}, "clp", "right hand side after adding row");
497        OSIUNITTEST_ASSERT_ERROR(eq(siC1rhs[1],2.1), {}, "clp", "right hand side after adding row");
498        OSIUNITTEST_ASSERT_ERROR(eq(siC1rhs[2],4.0), {}, "clp", "right hand side after adding row");
499        OSIUNITTEST_ASSERT_ERROR(eq(siC1rhs[3],5.0), {}, "clp", "right hand side after adding row");
500        OSIUNITTEST_ASSERT_ERROR(eq(siC1rhs[4],15.), {}, "clp", "right hand side after adding row");
501        OSIUNITTEST_ASSERT_ERROR(eq(siC1rhs[5],0.0), {}, "clp", "right hand side after adding row");
502
503        siC1rr  = siC1.getRowRange();
504        OSIUNITTEST_ASSERT_ERROR(eq(siC1rr[0],0.0), {}, "clp", "row range after adding row");
505        OSIUNITTEST_ASSERT_ERROR(eq(siC1rr[1],0.0), {}, "clp", "row range after adding row");
506        OSIUNITTEST_ASSERT_ERROR(eq(siC1rr[2],0.0), {}, "clp", "row range after adding row");
507        OSIUNITTEST_ASSERT_ERROR(eq(siC1rr[3],5.0-1.8), {}, "clp", "row range after adding row");
508        OSIUNITTEST_ASSERT_ERROR(eq(siC1rr[4],15.0-3.0), {}, "clp", "row range after adding row");
509        OSIUNITTEST_ASSERT_ERROR(eq(siC1rr[5],0.0), {}, "clp", "row range after adding row");
510
511        lhs = siC1;
512      }
513      // Test that lhs has correct values even though siC1 has gone out of scope
514      OSIUNITTEST_ASSERT_ERROR(lhs.rowrange_ == NULL, {}, "clp", "freed origin after assignment");
515      OSIUNITTEST_ASSERT_ERROR(lhs.rowsense_ == NULL, {}, "clp", "freed origin after assignment");
516      OSIUNITTEST_ASSERT_ERROR(lhs.rhs_ == NULL, {}, "clp", "freed origin after assignment");
517      OSIUNITTEST_ASSERT_ERROR(lhs.matrixByRow_ == NULL, {}, "clp", "freed origin after assignment");
518
519      const char * lhsrs  = lhs.getRowSense();
520      OSIUNITTEST_ASSERT_ERROR(lhsrs[0] == 'G', {}, "clp", "row sense after assignment");
521      OSIUNITTEST_ASSERT_ERROR(lhsrs[1] == 'L', {}, "clp", "row sense after assignment");
522      OSIUNITTEST_ASSERT_ERROR(lhsrs[2] == 'E', {}, "clp", "row sense after assignment");
523      OSIUNITTEST_ASSERT_ERROR(lhsrs[3] == 'R', {}, "clp", "row sense after assignment");
524      OSIUNITTEST_ASSERT_ERROR(lhsrs[4] == 'R', {}, "clp", "row sense after assignment");
525      OSIUNITTEST_ASSERT_ERROR(lhsrs[5] == 'N', {}, "clp", "row sense after assignment");
526     
527      const double * lhsrhs = lhs.getRightHandSide();
528      OSIUNITTEST_ASSERT_ERROR(eq(lhsrhs[0],2.5), {}, "clp", "right hand side after assignment");
529      OSIUNITTEST_ASSERT_ERROR(eq(lhsrhs[1],2.1), {}, "clp", "right hand side after assignment");
530      OSIUNITTEST_ASSERT_ERROR(eq(lhsrhs[2],4.0), {}, "clp", "right hand side after assignment");
531      OSIUNITTEST_ASSERT_ERROR(eq(lhsrhs[3],5.0), {}, "clp", "right hand side after assignment");
532      OSIUNITTEST_ASSERT_ERROR(eq(lhsrhs[4],15.), {}, "clp", "right hand side after assignment");
533      OSIUNITTEST_ASSERT_ERROR(eq(lhsrhs[5],0.0), {}, "clp", "right hand side after assignment");
534     
535      const double *lhsrr = lhs.getRowRange();
536      OSIUNITTEST_ASSERT_ERROR(eq(lhsrr[0],0.0), {}, "clp", "row range after assignment");
537      OSIUNITTEST_ASSERT_ERROR(eq(lhsrr[1],0.0), {}, "clp", "row range after assignment");
538      OSIUNITTEST_ASSERT_ERROR(eq(lhsrr[2],0.0), {}, "clp", "row range after assignment");
539      OSIUNITTEST_ASSERT_ERROR(eq(lhsrr[3],5.0-1.8), {}, "clp", "row range after assignment");
540      OSIUNITTEST_ASSERT_ERROR(eq(lhsrr[4],15.0-3.0), {}, "clp", "row range after assignment");
541      OSIUNITTEST_ASSERT_ERROR(eq(lhsrr[5],0.0), {}, "clp", "row range after assignment");
542     
543      const CoinPackedMatrix * lhsmbr = lhs.getMatrixByRow();
544      OSIUNITTEST_ASSERT_ERROR(lhsmbr != NULL, {}, "clp", "matrix by row after assignment");
545      OSIUNITTEST_ASSERT_ERROR(lhsmbr->getMajorDim()    ==  6, return, "clp", "matrix by row after assignment: major dim");
546      OSIUNITTEST_ASSERT_ERROR(lhsmbr->getNumElements() == 14, return, "clp", "matrix by row after assignment: num elements");
547
548      const double * ev = lhsmbr->getElements();
549      OSIUNITTEST_ASSERT_ERROR(eq(ev[ 0], 3.0), {}, "clp", "matrix by row after assignment: elements");
550      OSIUNITTEST_ASSERT_ERROR(eq(ev[ 1], 1.0), {}, "clp", "matrix by row after assignment: elements");
551      OSIUNITTEST_ASSERT_ERROR(eq(ev[ 2],-2.0), {}, "clp", "matrix by row after assignment: elements");
552      OSIUNITTEST_ASSERT_ERROR(eq(ev[ 3],-1.0), {}, "clp", "matrix by row after assignment: elements");
553      OSIUNITTEST_ASSERT_ERROR(eq(ev[ 4],-1.0), {}, "clp", "matrix by row after assignment: elements");
554      OSIUNITTEST_ASSERT_ERROR(eq(ev[ 5], 2.0), {}, "clp", "matrix by row after assignment: elements");
555      OSIUNITTEST_ASSERT_ERROR(eq(ev[ 6], 1.1), {}, "clp", "matrix by row after assignment: elements");
556      OSIUNITTEST_ASSERT_ERROR(eq(ev[ 7], 1.0), {}, "clp", "matrix by row after assignment: elements");
557      OSIUNITTEST_ASSERT_ERROR(eq(ev[ 8], 1.0), {}, "clp", "matrix by row after assignment: elements");
558      OSIUNITTEST_ASSERT_ERROR(eq(ev[ 9], 2.8), {}, "clp", "matrix by row after assignment: elements");
559      OSIUNITTEST_ASSERT_ERROR(eq(ev[10],-1.2), {}, "clp", "matrix by row after assignment: elements");
560      OSIUNITTEST_ASSERT_ERROR(eq(ev[11], 5.6), {}, "clp", "matrix by row after assignment: elements");
561      OSIUNITTEST_ASSERT_ERROR(eq(ev[12], 1.0), {}, "clp", "matrix by row after assignment: elements");
562      OSIUNITTEST_ASSERT_ERROR(eq(ev[13], 1.9), {}, "clp", "matrix by row after assignment: elements");
563     
564      const CoinBigIndex * mi = lhsmbr->getVectorStarts();
565      OSIUNITTEST_ASSERT_ERROR(mi[0] ==  0, {}, "clp", "matrix by row after assignment: vector starts");
566      OSIUNITTEST_ASSERT_ERROR(mi[1] ==  5, {}, "clp", "matrix by row after assignment: vector starts");
567      OSIUNITTEST_ASSERT_ERROR(mi[2] ==  7, {}, "clp", "matrix by row after assignment: vector starts");
568      OSIUNITTEST_ASSERT_ERROR(mi[3] ==  9, {}, "clp", "matrix by row after assignment: vector starts");
569      OSIUNITTEST_ASSERT_ERROR(mi[4] == 11, {}, "clp", "matrix by row after assignment: vector starts");
570      OSIUNITTEST_ASSERT_ERROR(mi[5] == 14, {}, "clp", "matrix by row after assignment: vector starts");
571     
572      const int * ei = lhsmbr->getIndices();
573      OSIUNITTEST_ASSERT_ERROR(ei[ 0] == 0, {}, "clp", "matrix by row after assignment: indices");
574      OSIUNITTEST_ASSERT_ERROR(ei[ 1] == 1, {}, "clp", "matrix by row after assignment: indices");
575      OSIUNITTEST_ASSERT_ERROR(ei[ 2] == 3, {}, "clp", "matrix by row after assignment: indices");
576      OSIUNITTEST_ASSERT_ERROR(ei[ 3] == 4, {}, "clp", "matrix by row after assignment: indices");
577      OSIUNITTEST_ASSERT_ERROR(ei[ 4] == 7, {}, "clp", "matrix by row after assignment: indices");
578      OSIUNITTEST_ASSERT_ERROR(ei[ 5] == 1, {}, "clp", "matrix by row after assignment: indices");
579      OSIUNITTEST_ASSERT_ERROR(ei[ 6] == 2, {}, "clp", "matrix by row after assignment: indices");
580      OSIUNITTEST_ASSERT_ERROR(ei[ 7] == 2, {}, "clp", "matrix by row after assignment: indices");
581      OSIUNITTEST_ASSERT_ERROR(ei[ 8] == 5, {}, "clp", "matrix by row after assignment: indices");
582      OSIUNITTEST_ASSERT_ERROR(ei[ 9] == 3, {}, "clp", "matrix by row after assignment: indices");
583      OSIUNITTEST_ASSERT_ERROR(ei[10] == 6, {}, "clp", "matrix by row after assignment: indices");
584      OSIUNITTEST_ASSERT_ERROR(ei[11] == 0, {}, "clp", "matrix by row after assignment: indices");
585      OSIUNITTEST_ASSERT_ERROR(ei[12] == 4, {}, "clp", "matrix by row after assignment: indices");
586      OSIUNITTEST_ASSERT_ERROR(ei[13] == 7, {}, "clp", "matrix by row after assignment: indices");
587    }
588  }
589
590  // Test ClpPlusMinusOneMatrix by way of loadProblem(ClpMatrixBase, ... )
591  { int pos_start[4] = {0,5,9,12};
592    int neg_start[4] = {3,7,11,12};
593    int col[12] = {0,1,2,3,4,5,6,7,0,1,2,3};
594    double rhs[3] = {0.0,0.0,0.0};
595    double cost[8];
596    double var_lb[8];
597    double var_ub[8];
598    for (int i = 0 ; i < 8 ; i++) {
599      cost[i] = 1.0;
600      var_lb[i] = 0.0;
601      var_ub[i] = 1.0;
602    }
603    ClpPlusMinusOneMatrix pmone_matrix(3,8,false,col,pos_start,neg_start);
604    OsiClpSolverInterface clpSi;
605    OSIUNITTEST_CATCH_ERROR(
606        {clpSi.loadProblem(pmone_matrix,var_lb,var_ub,cost,rhs,rhs);
607         clpSi.initialSolve();},
608        {},"clp","loadProblem(ClpMatrixBase, ...)")
609  }
610
611
612  // Test add/delete columns
613  {   
614    OsiClpSolverInterface m;
615    std::string fn = mpsDir+"p0033";
616    m.readMps(fn.c_str(),"mps");
617    double inf = m.getInfinity();
618
619    CoinPackedVector c0;
620    c0.insert(0, 4);
621    c0.insert(1, 1);
622    m.addCol(c0, 0, inf, 3);
623    m.initialSolve();
624    double objValue = m.getObjValue();
625    CoinRelFltEq eq(1.0e-2);
626    OSIUNITTEST_ASSERT_ERROR(eq(objValue,2520.57), {}, "clp", "objvalue after adding col");
627
628    // Try deleting first column
629    int * d = new int[1];
630    d[0]=0;
631    m.deleteCols(1,d);
632    delete [] d;
633    d=NULL;
634    m.resolve();
635    objValue = m.getObjValue();
636    OSIUNITTEST_ASSERT_ERROR(eq(objValue,2520.57), {}, "clp", "objvalue after deleting first col");
637
638    // Try deleting column we added
639    int iCol = m.getNumCols()-1;
640    m.deleteCols(1,&iCol);
641    m.resolve();
642    objValue = m.getObjValue();
643    OSIUNITTEST_ASSERT_ERROR(eq(objValue,2520.57), {}, "clp", "objvalue after deleting added col");
644  }
645
646  // Test branch and bound
647  {   
648    OsiClpSolverInterface m;
649    std::string fn = mpsDir+"p0033";
650    m.readMps(fn.c_str(),"mps");
651    // reduce printout
652    m.setHintParam(OsiDoReducePrint,true,OsiHintTry);
653    // test maximization
654    int n  = m.getNumCols();
655    int i;
656    double * obj2 = new double[n];
657    const double * obj = m.getObjCoefficients();
658    for (i=0;i<n;i++) {
659      obj2[i]=-obj[i];
660    }
661    m.setObjective(obj2);
662    delete [] obj2;
663    m.setObjSense(-1.0);
664    // Save bounds
665    double * saveUpper = CoinCopyOfArray(m.getColUpper(),n);
666    double * saveLower = CoinCopyOfArray(m.getColLower(),n);
667    m.branchAndBound();
668    // reset bounds
669    m.setColUpper(saveUpper);
670    m.setColLower(saveLower);
671    delete [] saveUpper;
672    delete [] saveLower;
673    // reset cutoff - otherwise we won't get solution
674    m.setDblParam(OsiDualObjectiveLimit,-COIN_DBL_MAX);
675    m.branchAndBound();
676    double objValue = m.getObjValue();
677    CoinRelFltEq eq(1.0e-2);
678    OSIUNITTEST_ASSERT_ERROR(eq(objValue,-3089), {}, "clp", "branch and bound");
679    const double * cs = m.getColSolution();
680    if( OsiUnitTest::verbosity >= 2 )
681      for ( i=0;i<n;i++) {
682        if (cs[i]>1.0e-7)
683          printf("%d has value %g\n",i,cs[i]);
684    }
685  }
686
687  // Test infeasible bounds
688  {
689    OsiClpSolverInterface solver;
690    std::string fn = mpsDir+"exmip1";
691    solver.readMps(fn.c_str(),"mps");
692    int index[]={0};
693    double value[]={1.0,0.0};
694    solver.setColSetBounds(index,index+1,value);
695    solver.setHintParam(OsiDoPresolveInInitial, false, OsiHintDo);
696    solver.initialSolve();
697    OSIUNITTEST_ASSERT_ERROR(!solver.isProvenOptimal(), {}, "clp", "infeasible bounds");
698  }
699
700  // Build a model
701  {   
702    OsiClpSolverInterface model;
703    std::string fn = mpsDir+"p0033";
704    model.readMps(fn.c_str(),"mps");
705    // Point to data
706    int numberRows = model.getNumRows();
707    const double * rowLower = model.getRowLower();
708    const double * rowUpper = model.getRowUpper();
709    int numberColumns = model.getNumCols();
710    const double * columnLower = model.getColLower();
711    const double * columnUpper = model.getColUpper();
712    const double * columnObjective = model.getObjCoefficients();
713    // get row copy
714    CoinPackedMatrix rowCopy = *model.getMatrixByRow();
715    const int * column = rowCopy.getIndices();
716    const int * rowLength = rowCopy.getVectorLengths();
717    const CoinBigIndex * rowStart = rowCopy.getVectorStarts();
718    const double * element = rowCopy.getElements();
719   
720    // solve
721    model.initialSolve();
722    // Now build new model
723    CoinModel build;
724    // Row bounds
725    int iRow;
726    for (iRow=0;iRow<numberRows;iRow++) {
727      build.setRowBounds(iRow,rowLower[iRow],rowUpper[iRow]);
728      build.setRowName(iRow,model.getRowName(iRow).c_str());
729    }
730    // Column bounds and objective
731    int iColumn;
732    for (iColumn=0;iColumn<numberColumns;iColumn++) {
733      build.setColumnLower(iColumn,columnLower[iColumn]);
734      build.setColumnUpper(iColumn,columnUpper[iColumn]);
735      build.setObjective(iColumn,columnObjective[iColumn]);
736      build.setColumnName(iColumn,model.getColName(iColumn).c_str());
737    }
738    // Adds elements one by one by row (backwards by row)
739    for (iRow=numberRows-1;iRow>=0;iRow--) {
740      int start = rowStart[iRow];
741      for (int j=start;j<start+rowLength[iRow];j++) 
742        build(iRow,column[j],element[j]);
743    }
744    // Now create Model
745    OsiClpSolverInterface model2;
746    model2.loadFromCoinModel(build);
747    model2.initialSolve();
748    // Save - should be continuous
749    model2.writeMps("continuous");
750    int * whichInteger = new int[numberColumns];
751    for (iColumn=0;iColumn<numberColumns;iColumn++) 
752      whichInteger[iColumn]=iColumn;
753    // mark as integer
754    model2.setInteger(whichInteger,numberColumns);
755    delete [] whichInteger;
756    // save - should be integer
757    model2.writeMps("integer");
758    // check names are there
759    if (OsiUnitTest::verbosity >= 1)
760      for (iColumn=0;iColumn<numberColumns;iColumn++)
761        printf("%d name %s\n",iColumn,model2.getColName(iColumn).c_str());
762   
763    // Now do with strings attached
764    // Save build to show how to go over rows
765    CoinModel saveBuild = build;
766    build = CoinModel();
767    // Column bounds
768    for (iColumn=0;iColumn<numberColumns;iColumn++) {
769      build.setColumnLower(iColumn,columnLower[iColumn]);
770      build.setColumnUpper(iColumn,columnUpper[iColumn]);
771    }
772    // Objective - half the columns as is and half with multiplier of "1.0+multiplier"
773    // Pick up from saveBuild (for no reason at all)
774    for (iColumn=0;iColumn<numberColumns;iColumn++) {
775      double value = saveBuild.objective(iColumn);
776      if (iColumn*2<numberColumns) {
777        build.setObjective(iColumn,columnObjective[iColumn]);
778      } else {
779        // create as string
780        char temp[100];
781        sprintf(temp,"%g + abs(%g*multiplier)",value,value);
782        build.setObjective(iColumn,temp);
783      }
784    }
785    // It then adds rows one by one but for half the rows sets their values
786    //      with multiplier of "1.0+1.5*multiplier"
787    for (iRow=0;iRow<numberRows;iRow++) {
788      if (iRow*2<numberRows) {
789        // add row in simple way
790        int start = rowStart[iRow];
791        build.addRow(rowLength[iRow],column+start,element+start,
792                     rowLower[iRow],rowUpper[iRow]);
793      } else {
794        // As we have to add one by one let's get from saveBuild
795        CoinModelLink triple=saveBuild.firstInRow(iRow);
796        while (triple.column()>=0) {
797          int iColumn = triple.column();
798          if (iColumn*2<numberColumns) {
799            // just value as normal
800            build(iRow,triple.column(),triple.value());
801          } else {
802            // create as string
803            char temp[100];
804            sprintf(temp,"%g + (1.5*%g*multiplier)",triple.value(), triple.value());
805            build(iRow,iColumn,temp);
806          }
807          triple=saveBuild.next(triple);
808        }
809        // but remember to do rhs
810        build.setRowLower(iRow,rowLower[iRow]);
811        build.setRowUpper(iRow,rowUpper[iRow]);
812      }
813    }
814    // If small switch on error printing
815    if (numberColumns<50)
816      build.setLogLevel(1);
817    // should fail as we never set multiplier
818    OSIUNITTEST_ASSERT_ERROR(model2.loadFromCoinModel(build) != 0, {}, "clp", "build model with missing multipliers");
819    build.associateElement("multiplier",0.0);
820    OSIUNITTEST_ASSERT_ERROR(model2.loadFromCoinModel(build) == 0, {}, "clp", "build model");
821    model2.initialSolve();
822    // It then loops with multiplier going from 0.0 to 2.0 in increments of 0.1
823    for (double multiplier=0.0;multiplier<2.0;multiplier+= 0.1) {
824      build.associateElement("multiplier",multiplier);
825      OSIUNITTEST_ASSERT_ERROR(model2.loadFromCoinModel(build,true) == 0, {}, "clp", "build model with increasing multiplier");
826      model2.resolve();
827    }
828  }
829
830  // Solve an lp by hand
831  {   
832    OsiClpSolverInterface m;
833    std::string fn = mpsDir+"p0033";
834    m.readMps(fn.c_str(),"mps");
835    m.setObjSense(-1.0);
836    m.getModelPtr()->messageHandler()->setLogLevel(4);
837    m.initialSolve();
838    m.getModelPtr()->factorization()->maximumPivots(5);
839    m.setObjSense(1.0);
840    // clone to test pivot as well as primalPivot
841    OsiSolverInterface * mm =m.clone();
842    // enable special mode
843    m.enableSimplexInterface(true);
844    // need to do after clone
845    mm->enableSimplexInterface(true);
846    // we happen to know that variables are 0-1 and rows are L
847    int numberIterations=0;
848    int numberColumns = m.getNumCols();
849    int numberRows = m.getNumRows();
850    double * fakeCost = new double[numberColumns];
851    double * duals = new double [numberRows];
852    double * djs = new double [numberColumns];
853    const double * solution = m.getColSolution();
854    memcpy(fakeCost,m.getObjCoefficients(),numberColumns*sizeof(double));
855    while (1) {
856      const double * dj;
857      const double * dual;
858      if ((numberIterations&1)==0) {
859        // use given ones
860        dj = m.getReducedCost();
861        dual = m.getRowPrice();
862      } else {
863        // create
864        dj = djs;
865        dual = duals;
866        m.getReducedGradient(djs,duals,fakeCost);
867      }
868      int i;
869      int colIn=9999;
870      int direction=1;
871      double best=1.0e-6;
872      // find most negative reduced cost
873      // Should check basic - but should be okay on this problem
874      for (i=0;i<numberRows;i++) {
875        double value=dual[i];
876        if (value>best) {
877          direction=-1;
878          best=value;
879          colIn=-i-1;
880        }
881      }
882      for (i=0;i<numberColumns;i++) {
883        double value=dj[i];
884        if (value<-best&&solution[i]<1.0e-6) {
885          direction=1;
886          best=-value;
887          colIn=i;
888        } else if (value>best&&solution[i]>1.0-1.0e-6) {
889          direction=-1;
890          best=value;
891          colIn=i;
892        }
893      }
894      if (colIn==9999)
895        break; // should be optimal
896      int colOut;
897      int outStatus;
898      double theta;
899      OSIUNITTEST_ASSERT_ERROR(m.primalPivotResult(colIn,direction,colOut,outStatus,theta,NULL) == 0, {}, "clp", "solve model by hand");
900      if (OsiUnitTest::verbosity >= 1)
901        printf("out %d, direction %d theta %g\n", colOut,outStatus,theta);
902      if (colIn!=colOut) {
903        OSIUNITTEST_ASSERT_ERROR(mm->pivot(colIn,colOut,outStatus) >= 0, {}, "clp", "solve model by hand");
904      } else {
905        // bound flip (so pivot does not make sense)
906        OSIUNITTEST_ASSERT_ERROR(mm->primalPivotResult(colIn,direction,colOut,outStatus,theta,NULL) == 0, {}, "clp", "solve model by hand");
907      }
908      numberIterations++;
909    }
910    delete mm;
911    delete [] fakeCost;
912    delete [] duals;
913    delete [] djs;
914    // exit special mode
915    m.disableSimplexInterface();
916    m.getModelPtr()->messageHandler()->setLogLevel(4);
917    m.messageHandler()->setLogLevel(0);
918    m.resolve();
919    OSIUNITTEST_ASSERT_ERROR(m.getIterationCount() == 0, {}, "clp", "resolve after solving model by hand");
920    m.setObjSense(-1.0);
921    m.initialSolve();
922  }
923
924# if 0
925/*
926  This section stops working without setObjectiveAndRefresh. Assertion failure
927  down in the guts of clp, likely due to reduced costs not properly updated.
928  Leave the code in for a bit so it's easily recoverable if anyone actually
929  yells about the loss. There was no response to a public announcement
930  of intent to delete, but sometimes it takes a whack on the head to get
931  peoples' attention. At some point, it'd be good to come back through and
932  make this work again. -- lh, 100828 --
933*/
934  // Do parametrics on the objective by hand
935  {
936    std::cout << " Beginning Osi Simplex mode 2 ... " << std::endl ;
937    OsiClpSolverInterface m;
938    std::string fn = mpsDir+"p0033";
939    m.readMps(fn.c_str(),"mps");
940    ClpSimplex * simplex = m.getModelPtr();
941    simplex->messageHandler()->setLogLevel(4);
942    m.initialSolve();
943    simplex->factorization()->maximumPivots(5);
944    simplex->messageHandler()->setLogLevel(63);
945    m.setObjSense(1.0);
946    // enable special mode
947    m.enableSimplexInterface(true);
948    int numberIterations=0;
949    int numberColumns = m.getNumCols();
950    int numberRows = m.getNumRows();
951    double * changeCost = new double[numberColumns];
952    double * duals = new double [numberRows];
953    double * djs = new double [numberColumns];
954    // As getReducedGradient mucks about with innards of Clp get arrays to save
955    double * dualsNow = new double [numberRows];
956    double * djsNow = new double [numberColumns];
957   
958    const double * solution = m.getColSolution();
959    int i;
960    // Set up change cost
961    for (i=0;i<numberColumns;i++)
962      changeCost[i]=1.0+0.1*i;;
963    // Range of investigation
964    double totalChange=100.0;
965    double totalDone=0.0;
966    while (true) {
967      std::cout << " Starting iterations ... " << std::endl ;
968      // Save current
969      // (would be more accurate to start from scratch)
970      memcpy(djsNow, m.getReducedCost(),numberColumns*sizeof(double));
971      memcpy(dualsNow,m.getRowPrice(),numberRows*sizeof(double));
972      // Get reduced gradient of changeCost
973      m.getReducedGradient(djs,duals,changeCost);
974      int colIn=9999;
975      int direction=1;
976      // We are going up to totalChange but we have done some
977      double best=totalChange-totalDone;
978      // find best ratio
979      // Should check basic - but should be okay on this problem
980      // We are cheating here as we know L rows
981      // Really should be using L and U status but I modified from previous example
982      for (i=0;i<numberRows;i++) {
983        if (simplex->getRowStatus(i)==ClpSimplex::basic) {
984          assert (fabs(dualsNow[i])<1.0e-4&&fabs(duals[i])<1.0e-4);
985        } else {
986          assert (dualsNow[i]<1.0e-4);
987          if (duals[i]>1.0e-8) {
988            if (dualsNow[i]+best*duals[i]>0.0) {
989              best = CoinMax(-dualsNow[i]/duals[i],0.0);
990              direction=-1;
991              colIn=-i-1;
992            }
993          }
994        }
995      }
996      for (i=0;i<numberColumns;i++) {
997        if (simplex->getColumnStatus(i)==ClpSimplex::basic) {
998          assert (fabs(djsNow[i])<1.0e-4&&fabs(djs[i])<1.0e-4);
999        } else {
1000          if (solution[i]<1.0e-6) {
1001            assert (djsNow[i]>-1.0e-4);
1002            if (djs[i]<-1.0e-8) {
1003              if (djsNow[i]+best*djs[i]<0.0) {
1004                best = CoinMax(-djsNow[i]/djs[i],0.0);
1005                direction=1;
1006                colIn=i;
1007              }
1008            }
1009          } else if (solution[i]>1.0-1.0e-6) {
1010            assert (djsNow[i]<1.0e-4);
1011            if (djs[i]>1.0e-8) {
1012              if (djsNow[i]+best*djs[i]>0.0) {
1013                best = CoinMax(-djsNow[i]/djs[i],0.0);
1014                direction=-1;
1015                colIn=i;
1016              }
1017            }
1018          }
1019        }
1020      }
1021      if (colIn==9999)
1022        break; // should be optimal
1023      // update objective - djs is spare array
1024      const double * obj = m.getObjCoefficients();
1025      for (i=0;i<numberColumns;i++) {
1026        djs[i]=obj[i]+best*changeCost[i];
1027      }
1028      totalDone += best;
1029      printf("Best change %g, total %g\n",best,totalDone);
1030      m.setObjectiveAndRefresh(djs);
1031      int colOut;
1032      int outStatus;
1033      double theta;
1034      int returnCode=m.primalPivotResult(colIn,direction,colOut,outStatus,theta,NULL);
1035      assert (!returnCode);
1036      double objValue = m.getObjValue(); // printed one may not be accurate
1037      printf("in %d out %d, direction %d theta %g, objvalue %g\n",
1038             colIn,colOut,outStatus,theta,objValue);
1039      numberIterations++;
1040    }
1041    // update objective to totalChange- djs is spare array
1042    const double * obj = m.getObjCoefficients();
1043    double best = totalChange-totalDone;
1044    for (i=0;i<numberColumns;i++) {
1045      djs[i]=obj[i]+best*changeCost[i];
1046    }
1047    delete [] changeCost;
1048    delete [] duals;
1049    delete [] djs;
1050    delete [] dualsNow;
1051    delete [] djsNow;
1052    // exit special mode
1053    std::cout << " Finished simplex mode 2 ; checking result." << std::endl ;
1054    m.disableSimplexInterface();
1055    simplex->messageHandler()->setLogLevel(4);
1056    m.resolve();
1057    assert (!m.getIterationCount());
1058  }
1059# endif
1060
1061  // Solve an lp when interface is on
1062  {   
1063    OsiClpSolverInterface m;
1064    std::string fn = mpsDir+"p0033";
1065    m.readMps(fn.c_str(),"mps");
1066    // enable special mode
1067    m.setHintParam(OsiDoScale,false,OsiHintDo);
1068    m.setHintParam(OsiDoPresolveInInitial,false,OsiHintDo);
1069    m.setHintParam(OsiDoDualInInitial,false,OsiHintDo);
1070    m.setHintParam(OsiDoPresolveInResolve,false,OsiHintDo);
1071    m.setHintParam(OsiDoDualInResolve,false,OsiHintDo);
1072    m.enableSimplexInterface(true);
1073    m.initialSolve();
1074  }
1075
1076  // Check tableau stuff when simplex interface is on
1077  {   
1078    OsiClpSolverInterface m;
1079    /*
1080       Wolsey : Page 130
1081       max 4x1 -  x2
1082       7x1 - 2x2    <= 14
1083       x2    <= 3
1084       2x1 - 2x2    <= 3
1085       x1 in Z+, x2 >= 0
1086    */
1087   
1088    double inf_ = m.getInfinity();
1089    int n_cols = 2;
1090    int n_rows = 3;
1091   
1092    double obj[2] = {-4.0, 1.0};
1093    double collb[2] = {0.0, 0.0};
1094    double colub[2] = {inf_, inf_};
1095    double rowlb[3] = {-inf_, -inf_, -inf_};
1096    double rowub[3] = {14.0, 3.0, 3.0};
1097   
1098    int rowIndices[5] =  {0,     2,    0,    1,    2};
1099    int colIndices[5] =  {0,     0,    1,    1,    1};
1100    double elements[5] = {7.0, 2.0, -2.0,  1.0, -2.0};
1101    CoinPackedMatrix M(true, rowIndices, colIndices, elements, 5);
1102   
1103    m.loadProblem(M, collb, colub, obj, rowlb, rowub);
1104    m.enableSimplexInterface(true);
1105   
1106    m.initialSolve();
1107   
1108    //check that the tableau matches wolsey (B-1 A)
1109    // slacks in second part of binvA
1110    double * binvA = (double*) malloc((n_cols+n_rows) * sizeof(double));
1111   
1112    if (OsiUnitTest::verbosity >= 2)
1113      printf("B-1 A");
1114    int i;
1115    for( i = 0; i < n_rows; i++){
1116      m.getBInvARow(i, binvA,binvA+n_cols);
1117      if (OsiUnitTest::verbosity >= 2) {
1118        printf("\nrow: %d -> ",i);
1119        for(int j=0; j < n_cols+n_rows; j++)
1120          printf("%g, ", binvA[j]);
1121      }
1122    }
1123    if (OsiUnitTest::verbosity >= 2) {
1124      printf("\n");
1125      printf("And by column");
1126    }
1127    for( i = 0; i < n_cols+n_rows; i++){
1128      m.getBInvACol(i, binvA);
1129      if (OsiUnitTest::verbosity >= 2) {
1130        printf("\ncolumn: %d -> ",i);
1131        for(int j=0; j < n_rows; j++)
1132          printf("%g, ", binvA[j]);
1133      }
1134    }
1135    if (OsiUnitTest::verbosity >= 2)
1136      printf("\n");
1137
1138    // and when doing as expert
1139    m.setSpecialOptions(m.specialOptions()|512);
1140    ClpSimplex * clp = m.getModelPtr();
1141    /* Do twice -
1142       first time with enableSimplexInterface still set
1143       then without and with scaling
1144    */
1145    for (int iPass=0;iPass<2;iPass++) {
1146      const double * rowScale = clp->rowScale();
1147      const double * columnScale = clp->columnScale();
1148#     if 0
1149      if (!iPass)
1150        assert (!rowScale);
1151      else
1152        assert (rowScale); // only true for this example
1153#     endif
1154      /* has to be exactly correct as in OsiClpsolverInterface.cpp
1155         (also redo each pass as may change
1156      */
1157      CoinIndexedVector * rowArray = clp->rowArray(1);
1158      CoinIndexedVector * columnArray = clp->columnArray(0);
1159      int n;
1160      int * which;
1161      double * array;
1162      if (OsiUnitTest::verbosity >= 2)
1163        printf("B-1 A");
1164      for( i = 0; i < n_rows; i++){
1165        m.getBInvARow(i, binvA,binvA+n_cols);
1166        if (OsiUnitTest::verbosity >= 2)
1167          printf("\nrow: %d -> ",i);
1168        int j;
1169        // First columns
1170        n = columnArray->getNumElements();
1171        which = columnArray->getIndices();
1172        array = columnArray->denseVector();
1173        for(j=0; j < n; j++){
1174          int k=which[j];
1175          if (OsiUnitTest::verbosity >= 2) {
1176            if (!columnScale) {
1177              printf("(%d %g), ", k, array[k]);
1178            } else {
1179              printf("(%d %g), ", k, array[k]/columnScale[k]);
1180            }
1181          }
1182          // zero out
1183          array[k]=0.0;
1184        }
1185        // say empty
1186        columnArray->setNumElements(0);
1187        // and check (would not be in any production code)
1188        columnArray->checkClear();
1189        // now rows
1190        n = rowArray->getNumElements();
1191        which = rowArray->getIndices();
1192        array = rowArray->denseVector();
1193        for(j=0; j < n; j++){
1194          int k=which[j];
1195          if (OsiUnitTest::verbosity >= 2) {
1196            if (!rowScale) {
1197              printf("(%d %g), ", k+n_cols, array[k]);
1198            } else {
1199              printf("(%d %g), ", k+n_cols, array[k]*rowScale[k]);
1200            }
1201          }
1202          // zero out
1203          array[k]=0.0;
1204        }
1205        // say empty
1206        rowArray->setNumElements(0);
1207        // and check (would not be in any production code)
1208        rowArray->checkClear();
1209      }
1210      if (OsiUnitTest::verbosity >= 2) {
1211        printf("\n");
1212        printf("And by column (trickier)");
1213      }
1214      const int * pivotVariable = clp->pivotVariable();
1215      for( i = 0; i < n_cols+n_rows; i++){
1216        m.getBInvACol(i, binvA);
1217        if (OsiUnitTest::verbosity >= 2)
1218          printf("\ncolumn: %d -> ",i);
1219        n = rowArray->getNumElements();
1220        which = rowArray->getIndices();
1221        array = rowArray->denseVector();
1222        for(int j=0; j < n; j++){
1223          int k=which[j];
1224          // need to know pivot variable for +1/-1 (slack) and row/column scaling
1225          int pivot = pivotVariable[k];
1226          if (OsiUnitTest::verbosity >= 2) {
1227            if (pivot<n_cols) {
1228              // scaled coding is in just in case
1229              if (!columnScale) {
1230                printf("(%d %g), ", k, array[k]);
1231              } else {
1232                printf("(%d %g), ", k, array[k]*columnScale[pivot]);
1233              }
1234            } else {
1235              if (!rowScale) {
1236                printf("(%d %g), ", k, -array[k]);
1237              } else {
1238                printf("(%d %g), ", k, -array[k]/rowScale[pivot-n_cols]);
1239              }
1240            }
1241          }
1242          // zero out
1243          array[k]=0.0;
1244        }
1245        // say empty
1246        rowArray->setNumElements(0);
1247        // and check (would not be in any production code)
1248        rowArray->checkClear();
1249      }
1250      if (OsiUnitTest::verbosity >= 2)
1251        printf("\n");
1252      // now deal with next pass
1253      if (!iPass) {
1254        m.disableSimplexInterface();
1255        // see if we can get scaling for testing
1256        clp->scaling(1);
1257        m.enableFactorization();
1258      } else {
1259        // may not be needed - but cleaner
1260        m.disableFactorization();
1261      }
1262    }
1263    m.setSpecialOptions(m.specialOptions()&~512);
1264    free(binvA);
1265    // Do using newer interface
1266    m.enableFactorization();
1267    {
1268      CoinIndexedVector * rowArray = new CoinIndexedVector(n_rows);
1269      CoinIndexedVector * columnArray = new CoinIndexedVector(n_cols);
1270      int n;
1271      int * which;
1272      double * array;
1273      if (OsiUnitTest::verbosity >= 2)
1274        printf("B-1 A");
1275      for( i = 0; i < n_rows; i++){
1276        m.getBInvARow(i, columnArray,rowArray);
1277        if (OsiUnitTest::verbosity >= 2)
1278          printf("\nrow: %d -> ",i);
1279        int j;
1280        // First columns
1281        n = columnArray->getNumElements();
1282        which = columnArray->getIndices();
1283        array = columnArray->denseVector();
1284        for(j=0; j < n; j++){
1285          int k=which[j];
1286          if (OsiUnitTest::verbosity >= 2)
1287            printf("(%d %g), ", k, array[k]);
1288          // zero out
1289          array[k]=0.0;
1290        }
1291        // say empty (if I had not zeroed array[k] I would use ->clear())
1292        columnArray->setNumElements(0);
1293        // and check (would not be in any production code)
1294        columnArray->checkClear();
1295        // now rows
1296        n = rowArray->getNumElements();
1297        which = rowArray->getIndices();
1298        array = rowArray->denseVector();
1299        for(j=0; j < n; j++){
1300          int k=which[j];
1301          if (OsiUnitTest::verbosity >= 2)
1302            printf("(%d %g), ", k+n_cols, array[k]);
1303          // zero out
1304          array[k]=0.0;
1305        }
1306        // say empty
1307        rowArray->setNumElements(0);
1308        // and check (would not be in any production code)
1309        rowArray->checkClear();
1310      }
1311      if (OsiUnitTest::verbosity >= 2) {
1312        printf("\n");
1313        printf("And by column");
1314      }
1315      for( i = 0; i < n_cols+n_rows; i++){
1316        m.getBInvACol(i, rowArray);
1317        if (OsiUnitTest::verbosity >= 2)
1318          printf("\ncolumn: %d -> ",i);
1319        n = rowArray->getNumElements();
1320        which = rowArray->getIndices();
1321        array = rowArray->denseVector();
1322        for(int j=0; j < n; j++){
1323          int k=which[j];
1324          if (OsiUnitTest::verbosity >= 2)
1325            printf("(%d %g), ", k, array[k]);
1326          // zero out
1327          array[k]=0.0;
1328        }
1329        // say empty
1330        rowArray->setNumElements(0);
1331        // and check (would not be in any production code)
1332        rowArray->checkClear();
1333      }
1334      if (OsiUnitTest::verbosity >= 2)
1335        printf("\n");
1336      delete rowArray;
1337      delete columnArray;
1338    }
1339    // may not be needed - but cleaner
1340    m.disableFactorization();
1341  }
1342
1343  // Check tableau stuff when simplex interface is off
1344  {   
1345    OsiClpSolverInterface m;
1346    /*
1347       Wolsey : Page 130
1348       max 4x1 -  x2
1349       7x1 - 2x2    <= 14
1350       x2    <= 3
1351       2x1 - 2x2    <= 3
1352       x1 in Z+, x2 >= 0
1353    */
1354   
1355    double inf_ = m.getInfinity();
1356    int n_cols = 2;
1357    int n_rows = 3;
1358   
1359    double obj[2] = {-4.0, 1.0};
1360    double collb[2] = {0.0, 0.0};
1361    double colub[2] = {inf_, inf_};
1362    double rowlb[3] = {-inf_, -inf_, -inf_};
1363    double rowub[3] = {14.0, 3.0, 3.0};
1364   
1365    int rowIndices[5] =  {0,     2,    0,    1,    2};
1366    int colIndices[5] =  {0,     0,    1,    1,    1};
1367    double elements[5] = {7.0, 2.0, -2.0,  1.0, -2.0};
1368    CoinPackedMatrix M(true, rowIndices, colIndices, elements, 5);
1369   
1370    m.loadProblem(M, collb, colub, obj, rowlb, rowub);
1371    // Test with scaling
1372    m.setHintParam(OsiDoScale,true);
1373    // Tell code to keep factorization
1374    m.setupForRepeatedUse(3,0);
1375   
1376    m.initialSolve();
1377    // Repeated stuff only in resolve
1378    // m.resolve(); (now in both (just for (3,0))
1379   
1380    //check that the tableau matches wolsey (B-1 A)
1381    // slacks in second part of binvA
1382    double * binvA = (double*) malloc((n_cols+n_rows) * sizeof(double));
1383    // and getBasics
1384    int * pivots = new int[n_rows];
1385    m.getBasics(pivots);
1386   
1387    if (OsiUnitTest::verbosity >= 2)
1388      printf("B-1 A");
1389    int i;
1390    for( i = 0; i < n_rows; i++){
1391      m.getBInvARow(i, binvA,binvA+n_cols);
1392      if (OsiUnitTest::verbosity >= 2) {
1393        printf("\nrow: %d (pivot %d) -> ",i,pivots[i]);
1394        for(int j=0; j < n_cols+n_rows; j++)
1395          printf("%g, ", binvA[j]);
1396      }
1397    }
1398    if (OsiUnitTest::verbosity >= 2) {
1399      printf("\n");
1400      printf("And by column");
1401    }
1402    for( i = 0; i < n_cols+n_rows; i++){
1403      m.getBInvACol(i, binvA);
1404      if (OsiUnitTest::verbosity >= 2) {
1405        printf("\ncolumn: %d -> ",i);
1406        for(int j=0; j < n_rows; j++)
1407          printf("%g, ", binvA[j]);
1408      }
1409    }
1410    if (OsiUnitTest::verbosity >= 2)
1411      printf("\n");
1412    free(binvA);
1413    delete [] pivots;
1414  }
1415
1416  // Do common solverInterface testing
1417  {
1418    OsiClpSolverInterface m;
1419    OsiSolverInterfaceCommonUnitTest(&m, mpsDir, netlibDir);
1420  }
1421}
Note: See TracBrowser for help on using the repository browser.