source: trunk/Cbc/src/CbcGenSolution.cpp @ 1424

Last change on this file since 1424 was 1374, checked in by bjarni, 10 years ago

Renamed all CbcGenXXXX.cpp_lou back to CbcGenXXXX.cpp (same for .hpp) after reorganizing the CbcParam? and CbcSolver? files

File size: 18.4 KB
Line 
1/*
2  Copyright (C) 2007, Lou Hafer, International Business Machines Corporation
3  and others.  All Rights Reserved.
4
5  This file is part of cbc-generic.
6*/
7
8
9#include <string>
10#include <cstdio>
11#include <cassert>
12
13#include "CoinHelperFunctions.hpp"
14#include "CoinSort.hpp"
15#include "CoinFileIO.hpp"
16
17#include "CbcGenCtlBlk.hpp"
18#include "CbcGenParam.hpp"
19
20namespace {
21
22char svnid[] = "$Id: CbcGenSolution.cpp 1173 2009-06-04 09:44:10Z forrest $" ;
23
24}
25
26namespace {
27
28/*
29  Helper routine to generate masks for selecting names to print. Returns true
30  if masks are generated without error, false otherwise.
31
32  This is John Forrest's code, shamelessly stolen from CoinSolve and tweaked
33  just enough to allow it to be yanked out of the main body of CoinSolve.
34
35  Returns the number of generated masks, -1 on error.
36*/
37
38int generateMasks (std::string proto, int longestName,
39                   int *&maskStarts, char **&finalMasks)
40
41{
42    int nAst = 0 ;
43    const char *pMask2 = proto.c_str() ;
44    char pMask[100] ;
45    int iChar ;
46    int lengthMask = strlen(pMask2) ;
47    assert (lengthMask < 100) ;
48
49    int maxMasks = -1 ;
50
51    maskStarts = 0 ;
52    finalMasks = 0 ;
53    /*
54      Remove surrounding matched quotes if present. Abort if unmatched.
55    */
56    if (*pMask2 == '"') {
57        if (pMask2[lengthMask-1] != '"') {
58            printf("generateMasks: Mismatched \" in mask %s\n", pMask2) ;
59            return (-1) ;
60        } else {
61            strcpy(pMask, pMask2 + 1) ;
62            *strchr(pMask, '"') = '\0' ;
63        }
64    } else if (*pMask2 == '\'') {
65        if (pMask2[lengthMask-1] != '\'') {
66            printf("mismatched ' in mask %s\n", pMask2) ;
67            return (maxMasks) ;
68        } else {
69            strcpy(pMask, pMask2 + 1) ;
70            *strchr(pMask, '\'') = '\0' ;
71        }
72    } else {
73        strcpy(pMask, pMask2) ;
74    }
75    /*
76      Mask should not be longer than longest name.
77    */
78    if (lengthMask > longestName) {
79        printf("mask %s too long - skipping\n", pMask) ;
80        return (maxMasks) ;
81    }
82    /*
83      Expand `*' to multiple masks with varying number of `?' characters.
84    */
85    maxMasks = 1 ;
86    for (iChar = 0; iChar < lengthMask; iChar++) {
87        if (pMask[iChar] == '*') {
88            nAst++ ;
89            maxMasks *= (longestName + 1) ;
90        }
91    }
92    int nEntries = 1 ;
93    maskStarts = new int[longestName+2] ;
94    char ** masks = new char * [maxMasks] ;
95    char ** newMasks = new char * [maxMasks] ;
96    int i ;
97    for (i = 0; i < maxMasks; i++) {
98        masks[i] = new char[longestName+1] ;
99        newMasks[i] = new char[longestName+1] ;
100    }
101    strcpy(masks[0], pMask) ;
102    for (int iAst = 0; iAst < nAst; iAst++) {
103        int nOldEntries = nEntries ;
104        nEntries = 0 ;
105        for (int iEntry = 0; iEntry < nOldEntries; iEntry++) {
106            char * oldMask = masks[iEntry] ;
107            char * ast = strchr(oldMask, '*') ;
108            assert (ast) ;
109            int length = strlen(oldMask) - 1 ;
110            int nBefore = ast - oldMask ;
111            int nAfter = length - nBefore ;
112            // and add null
113            nAfter++ ;
114            for (int i = 0; i <= longestName - length; i++) {
115                char * maskOut = newMasks[nEntries] ;
116                memcpy(maskOut, oldMask, nBefore) ;
117                for (int k = 0; k < i; k++)
118                    maskOut[k+nBefore] = '?' ;
119                memcpy(maskOut + nBefore + i, ast + 1, nAfter) ;
120                nEntries++ ;
121                assert (nEntries <= maxMasks) ;
122            }
123        }
124        char ** temp = masks ;
125        masks = newMasks ;
126        newMasks = temp ;
127    }
128    /*
129      Trim trailing blanks and record final length.
130    */
131    int * sort = new int[nEntries] ;
132    for (i = 0; i < nEntries; i++) {
133        char * maskThis = masks[i] ;
134        int length = strlen(maskThis) ;
135        while (maskThis[length-1] == ' ')
136            length-- ;
137        maskThis[length] = '\0' ;
138        sort[i] = length ;
139    }
140    /*
141      Sort by length.
142    */
143    CoinSort_2(sort, sort + nEntries, masks) ;
144    int lastLength = -1 ;
145    for (i = 0; i < nEntries; i++) {
146        int length = sort[i] ;
147        while (length > lastLength)
148            maskStarts[++lastLength] = i ;
149    }
150    maskStarts[++lastLength] = nEntries ;
151    delete [] sort ;
152    for (i = 0; i < maxMasks; i++)
153        delete [] newMasks[i] ;
154    delete [] newMasks ;
155
156    finalMasks = masks ;
157
158    return (maxMasks) ;
159}
160
161/*
162  Utility routine to check a string against the array of masks. Borrowed
163  from CoinSolve.
164*/
165
166bool maskMatches (const int *starts, char **masks, const char *checkC)
167
168{
169    int length = strlen(checkC);
170    while (checkC[length-1] == ' ')
171        length--;
172    for (int i = starts[length]; i < starts[length+1]; i++) {
173        char * thisMask = masks[i];
174        int k;
175        for ( k = 0; k < length; k++) {
176            if (thisMask[k] != '?' && thisMask[k] != checkC[k])
177                break;
178        }
179        if (k == length)
180            return true;
181    }
182    return (false) ;
183}
184
185}  // end unnamed namespace
186
187/*
188  Routine to write out the solution. Minimally adapted from John's code in
189  CoinSolve to break out a few subroutines and use generic OSI functions.
190
191  The print mode is established by the printingOptions command, and the integer
192  coding is established by the order of declaration of the keyword parameters.
193  As of 060920, known modes are
194
195    normal (0)    print nonzero primal variables
196    integer (1)   print nonzero integer primal variables
197    special (2)   print in a format suitable for OsiRowCutDebugger
198    rows (3)      `normal', plus rows with nonzero row activity
199    all (4)       all primal variables and row activities
200*/
201
202int CbcGenParamUtils::doSolutionParam (CoinParam *param)
203
204{
205    assert (param != 0) ;
206    CbcGenParam *genParam = dynamic_cast<CbcGenParam *>(param) ;
207    assert (genParam != 0) ;
208    CbcGenCtlBlk *ctlBlk = genParam->obj() ;
209    assert (ctlBlk != 0) ;
210    CbcModel *model = ctlBlk->model_ ;
211    assert (model != 0) ;
212    /*
213      Setup to return nonfatal/fatal error (1/-1) by default.
214    */
215    int retval ;
216    if (CoinParamUtils::isInteractive()) {
217        retval = 1 ;
218    } else {
219        retval = -1 ;
220    }
221    /*
222      It's hard to print a solution we don't have.
223    */
224    if (ctlBlk->bab_.haveAnswer_ == false) {
225        std::cout << "There is no solution available to print." << std::endl ;
226        return (retval) ;
227    }
228    OsiSolverInterface *osi = ctlBlk->bab_.answerSolver_ ;
229    assert (osi != 0) ;
230    /*
231      Figure out where we're going to write the solution. As special cases,
232      `$' says `use the previous output file' and `-' says `use stdout'.
233
234      cbc will also accept `no string value' as stdout, but that'd be a real pain
235      in this architecture.
236    */
237    std::string field = genParam->strVal() ;
238    std::string fileName ;
239    if (field == "$") {
240        fileName = ctlBlk->lastSolnOut_ ;
241        field = fileName ;
242    }
243    if (field == "-") {
244        fileName = "stdout" ;
245        field = fileName ;
246    } else {
247        fileName = field ;
248    }
249    if (!(fileName == "stdout" || fileName == "stderr")) {
250        if (fileName[0] == '~') {
251            char dirSep = CoinFindDirSeparator() ;
252            if (fileName[1] == dirSep) {
253                char *environVar = getenv("HOME") ;
254                if (environVar) {
255                    std::string home(environVar) ;
256                    fileName = home + fileName.substr(1) ;
257                }
258            }
259        }
260        if (!(fileAbsPath(fileName) || fileName.substr(0, 2) == "./")) {
261            fileName = ctlBlk->dfltDirectory_ + fileName ;
262        }
263    }
264    /*
265      See if we can open the file. Bail out if the open fails.
266    */
267    FILE *fp ;
268    if (fileName == "stdout") {
269        fp = stdout ;
270    } else if (fileName == "stderr") {
271        fp = stderr ;
272    } else {
273        fp = fopen(fileName.c_str(), "w") ;
274        fp = fopen(fileName.c_str(), "w") ;
275    }
276    if (!fp) {
277        std::cout
278            << "Unable to open file `" << fileName
279            << "', original name '" << genParam->strVal() << "'." << std::endl ;
280        return (retval) ;
281    } else {
282        std::cout
283            << "Writing solution to `" << fileName << "'." << std::endl ;
284    }
285
286    int m = osi->getNumRows() ;
287    int n = osi->getNumCols() ;
288    int iColumn ;
289    const double *primalColSolution = osi->getColSolution() ;
290    /*
291      Separate printing a solution for humans from printing a solution for the
292      row cut debugger. For the row cut debugger, we want to produce C++ code
293      that can be pasted into the debugger's set of known problems.
294    */
295    if (ctlBlk->printMode_ == 2) {
296        int k = 0 ;
297        bool newLine = true ;
298        bool comma = false ;
299        fprintf(fp, "  int intIndicesV[] = {") ;
300        for (iColumn = 0 ; iColumn < n ; iColumn++ ) {
301            if (fabs(primalColSolution[iColumn]) > 0.5 && osi->isInteger(iColumn)) {
302                if (newLine) {
303                    fprintf(fp, "\n\t") ;
304                    newLine = false ;
305                } else {
306                    fprintf(fp, ", ") ;
307                }
308                fprintf(fp, "%d", iColumn) ;
309                if (++k == 10) {
310                    k = 0 ;
311                    newLine = true ;
312                }
313            }
314        }
315        fprintf(fp, "\n      } ;\n") ;
316        k = 0 ;
317        newLine = true ;
318        fprintf(fp, "  double intSolnV[] = {") ;
319        for (iColumn = 0 ; iColumn < n ; iColumn++) {
320            double value = primalColSolution[iColumn] ;
321            if (fabs(value) > 0.5 && osi->isInteger(iColumn)) {
322                if (newLine) {
323                    fprintf(fp, "\n\t") ;
324                    newLine = false ;
325                } else {
326                    fprintf(fp, ", ") ;
327                }
328                if (value > 0) {
329                    value = floor(value + .5)  ;
330                } else {
331                    value = ceil(value - .5) ;
332                }
333                int ivalue = static_cast<int>(value) ;
334                fprintf(fp, "%d.0", ivalue) ;
335                if (++k == 10) {
336                    k = 0 ;
337                    newLine = true ;
338                }
339            }
340        }
341        fprintf(fp, "\n      } ;\n") ;
342        return (0) ;
343    }
344    /*
345      Begin the code to generate output meant for a human.  What's our longest
346      name? Scan the names we're going to print. printMode_ of 3 or 4 requires we
347      scan the row names too. Force between 8 and 20 characters in any event.
348    */
349    int longestName = 0 ;
350    for (int j = 0 ; j < n ; j++) {
351        int len = osi->getColName(j).length() ;
352        longestName = CoinMax(longestName, len) ;
353    }
354    if (ctlBlk->printMode_ >= 3) {
355        for (int i = 0 ; i < m ; i++) {
356            int len = osi->getRowName(i).length() ;
357            longestName = CoinMax(longestName, len) ;
358        }
359    }
360    /*
361      Generate masks if we need to do so.
362    */
363    bool doMask = ctlBlk->printMask_ != "" ;
364    int *maskStarts = NULL ;
365    int maxMasks = 0 ;
366    char **masks = NULL ;
367    if (doMask) {
368        maxMasks = generateMasks(ctlBlk->printMask_, longestName, maskStarts, masks) ;
369        if (maxMasks < 0) {
370            return (retval) ;
371        }
372    }
373    /*
374      Force the space allocated to names to be between 8 and 20 characters.
375    */
376    if (longestName < 8) {
377        longestName = 8 ;
378    } else if (longestName > 20) {
379        longestName = 20 ;
380    }
381    /*
382      Print requested components of the row solution. Only modes 3 (rows) and 4
383      (all) will print row information. For the rows that we print, print both
384      the row activity and the value of the associated dual.
385
386      Which to print? Violated constraints will always be flagged to print.
387      Otherwise, if m < 50 or all rows are requested, print all rows.  Otherwise,
388      print tight constraints (non-zero dual).
389
390      All of this is filtered through printMask, if specified.
391    */
392    double primalTolerance ;
393    osi->getDblParam(OsiPrimalTolerance, primalTolerance) ;
394
395    int iRow ;
396    if (ctlBlk->printMode_ >= 3) {
397        const double *dualRowSolution = osi->getRowPrice() ;
398        const double *primalRowSolution = osi->getRowActivity() ;
399        const double *rowLower = osi->getRowLower() ;
400        const double *rowUpper = osi->getRowUpper() ;
401
402        fprintf(fp, "\n   %7s %-*s%15s%15s\n\n",
403                "Index", longestName, "Row", "Activity", "Dual") ;
404
405        for (iRow = 0 ; iRow < m ; iRow++) {
406            bool violated = false ;
407            bool print = false ;
408            if (primalRowSolution[iRow] > rowUpper[iRow] + primalTolerance ||
409                    primalRowSolution[iRow] < rowLower[iRow] - primalTolerance) {
410                violated = true ;
411                print = true ;
412            } else {
413                if (m < 50 || ctlBlk->printMode_ >= 4) {
414                    print = true ;
415                } else if (fabs(dualRowSolution[iRow]) > 1.0e-8) {
416                    print = true ;
417                }
418            }
419            const char *name = osi->getRowName(iRow).c_str() ;
420            if (doMask && !maskMatches(maskStarts, masks, name)) {
421                print = false ;
422            }
423            if (print) {
424                if (violated) {
425                    fprintf(fp, "** ") ;
426                } else {
427                    fprintf(fp, "%3s", " ") ;
428                }
429                fprintf(fp, "%7d %-*s%15.8g%15.8g\n", iRow, longestName, name,
430                        primalRowSolution[iRow], dualRowSolution[iRow]) ;
431            }
432        }
433        fprintf(fp, "\n") ;
434    }
435    /*
436      Now do the columns. This first block handles all modes except 2 (special).
437      Out-of-bounds variables are flagged with `**'. If there are less than 50
438      variables, all are printed. All of this is filtered through `integer only'
439      and can be further filtered using printMask.
440    */
441    if (ctlBlk->printMode_ != 2) {
442        const double *columnLower = osi->getColLower() ;
443        const double *columnUpper = osi->getColUpper() ;
444        const double *dualColSolution = osi->getReducedCost() ;
445
446        fprintf(fp, "\n   %7s %-*s%15s%15s\n\n",
447                "Index", longestName, "Column", "Value", "Reduced Cost") ;
448
449        for (iColumn = 0 ; iColumn < n ; iColumn++) {
450            bool violated = false ;
451            bool print = false ;
452            if (primalColSolution[iColumn] > columnUpper[iColumn] + primalTolerance ||
453                    primalColSolution[iColumn] < columnLower[iColumn] - primalTolerance) {
454                violated = true ;
455                print = true ;
456            } else {
457                if (n < 50 || ctlBlk->printMode_ == 4) {
458                    print = true ;
459                } else if (fabs(primalColSolution[iColumn]) > 1.0e-8) {
460                    if (ctlBlk->printMode_ == 1) {
461                        print = osi->isInteger(iColumn) ;
462                    } else {
463                        print = true ;
464                    }
465                }
466            }
467            const char *name = osi->getColName(iColumn).c_str() ;
468            if (doMask && !maskMatches(maskStarts, masks, name)) {
469                print = false ;
470            }
471            if (print) {
472                if (violated) {
473                    fprintf(fp, "** ") ;
474                } else {
475                    fprintf(fp, "%3s", " ") ;
476                }
477                fprintf(fp, "%7d %-*s%15.8g%15.8g\n", iColumn, longestName, name,
478                        primalColSolution[iColumn], dualColSolution[iColumn]) ;
479            }
480        }
481    }
482    /*
483      Close out the file, but don't close stdout. Delete any masks.
484    */
485    if (fp != stdout) {
486        fclose(fp) ;
487    }
488
489    if (masks) {
490        delete [] maskStarts ;
491        for (int i = 0 ; i < maxMasks ; i++)
492            delete [] masks[i] ;
493        delete [] masks ;
494    }
495
496    return (0) ;
497}
498
499
500/*
501  Routine to do initial verification of a print mask. We don't generate the
502  full set of print masks here, but we do some verification checks to make sure
503  it's valid.
504*/
505
506int CbcGenParamUtils::doPrintMaskParam (CoinParam *param)
507
508{
509    assert (param != 0) ;
510    CbcGenParam *genParam = dynamic_cast<CbcGenParam *>(param) ;
511    assert (genParam != 0) ;
512    CbcGenCtlBlk *ctlBlk = genParam->obj() ;
513    assert (ctlBlk != 0) ;
514    /*
515      Setup to return nonfatal/fatal error (1/-1) by default.
516    */
517    int retval ;
518    if (CoinParamUtils::isInteractive()) {
519        retval = 1 ;
520    } else {
521        retval = -1 ;
522    }
523    /*
524      Now do a bit of verification of the mask. It should be non-null and, if
525      quoted, the quotes should be matched. Aribtrarily put the absolute maximum
526      length at 50 characters. If we have a model loaded, that'll be tightened to
527      the length of the longest name.
528    */
529    std::string maskProto = param->strVal() ;
530    int maskLen = maskProto.length() ;
531    if (maskLen <= 0 || maskLen > 50) {
532        std::cerr
533            << "Mask |" << maskProto
534            << "| is " << maskLen << " characters; should be between "
535            << 0 << " and " << 50 << "." << std::endl ;
536        return (retval) ;
537    }
538    /*
539      Remove surrounding matched quotes if present. Abort if unmatched.
540    */
541    if (maskProto[0] == '"' || maskProto[0] == '\'') {
542        char quoteChar = maskProto[0] ;
543        if (maskProto[maskLen-1] != quoteChar) {
544            std::cerr
545                << "Mismatched quotes around mask |" << maskProto
546                << "|." << std::endl ;
547            return (retval) ;
548        } else {
549            maskProto = maskProto.substr(1, maskLen - 2) ;
550        }
551    }
552    /*
553      Mask should not be longer than longest name. Of course, if we don't have a
554      model, we can't do this check.
555    */
556    if (ctlBlk->goodModel_) {
557        CbcModel *model = ctlBlk->model_ ;
558        assert (model != 0) ;
559        OsiSolverInterface *osi = model->solver() ;
560        assert (osi != 0) ;
561        int longestName = 0 ;
562        int n = osi->getNumCols() ;
563        for (int j = 0 ; j < n ; j++) {
564            int len = osi->getColName(j).length() ;
565            longestName = CoinMax(longestName, len) ;
566        }
567        int m = osi->getNumRows() ;
568        for (int i = 0 ; i < m ; i++) {
569            int len = osi->getRowName(i).length() ;
570            longestName = CoinMax(longestName, len) ;
571        }
572        if (maskLen > longestName) {
573            std::cerr
574                << "Mask |" << maskProto << "| has " << maskLen << " chars; this"
575                << " is longer than the longest name (" << longestName
576                << " chars)." << std::endl ;
577            return (retval) ;
578        }
579    }
580
581    ctlBlk->printMask_ = maskProto ;
582
583    return (0) ;
584}
Note: See TracBrowser for help on using the repository browser.