Ignore:
Timestamp:
Jun 5, 2010 3:42:36 PM (10 years ago)
Author:
stefan
Message:

merge split branch into trunk

Location:
trunk
Files:
3 edited

Legend:

Unmodified
Added
Removed
  • trunk

    • Property svn:externals
      •  

        old new  
        1 BuildTools    https://projects.coin-or.org/svn/BuildTools/stable/0.6
        2 ThirdParty/Blas https://projects.coin-or.org/svn/BuildTools/ThirdParty/Blas/stable/1.0
        3 ThirdParty/Lapack https://projects.coin-or.org/svn/BuildTools/ThirdParty/Lapack/stable/1.0
        4 Data/Netlib   https://projects.coin-or.org/svn/Data/stable/1.0/Netlib
        5 Data/Sample   https://projects.coin-or.org/svn/Data/stable/1.0/Sample
        6 CoinUtils     https://projects.coin-or.org/svn/CoinUtils/stable/2.6/CoinUtils
        7 
         1BuildTools    https://projects.coin-or.org/svn/BuildTools/trunk
         2ThirdParty/Blas https://projects.coin-or.org/svn/BuildTools/ThirdParty/Blas/trunk
         3ThirdParty/Lapack https://projects.coin-or.org/svn/BuildTools/ThirdParty/Lapack/trunk
         4#Data/Netlib   https://projects.coin-or.org/svn/Data/trunk/Netlib
         5Data/Sample   https://projects.coin-or.org/svn/Data/trunk/Sample
         6CoinUtils     https://projects.coin-or.org/svn/CoinUtils/trunk/CoinUtils
    • Property svn:mergeinfo changed
      /branches/split (added)merged: 1522,​1527-1531,​1542-1550,​1558
  • trunk/Clp/src

  • trunk/Clp/src/unitTest.cpp

    r1525 r1559  
    13131313          CoinMpsIO m;
    13141314          std::string fn = dirSample + "exmip1";
    1315           m.readMps(fn.c_str(), "mps");
    1316           ClpSimplex solution;
    1317           solution.loadProblem(*m.getMatrixByCol(), m.getColLower(), m.getColUpper(),
    1318                                m.getObjCoefficients(),
    1319                                m.getRowLower(), m.getRowUpper());
    1320           solution.dual();
    1321           // Test event handling
    1322           MyEventHandler handler;
    1323           solution.passInEventHandler(&handler);
    1324           int numberRows = solution.numberRows();
    1325           // make sure values pass has something to do
    1326           for (int i = 0; i < numberRows; i++)
    1327                solution.setRowStatus(i, ClpSimplex::basic);
    1328           solution.primal(1);
    1329           assert (solution.secondaryStatus() == 102); // Came out at end of pass
     1315          if (m.readMps(fn.c_str(), "mps") == 0) {
     1316               ClpSimplex solution;
     1317               solution.loadProblem(*m.getMatrixByCol(), m.getColLower(), m.getColUpper(),
     1318                                    m.getObjCoefficients(),
     1319                                    m.getRowLower(), m.getRowUpper());
     1320               solution.dual();
     1321               // Test event handling
     1322               MyEventHandler handler;
     1323               solution.passInEventHandler(&handler);
     1324               int numberRows = solution.numberRows();
     1325               // make sure values pass has something to do
     1326               for (int i = 0; i < numberRows; i++)
     1327                    solution.setRowStatus(i, ClpSimplex::basic);
     1328               solution.primal(1);
     1329               assert (solution.secondaryStatus() == 102); // Came out at end of pass
     1330          } else {
     1331               std::cerr << "Error reading exmip1 from sample data. Skipping test." << std::endl;
     1332          }
    13301333     }
    13311334     // Test Message handler
     
    13341337          std::string fn = dirSample + "exmip1";
    13351338          //fn = "Test/subGams4";
    1336           m.readMps(fn.c_str(), "mps");
    1337           ClpSimplex model;
    1338           model.loadProblem(*m.getMatrixByCol(), m.getColLower(), m.getColUpper(),
    1339                             m.getObjCoefficients(),
    1340                             m.getRowLower(), m.getRowUpper());
    1341           // Message handler
    1342           MyMessageHandler messageHandler(&model);
    1343           std::cout << "Testing derived message handler" << std::endl;
    1344           model.passInMessageHandler(&messageHandler);
    1345           model.messagesPointer()->setDetailMessage(1, 102);
    1346           model.setFactorizationFrequency(10);
    1347           model.primal();
    1348           model.primal(0, 3);
    1349           model.setObjCoeff(3, -2.9473684210526314);
    1350           model.primal(0, 3);
    1351           // Write saved solutions
    1352           int nc = model.getNumCols();
    1353           int s;
    1354           std::deque<StdVectorDouble> fep = messageHandler.getFeasibleExtremePoints();
    1355           int numSavedSolutions = fep.size();
    1356           for ( s = 0; s < numSavedSolutions; ++s ) {
    1357                const StdVectorDouble & solnVec = fep[s];
    1358                for ( int c = 0; c < nc; ++c ) {
    1359                     if (fabs(solnVec[c]) > 1.0e-8)
    1360                          std::cout << "Saved Solution: " << s << " ColNum: " << c << " Value: " << solnVec[c] << std::endl;
    1361                }
    1362           }
    1363           // Solve again without scaling
    1364           // and maximize then minimize
    1365           messageHandler.clearFeasibleExtremePoints();
    1366           model.scaling(0);
    1367           model.setOptimizationDirection(-1);
    1368           model.primal();
    1369           model.setOptimizationDirection(1);
    1370           model.primal();
    1371           fep = messageHandler.getFeasibleExtremePoints();
    1372           numSavedSolutions = fep.size();
    1373           for ( s = 0; s < numSavedSolutions; ++s ) {
    1374                const StdVectorDouble & solnVec = fep[s];
    1375                for ( int c = 0; c < nc; ++c ) {
    1376                     if (fabs(solnVec[c]) > 1.0e-8)
    1377                          std::cout << "Saved Solution: " << s << " ColNum: " << c << " Value: " << solnVec[c] << std::endl;
    1378                }
     1339          if (m.readMps(fn.c_str(), "mps") == 0) {
     1340               ClpSimplex model;
     1341               model.loadProblem(*m.getMatrixByCol(), m.getColLower(), m.getColUpper(),
     1342                                 m.getObjCoefficients(),
     1343                                 m.getRowLower(), m.getRowUpper());
     1344               // Message handler
     1345               MyMessageHandler messageHandler(&model);
     1346               std::cout << "Testing derived message handler" << std::endl;
     1347               model.passInMessageHandler(&messageHandler);
     1348               model.messagesPointer()->setDetailMessage(1, 102);
     1349               model.setFactorizationFrequency(10);
     1350               model.primal();
     1351               model.primal(0, 3);
     1352               model.setObjCoeff(3, -2.9473684210526314);
     1353               model.primal(0, 3);
     1354               // Write saved solutions
     1355               int nc = model.getNumCols();
     1356               int s;
     1357               std::deque<StdVectorDouble> fep = messageHandler.getFeasibleExtremePoints();
     1358               int numSavedSolutions = fep.size();
     1359               for ( s = 0; s < numSavedSolutions; ++s ) {
     1360                    const StdVectorDouble & solnVec = fep[s];
     1361                    for ( int c = 0; c < nc; ++c ) {
     1362                         if (fabs(solnVec[c]) > 1.0e-8)
     1363                              std::cout << "Saved Solution: " << s << " ColNum: " << c << " Value: " << solnVec[c] << std::endl;
     1364                    }
     1365               }
     1366               // Solve again without scaling
     1367               // and maximize then minimize
     1368               messageHandler.clearFeasibleExtremePoints();
     1369               model.scaling(0);
     1370               model.setOptimizationDirection(-1);
     1371               model.primal();
     1372               model.setOptimizationDirection(1);
     1373               model.primal();
     1374               fep = messageHandler.getFeasibleExtremePoints();
     1375               numSavedSolutions = fep.size();
     1376               for ( s = 0; s < numSavedSolutions; ++s ) {
     1377                    const StdVectorDouble & solnVec = fep[s];
     1378                    for ( int c = 0; c < nc; ++c ) {
     1379                         if (fabs(solnVec[c]) > 1.0e-8)
     1380                              std::cout << "Saved Solution: " << s << " ColNum: " << c << " Value: " << solnVec[c] << std::endl;
     1381                    }
     1382               }
     1383          } else {
     1384               std::cerr << "Error reading exmip1 from sample data. Skipping test." << std::endl;
    13791385          }
    13801386     }
     
    13841390          CoinMpsIO m;
    13851391          std::string fn = dirSample + "exmip1";
    1386           m.readMps(fn.c_str(), "mps");
    1387           ClpSimplex model;
    1388           model.loadProblem(*m.getMatrixByCol(), m.getColLower(), m.getColUpper(),
    1389                             m.getObjCoefficients(),
    1390                             m.getRowLower(), m.getRowUpper());
    1391           model.primal();
    1392           int which[13] = {0, 1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12};
    1393           double costIncrease[13];
    1394           int sequenceIncrease[13];
    1395           double costDecrease[13];
    1396           int sequenceDecrease[13];
    1397           // ranging
    1398           model.dualRanging(13, which, costIncrease, sequenceIncrease,
    1399                             costDecrease, sequenceDecrease);
    1400           int i;
    1401           for ( i = 0; i < 13; i++)
    1402                printf("%d increase %g %d, decrease %g %d\n",
    1403                       i, costIncrease[i], sequenceIncrease[i],
    1404                       costDecrease[i], sequenceDecrease[i]);
    1405           assert (fabs(costDecrease[3]) < 1.0e-4);
    1406           assert (fabs(costIncrease[7] - 1.0) < 1.0e-4);
    1407           model.setOptimizationDirection(-1);
    1408           {
    1409                int j;
    1410                double * obj = model.objective();
    1411                int n = model.numberColumns();
    1412                for (j = 0; j < n; j++)
    1413                     obj[j] *= -1.0;
    1414           }
    1415           double costIncrease2[13];
    1416           int sequenceIncrease2[13];
    1417           double costDecrease2[13];
    1418           int sequenceDecrease2[13];
    1419           // ranging
    1420           model.dualRanging(13, which, costIncrease2, sequenceIncrease2,
    1421                             costDecrease2, sequenceDecrease2);
    1422           for (i = 0; i < 13; i++) {
    1423                assert (fabs(costIncrease[i] - costDecrease2[i]) < 1.0e-6);
    1424                assert (fabs(costDecrease[i] - costIncrease2[i]) < 1.0e-6);
    1425                assert (sequenceIncrease[i] == sequenceDecrease2[i]);
    1426                assert (sequenceDecrease[i] == sequenceIncrease2[i]);
    1427           }
    1428           // Now delete all rows and see what happens
    1429           model.deleteRows(model.numberRows(), which);
    1430           model.primal();
    1431           // ranging
    1432           if (!model.dualRanging(8, which, costIncrease, sequenceIncrease,
    1433                                  costDecrease, sequenceDecrease)) {
    1434                for (i = 0; i < 8; i++) {
     1392          if (m.readMps(fn.c_str(), "mps") == 0) {
     1393               ClpSimplex model;
     1394               model.loadProblem(*m.getMatrixByCol(), m.getColLower(), m.getColUpper(),
     1395                                 m.getObjCoefficients(),
     1396                                 m.getRowLower(), m.getRowUpper());
     1397               model.primal();
     1398               int which[13] = {0, 1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12};
     1399               double costIncrease[13];
     1400               int sequenceIncrease[13];
     1401               double costDecrease[13];
     1402               int sequenceDecrease[13];
     1403               // ranging
     1404               model.dualRanging(13, which, costIncrease, sequenceIncrease,
     1405                                 costDecrease, sequenceDecrease);
     1406               int i;
     1407               for ( i = 0; i < 13; i++)
    14351408                    printf("%d increase %g %d, decrease %g %d\n",
    14361409                           i, costIncrease[i], sequenceIncrease[i],
    14371410                           costDecrease[i], sequenceDecrease[i]);
    1438                }
     1411               assert (fabs(costDecrease[3]) < 1.0e-4);
     1412               assert (fabs(costIncrease[7] - 1.0) < 1.0e-4);
     1413               model.setOptimizationDirection(-1);
     1414               {
     1415                    int j;
     1416                    double * obj = model.objective();
     1417                    int n = model.numberColumns();
     1418                    for (j = 0; j < n; j++)
     1419                         obj[j] *= -1.0;
     1420               }
     1421               double costIncrease2[13];
     1422               int sequenceIncrease2[13];
     1423               double costDecrease2[13];
     1424               int sequenceDecrease2[13];
     1425               // ranging
     1426               model.dualRanging(13, which, costIncrease2, sequenceIncrease2,
     1427                                 costDecrease2, sequenceDecrease2);
     1428               for (i = 0; i < 13; i++) {
     1429                    assert (fabs(costIncrease[i] - costDecrease2[i]) < 1.0e-6);
     1430                    assert (fabs(costDecrease[i] - costIncrease2[i]) < 1.0e-6);
     1431                    assert (sequenceIncrease[i] == sequenceDecrease2[i]);
     1432                    assert (sequenceDecrease[i] == sequenceIncrease2[i]);
     1433               }
     1434               // Now delete all rows and see what happens
     1435               model.deleteRows(model.numberRows(), which);
     1436               model.primal();
     1437               // ranging
     1438               if (!model.dualRanging(8, which, costIncrease, sequenceIncrease,
     1439                                      costDecrease, sequenceDecrease)) {
     1440                    for (i = 0; i < 8; i++) {
     1441                         printf("%d increase %g %d, decrease %g %d\n",
     1442                                i, costIncrease[i], sequenceIncrease[i],
     1443                                costDecrease[i], sequenceDecrease[i]);
     1444                    }
     1445               }
     1446          } else {
     1447               std::cerr << "Error reading exmip1 from sample data. Skipping test." << std::endl;
    14391448          }
    14401449     }
     
    14431452          CoinMpsIO m;
    14441453          std::string fn = dirSample + "exmip1";
    1445           m.readMps(fn.c_str(), "mps");
    1446           ClpSimplex model;
    1447           model.loadProblem(*m.getMatrixByCol(), m.getColLower(), m.getColUpper(),
    1448                             m.getObjCoefficients(),
    1449                             m.getRowLower(), m.getRowUpper());
    1450           model.primal();
    1451           int which[13] = {0, 1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12};
    1452           double valueIncrease[13];
    1453           int sequenceIncrease[13];
    1454           double valueDecrease[13];
    1455           int sequenceDecrease[13];
    1456           // ranging
    1457           model.primalRanging(13, which, valueIncrease, sequenceIncrease,
    1458                               valueDecrease, sequenceDecrease);
    1459           int i;
    1460           for ( i = 0; i < 13; i++)
    1461                printf("%d increase %g %d, decrease %g %d\n",
    1462                       i, valueIncrease[i], sequenceIncrease[i],
    1463                       valueDecrease[i], sequenceDecrease[i]);
    1464           assert (fabs(valueIncrease[3] - 0.642857) < 1.0e-4);
    1465           assert (fabs(valueIncrease[8] - 2.95113) < 1.0e-4);
     1454          if (m.readMps(fn.c_str(), "mps") == 0) {
     1455               ClpSimplex model;
     1456               model.loadProblem(*m.getMatrixByCol(), m.getColLower(), m.getColUpper(),
     1457                                 m.getObjCoefficients(),
     1458                                 m.getRowLower(), m.getRowUpper());
     1459               model.primal();
     1460               int which[13] = {0, 1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12};
     1461               double valueIncrease[13];
     1462               int sequenceIncrease[13];
     1463               double valueDecrease[13];
     1464               int sequenceDecrease[13];
     1465               // ranging
     1466               model.primalRanging(13, which, valueIncrease, sequenceIncrease,
     1467                                   valueDecrease, sequenceDecrease);
     1468               int i;
     1469               for ( i = 0; i < 13; i++)
     1470                    printf("%d increase %g %d, decrease %g %d\n",
     1471                           i, valueIncrease[i], sequenceIncrease[i],
     1472                           valueDecrease[i], sequenceDecrease[i]);
     1473               assert (fabs(valueIncrease[3] - 0.642857) < 1.0e-4);
     1474               assert (fabs(valueIncrease[8] - 2.95113) < 1.0e-4);
     1475          } else {
     1476               std::cerr << "Error reading exmip1 from sample data. Skipping test." << std::endl;
     1477          }
    14661478#if 0
    14671479          // out until I find optimization bug
     
    16401652               // probable cause is that gz not there
    16411653               fprintf(stderr, "Unable to open finnis.mps in %s\n", dirSample.c_str());
    1642                fprintf(stderr, "Most probable cause is finnis.mps is gzipped i.e. finnis.mps.gz and libz has not been activated\n");
     1654               fprintf(stderr, "Most probable cause is that sample data is not available, or finnis.mps is gzipped i.e. finnis.mps.gz and libz has not been activated\n");
    16431655               fprintf(stderr, "Either gunzip files or edit Makefiles/Makefile.location to get libz\n");
    1644                exit(999);
    1645           }
    1646           ClpModel model;
    1647           model.loadProblem(*m.getMatrixByCol(), m.getColLower(),
    1648                             m.getColUpper(),
    1649                             m.getObjCoefficients(),
    1650                             m.getRowLower(), m.getRowUpper());
    1651           ClpSimplex solution(model);
    1652 
    1653           solution.scaling(1);
    1654           solution.setDualBound(1.0e8);
    1655           //solution.factorization()->maximumPivots(1);
    1656           //solution.setLogLevel(3);
    1657           solution.setDualTolerance(1.0e-7);
    1658           // set objective sense,
    1659           ClpDualRowSteepest steep;
    1660           solution.setDualRowPivotAlgorithm(steep);
    1661           solution.setDblParam(ClpObjOffset, m.objectiveOffset());
    1662           solution.dual();
     1656          } else {
     1657               ClpModel model;
     1658               model.loadProblem(*m.getMatrixByCol(), m.getColLower(),
     1659                                 m.getColUpper(),
     1660                                 m.getObjCoefficients(),
     1661                                 m.getRowLower(), m.getRowUpper());
     1662               ClpSimplex solution(model);
     1663
     1664               solution.scaling(1);
     1665               solution.setDualBound(1.0e8);
     1666               //solution.factorization()->maximumPivots(1);
     1667               //solution.setLogLevel(3);
     1668               solution.setDualTolerance(1.0e-7);
     1669               // set objective sense,
     1670               ClpDualRowSteepest steep;
     1671               solution.setDualRowPivotAlgorithm(steep);
     1672               solution.setDblParam(ClpObjOffset, m.objectiveOffset());
     1673               solution.dual();
     1674          }
    16631675     }
    16641676     // test normal solution
     
    16661678          CoinMpsIO m;
    16671679          std::string fn = dirSample + "afiro";
    1668           m.readMps(fn.c_str(), "mps");
    1669           ClpSimplex solution;
    1670           ClpModel model;
    1671           // do twice - without and with scaling
    1672           int iPass;
    1673           for (iPass = 0; iPass < 2; iPass++) {
    1674                // explicit row objective for testing
    1675                int nr = m.getNumRows();
    1676                double * rowObj = new double[nr];
    1677                CoinFillN(rowObj, nr, 0.0);
    1678                model.loadProblem(*m.getMatrixByCol(), m.getColLower(), m.getColUpper(),
    1679                                  m.getObjCoefficients(),
    1680                                  m.getRowLower(), m.getRowUpper(), rowObj);
    1681                delete [] rowObj;
    1682                solution = ClpSimplex(model);
    1683                if (iPass) {
    1684                     solution.scaling();
    1685                }
    1686                solution.dual();
    1687                solution.dual();
    1688                // test optimal
    1689                assert (solution.status() == 0);
    1690                int numberColumns = solution.numberColumns();
    1691                int numberRows = solution.numberRows();
    1692                CoinPackedVector colsol(numberColumns, solution.primalColumnSolution());
    1693                double * objective = solution.objective();
     1680          if (m.readMps(fn.c_str(), "mps") == 0) {
     1681               ClpSimplex solution;
     1682               ClpModel model;
     1683               // do twice - without and with scaling
     1684               int iPass;
     1685               for (iPass = 0; iPass < 2; iPass++) {
     1686                    // explicit row objective for testing
     1687                    int nr = m.getNumRows();
     1688                    double * rowObj = new double[nr];
     1689                    CoinFillN(rowObj, nr, 0.0);
     1690                    model.loadProblem(*m.getMatrixByCol(), m.getColLower(), m.getColUpper(),
     1691                                      m.getObjCoefficients(),
     1692                                      m.getRowLower(), m.getRowUpper(), rowObj);
     1693                    delete [] rowObj;
     1694                    solution = ClpSimplex(model);
     1695                    if (iPass) {
     1696                         solution.scaling();
     1697                    }
     1698                    solution.dual();
     1699                    solution.dual();
     1700                    // test optimal
     1701                    assert (solution.status() == 0);
     1702                    int numberColumns = solution.numberColumns();
     1703                    int numberRows = solution.numberRows();
     1704                    CoinPackedVector colsol(numberColumns, solution.primalColumnSolution());
     1705                    double * objective = solution.objective();
    16941706#ifndef NDEBUG
    1695                double objValue = colsol.dotProduct(objective);
     1707                    double objValue = colsol.dotProduct(objective);
    16961708#endif
    1697                CoinRelFltEq eq(1.0e-8);
    1698                assert(eq(objValue, -4.6475314286e+02));
    1699                solution.dual();
    1700                assert(eq(solution.objectiveValue(), -4.6475314286e+02));
    1701                double * lower = solution.columnLower();
    1702                double * upper = solution.columnUpper();
    1703                double * sol = solution.primalColumnSolution();
    1704                double * result = new double[numberColumns];
    1705                CoinFillN ( result, numberColumns, 0.0);
    1706                solution.matrix()->transposeTimes(solution.dualRowSolution(), result);
    1707                int iRow , iColumn;
    1708                // see if feasible and dual feasible
    1709                for (iColumn = 0; iColumn < numberColumns; iColumn++) {
    1710                     double value = sol[iColumn];
    1711                     assert(value < upper[iColumn] + 1.0e-8);
    1712                     assert(value > lower[iColumn] - 1.0e-8);
    1713                     value = objective[iColumn] - result[iColumn];
    1714                     assert (value > -1.0e-5);
    1715                     if (sol[iColumn] > 1.0e-5)
    1716                          assert (value < 1.0e-5);
    1717                }
    1718                delete [] result;
    1719                result = new double[numberRows];
    1720                CoinFillN ( result, numberRows, 0.0);
    1721                solution.matrix()->times(colsol, result);
    1722                lower = solution.rowLower();
    1723                upper = solution.rowUpper();
    1724                sol = solution.primalRowSolution();
     1709                    CoinRelFltEq eq(1.0e-8);
     1710                    assert(eq(objValue, -4.6475314286e+02));
     1711                    solution.dual();
     1712                    assert(eq(solution.objectiveValue(), -4.6475314286e+02));
     1713                    double * lower = solution.columnLower();
     1714                    double * upper = solution.columnUpper();
     1715                    double * sol = solution.primalColumnSolution();
     1716                    double * result = new double[numberColumns];
     1717                    CoinFillN ( result, numberColumns, 0.0);
     1718                    solution.matrix()->transposeTimes(solution.dualRowSolution(), result);
     1719                    int iRow , iColumn;
     1720                    // see if feasible and dual feasible
     1721                    for (iColumn = 0; iColumn < numberColumns; iColumn++) {
     1722                         double value = sol[iColumn];
     1723                         assert(value < upper[iColumn] + 1.0e-8);
     1724                         assert(value > lower[iColumn] - 1.0e-8);
     1725                         value = objective[iColumn] - result[iColumn];
     1726                         assert (value > -1.0e-5);
     1727                         if (sol[iColumn] > 1.0e-5)
     1728                              assert (value < 1.0e-5);
     1729                    }
     1730                    delete [] result;
     1731                    result = new double[numberRows];
     1732                    CoinFillN ( result, numberRows, 0.0);
     1733                    solution.matrix()->times(colsol, result);
     1734                    lower = solution.rowLower();
     1735                    upper = solution.rowUpper();
     1736                    sol = solution.primalRowSolution();
    17251737#ifndef NDEBUG
    1726                for (iRow = 0; iRow < numberRows; iRow++) {
    1727                     double value = result[iRow];
    1728                     assert(eq(value, sol[iRow]));
    1729                     assert(value < upper[iRow] + 1.0e-8);
    1730                     assert(value > lower[iRow] - 1.0e-8);
    1731                }
     1738                    for (iRow = 0; iRow < numberRows; iRow++) {
     1739                         double value = result[iRow];
     1740                         assert(eq(value, sol[iRow]));
     1741                         assert(value < upper[iRow] + 1.0e-8);
     1742                         assert(value > lower[iRow] - 1.0e-8);
     1743                    }
    17321744#endif
    1733                delete [] result;
    1734                // test row objective
    1735                double * rowObjective = solution.rowObjective();
    1736                CoinDisjointCopyN(solution.dualRowSolution(), numberRows, rowObjective);
    1737                CoinDisjointCopyN(solution.dualColumnSolution(), numberColumns, objective);
    1738                // this sets up all slack basis
    1739                solution.createStatus();
    1740                solution.dual();
    1741                CoinFillN(rowObjective, numberRows, 0.0);
    1742                CoinDisjointCopyN(m.getObjCoefficients(), numberColumns, objective);
    1743                solution.dual();
     1745                    delete [] result;
     1746                    // test row objective
     1747                    double * rowObjective = solution.rowObjective();
     1748                    CoinDisjointCopyN(solution.dualRowSolution(), numberRows, rowObjective);
     1749                    CoinDisjointCopyN(solution.dualColumnSolution(), numberColumns, objective);
     1750                    // this sets up all slack basis
     1751                    solution.createStatus();
     1752                    solution.dual();
     1753                    CoinFillN(rowObjective, numberRows, 0.0);
     1754                    CoinDisjointCopyN(m.getObjCoefficients(), numberColumns, objective);
     1755                    solution.dual();
     1756               }
     1757          } else {
     1758               std::cerr << "Error reading afiro from sample data. Skipping test." << std::endl;
    17441759          }
    17451760     }
     
    17481763          CoinMpsIO m;
    17491764          std::string fn = dirSample + "brandy";
    1750           m.readMps(fn.c_str(), "mps");
    1751           ClpSimplex solution;
    1752           // do twice - without and with scaling
    1753           int iPass;
    1754           for (iPass = 0; iPass < 2; iPass++) {
    1755                solution.loadProblem(*m.getMatrixByCol(), m.getColLower(), m.getColUpper(),
    1756                                     m.getObjCoefficients(),
    1757                                     m.getRowLower(), m.getRowUpper());
    1758                if (iPass)
    1759                     solution.scaling();
    1760                solution.setOptimizationDirection(-1);
    1761                // test unbounded and ray
     1765          if (m.readMps(fn.c_str(), "mps") == 0) {
     1766               ClpSimplex solution;
     1767               // do twice - without and with scaling
     1768               int iPass;
     1769               for (iPass = 0; iPass < 2; iPass++) {
     1770                    solution.loadProblem(*m.getMatrixByCol(), m.getColLower(), m.getColUpper(),
     1771                                         m.getObjCoefficients(),
     1772                                         m.getRowLower(), m.getRowUpper());
     1773                    if (iPass)
     1774                         solution.scaling();
     1775                    solution.setOptimizationDirection(-1);
     1776                    // test unbounded and ray
    17621777#ifdef DUAL
    1763                solution.setDualBound(100.0);
    1764                solution.dual();
     1778                    solution.setDualBound(100.0);
     1779                    solution.dual();
    17651780#else
    1766                solution.primal();
     1781                    solution.primal();
    17671782#endif
    1768                assert (solution.status() == 2);
    1769                int numberColumns = solution.numberColumns();
    1770                int numberRows = solution.numberRows();
    1771                double * lower = solution.columnLower();
    1772                double * upper = solution.columnUpper();
    1773                double * sol = solution.primalColumnSolution();
    1774                double * ray = solution.unboundedRay();
    1775                double * objective = solution.objective();
    1776                double objChange = 0.0;
    1777                int iRow , iColumn;
    1778                // make sure feasible and columns form ray
    1779                for (iColumn = 0; iColumn < numberColumns; iColumn++) {
    1780                     double value = sol[iColumn];
    1781                     assert(value < upper[iColumn] + 1.0e-8);
    1782                     assert(value > lower[iColumn] - 1.0e-8);
    1783                     value = ray[iColumn];
    1784                     if (value > 0.0)
    1785                          assert(upper[iColumn] > 1.0e30);
    1786                     else if (value < 0.0)
    1787                          assert(lower[iColumn] < -1.0e30);
    1788                     objChange += value * objective[iColumn];
    1789                }
    1790                // make sure increasing objective
    1791                assert(objChange > 0.0);
    1792                double * result = new double[numberRows];
    1793                CoinFillN ( result, numberRows, 0.0);
    1794                solution.matrix()->times(sol, result);
    1795                lower = solution.rowLower();
    1796                upper = solution.rowUpper();
    1797                sol = solution.primalRowSolution();
     1783                    assert (solution.status() == 2);
     1784                    int numberColumns = solution.numberColumns();
     1785                    int numberRows = solution.numberRows();
     1786                    double * lower = solution.columnLower();
     1787                    double * upper = solution.columnUpper();
     1788                    double * sol = solution.primalColumnSolution();
     1789                    double * ray = solution.unboundedRay();
     1790                    double * objective = solution.objective();
     1791                    double objChange = 0.0;
     1792                    int iRow , iColumn;
     1793                    // make sure feasible and columns form ray
     1794                    for (iColumn = 0; iColumn < numberColumns; iColumn++) {
     1795                         double value = sol[iColumn];
     1796                         assert(value < upper[iColumn] + 1.0e-8);
     1797                         assert(value > lower[iColumn] - 1.0e-8);
     1798                         value = ray[iColumn];
     1799                         if (value > 0.0)
     1800                              assert(upper[iColumn] > 1.0e30);
     1801                         else if (value < 0.0)
     1802                              assert(lower[iColumn] < -1.0e30);
     1803                         objChange += value * objective[iColumn];
     1804                    }
     1805                    // make sure increasing objective
     1806                    assert(objChange > 0.0);
     1807                    double * result = new double[numberRows];
     1808                    CoinFillN ( result, numberRows, 0.0);
     1809                    solution.matrix()->times(sol, result);
     1810                    lower = solution.rowLower();
     1811                    upper = solution.rowUpper();
     1812                    sol = solution.primalRowSolution();
    17981813#ifndef NDEBUG
    1799                for (iRow = 0; iRow < numberRows; iRow++) {
    1800                     double value = result[iRow];
    1801                     assert(eq(value, sol[iRow]));
    1802                     assert(value < upper[iRow] + 2.0e-8);
    1803                     assert(value > lower[iRow] - 2.0e-8);
    1804                }
     1814                    for (iRow = 0; iRow < numberRows; iRow++) {
     1815                         double value = result[iRow];
     1816                         assert(eq(value, sol[iRow]));
     1817                         assert(value < upper[iRow] + 2.0e-8);
     1818                         assert(value > lower[iRow] - 2.0e-8);
     1819                    }
    18051820#endif
    1806                CoinFillN ( result, numberRows, 0.0);
    1807                solution.matrix()->times(ray, result);
    1808                // there may be small differences (especially if scaled)
    1809                for (iRow = 0; iRow < numberRows; iRow++) {
    1810                     double value = result[iRow];
    1811                     if (value > 1.0e-8)
    1812                          assert(upper[iRow] > 1.0e30);
    1813                     else if (value < -1.0e-8)
    1814                          assert(lower[iRow] < -1.0e30);
    1815                }
    1816                delete [] result;
    1817                delete [] ray;
     1821                    CoinFillN ( result, numberRows, 0.0);
     1822                    solution.matrix()->times(ray, result);
     1823                    // there may be small differences (especially if scaled)
     1824                    for (iRow = 0; iRow < numberRows; iRow++) {
     1825                         double value = result[iRow];
     1826                         if (value > 1.0e-8)
     1827                              assert(upper[iRow] > 1.0e30);
     1828                         else if (value < -1.0e-8)
     1829                              assert(lower[iRow] < -1.0e30);
     1830                    }
     1831                    delete [] result;
     1832                    delete [] ray;
     1833               }
     1834          } else {
     1835               std::cerr << "Error reading brandy from sample data. Skipping test." << std::endl;
    18181836          }
    18191837     }
     
    18221840          CoinMpsIO m;
    18231841          std::string fn = dirSample + "brandy";
    1824           m.readMps(fn.c_str(), "mps");
    1825           ClpSimplex solution;
    1826           // do twice - without and with scaling
    1827           int iPass;
    1828           for (iPass = 0; iPass < 2; iPass++) {
    1829                solution.loadProblem(*m.getMatrixByCol(), m.getColLower(), m.getColUpper(),
    1830                                     m.getObjCoefficients(),
    1831                                     m.getRowLower(), m.getRowUpper());
    1832                if (iPass)
    1833                     solution.scaling();
    1834                // test infeasible and ray
    1835                solution.columnUpper()[0] = 0.0;
     1842          if (m.readMps(fn.c_str(), "mps") == 0) {
     1843               ClpSimplex solution;
     1844               // do twice - without and with scaling
     1845               int iPass;
     1846               for (iPass = 0; iPass < 2; iPass++) {
     1847                    solution.loadProblem(*m.getMatrixByCol(), m.getColLower(), m.getColUpper(),
     1848                                         m.getObjCoefficients(),
     1849                                         m.getRowLower(), m.getRowUpper());
     1850                    if (iPass)
     1851                         solution.scaling();
     1852                    // test infeasible and ray
     1853                    solution.columnUpper()[0] = 0.0;
    18361854#ifdef DUAL
    1837                solution.setDualBound(100.0);
    1838                solution.dual();
     1855                    solution.setDualBound(100.0);
     1856                    solution.dual();
    18391857#else
    1840                solution.primal();
     1858                    solution.primal();
    18411859#endif
    1842                assert (solution.status() == 1);
    1843                int numberColumns = solution.numberColumns();
    1844                int numberRows = solution.numberRows();
    1845                double * lower = solution.rowLower();
    1846                double * upper = solution.rowUpper();
    1847                double * ray = solution.infeasibilityRay();
    1848                assert(ray);
    1849                // construct proof of infeasibility
    1850                int iRow , iColumn;
    1851                double lo = 0.0, up = 0.0;
    1852                int nl = 0, nu = 0;
    1853                for (iRow = 0; iRow < numberRows; iRow++) {
    1854                     if (lower[iRow] > -1.0e20) {
    1855                          lo += ray[iRow] * lower[iRow];
    1856                     } else {
    1857                          if (ray[iRow] > 1.0e-8)
    1858                               nl++;
    1859                     }
    1860                     if (upper[iRow] < 1.0e20) {
    1861                          up += ray[iRow] * upper[iRow];
    1862                     } else {
    1863                          if (ray[iRow] > 1.0e-8)
    1864                               nu++;
    1865                     }
    1866                }
    1867                if (nl)
    1868                     lo = -1.0e100;
    1869                if (nu)
    1870                     up = 1.0e100;
    1871                double * result = new double[numberColumns];
    1872                double lo2 = 0.0, up2 = 0.0;
    1873                CoinFillN ( result, numberColumns, 0.0);
    1874                solution.matrix()->transposeTimes(ray, result);
    1875                lower = solution.columnLower();
    1876                upper = solution.columnUpper();
    1877                nl = nu = 0;
    1878                for (iColumn = 0; iColumn < numberColumns; iColumn++) {
    1879                     if (result[iColumn] > 1.0e-8) {
    1880                          if (lower[iColumn] > -1.0e20)
    1881                               lo2 += result[iColumn] * lower[iColumn];
    1882                          else
    1883                               nl++;
    1884                          if (upper[iColumn] < 1.0e20)
    1885                               up2 += result[iColumn] * upper[iColumn];
    1886                          else
    1887                               nu++;
    1888                     } else if (result[iColumn] < -1.0e-8) {
    1889                          if (lower[iColumn] > -1.0e20)
    1890                               up2 += result[iColumn] * lower[iColumn];
    1891                          else
    1892                               nu++;
    1893                          if (upper[iColumn] < 1.0e20)
    1894                               lo2 += result[iColumn] * upper[iColumn];
    1895                          else
    1896                               nl++;
    1897                     }
    1898                }
    1899                if (nl)
    1900                     lo2 = -1.0e100;
    1901                if (nu)
    1902                     up2 = 1.0e100;
    1903                // make sure inconsistency
    1904                assert(lo2 > up || up2 < lo);
    1905                delete [] result;
    1906                delete [] ray;
     1860                    assert (solution.status() == 1);
     1861                    int numberColumns = solution.numberColumns();
     1862                    int numberRows = solution.numberRows();
     1863                    double * lower = solution.rowLower();
     1864                    double * upper = solution.rowUpper();
     1865                    double * ray = solution.infeasibilityRay();
     1866                    assert(ray);
     1867                    // construct proof of infeasibility
     1868                    int iRow , iColumn;
     1869                    double lo = 0.0, up = 0.0;
     1870                    int nl = 0, nu = 0;
     1871                    for (iRow = 0; iRow < numberRows; iRow++) {
     1872                         if (lower[iRow] > -1.0e20) {
     1873                              lo += ray[iRow] * lower[iRow];
     1874                         } else {
     1875                              if (ray[iRow] > 1.0e-8)
     1876                                   nl++;
     1877                         }
     1878                         if (upper[iRow] < 1.0e20) {
     1879                              up += ray[iRow] * upper[iRow];
     1880                         } else {
     1881                              if (ray[iRow] > 1.0e-8)
     1882                                   nu++;
     1883                         }
     1884                    }
     1885                    if (nl)
     1886                         lo = -1.0e100;
     1887                    if (nu)
     1888                         up = 1.0e100;
     1889                    double * result = new double[numberColumns];
     1890                    double lo2 = 0.0, up2 = 0.0;
     1891                    CoinFillN ( result, numberColumns, 0.0);
     1892                    solution.matrix()->transposeTimes(ray, result);
     1893                    lower = solution.columnLower();
     1894                    upper = solution.columnUpper();
     1895                    nl = nu = 0;
     1896                    for (iColumn = 0; iColumn < numberColumns; iColumn++) {
     1897                         if (result[iColumn] > 1.0e-8) {
     1898                              if (lower[iColumn] > -1.0e20)
     1899                                   lo2 += result[iColumn] * lower[iColumn];
     1900                              else
     1901                                   nl++;
     1902                              if (upper[iColumn] < 1.0e20)
     1903                                   up2 += result[iColumn] * upper[iColumn];
     1904                              else
     1905                                   nu++;
     1906                         } else if (result[iColumn] < -1.0e-8) {
     1907                              if (lower[iColumn] > -1.0e20)
     1908                                   up2 += result[iColumn] * lower[iColumn];
     1909                              else
     1910                                   nu++;
     1911                              if (upper[iColumn] < 1.0e20)
     1912                                   lo2 += result[iColumn] * upper[iColumn];
     1913                              else
     1914                                   nl++;
     1915                         }
     1916                    }
     1917                    if (nl)
     1918                         lo2 = -1.0e100;
     1919                    if (nu)
     1920                         up2 = 1.0e100;
     1921                    // make sure inconsistency
     1922                    assert(lo2 > up || up2 < lo);
     1923                    delete [] result;
     1924                    delete [] ray;
     1925               }
     1926          } else {
     1927               std::cerr << "Error reading brandy from sample data. Skipping test." << std::endl;
    19071928          }
    19081929     }
     
    19111932          CoinMpsIO m;
    19121933          std::string fn = dirSample + "brandy";
    1913           m.readMps(fn.c_str(), "mps");
    1914           ClpSimplex solution;
    1915           solution.loadProblem(*m.getMatrixByCol(), m.getColLower(), m.getColUpper(),
    1916                                m.getObjCoefficients(),
    1917                                m.getRowLower(), m.getRowUpper());
    1918           solution.dual();
    1919           CoinRelFltEq eq(1.0e-8);
    1920           assert(eq(solution.objectiveValue(), 1.5185098965e+03));
    1921 
    1922           int numberColumns = solution.numberColumns();
    1923           int numberRows = solution.numberRows();
    1924           double * saveObj = new double [numberColumns];
    1925           double * saveLower = new double[numberRows+numberColumns];
    1926           double * saveUpper = new double[numberRows+numberColumns];
    1927           int * which = new int [numberRows+numberColumns];
    1928 
    1929           int numberElements = m.getMatrixByCol()->getNumElements();
    1930           int * starts = new int[numberRows+numberColumns];
    1931           int * index = new int[numberElements];
    1932           double * element = new double[numberElements];
    1933 
    1934           const CoinBigIndex * startM;
    1935           const int * lengthM;
    1936           const int * indexM;
    1937           const double * elementM;
    1938 
    1939           int n, nel;
    1940 
    1941           // delete non basic columns
    1942           n = 0;
    1943           nel = 0;
    1944           int iRow , iColumn;
    1945           const double * lower = m.getColLower();
    1946           const double * upper = m.getColUpper();
    1947           const double * objective = m.getObjCoefficients();
    1948           startM = m.getMatrixByCol()->getVectorStarts();
    1949           lengthM = m.getMatrixByCol()->getVectorLengths();
    1950           indexM = m.getMatrixByCol()->getIndices();
    1951           elementM = m.getMatrixByCol()->getElements();
    1952           starts[0] = 0;
    1953           for (iColumn = 0; iColumn < numberColumns; iColumn++) {
    1954                if (solution.getColumnStatus(iColumn) != ClpSimplex::basic) {
     1934          if (m.readMps(fn.c_str(), "mps") == 0) {
     1935               ClpSimplex solution;
     1936               solution.loadProblem(*m.getMatrixByCol(), m.getColLower(), m.getColUpper(),
     1937                                    m.getObjCoefficients(),
     1938                                    m.getRowLower(), m.getRowUpper());
     1939               solution.dual();
     1940               CoinRelFltEq eq(1.0e-8);
     1941               assert(eq(solution.objectiveValue(), 1.5185098965e+03));
     1942
     1943               int numberColumns = solution.numberColumns();
     1944               int numberRows = solution.numberRows();
     1945               double * saveObj = new double [numberColumns];
     1946               double * saveLower = new double[numberRows+numberColumns];
     1947               double * saveUpper = new double[numberRows+numberColumns];
     1948               int * which = new int [numberRows+numberColumns];
     1949
     1950               int numberElements = m.getMatrixByCol()->getNumElements();
     1951               int * starts = new int[numberRows+numberColumns];
     1952               int * index = new int[numberElements];
     1953               double * element = new double[numberElements];
     1954
     1955               const CoinBigIndex * startM;
     1956               const int * lengthM;
     1957               const int * indexM;
     1958               const double * elementM;
     1959
     1960               int n, nel;
     1961
     1962               // delete non basic columns
     1963               n = 0;
     1964               nel = 0;
     1965               int iRow , iColumn;
     1966               const double * lower = m.getColLower();
     1967               const double * upper = m.getColUpper();
     1968               const double * objective = m.getObjCoefficients();
     1969               startM = m.getMatrixByCol()->getVectorStarts();
     1970               lengthM = m.getMatrixByCol()->getVectorLengths();
     1971               indexM = m.getMatrixByCol()->getIndices();
     1972               elementM = m.getMatrixByCol()->getElements();
     1973               starts[0] = 0;
     1974               for (iColumn = 0; iColumn < numberColumns; iColumn++) {
     1975                    if (solution.getColumnStatus(iColumn) != ClpSimplex::basic) {
     1976                         saveObj[n] = objective[iColumn];
     1977                         saveLower[n] = lower[iColumn];
     1978                         saveUpper[n] = upper[iColumn];
     1979                         int j;
     1980                         for (j = startM[iColumn]; j < startM[iColumn] + lengthM[iColumn]; j++) {
     1981                              index[nel] = indexM[j];
     1982                              element[nel++] = elementM[j];
     1983                         }
     1984                         which[n++] = iColumn;
     1985                         starts[n] = nel;
     1986                    }
     1987               }
     1988               solution.deleteColumns(n, which);
     1989               solution.dual();
     1990               // Put back
     1991               solution.addColumns(n, saveLower, saveUpper, saveObj,
     1992                                   starts, index, element);
     1993               solution.dual();
     1994               assert(eq(solution.objectiveValue(), 1.5185098965e+03));
     1995               // Delete all columns and add back
     1996               n = 0;
     1997               nel = 0;
     1998               starts[0] = 0;
     1999               lower = m.getColLower();
     2000               upper = m.getColUpper();
     2001               objective = m.getObjCoefficients();
     2002               for (iColumn = 0; iColumn < numberColumns; iColumn++) {
    19552003                    saveObj[n] = objective[iColumn];
    19562004                    saveLower[n] = lower[iColumn];
     
    19642012                    starts[n] = nel;
    19652013               }
    1966           }
    1967           solution.deleteColumns(n, which);
    1968           solution.dual();
    1969           // Put back
    1970           solution.addColumns(n, saveLower, saveUpper, saveObj,
    1971                               starts, index, element);
    1972           solution.dual();
    1973           assert(eq(solution.objectiveValue(), 1.5185098965e+03));
    1974           // Delete all columns and add back
    1975           n = 0;
    1976           nel = 0;
    1977           starts[0] = 0;
    1978           lower = m.getColLower();
    1979           upper = m.getColUpper();
    1980           objective = m.getObjCoefficients();
    1981           for (iColumn = 0; iColumn < numberColumns; iColumn++) {
    1982                saveObj[n] = objective[iColumn];
    1983                saveLower[n] = lower[iColumn];
    1984                saveUpper[n] = upper[iColumn];
    1985                int j;
    1986                for (j = startM[iColumn]; j < startM[iColumn] + lengthM[iColumn]; j++) {
    1987                     index[nel] = indexM[j];
    1988                     element[nel++] = elementM[j];
    1989                }
    1990                which[n++] = iColumn;
    1991                starts[n] = nel;
    1992           }
    1993           solution.deleteColumns(n, which);
    1994           solution.dual();
    1995           // Put back
    1996           solution.addColumns(n, saveLower, saveUpper, saveObj,
    1997                               starts, index, element);
    1998           solution.dual();
    1999           assert(eq(solution.objectiveValue(), 1.5185098965e+03));
    2000 
    2001           // reload with original
    2002           solution.loadProblem(*m.getMatrixByCol(), m.getColLower(), m.getColUpper(),
    2003                                m.getObjCoefficients(),
    2004                                m.getRowLower(), m.getRowUpper());
    2005           // delete half rows
    2006           n = 0;
    2007           nel = 0;
    2008           lower = m.getRowLower();
    2009           upper = m.getRowUpper();
    2010           startM = m.getMatrixByRow()->getVectorStarts();
    2011           lengthM = m.getMatrixByRow()->getVectorLengths();
    2012           indexM = m.getMatrixByRow()->getIndices();
    2013           elementM = m.getMatrixByRow()->getElements();
    2014           starts[0] = 0;
    2015           for (iRow = 0; iRow < numberRows; iRow++) {
    2016                if ((iRow & 1) == 0) {
     2014               solution.deleteColumns(n, which);
     2015               solution.dual();
     2016               // Put back
     2017               solution.addColumns(n, saveLower, saveUpper, saveObj,
     2018                                   starts, index, element);
     2019               solution.dual();
     2020               assert(eq(solution.objectiveValue(), 1.5185098965e+03));
     2021
     2022               // reload with original
     2023               solution.loadProblem(*m.getMatrixByCol(), m.getColLower(), m.getColUpper(),
     2024                                    m.getObjCoefficients(),
     2025                                    m.getRowLower(), m.getRowUpper());
     2026               // delete half rows
     2027               n = 0;
     2028               nel = 0;
     2029               lower = m.getRowLower();
     2030               upper = m.getRowUpper();
     2031               startM = m.getMatrixByRow()->getVectorStarts();
     2032               lengthM = m.getMatrixByRow()->getVectorLengths();
     2033               indexM = m.getMatrixByRow()->getIndices();
     2034               elementM = m.getMatrixByRow()->getElements();
     2035               starts[0] = 0;
     2036               for (iRow = 0; iRow < numberRows; iRow++) {
     2037                    if ((iRow & 1) == 0) {
     2038                         saveLower[n] = lower[iRow];
     2039                         saveUpper[n] = upper[iRow];
     2040                         int j;
     2041                         for (j = startM[iRow]; j < startM[iRow] + lengthM[iRow]; j++) {
     2042                              index[nel] = indexM[j];
     2043                              element[nel++] = elementM[j];
     2044                         }
     2045                         which[n++] = iRow;
     2046                         starts[n] = nel;
     2047                    }
     2048               }
     2049               solution.deleteRows(n, which);
     2050               solution.dual();
     2051               // Put back
     2052               solution.addRows(n, saveLower, saveUpper,
     2053                                starts, index, element);
     2054               solution.dual();
     2055               assert(eq(solution.objectiveValue(), 1.5185098965e+03));
     2056               solution.writeMps("yy.mps");
     2057               // Delete all rows
     2058               n = 0;
     2059               nel = 0;
     2060               lower = m.getRowLower();
     2061               upper = m.getRowUpper();
     2062               starts[0] = 0;
     2063               for (iRow = 0; iRow < numberRows; iRow++) {
    20172064                    saveLower[n] = lower[iRow];
    20182065                    saveUpper[n] = upper[iRow];
     
    20252072                    starts[n] = nel;
    20262073               }
    2027           }
    2028           solution.deleteRows(n, which);
    2029           solution.dual();
    2030           // Put back
    2031           solution.addRows(n, saveLower, saveUpper,
    2032                            starts, index, element);
    2033           solution.dual();
    2034           assert(eq(solution.objectiveValue(), 1.5185098965e+03));
    2035           solution.writeMps("yy.mps");
    2036           // Delete all rows
    2037           n = 0;
    2038           nel = 0;
    2039           lower = m.getRowLower();
    2040           upper = m.getRowUpper();
    2041           starts[0] = 0;
    2042           for (iRow = 0; iRow < numberRows; iRow++) {
    2043                saveLower[n] = lower[iRow];
    2044                saveUpper[n] = upper[iRow];
    2045                int j;
    2046                for (j = startM[iRow]; j < startM[iRow] + lengthM[iRow]; j++) {
    2047                     index[nel] = indexM[j];
    2048                     element[nel++] = elementM[j];
    2049                }
    2050                which[n++] = iRow;
    2051                starts[n] = nel;
    2052           }
    2053           solution.deleteRows(n, which);
    2054           solution.dual();
    2055           // Put back
    2056           solution.addRows(n, saveLower, saveUpper,
    2057                            starts, index, element);
    2058           solution.dual();
    2059           solution.writeMps("xx.mps");
    2060           assert(eq(solution.objectiveValue(), 1.5185098965e+03));
    2061           // Zero out status array to give some interest
    2062           memset(solution.statusArray() + numberColumns, 0, numberRows);
    2063           solution.primal(1);
    2064           assert(eq(solution.objectiveValue(), 1.5185098965e+03));
    2065           // Delete all columns and rows
    2066           n = 0;
    2067           for (iColumn = 0; iColumn < numberColumns; iColumn++) {
    2068                which[n++] = iColumn;
    2069                starts[n] = nel;
    2070           }
    2071           solution.deleteColumns(n, which);
    2072           n = 0;
    2073           for (iRow = 0; iRow < numberRows; iRow++) {
    2074                which[n++] = iRow;
    2075                starts[n] = nel;
    2076           }
    2077           solution.deleteRows(n, which);
    2078 
    2079           delete [] saveObj;
    2080           delete [] saveLower;
    2081           delete [] saveUpper;
    2082           delete [] which;
    2083           delete [] starts;
    2084           delete [] index;
    2085           delete [] element;
     2074               solution.deleteRows(n, which);
     2075               solution.dual();
     2076               // Put back
     2077               solution.addRows(n, saveLower, saveUpper,
     2078                                starts, index, element);
     2079               solution.dual();
     2080               solution.writeMps("xx.mps");
     2081               assert(eq(solution.objectiveValue(), 1.5185098965e+03));
     2082               // Zero out status array to give some interest
     2083               memset(solution.statusArray() + numberColumns, 0, numberRows);
     2084               solution.primal(1);
     2085               assert(eq(solution.objectiveValue(), 1.5185098965e+03));
     2086               // Delete all columns and rows
     2087               n = 0;
     2088               for (iColumn = 0; iColumn < numberColumns; iColumn++) {
     2089                    which[n++] = iColumn;
     2090                    starts[n] = nel;
     2091               }
     2092               solution.deleteColumns(n, which);
     2093               n = 0;
     2094               for (iRow = 0; iRow < numberRows; iRow++) {
     2095                    which[n++] = iRow;
     2096                    starts[n] = nel;
     2097               }
     2098               solution.deleteRows(n, which);
     2099
     2100               delete [] saveObj;
     2101               delete [] saveLower;
     2102               delete [] saveUpper;
     2103               delete [] which;
     2104               delete [] starts;
     2105               delete [] index;
     2106               delete [] element;
     2107          } else {
     2108               std::cerr << "Error reading brandy from sample data. Skipping test." << std::endl;
     2109          }
    20862110     }
    20872111#if 1
     
    20902114          CoinMpsIO m;
    20912115          std::string fn = dirSample + "exmip1";
    2092           m.readMps(fn.c_str(), "mps");
    2093           ClpInterior solution;
    2094           solution.loadProblem(*m.getMatrixByCol(), m.getColLower(), m.getColUpper(),
    2095                                m.getObjCoefficients(),
    2096                                m.getRowLower(), m.getRowUpper());
    2097           solution.primalDual();
     2116          if (m.readMps(fn.c_str(), "mps") == 0) {
     2117               ClpInterior solution;
     2118               solution.loadProblem(*m.getMatrixByCol(), m.getColLower(), m.getColUpper(),
     2119                                    m.getObjCoefficients(),
     2120                                    m.getRowLower(), m.getRowUpper());
     2121               solution.primalDual();
     2122          } else {
     2123               std::cerr << "Error reading exmip1 from sample data. Skipping test." << std::endl;
     2124          }
    20982125     }
    20992126#endif
     
    22832310          CoinMpsIO m;
    22842311          std::string fn = dirSample + "exmip1";
    2285           m.readMps(fn.c_str(), "mps");
    2286           ClpSimplex solution;
    2287           solution.loadProblem(*m.getMatrixByCol(), m.getColLower(), m.getColUpper(),
    2288                                m.getObjCoefficients(),
    2289                                m.getRowLower(), m.getRowUpper());
    2290           //solution.dual();
    2291           // get quadratic part
    2292           int numberColumns = solution.numberColumns();
    2293           int * start = new int [numberColumns+1];
    2294           int * column = new int[numberColumns];
    2295           double * element = new double[numberColumns];
    2296           int i;
    2297           start[0] = 0;
    2298           int n = 0;
    2299           int kk = numberColumns - 1;
    2300           int kk2 = numberColumns - 1;
    2301           for (i = 0; i < numberColumns; i++) {
    2302                if (i >= kk) {
    2303                     column[n] = i;
    2304                     if (i >= kk2)
    2305                          element[n] = 1.0e-1;
    2306                     else
    2307                          element[n] = 0.0;
    2308                     n++;
    2309                }
    2310                start[i+1] = n;
    2311           }
    2312           // Load up objective
    2313           solution.loadQuadraticObjective(numberColumns, start, column, element);
    2314           delete [] start;
    2315           delete [] column;
    2316           delete [] element;
    2317           //solution.quadraticSLP(50,1.0e-4);
    2318           CoinRelFltEq eq(1.0e-4);
    2319           //assert(eq(objValue,820.0));
    2320           //solution.setLogLevel(63);
    2321           solution.primal();
    2322           printSol(solution);
    2323           //assert(eq(objValue,3.2368421));
    2324           //exit(77);
     2312          if (m.readMps(fn.c_str(), "mps") == 0) {
     2313               ClpSimplex solution;
     2314               solution.loadProblem(*m.getMatrixByCol(), m.getColLower(), m.getColUpper(),
     2315                                    m.getObjCoefficients(),
     2316                                    m.getRowLower(), m.getRowUpper());
     2317               //solution.dual();
     2318               // get quadratic part
     2319               int numberColumns = solution.numberColumns();
     2320               int * start = new int [numberColumns+1];
     2321               int * column = new int[numberColumns];
     2322               double * element = new double[numberColumns];
     2323               int i;
     2324               start[0] = 0;
     2325               int n = 0;
     2326               int kk = numberColumns - 1;
     2327               int kk2 = numberColumns - 1;
     2328               for (i = 0; i < numberColumns; i++) {
     2329                    if (i >= kk) {
     2330                         column[n] = i;
     2331                         if (i >= kk2)
     2332                              element[n] = 1.0e-1;
     2333                         else
     2334                              element[n] = 0.0;
     2335                         n++;
     2336                    }
     2337                    start[i+1] = n;
     2338               }
     2339               // Load up objective
     2340               solution.loadQuadraticObjective(numberColumns, start, column, element);
     2341               delete [] start;
     2342               delete [] column;
     2343               delete [] element;
     2344               //solution.quadraticSLP(50,1.0e-4);
     2345               CoinRelFltEq eq(1.0e-4);
     2346               //assert(eq(objValue,820.0));
     2347               //solution.setLogLevel(63);
     2348               solution.primal();
     2349               printSol(solution);
     2350               //assert(eq(objValue,3.2368421));
     2351               //exit(77);
     2352          } else {
     2353               std::cerr << "Error reading exmip1 from sample data. Skipping test." << std::endl;
     2354          }
    23252355     }
    23262356     // Test quadratic
     
    23292359          std::string fn = dirSample + "share2qp";
    23302360          //fn = "share2qpb";
    2331           m.readMps(fn.c_str(), "mps");
    2332           ClpSimplex model;
    2333           model.loadProblem(*m.getMatrixByCol(), m.getColLower(), m.getColUpper(),
    2334                             m.getObjCoefficients(),
    2335                             m.getRowLower(), m.getRowUpper());
    2336           model.dual();
    2337           // get quadratic part
    2338           int * start = NULL;
    2339           int * column = NULL;
    2340           double * element = NULL;
    2341           m.readQuadraticMps(NULL, start, column, element, 2);
    2342           int column2[200];
    2343           double element2[200];
    2344           int start2[80];
    2345           int j;
    2346           start2[0] = 0;
    2347           int nel = 0;
    2348           bool good = false;
    2349           for (j = 0; j < 79; j++) {
    2350                if (start[j] == start[j+1]) {
    2351                     column2[nel] = j;
    2352                     element2[nel] = 0.0;
    2353                     nel++;
    2354                } else {
    2355                     int i;
    2356                     for (i = start[j]; i < start[j+1]; i++) {
    2357                          column2[nel] = column[i];
    2358                          element2[nel++] = element[i];
    2359                     }
    2360                }
    2361                start2[j+1] = nel;
    2362           }
    2363           // Load up objective
    2364           if (good)
    2365                model.loadQuadraticObjective(model.numberColumns(), start2, column2, element2);
    2366           else
    2367                model.loadQuadraticObjective(model.numberColumns(), start, column, element);
    2368           delete [] start;
    2369           delete [] column;
    2370           delete [] element;
    2371           int numberColumns = model.numberColumns();
    2372           model.scaling(0);
     2361          if (m.readMps(fn.c_str(), "mps") == 0) {
     2362               ClpSimplex model;
     2363               model.loadProblem(*m.getMatrixByCol(), m.getColLower(), m.getColUpper(),
     2364                                 m.getObjCoefficients(),
     2365                                 m.getRowLower(), m.getRowUpper());
     2366               model.dual();
     2367               // get quadratic part
     2368               int * start = NULL;
     2369               int * column = NULL;
     2370               double * element = NULL;
     2371               m.readQuadraticMps(NULL, start, column, element, 2);
     2372               int column2[200];
     2373               double element2[200];
     2374               int start2[80];
     2375               int j;
     2376               start2[0] = 0;
     2377               int nel = 0;
     2378               bool good = false;
     2379               for (j = 0; j < 79; j++) {
     2380                    if (start[j] == start[j+1]) {
     2381                         column2[nel] = j;
     2382                         element2[nel] = 0.0;
     2383                         nel++;
     2384                    } else {
     2385                         int i;
     2386                         for (i = start[j]; i < start[j+1]; i++) {
     2387                              column2[nel] = column[i];
     2388                              element2[nel++] = element[i];
     2389                         }
     2390                    }
     2391                    start2[j+1] = nel;
     2392               }
     2393               // Load up objective
     2394               if (good)
     2395                    model.loadQuadraticObjective(model.numberColumns(), start2, column2, element2);
     2396               else
     2397                    model.loadQuadraticObjective(model.numberColumns(), start, column, element);
     2398               delete [] start;
     2399               delete [] column;
     2400               delete [] element;
     2401               int numberColumns = model.numberColumns();
     2402               model.scaling(0);
    23732403#if 0
    2374           model.nonlinearSLP(50, 1.0e-4);
     2404               model.nonlinearSLP(50, 1.0e-4);
    23752405#else
    2376           // Get feasible
    2377           ClpObjective * saveObjective = model.objectiveAsObject()->clone();
    2378           ClpLinearObjective zeroObjective(NULL, numberColumns);
    2379           model.setObjective(&zeroObjective);
    2380           model.dual();
    2381           model.setObjective(saveObjective);
    2382           delete saveObjective;
     2406               // Get feasible
     2407               ClpObjective * saveObjective = model.objectiveAsObject()->clone();
     2408               ClpLinearObjective zeroObjective(NULL, numberColumns);
     2409               model.setObjective(&zeroObjective);
     2410               model.dual();
     2411               model.setObjective(saveObjective);
     2412               delete saveObjective;
    23832413#endif
    2384           //model.setLogLevel(63);
    2385           //exit(77);
    2386           model.setFactorizationFrequency(10);
    2387           model.primal();
    2388           printSol(model);
     2414               //model.setLogLevel(63);
     2415               //exit(77);
     2416               model.setFactorizationFrequency(10);
     2417               model.primal();
     2418               printSol(model);
    23892419#ifndef NDEBUG
    2390           double objValue = model.getObjValue();
     2420               double objValue = model.getObjValue();
    23912421#endif
    2392           CoinRelFltEq eq(1.0e-4);
    2393           assert(eq(objValue, -400.92));
    2394           // and again for barrier
    2395           model.barrier(false);
    2396           //printSol(model);
    2397           model.allSlackBasis();
    2398           model.primal();
    2399           //printSol(model);
     2422               CoinRelFltEq eq(1.0e-4);
     2423               assert(eq(objValue, -400.92));
     2424               // and again for barrier
     2425               model.barrier(false);
     2426               //printSol(model);
     2427               model.allSlackBasis();
     2428               model.primal();
     2429               //printSol(model);
     2430          } else {
     2431               std::cerr << "Error reading share2qp from sample data. Skipping test." << std::endl;
     2432          }
    24002433     }
    24012434     if (0) {
     
    24032436          std::string fn = "./beale";
    24042437          //fn = "./jensen";
    2405           m.readMps(fn.c_str(), "mps");
    2406           ClpSimplex solution;
    2407           solution.loadProblem(*m.getMatrixByCol(), m.getColLower(), m.getColUpper(),
    2408                                m.getObjCoefficients(),
    2409                                m.getRowLower(), m.getRowUpper());
    2410           solution.setDblParam(ClpObjOffset, m.objectiveOffset());
    2411           solution.dual();
    2412           // get quadratic part
    2413           int * start = NULL;
    2414           int * column = NULL;
    2415           double * element = NULL;
    2416           m.readQuadraticMps(NULL, start, column, element, 2);
    2417           // Load up objective
    2418           solution.loadQuadraticObjective(solution.numberColumns(), start, column, element);
    2419           delete [] start;
    2420           delete [] column;
    2421           delete [] element;
    2422           solution.primal(1);
    2423           solution.nonlinearSLP(50, 1.0e-4);
    2424           double objValue = solution.getObjValue();
    2425           CoinRelFltEq eq(1.0e-4);
    2426           assert(eq(objValue, 0.5));
    2427           solution.primal();
    2428           objValue = solution.getObjValue();
    2429           assert(eq(objValue, 0.5));
     2438          if (m.readMps(fn.c_str(), "mps") == 0) {
     2439               ClpSimplex solution;
     2440               solution.loadProblem(*m.getMatrixByCol(), m.getColLower(), m.getColUpper(),
     2441                                    m.getObjCoefficients(),
     2442                                    m.getRowLower(), m.getRowUpper());
     2443               solution.setDblParam(ClpObjOffset, m.objectiveOffset());
     2444               solution.dual();
     2445               // get quadratic part
     2446               int * start = NULL;
     2447               int * column = NULL;
     2448               double * element = NULL;
     2449               m.readQuadraticMps(NULL, start, column, element, 2);
     2450               // Load up objective
     2451               solution.loadQuadraticObjective(solution.numberColumns(), start, column, element);
     2452               delete [] start;
     2453               delete [] column;
     2454               delete [] element;
     2455               solution.primal(1);
     2456               solution.nonlinearSLP(50, 1.0e-4);
     2457               double objValue = solution.getObjValue();
     2458               CoinRelFltEq eq(1.0e-4);
     2459               assert(eq(objValue, 0.5));
     2460               solution.primal();
     2461               objValue = solution.getObjValue();
     2462               assert(eq(objValue, 0.5));
     2463          } else {
     2464               std::cerr << "Error reading beale.mps. Skipping test." << std::endl;
     2465          }
    24302466     }
    24312467#endif
Note: See TracChangeset for help on using the changeset viewer.