source: branches/pre/ClpNetworkMatrix.cpp @ 180

Last change on this file since 180 was 180, checked in by forrest, 18 years ago

new stuff

  • Property svn:eol-style set to native
  • Property svn:keywords set to Author Date Id Revision
File size: 21.1 KB
Line 
1// Copyright (C) 2003, International Business Machines
2// Corporation and others.  All Rights Reserved.
3
4
5#include <cstdio>
6
7#include "CoinPragma.hpp"
8#include "CoinIndexedVector.hpp"
9#include "CoinHelperFunctions.hpp"
10
11#include "ClpSimplex.hpp"
12#include "ClpFactorization.hpp"
13// at end to get min/max!
14#include "ClpNetworkMatrix.hpp"
15#include "ClpPlusMinusOneMatrix.hpp"
16#include "ClpMessage.hpp"
17
18//#############################################################################
19// Constructors / Destructor / Assignment
20//#############################################################################
21
22//-------------------------------------------------------------------
23// Default Constructor
24//-------------------------------------------------------------------
25ClpNetworkMatrix::ClpNetworkMatrix () 
26  : ClpMatrixBase()
27{
28  setType(11);
29  elements_ = NULL;
30  starts_ = NULL;
31  lengths_=NULL;
32  indices_=NULL;
33  numberRows_=0;
34  numberColumns_=0;
35  trueNetwork_=false;
36}
37
38/* Constructor from two arrays */
39ClpNetworkMatrix::ClpNetworkMatrix(int numberColumns, const int * head,
40                                   const int * tail)
41  : ClpMatrixBase()
42{
43  setType(11);
44  elements_ = NULL;
45  starts_ = NULL;
46  lengths_=NULL;
47  indices_=new int[2*numberColumns];;
48  numberRows_=-1;
49  numberColumns_=numberColumns;
50  trueNetwork_=true;
51  int iColumn;
52  CoinBigIndex j=0;
53  for (iColumn=0;iColumn<numberColumns_;iColumn++, j+=2) {
54    int iRow = head[iColumn];
55    numberRows_ = max(numberRows_,iRow);
56    indices_[j]=iRow;
57    iRow = tail[iColumn];
58    numberRows_ = max(numberRows_,iRow);
59    indices_[j+1]=iRow;
60  }
61  numberRows_++;
62}
63//-------------------------------------------------------------------
64// Copy constructor
65//-------------------------------------------------------------------
66ClpNetworkMatrix::ClpNetworkMatrix (const ClpNetworkMatrix & rhs) 
67: ClpMatrixBase(rhs)
68{ 
69  elements_ = NULL;
70  starts_ = NULL;
71  lengths_=NULL;
72  indices_=NULL;
73  numberRows_=rhs.numberRows_;
74  numberColumns_=rhs.numberColumns_;
75  trueNetwork_=rhs.trueNetwork_;
76  if (numberColumns_) {
77    indices_ = new int [ 2*numberColumns_];
78    memcpy(indices_,rhs.indices_,2*numberColumns_*sizeof(int));
79  }
80}
81
82ClpNetworkMatrix::ClpNetworkMatrix (const CoinPackedMatrix & rhs) 
83  : ClpMatrixBase()
84{ 
85  setType(11);
86  elements_ = NULL;
87  starts_ = NULL;
88  lengths_=NULL;
89  indices_=NULL;
90  int iColumn;
91  assert (rhs.isColOrdered());
92  // get matrix data pointers
93  const int * row = rhs.getIndices();
94  const CoinBigIndex * columnStart = rhs.getVectorStarts();
95  const int * columnLength = rhs.getVectorLengths(); 
96  const double * elementByColumn = rhs.getElements();
97  numberColumns_ = rhs.getNumCols();
98  int goodNetwork=1;
99  numberRows_=-1;
100  indices_ = new int[2*numberColumns_];
101  CoinBigIndex j=0;
102  for (iColumn=0;iColumn<numberColumns_;iColumn++, j+=2) {
103    CoinBigIndex k=columnStart[iColumn];
104    int iRow;
105    switch (columnLength[iColumn]) {
106    case 0:
107      goodNetwork=-1; // not classic network
108      indices_[j]=-1;
109      indices_[j+1]=-1;
110      break;
111     
112    case 1:
113      goodNetwork=-1; // not classic network
114      if (fabs(elementByColumn[k]-1.0)<1.0e-10) {
115        indices_[j] = -1;
116        iRow = row[k];
117        numberRows_ = max(numberRows_,iRow);
118        indices_[j+1]=iRow;
119      } else if (fabs(elementByColumn[k]+1.0)<1.0e-10) {
120        indices_[j+1] = -1;
121        iRow = row[k];
122        numberRows_ = max(numberRows_,iRow);
123        indices_[j]=iRow;
124      } else {
125        goodNetwork = 0; // not a network
126      }
127      break;
128     
129    case 2:
130      if (fabs(elementByColumn[k]-1.0)<1.0e-10) {
131        if (fabs(elementByColumn[k+1]+1.0)<1.0e-10) {
132          iRow = row[k];
133          numberRows_ = max(numberRows_,iRow);
134          indices_[j+1]=iRow;
135          iRow = row[k+1];
136          numberRows_ = max(numberRows_,iRow);
137          indices_[j] = iRow;
138        } else {
139          goodNetwork = 0; // not a network
140        }
141      } else if (fabs(elementByColumn[k]+1.0)<1.0e-10) {
142        if (fabs(elementByColumn[k+1]-1.0)<1.0e-10) {
143          iRow = row[k];
144          numberRows_ = max(numberRows_,iRow);
145          indices_[j]=iRow;
146          iRow = row[k+1];
147          numberRows_ = max(numberRows_,iRow);
148          indices_[j+1] = iRow;
149        } else {
150          goodNetwork = 0; // not a network
151        }
152      } else {
153        goodNetwork = 0; // not a network
154      }
155      break;
156
157    default:
158      goodNetwork = 0; // not a network
159      break;
160    }
161    if (!goodNetwork)
162      break;
163  }
164  if (!goodNetwork) {
165    delete [] indices_;
166    // put in message
167    printf("Not a network - can test if indices_ null\n");
168    indices_=NULL;
169    numberRows_=0;
170    numberColumns_=0;
171  } else {
172    numberRows_ ++; //  correct
173    trueNetwork_ = goodNetwork>0;
174  }
175}
176
177//-------------------------------------------------------------------
178// Destructor
179//-------------------------------------------------------------------
180ClpNetworkMatrix::~ClpNetworkMatrix ()
181{
182  delete [] elements_;
183  delete [] starts_;
184  delete [] lengths_;
185  delete [] indices_;
186}
187
188//----------------------------------------------------------------
189// Assignment operator
190//-------------------------------------------------------------------
191ClpNetworkMatrix &
192ClpNetworkMatrix::operator=(const ClpNetworkMatrix& rhs)
193{
194  if (this != &rhs) {
195    ClpMatrixBase::operator=(rhs);
196    delete [] elements_;
197    delete [] starts_;
198    delete [] lengths_;
199    delete [] indices_;
200    elements_ = NULL;
201    starts_ = NULL;
202    lengths_=NULL;
203    indices_=NULL;
204    numberRows_=rhs.numberRows_;
205    numberColumns_=rhs.numberColumns_;
206    trueNetwork_=rhs.trueNetwork_;
207    if (numberColumns_) {
208      indices_ = new int [ 2*numberColumns_];
209      memcpy(indices_,rhs.indices_,2*numberColumns_*sizeof(int));
210    }
211  }
212  return *this;
213}
214//-------------------------------------------------------------------
215// Clone
216//-------------------------------------------------------------------
217ClpMatrixBase * ClpNetworkMatrix::clone() const
218{
219  return new ClpNetworkMatrix(*this);
220}
221
222/* Returns a new matrix in reverse order without gaps */
223ClpMatrixBase * 
224ClpNetworkMatrix::reverseOrderedCopy() const
225{
226  // count number in each row
227  int * tempP = new int [numberRows_];
228  int * tempN = new int [numberRows_];
229  memset(tempP,0,numberRows_*sizeof(int));
230  memset(tempN,0,numberRows_*sizeof(int));
231  CoinBigIndex j=0;
232  int i;
233  for (i=0;i<numberColumns_;i++,j+=2) {
234    int iRow = indices_[j];
235    tempN[iRow]++;
236    iRow = indices_[j+1];
237    tempP[iRow]++;
238  }
239  int * newIndices = new int [2*numberColumns_];
240  int * newP = new int [numberRows_+1];
241  int * newN = new int[numberRows_];
242  int iRow;
243  j=0;
244  // do starts
245  for (iRow=0;iRow<numberRows_;iRow++) {
246    newP[iRow]=j;
247    j += tempP[iRow];
248    tempP[iRow]=newP[iRow];
249    newN[iRow] = j;
250    j += tempN[iRow];
251    tempN[iRow]=newN[iRow];
252  }
253  newP[numberRows_]=j;
254  j=0;
255  for (i=0;i<numberColumns_;i++,j+=2) {
256    int iRow = indices_[j];
257    int put = tempN[iRow];
258    newIndices[put++] = i;
259    tempN[iRow] = put;
260    iRow = indices_[j+1];
261    put = tempP[iRow];
262    newIndices[put++] = i;
263    tempP[iRow] = put;
264  }
265  delete [] tempP;
266  delete [] tempN;
267  ClpPlusMinusOneMatrix * newCopy = new ClpPlusMinusOneMatrix();
268  newCopy->passInCopy(numberRows_, numberColumns_,
269                      false,  newIndices, newP, newN);
270  return newCopy;
271}
272//unscaled versions
273void 
274ClpNetworkMatrix::times(double scalar,
275                   const double * x, double * y) const
276{
277  int iColumn;
278  CoinBigIndex j=0;
279  if (trueNetwork_) {
280    for (iColumn=0;iColumn<numberColumns_;iColumn++,j+=2) {
281      double value = scalar*x[iColumn];
282      if (value) {
283        int iRowM = indices_[j];
284        int iRowP = indices_[j+1];
285        y[iRowM] -= value;
286        y[iRowP] += value;
287      }
288    }
289  } else {
290    // skip negative rows
291    for (iColumn=0;iColumn<numberColumns_;iColumn++,j+=2) {
292      double value = scalar*x[iColumn];
293      if (value) {
294        int iRowM = indices_[j];
295        int iRowP = indices_[j+1];
296        if (iRowM>=0)
297          y[iRowM] -= value;
298        if (iRowP>=0)
299          y[iRowP] += value;
300      }
301    }
302  }
303}
304void 
305ClpNetworkMatrix::transposeTimes(double scalar,
306                                const double * x, double * y) const
307{
308  int iColumn;
309  CoinBigIndex j=0;
310  if (trueNetwork_) {
311    for (iColumn=0;iColumn<numberColumns_;iColumn++,j+=2) {
312      double value = y[iColumn];
313      int iRowM = indices_[j];
314      int iRowP = indices_[j+1];
315      value -= scalar*x[iRowM];
316      value += scalar*x[iRowP];
317      y[iColumn] = value;
318    }
319  } else {
320    // skip negative rows
321    for (iColumn=0;iColumn<numberColumns_;iColumn++,j+=2) {
322      double value = y[iColumn];
323      int iRowM = indices_[j];
324      int iRowP = indices_[j+1];
325      if (iRowM>=0)
326        value -= scalar*x[iRowM];
327      if (iRowP>=0)
328        value += scalar*x[iRowP];
329      y[iColumn] = value;
330    }
331  }
332}
333void 
334ClpNetworkMatrix::times(double scalar,
335                       const double * x, double * y,
336                       const double * rowScale, 
337                       const double * columnScale) const
338{
339  // we know it is not scaled
340  times(scalar, x, y);
341}
342void 
343ClpNetworkMatrix::transposeTimes( double scalar,
344                                 const double * x, double * y,
345                                 const double * rowScale, 
346                                 const double * columnScale) const
347{
348  // we know it is not scaled
349  transposeTimes(scalar, x, y);
350}
351/* Return <code>x * A + y</code> in <code>z</code>.
352        Squashes small elements and knows about ClpSimplex */
353void 
354ClpNetworkMatrix::transposeTimes(const ClpSimplex * model, double scalar,
355                              const CoinIndexedVector * rowArray,
356                              CoinIndexedVector * y,
357                              CoinIndexedVector * columnArray) const
358{
359  // we know it is not scaled
360  columnArray->clear();
361  double * pi = rowArray->denseVector();
362  int numberNonZero=0;
363  int * index = columnArray->getIndices();
364  double * array = columnArray->denseVector();
365  int numberInRowArray = rowArray->getNumElements();
366  // maybe I need one in OsiSimplex
367  double zeroTolerance = model->factorization()->zeroTolerance();
368  int numberRows = model->numberRows();
369  ClpPlusMinusOneMatrix* rowCopy =
370    dynamic_cast< ClpPlusMinusOneMatrix*>(model->rowCopy());
371  bool packed = rowArray->packedMode();
372  if (numberInRowArray>0.3*numberRows||!rowCopy) {
373    // do by column
374    int iColumn;
375    assert (!y->getNumElements());
376    CoinBigIndex j=0;
377    if (packed) {
378      // need to expand pi into y
379      assert(y->capacity()>=numberRows);
380      double * piOld = pi;
381      pi = y->denseVector();
382      const int * whichRow = rowArray->getIndices();
383      int i;
384      // modify pi so can collapse to one loop
385      for (i=0;i<numberInRowArray;i++) {
386        int iRow = whichRow[i];
387        pi[iRow]=scalar*piOld[i];
388      }
389      if (trueNetwork_) {
390        for (iColumn=0;iColumn<numberColumns_;iColumn++,j+=2) {
391          double value = 0.0;
392          int iRowM = indices_[j];
393          int iRowP = indices_[j+1];
394          value -= pi[iRowM];
395          value += pi[iRowP];
396          if (fabs(value)>zeroTolerance) {
397            array[numberNonZero]=value;
398            index[numberNonZero++]=iColumn;
399          }
400        }
401      } else {
402        // skip negative rows
403        for (iColumn=0;iColumn<numberColumns_;iColumn++,j+=2) {
404          double value = 0.0;
405          int iRowM = indices_[j];
406          int iRowP = indices_[j+1];
407          if (iRowM>=0)
408            value -= pi[iRowM];
409          if (iRowP>=0)
410            value += pi[iRowP];
411          if (fabs(value)>zeroTolerance) {
412            array[numberNonZero]=value;
413            index[numberNonZero++]=iColumn;
414          }
415        }
416      }
417      for (i=0;i<numberInRowArray;i++) {
418        int iRow = whichRow[i];
419        pi[iRow]=0.0;
420      }
421    } else {
422      if (trueNetwork_) {
423        for (iColumn=0;iColumn<numberColumns_;iColumn++,j+=2) {
424          double value = 0.0;
425          int iRowM = indices_[j];
426          int iRowP = indices_[j+1];
427          value -= scalar*pi[iRowM];
428          value += scalar*pi[iRowP];
429          if (fabs(value)>zeroTolerance) {
430            index[numberNonZero++]=iColumn;
431            array[iColumn]=value;
432          }
433        }
434      } else {
435        // skip negative rows
436        for (iColumn=0;iColumn<numberColumns_;iColumn++,j+=2) {
437          double value = 0.0;
438          int iRowM = indices_[j];
439          int iRowP = indices_[j+1];
440          if (iRowM>=0)
441            value -= scalar*pi[iRowM];
442          if (iRowP>=0)
443            value += scalar*pi[iRowP];
444          if (fabs(value)>zeroTolerance) {
445            index[numberNonZero++]=iColumn;
446            array[iColumn]=value;
447          }
448        }
449      }
450    }
451    columnArray->setNumElements(numberNonZero);
452  } else {
453    // do by row
454    rowCopy->transposeTimesByRow(model, scalar, rowArray, y, columnArray);
455  }
456}
457/* Return <code>x *A in <code>z</code> but
458   just for indices in y.
459   Squashes small elements and knows about ClpSimplex */
460void 
461ClpNetworkMatrix::subsetTransposeTimes(const ClpSimplex * model,
462                              const CoinIndexedVector * rowArray,
463                              const CoinIndexedVector * y,
464                              CoinIndexedVector * columnArray) const
465{
466  columnArray->clear();
467  double * pi = rowArray->denseVector();
468  int numberNonZero=0;
469  int * index = columnArray->getIndices();
470  double * array = columnArray->denseVector();
471  // maybe I need one in OsiSimplex
472  double zeroTolerance = model->factorization()->zeroTolerance();
473  int jColumn;
474  int numberToDo = y->getNumElements();
475  const int * which = y->getIndices();
476  bool packed = rowArray->packedMode();
477  if (packed) {
478    // Must line up with y
479    // need to expand pi into y
480    int numberInRowArray = rowArray->getNumElements();
481    assert(y->capacity()>=model->numberRows());
482    double * piOld = pi;
483    pi = y->denseVector();
484    const int * whichRow = rowArray->getIndices();
485    int i;
486    for (i=0;i<numberInRowArray;i++) {
487      int iRow = whichRow[i];
488      pi[iRow]=piOld[i];
489    }
490    if (trueNetwork_) {
491      for (jColumn=0;jColumn<numberToDo;jColumn++) {
492        int iColumn = which[jColumn];
493        double value = 0.0;
494        CoinBigIndex j=iColumn<<1;
495        int iRowM = indices_[j];
496        int iRowP = indices_[j+1];
497        value -= pi[iRowM];
498        value += pi[iRowP];
499        array[jColumn]=value;
500      }
501    } else {
502      // skip negative rows
503      for (jColumn=0;jColumn<numberToDo;jColumn++) {
504        int iColumn = which[jColumn];
505        double value = 0.0;
506        CoinBigIndex j=iColumn<<1;
507        int iRowM = indices_[j];
508        int iRowP = indices_[j+1];
509        if (iRowM>=0)
510          value -= pi[iRowM];
511        if (iRowP>=0)
512          value += pi[iRowP];
513        array[jColumn]=value;
514      }
515    }
516    for (i=0;i<numberInRowArray;i++) {
517      int iRow = whichRow[i];
518      pi[iRow]=0.0;
519    }
520  } else {
521    if (trueNetwork_) {
522      for (jColumn=0;jColumn<numberToDo;jColumn++) {
523        int iColumn = which[jColumn];
524        double value = 0.0;
525        CoinBigIndex j=iColumn<<1;
526        int iRowM = indices_[j];
527        int iRowP = indices_[j+1];
528        value -= pi[iRowM];
529        value += pi[iRowP];
530        if (fabs(value)>zeroTolerance) {
531          index[numberNonZero++]=iColumn;
532          array[iColumn]=value;
533        }
534      }
535    } else {
536      // skip negative rows
537      for (jColumn=0;jColumn<numberToDo;jColumn++) {
538        int iColumn = which[jColumn];
539        double value = 0.0;
540        CoinBigIndex j=iColumn<<1;
541        int iRowM = indices_[j];
542        int iRowP = indices_[j+1];
543        if (iRowM>=0)
544          value -= pi[iRowM];
545        if (iRowP>=0)
546          value += pi[iRowP];
547        if (fabs(value)>zeroTolerance) {
548          index[numberNonZero++]=iColumn;
549          array[iColumn]=value;
550        }
551      }
552    }
553  }
554}
555/* Returns number of elements in basis
556   column is basic if entry >=0 */
557CoinBigIndex
558ClpNetworkMatrix::numberInBasis(const int * columnIsBasic) const 
559{
560  int i;
561  CoinBigIndex numberElements=0;
562  if (trueNetwork_) {
563    for (i=0;i<numberColumns_;i++) {
564      if (columnIsBasic[i]>=0) 
565        numberElements ++;
566    }
567    numberElements *= 2;
568  } else {
569    for (i=0;i<numberColumns_;i++) {
570      if (columnIsBasic[i]>=0) {
571        CoinBigIndex j=i<<1;
572        int iRowM = indices_[j];
573        int iRowP = indices_[j+1];
574        if (iRowM>=0)
575          numberElements ++;
576        if (iRowP>=0)
577          numberElements ++;
578      }
579    }
580  }
581  return numberElements;
582}
583// Fills in basis (Returns number of elements and updates numberBasic)
584CoinBigIndex
585ClpNetworkMatrix::fillBasis(const ClpSimplex * model,
586                                const int * columnIsBasic, int & numberBasic,
587                                int * indexRowU, int * indexColumnU,
588                                double * elementU) const 
589{
590#ifdef CLPDEBUG
591  const double * rowScale = model->rowScale();
592  assert (!rowScale);
593#endif
594  int i;
595  CoinBigIndex numberElements=0;
596  if (trueNetwork_) {
597    for (i=0;i<numberColumns_;i++) {
598      if (columnIsBasic[i]>=0) {
599        CoinBigIndex j=i<<1;
600        int iRowM = indices_[j];
601        int iRowP = indices_[j+1];
602        indexRowU[numberElements]=iRowM;
603        indexColumnU[numberElements]=numberBasic;
604        elementU[numberElements]=-1.0;
605        indexRowU[numberElements+1]=iRowP;
606        indexColumnU[numberElements+1]=numberBasic;
607        elementU[numberElements+1]=1.0;
608        numberElements+=2;
609        numberBasic++;
610      }
611    }
612  } else {
613    for (i=0;i<numberColumns_;i++) {
614      if (columnIsBasic[i]>=0) {
615        CoinBigIndex j=i<<1;
616        int iRowM = indices_[j];
617        int iRowP = indices_[j+1];
618        if (iRowM>=0) {
619          indexRowU[numberElements]=iRowM;
620          indexColumnU[numberElements]=numberBasic;
621          elementU[numberElements++]=-1.0;
622        }
623        if (iRowP>=0) {
624          indexRowU[numberElements]=iRowP;
625          indexColumnU[numberElements]=numberBasic;
626          elementU[numberElements++]=1.0;
627        }
628        numberBasic++;
629      }
630    }
631  }
632  return numberElements;
633}
634/* If element NULL returns number of elements in column part of basis,
635   If not NULL fills in as well */
636CoinBigIndex
637ClpNetworkMatrix::fillBasis(const ClpSimplex * model,
638                                 const int * whichColumn, 
639                                 int numberBasic,
640                                 int numberColumnBasic,
641                                 int * indexRowU, int * indexColumnU,
642                                 double * elementU) const 
643{
644  int i;
645  CoinBigIndex numberElements=0;
646  if (elementU!=NULL) {
647    if (trueNetwork_) {
648      for (i=0;i<numberColumnBasic;i++) {
649        int iColumn = whichColumn[i];
650        CoinBigIndex j=iColumn<<1;
651        int iRowM = indices_[j];
652        int iRowP = indices_[j+1];
653        indexRowU[numberElements]=iRowM;
654        indexColumnU[numberElements]=numberBasic;
655        elementU[numberElements]=-1.0;
656        indexRowU[numberElements+1]=iRowP;
657        indexColumnU[numberElements+1]=numberBasic;
658        elementU[numberElements+1]=1.0;
659        numberElements+=2;
660        numberBasic++;
661      }
662    } else {
663      for (i=0;i<numberColumnBasic;i++) {
664        int iColumn = whichColumn[i];
665        CoinBigIndex j=iColumn<<1;
666        int iRowM = indices_[j];
667        int iRowP = indices_[j+1];
668        if (iRowM>=0) {
669          indexRowU[numberElements]=iRowM;
670          indexColumnU[numberElements]=numberBasic;
671          elementU[numberElements++]=-1.0;
672        }
673        if (iRowP>=0) {
674          indexRowU[numberElements]=iRowP;
675          indexColumnU[numberElements]=numberBasic;
676          elementU[numberElements++]=1.0;
677        }
678        numberBasic++;
679      }
680    }
681  } else {
682    if (trueNetwork_) {
683      numberElements = 2*numberColumnBasic;
684    } else {
685      for (i=0;i<numberColumnBasic;i++) {
686        int iColumn = whichColumn[i];
687        CoinBigIndex j=iColumn<<1;
688        int iRowM = indices_[j];
689        int iRowP = indices_[j+1];
690        if (iRowM>=0)
691          numberElements ++;
692        if (iRowP>=0)
693          numberElements ++;
694      }
695    }
696  }
697  return numberElements;
698}
699/* Unpacks a column into an CoinIndexedvector
700 */
701void 
702ClpNetworkMatrix::unpack(const ClpSimplex * model,CoinIndexedVector * rowArray,
703                   int iColumn) const 
704{
705  CoinBigIndex j=iColumn<<1;
706  int iRowM = indices_[j];
707  int iRowP = indices_[j+1];
708  if (iRowM>=0) 
709    rowArray->add(iRowM,-1.0);
710  if (iRowP>=0) 
711    rowArray->add(iRowP,1.0);
712}
713/* Unpacks a column into an CoinIndexedvector
714** in packed foramt
715Note that model is NOT const.  Bounds and objective could
716be modified if doing column generation (just for this variable) */
717void 
718ClpNetworkMatrix::unpackPacked(ClpSimplex * model,
719                            CoinIndexedVector * rowArray,
720                            int iColumn) const
721{
722  int * index = rowArray->getIndices();
723  double * array = rowArray->denseVector();
724  int number = 0;
725  CoinBigIndex j=iColumn<<1;
726  int iRowM = indices_[j];
727  int iRowP = indices_[j+1];
728  if (iRowM>=0) {
729    array[number]=-1.0;
730    index[number++]=iRowM;
731  }
732  if (iRowP>=0) {
733    array[number]=1.0;
734    index[number++]=iRowP;
735  }
736  rowArray->setNumElements(number);
737  rowArray->setPackedMode(true);
738}
739/* Adds multiple of a column into an CoinIndexedvector
740      You can use quickAdd to add to vector */
741void 
742ClpNetworkMatrix::add(const ClpSimplex * model,CoinIndexedVector * rowArray,
743                   int iColumn, double multiplier) const 
744{
745  CoinBigIndex j=iColumn<<1;
746  int iRowM = indices_[j];
747  int iRowP = indices_[j+1];
748  if (iRowM>=0) 
749    rowArray->quickAdd(iRowM,-multiplier);
750  if (iRowP>=0) 
751    rowArray->quickAdd(iRowP,multiplier);
752}
753
754// Return a complete CoinPackedMatrix
755CoinPackedMatrix * 
756ClpNetworkMatrix::getPackedMatrix() const 
757{
758  return new CoinPackedMatrix(true,numberRows_,numberColumns_,
759                              2*numberColumns_,
760                              getElements(),indices_,
761                              getVectorStarts(),getVectorLengths());
762
763}
764/* A vector containing the elements in the packed matrix. Note that there
765   might be gaps in this list, entries that do not belong to any
766   major-dimension vector. To get the actual elements one should look at
767   this vector together with vectorStarts and vectorLengths. */
768const double * 
769ClpNetworkMatrix::getElements() const 
770{
771  assert (trueNetwork_); // fix later
772  if (!elements_) {
773    elements_ = new double [2*numberColumns_];
774    int i;
775    for (i=0;i<2*numberColumns_;i+=2) {
776      elements_[i]=-1.0;
777      elements_[i+1]=1.0;
778    }
779  }
780  return elements_;
781}
782
783const CoinBigIndex * 
784ClpNetworkMatrix::getVectorStarts() const 
785{
786  assert (trueNetwork_); // fix later
787  if (!starts_) {
788    starts_ = new int [numberColumns_+1];
789    int i;
790    for (i=0;i<numberColumns_+1;i++) {
791      starts_[i]=i;
792    }
793  }
794  return starts_;
795}
796/* The lengths of the major-dimension vectors. */
797const int * 
798ClpNetworkMatrix::getVectorLengths() const
799{
800  assert (trueNetwork_); // fix later
801  if (!lengths_) {
802    lengths_ = new int [numberColumns_];
803    int i;
804    for (i=0;i<numberColumns_;i++) {
805      lengths_[i]=2;
806    }
807  }
808  return lengths_;
809}
810/* Delete the columns whose indices are listed in <code>indDel</code>. */
811void ClpNetworkMatrix::deleteCols(const int numDel, const int * indDel) 
812{
813  abort();
814}
815/* Delete the rows whose indices are listed in <code>indDel</code>. */
816void ClpNetworkMatrix::deleteRows(const int numDel, const int * indDel) 
817{
818  abort();
819}
Note: See TracBrowser for help on using the repository browser.