1 /* $Id: CoinIndexedVector.cpp 2156 2019-08-06 20:33:03Z stefan $ */
2 // Copyright (C) 2000, International Business Machines
3 // Corporation and others.  All Rights Reserved.
4 // This code is licensed under the terms of the Eclipse Public License (EPL).
5 
6 #if defined(_MSC_VER)
7 // Turn off compiler warning about long names
8 #pragma warning(disable : 4786)
9 #endif
10 
11 #include <cassert>
12 #include <cstdio>
13 
14 #include "CoinTypes.hpp"
15 #include "CoinFloatEqual.hpp"
16 #include "CoinHelperFunctions.hpp"
17 #include "CoinIndexedVector.hpp"
18 #include "CoinTypes.hpp"
19 //#############################################################################
20 #define WARN_USELESS 0
clear()21 void CoinIndexedVector::clear()
22 {
23 #if WARN_USELESS
24   int nNonZero = 0;
25   for (int i = 0; i < capacity_; i++) {
26     if (elements_[i])
27       nNonZero++;
28   }
29   assert(nNonZero <= nElements_);
30 #if WARN_USELESS > 1
31   if (nNonZero != nElements_)
32     printf("Vector said it had %d nonzeros - only had %d\n",
33       nElements_, nNonZero);
34 #endif
35   if (!nNonZero && nElements_)
36     printf("Vector said it had %d nonzeros - but is already empty\n",
37       nElements_);
38 #endif
39   assert(nElements_ <= capacity_);
40   if (!packedMode_) {
41 #ifndef NDEBUG
42     for (int i = 0; i < nElements_; i++)
43       assert(indices_[i] >= 0 && indices_[i] < capacity_);
44 #endif
45     if (3 * nElements_ < capacity_) {
46       int i = 0;
47       if ((nElements_ & 1) != 0) {
48         elements_[indices_[0]] = 0.0;
49         i = 1;
50       }
51       for (; i < nElements_; i += 2) {
52         int i0 = indices_[i];
53         int i1 = indices_[i + 1];
54         elements_[i0] = 0.0;
55         elements_[i1] = 0.0;
56       }
57     } else {
58       CoinZeroN(elements_, capacity_);
59     }
60   } else {
61     CoinZeroN(elements_, nElements_);
62   }
63   nElements_ = 0;
64   packedMode_ = false;
65   //checkClear();
66 }
reallyClear()67 void CoinIndexedVector::reallyClear()
68 {
69   CoinZeroN(elements_, capacity_);
70   nElements_ = 0;
71   packedMode_ = false;
72 }
73 //#############################################################################
74 
empty()75 void CoinIndexedVector::empty()
76 {
77   delete[] indices_;
78   indices_ = NULL;
79   if (elements_)
80     delete[](elements_ - offset_);
81   elements_ = NULL;
82   nElements_ = 0;
83   capacity_ = 0;
84   packedMode_ = false;
85 }
86 
87 //#############################################################################
88 
89 CoinIndexedVector &
operator =(const CoinIndexedVector & rhs)90 CoinIndexedVector::operator=(const CoinIndexedVector &rhs)
91 {
92   if (this != &rhs) {
93     clear();
94     packedMode_ = rhs.packedMode_;
95     if (!packedMode_)
96       gutsOfSetVector(rhs.capacity_, rhs.nElements_,
97         rhs.indices_, rhs.elements_);
98     else
99       gutsOfSetPackedVector(rhs.capacity_, rhs.nElements_,
100         rhs.indices_, rhs.elements_);
101   }
102   return *this;
103 }
104 /* Copy the contents of one vector into another.  If multiplier is 1
105    It is the equivalent of = but if vectors are same size does
106    not re-allocate memory just clears and copies */
copy(const CoinIndexedVector & rhs,double multiplier)107 void CoinIndexedVector::copy(const CoinIndexedVector &rhs, double multiplier)
108 {
109   if (capacity_ == rhs.capacity_) {
110     // can do fast
111     clear();
112     packedMode_ = rhs.packedMode_;
113     nElements_ = 0;
114     if (!packedMode_) {
115       for (int i = 0; i < rhs.nElements_; i++) {
116         int index = rhs.indices_[i];
117         double value = rhs.elements_[index] * multiplier;
118         if (fabs(value) < COIN_INDEXED_TINY_ELEMENT)
119           value = COIN_INDEXED_REALLY_TINY_ELEMENT;
120         elements_[index] = value;
121         indices_[nElements_++] = index;
122       }
123     } else {
124       for (int i = 0; i < rhs.nElements_; i++) {
125         int index = rhs.indices_[i];
126         double value = rhs.elements_[i] * multiplier;
127         if (fabs(value) < COIN_INDEXED_TINY_ELEMENT)
128           value = COIN_INDEXED_REALLY_TINY_ELEMENT;
129         elements_[nElements_] = value;
130         indices_[nElements_++] = index;
131       }
132     }
133   } else {
134     // do as two operations
135     *this = rhs;
136     (*this) *= multiplier;
137   }
138 }
139 
140 //#############################################################################
141 #ifndef CLP_NO_VECTOR
142 CoinIndexedVector &
operator =(const CoinPackedVectorBase & rhs)143 CoinIndexedVector::operator=(const CoinPackedVectorBase &rhs)
144 {
145   clear();
146   packedMode_ = false;
147   gutsOfSetVector(rhs.getNumElements(),
148     rhs.getIndices(), rhs.getElements());
149   return *this;
150 }
151 #endif
152 //#############################################################################
153 
borrowVector(int size,int numberIndices,int * inds,double * elems)154 void CoinIndexedVector::borrowVector(int size, int numberIndices, int *inds, double *elems)
155 {
156   empty();
157   capacity_ = size;
158   nElements_ = numberIndices;
159   indices_ = inds;
160   elements_ = elems;
161 
162   // whole point about borrowvector is that it is lightweight so no testing is done
163 }
164 
165 //#############################################################################
166 
returnVector()167 void CoinIndexedVector::returnVector()
168 {
169   indices_ = NULL;
170   elements_ = NULL;
171   nElements_ = 0;
172   capacity_ = 0;
173   packedMode_ = false;
174 }
175 
176 //#############################################################################
177 
setVector(int size,const int * inds,const double * elems)178 void CoinIndexedVector::setVector(int size, const int *inds, const double *elems)
179 {
180   clear();
181   gutsOfSetVector(size, inds, elems);
182 }
183 //#############################################################################
184 
setVector(int size,int numberIndices,const int * inds,const double * elems)185 void CoinIndexedVector::setVector(int size, int numberIndices, const int *inds, const double *elems)
186 {
187   clear();
188   gutsOfSetVector(size, numberIndices, inds, elems);
189 }
190 //#############################################################################
191 
setConstant(int size,const int * inds,double value)192 void CoinIndexedVector::setConstant(int size, const int *inds, double value)
193 {
194   clear();
195   gutsOfSetConstant(size, inds, value);
196 }
197 
198 //#############################################################################
199 
setFull(int size,const double * elems)200 void CoinIndexedVector::setFull(int size, const double *elems)
201 {
202   // Clear out any values presently stored
203   clear();
204 
205 #ifndef COIN_FAST_CODE
206   if (size < 0)
207     throw CoinError("negative number of indices", "setFull", "CoinIndexedVector");
208 #endif
209   reserve(size);
210   nElements_ = 0;
211   // elements_ array is all zero
212   int i;
213   for (i = 0; i < size; i++) {
214     int indexValue = i;
215     if (fabs(elems[i]) >= COIN_INDEXED_TINY_ELEMENT) {
216       elements_[indexValue] = elems[i];
217       indices_[nElements_++] = indexValue;
218     }
219   }
220 }
221 //#############################################################################
222 
223 /** Access the i'th element of the full storage vector.  */
224 double &
operator [](int index) const225   CoinIndexedVector::operator[](int index) const
226 {
227   assert(!packedMode_);
228 #ifndef COIN_FAST_CODE
229   if (index >= capacity_)
230     throw CoinError("index >= capacity()", "[]", "CoinIndexedVector");
231   if (index < 0)
232     throw CoinError("index < 0", "[]", "CoinIndexedVector");
233 #endif
234   double *where = elements_ + index;
235   return *where;
236 }
237 //#############################################################################
238 
setElement(int index,double element)239 void CoinIndexedVector::setElement(int index, double element)
240 {
241 #ifndef COIN_FAST_CODE
242   if (index >= nElements_)
243     throw CoinError("index >= size()", "setElement", "CoinIndexedVector");
244   if (index < 0)
245     throw CoinError("index < 0", "setElement", "CoinIndexedVector");
246 #endif
247   elements_[indices_[index]] = element;
248 }
249 
250 //#############################################################################
251 
insert(int index,double element)252 void CoinIndexedVector::insert(int index, double element)
253 {
254 #ifndef COIN_FAST_CODE
255   if (index < 0)
256     throw CoinError("index < 0", "setElement", "CoinIndexedVector");
257 #endif
258   if (index >= capacity_)
259     reserve(index + 1);
260 #ifndef COIN_FAST_CODE
261   if (elements_[index])
262     throw CoinError("Index already exists", "insert", "CoinIndexedVector");
263 #endif
264   indices_[nElements_++] = index;
265   elements_[index] = element;
266 }
267 
268 //#############################################################################
269 
add(int index,double element)270 void CoinIndexedVector::add(int index, double element)
271 {
272 #ifndef COIN_FAST_CODE
273   if (index < 0)
274     throw CoinError("index < 0", "setElement", "CoinIndexedVector");
275 #endif
276   if (index >= capacity_)
277     reserve(index + 1);
278   if (elements_[index]) {
279     element += elements_[index];
280     if (fabs(element) >= COIN_INDEXED_TINY_ELEMENT) {
281       elements_[index] = element;
282     } else {
283       elements_[index] = COIN_INDEXED_REALLY_TINY_ELEMENT;
284     }
285   } else if (fabs(element) >= COIN_INDEXED_TINY_ELEMENT) {
286     indices_[nElements_++] = index;
287     assert(nElements_ <= capacity_);
288     elements_[index] = element;
289   }
290 }
291 
292 //#############################################################################
293 
clean(double tolerance)294 int CoinIndexedVector::clean(double tolerance)
295 {
296   int number = nElements_;
297   int i;
298   nElements_ = 0;
299   assert(!packedMode_);
300   for (i = 0; i < number; i++) {
301     int indexValue = indices_[i];
302     if (fabs(elements_[indexValue]) >= tolerance) {
303       indices_[nElements_++] = indexValue;
304     } else {
305       elements_[indexValue] = 0.0;
306     }
307   }
308   return nElements_;
309 }
310 #ifndef NDEBUG
311 //#############################################################################
312 // For debug check vector is clear i.e. no elements
checkClear()313 void CoinIndexedVector::checkClear()
314 {
315 #ifndef NDEBUG
316   //printf("checkClear %p\n",this);
317   assert(!nElements_);
318   //assert(!packedMode_);
319   int i;
320   for (i = 0; i < capacity_; i++) {
321     assert(!elements_[i]);
322   }
323   // check mark array zeroed
324   char *mark = reinterpret_cast< char * >(indices_ + capacity_);
325   for (i = 0; i < capacity_; i++) {
326     assert(!mark[i]);
327   }
328 #else
329   if (nElements_) {
330     printf("%d nElements_ - checkClear\n", nElements_);
331     abort();
332   }
333   if (packedMode_) {
334     printf("packed mode when empty - checkClear\n");
335     abort();
336   }
337   int i;
338   int n = 0;
339   int k = -1;
340   for (i = 0; i < capacity_; i++) {
341     if (elements_[i]) {
342       n++;
343       if (k < 0)
344         k = i;
345     }
346   }
347   if (n) {
348     printf("%d elements, first %d - checkClear\n", n, k);
349     abort();
350   }
351 #endif
352 }
353 // For debug check vector is clean i.e. elements match indices
checkClean()354 void CoinIndexedVector::checkClean()
355 {
356   //printf("checkClean %p\n",this);
357   int i;
358   if (packedMode_) {
359     for (i = 0; i < nElements_; i++)
360       assert(elements_[i]);
361     for (; i < capacity_; i++)
362       assert(!elements_[i]);
363   } else {
364     double *copy = new double[capacity_];
365     CoinMemcpyN(elements_, capacity_, copy);
366     for (i = 0; i < nElements_; i++) {
367       int indexValue = indices_[i];
368       assert(copy[indexValue]);
369       copy[indexValue] = 0.0;
370     }
371     for (i = 0; i < capacity_; i++)
372       assert(!copy[i]);
373     delete[] copy;
374   }
375 #ifndef NDEBUG
376   // check mark array zeroed
377   char *mark = reinterpret_cast< char * >(indices_ + capacity_);
378   for (i = 0; i < capacity_; i++) {
379     assert(!mark[i]);
380   }
381 #endif
382 }
383 #endif
384 //#############################################################################
385 #ifndef CLP_NO_VECTOR
append(const CoinPackedVectorBase & caboose)386 void CoinIndexedVector::append(const CoinPackedVectorBase &caboose)
387 {
388   const int cs = caboose.getNumElements();
389 
390   const int *cind = caboose.getIndices();
391   const double *celem = caboose.getElements();
392   int maxIndex = -1;
393   int i;
394   for (i = 0; i < cs; i++) {
395     int indexValue = cind[i];
396     if (indexValue < 0)
397       throw CoinError("negative index", "append", "CoinIndexedVector");
398     if (maxIndex < indexValue)
399       maxIndex = indexValue;
400   }
401   reserve(maxIndex + 1);
402   bool needClean = false;
403   int numberDuplicates = 0;
404   for (i = 0; i < cs; i++) {
405     int indexValue = cind[i];
406     if (elements_[indexValue]) {
407       numberDuplicates++;
408       elements_[indexValue] += celem[i];
409       if (fabs(elements_[indexValue]) < COIN_INDEXED_TINY_ELEMENT)
410         needClean = true; // need to go through again
411     } else {
412       if (fabs(celem[i]) >= COIN_INDEXED_TINY_ELEMENT) {
413         elements_[indexValue] = celem[i];
414         indices_[nElements_++] = indexValue;
415       }
416     }
417   }
418   if (needClean) {
419     // go through again
420     int size = nElements_;
421     nElements_ = 0;
422     for (i = 0; i < size; i++) {
423       int indexValue = indices_[i];
424       double value = elements_[indexValue];
425       if (fabs(value) >= COIN_INDEXED_TINY_ELEMENT) {
426         indices_[nElements_++] = indexValue;
427       } else {
428         elements_[indexValue] = 0.0;
429       }
430     }
431   }
432   if (numberDuplicates)
433     throw CoinError("duplicate index", "append", "CoinIndexedVector");
434 }
435 #endif
436 //#############################################################################
437 
swap(int i,int j)438 void CoinIndexedVector::swap(int i, int j)
439 {
440   if (i >= nElements_)
441     throw CoinError("index i >= size()", "swap", "CoinIndexedVector");
442   if (i < 0)
443     throw CoinError("index i < 0", "swap", "CoinIndexedVector");
444   if (j >= nElements_)
445     throw CoinError("index j >= size()", "swap", "CoinIndexedVector");
446   if (j < 0)
447     throw CoinError("index j < 0", "swap", "CoinIndexedVector");
448 
449   // Swap positions i and j of the
450   // indices array
451 
452   int isave = indices_[i];
453   indices_[i] = indices_[j];
454   indices_[j] = isave;
455 }
456 
457 //#############################################################################
458 
truncate(int n)459 void CoinIndexedVector::truncate(int n)
460 {
461   reserve(n);
462 }
463 
464 //#############################################################################
465 
operator +=(double value)466 void CoinIndexedVector::operator+=(double value)
467 {
468   assert(!packedMode_);
469   int i, indexValue;
470   for (i = 0; i < nElements_; i++) {
471     indexValue = indices_[i];
472     double newValue = elements_[indexValue] + value;
473     if (fabs(newValue) >= COIN_INDEXED_TINY_ELEMENT)
474       elements_[indexValue] = newValue;
475     else
476       elements_[indexValue] = COIN_INDEXED_REALLY_TINY_ELEMENT;
477   }
478 }
479 
480 //-----------------------------------------------------------------------------
481 
operator -=(double value)482 void CoinIndexedVector::operator-=(double value)
483 {
484   assert(!packedMode_);
485   int i, indexValue;
486   for (i = 0; i < nElements_; i++) {
487     indexValue = indices_[i];
488     double newValue = elements_[indexValue] - value;
489     if (fabs(newValue) >= COIN_INDEXED_TINY_ELEMENT)
490       elements_[indexValue] = newValue;
491     else
492       elements_[indexValue] = COIN_INDEXED_REALLY_TINY_ELEMENT;
493   }
494 }
495 
496 //-----------------------------------------------------------------------------
497 
operator *=(double value)498 void CoinIndexedVector::operator*=(double value)
499 {
500   assert(!packedMode_);
501   int i, indexValue;
502   for (i = 0; i < nElements_; i++) {
503     indexValue = indices_[i];
504     double newValue = elements_[indexValue] * value;
505     if (fabs(newValue) >= COIN_INDEXED_TINY_ELEMENT)
506       elements_[indexValue] = newValue;
507     else
508       elements_[indexValue] = COIN_INDEXED_REALLY_TINY_ELEMENT;
509   }
510 }
511 
512 //-----------------------------------------------------------------------------
513 
operator /=(double value)514 void CoinIndexedVector::operator/=(double value)
515 {
516   assert(!packedMode_);
517   int i, indexValue;
518   for (i = 0; i < nElements_; i++) {
519     indexValue = indices_[i];
520     double newValue = elements_[indexValue] / value;
521     if (fabs(newValue) >= COIN_INDEXED_TINY_ELEMENT)
522       elements_[indexValue] = newValue;
523     else
524       elements_[indexValue] = COIN_INDEXED_REALLY_TINY_ELEMENT;
525   }
526 }
527 //#############################################################################
528 
reserve(int n)529 void CoinIndexedVector::reserve(int n)
530 {
531   int i;
532   int nPlus;
533   if (sizeof(int) == 4 * sizeof(char))
534     nPlus = (n + 3) >> 2;
535   else
536     nPlus = (n + 7) >> 4;
537 #ifdef COIN_AVX2
538   nPlus += 16;
539 #endif
540   // don't make allocated space smaller but do take off values
541   if (n + nPlus < capacity_) {
542 #ifndef COIN_FAST_CODE
543     if (n < 0)
544       throw CoinError("negative capacity", "reserve", "CoinIndexedVector");
545 #endif
546     int nNew = 0;
547     for (i = 0; i < nElements_; i++) {
548       int indexValue = indices_[i];
549       if (indexValue < n) {
550         indices_[nNew++] = indexValue;
551       } else {
552         elements_[indexValue] = 0.0;
553       }
554     }
555     nElements_ = nNew;
556   } else if (n > capacity_) {
557 
558     // save pointers to existing data
559     int *tempIndices = indices_;
560     double *tempElements = elements_;
561     double *delTemp = elements_ - offset_;
562 
563     // allocate new space
564     indices_ = new int[n + nPlus];
565     CoinZeroN(indices_ + n, nPlus);
566     // align on 64 byte boundary
567     double *temp = new double[n + 9 + nPlus];
568     offset_ = 0;
569     CoinInt64 xx = reinterpret_cast< CoinInt64 >(temp);
570     int iBottom = static_cast< int >(xx & 63);
571     //if (iBottom)
572     offset_ = (64 - iBottom) >> 3;
573     elements_ = temp + offset_;
574     ;
575 
576     // copy data to new space
577     // and zero out part of array
578     if (nElements_ > 0) {
579       CoinMemcpyN(tempIndices, nElements_, indices_);
580       CoinMemcpyN(tempElements, capacity_, elements_);
581       CoinZeroN(elements_ + capacity_, n - capacity_);
582     } else {
583       CoinZeroN(elements_, n);
584     }
585     capacity_ = n;
586 
587     // free old data
588     if (tempElements)
589       delete[] delTemp;
590     delete[] tempIndices;
591   }
592 }
593 
594 //#############################################################################
595 
CoinIndexedVector()596 CoinIndexedVector::CoinIndexedVector()
597   : indices_(NULL)
598   , elements_(NULL)
599   , nElements_(0)
600   , capacity_(0)
601   , offset_(0)
602   , packedMode_(false)
603 {
604 }
605 
CoinIndexedVector(int size)606 CoinIndexedVector::CoinIndexedVector(int size)
607   : indices_(NULL)
608   , elements_(NULL)
609   , nElements_(0)
610   , capacity_(0)
611   , offset_(0)
612   , packedMode_(false)
613 {
614   // Get space
615   reserve(size);
616 }
617 
618 //-----------------------------------------------------------------------------
619 
CoinIndexedVector(int size,const int * inds,const double * elems)620 CoinIndexedVector::CoinIndexedVector(int size,
621   const int *inds, const double *elems)
622   : indices_(NULL)
623   , elements_(NULL)
624   , nElements_(0)
625   , capacity_(0)
626   , offset_(0)
627   , packedMode_(false)
628 {
629   gutsOfSetVector(size, inds, elems);
630 }
631 
632 //-----------------------------------------------------------------------------
633 
CoinIndexedVector(int size,const int * inds,double value)634 CoinIndexedVector::CoinIndexedVector(int size,
635   const int *inds, double value)
636   : indices_(NULL)
637   , elements_(NULL)
638   , nElements_(0)
639   , capacity_(0)
640   , offset_(0)
641   , packedMode_(false)
642 {
643   gutsOfSetConstant(size, inds, value);
644 }
645 
646 //-----------------------------------------------------------------------------
647 
CoinIndexedVector(int size,const double * element)648 CoinIndexedVector::CoinIndexedVector(int size, const double *element)
649   : indices_(NULL)
650   , elements_(NULL)
651   , nElements_(0)
652   , capacity_(0)
653   , offset_(0)
654   , packedMode_(false)
655 {
656   setFull(size, element);
657 }
658 
659 //-----------------------------------------------------------------------------
660 #ifndef CLP_NO_VECTOR
CoinIndexedVector(const CoinPackedVectorBase & rhs)661 CoinIndexedVector::CoinIndexedVector(const CoinPackedVectorBase &rhs)
662   : indices_(NULL)
663   , elements_(NULL)
664   , nElements_(0)
665   , capacity_(0)
666   , offset_(0)
667   , packedMode_(false)
668 {
669   gutsOfSetVector(rhs.getNumElements(),
670     rhs.getIndices(), rhs.getElements());
671 }
672 #endif
673 //-----------------------------------------------------------------------------
674 
CoinIndexedVector(const CoinIndexedVector & rhs)675 CoinIndexedVector::CoinIndexedVector(const CoinIndexedVector &rhs)
676   : indices_(NULL)
677   , elements_(NULL)
678   , nElements_(0)
679   , capacity_(0)
680   , offset_(0)
681   , packedMode_(false)
682 {
683   if (!rhs.packedMode_)
684     gutsOfSetVector(rhs.capacity_, rhs.nElements_, rhs.indices_, rhs.elements_);
685   else
686     gutsOfSetPackedVector(rhs.capacity_, rhs.nElements_, rhs.indices_, rhs.elements_);
687 }
688 
689 //-----------------------------------------------------------------------------
690 
CoinIndexedVector(const CoinIndexedVector * rhs)691 CoinIndexedVector::CoinIndexedVector(const CoinIndexedVector *rhs)
692   : indices_(NULL)
693   , elements_(NULL)
694   , nElements_(0)
695   , capacity_(0)
696   , offset_(0)
697   , packedMode_(false)
698 {
699   if (!rhs->packedMode_)
700     gutsOfSetVector(rhs->capacity_, rhs->nElements_, rhs->indices_, rhs->elements_);
701   else
702     gutsOfSetPackedVector(rhs->capacity_, rhs->nElements_, rhs->indices_, rhs->elements_);
703 }
704 
705 //-----------------------------------------------------------------------------
706 
~CoinIndexedVector()707 CoinIndexedVector::~CoinIndexedVector()
708 {
709   delete[] indices_;
710   if (elements_)
711     delete[](elements_ - offset_);
712 }
713 //#############################################################################
714 //#############################################################################
715 
716 /// Return the sum of two indexed vectors
717 CoinIndexedVector
operator +(const CoinIndexedVector & op2)718 CoinIndexedVector::operator+(
719   const CoinIndexedVector &op2)
720 {
721   assert(!packedMode_);
722   int i;
723   int nElements = nElements_;
724   int capacity = CoinMax(capacity_, op2.capacity_);
725   CoinIndexedVector newOne(*this);
726   newOne.reserve(capacity);
727   bool needClean = false;
728   // new one now can hold everything so just modify old and add new
729   for (i = 0; i < op2.nElements_; i++) {
730     int indexValue = op2.indices_[i];
731     double value = op2.elements_[indexValue];
732     double oldValue = elements_[indexValue];
733     if (!oldValue) {
734       if (fabs(value) >= COIN_INDEXED_TINY_ELEMENT) {
735         newOne.elements_[indexValue] = value;
736         newOne.indices_[nElements++] = indexValue;
737       }
738     } else {
739       value += oldValue;
740       newOne.elements_[indexValue] = value;
741       if (fabs(value) < COIN_INDEXED_TINY_ELEMENT) {
742         needClean = true;
743       }
744     }
745   }
746   newOne.nElements_ = nElements;
747   if (needClean) {
748     // go through again
749     newOne.nElements_ = 0;
750     for (i = 0; i < nElements; i++) {
751       int indexValue = newOne.indices_[i];
752       double value = newOne.elements_[indexValue];
753       if (fabs(value) >= COIN_INDEXED_TINY_ELEMENT) {
754         newOne.indices_[newOne.nElements_++] = indexValue;
755       } else {
756         newOne.elements_[indexValue] = 0.0;
757       }
758     }
759   }
760   return newOne;
761 }
762 
763 /// Return the difference of two indexed vectors
764 CoinIndexedVector
operator -(const CoinIndexedVector & op2)765 CoinIndexedVector::operator-(
766   const CoinIndexedVector &op2)
767 {
768   assert(!packedMode_);
769   int i;
770   int nElements = nElements_;
771   int capacity = CoinMax(capacity_, op2.capacity_);
772   CoinIndexedVector newOne(*this);
773   newOne.reserve(capacity);
774   bool needClean = false;
775   // new one now can hold everything so just modify old and add new
776   for (i = 0; i < op2.nElements_; i++) {
777     int indexValue = op2.indices_[i];
778     double value = op2.elements_[indexValue];
779     double oldValue = elements_[indexValue];
780     if (!oldValue) {
781       if (fabs(value) >= COIN_INDEXED_TINY_ELEMENT) {
782         newOne.elements_[indexValue] = -value;
783         newOne.indices_[nElements++] = indexValue;
784       }
785     } else {
786       value = oldValue - value;
787       newOne.elements_[indexValue] = value;
788       if (fabs(value) < COIN_INDEXED_TINY_ELEMENT) {
789         needClean = true;
790       }
791     }
792   }
793   newOne.nElements_ = nElements;
794   if (needClean) {
795     // go through again
796     newOne.nElements_ = 0;
797     for (i = 0; i < nElements; i++) {
798       int indexValue = newOne.indices_[i];
799       double value = newOne.elements_[indexValue];
800       if (fabs(value) >= COIN_INDEXED_TINY_ELEMENT) {
801         newOne.indices_[newOne.nElements_++] = indexValue;
802       } else {
803         newOne.elements_[indexValue] = 0.0;
804       }
805     }
806   }
807   return newOne;
808 }
809 
810 /// Return the element-wise product of two indexed vectors
811 CoinIndexedVector
operator *(const CoinIndexedVector & op2)812   CoinIndexedVector::operator*(
813     const CoinIndexedVector &op2)
814 {
815   assert(!packedMode_);
816   int i;
817   int nElements = nElements_;
818   int capacity = CoinMax(capacity_, op2.capacity_);
819   CoinIndexedVector newOne(*this);
820   newOne.reserve(capacity);
821   bool needClean = false;
822   // new one now can hold everything so just modify old and add new
823   for (i = 0; i < op2.nElements_; i++) {
824     int indexValue = op2.indices_[i];
825     double value = op2.elements_[indexValue];
826     double oldValue = elements_[indexValue];
827     if (oldValue) {
828       value *= oldValue;
829       newOne.elements_[indexValue] = value;
830       if (fabs(value) < COIN_INDEXED_TINY_ELEMENT) {
831         needClean = true;
832       }
833     }
834   }
835 
836   newOne.nElements_ = nElements;
837 
838   if (needClean) {
839     // go through again
840     newOne.nElements_ = 0;
841     for (i = 0; i < nElements; i++) {
842       int indexValue = newOne.indices_[i];
843       double value = newOne.elements_[indexValue];
844       if (fabs(value) >= COIN_INDEXED_TINY_ELEMENT) {
845         newOne.indices_[newOne.nElements_++] = indexValue;
846       } else {
847         newOne.elements_[indexValue] = 0.0;
848       }
849     }
850   }
851   return newOne;
852 }
853 
854 /// Return the element-wise ratio of two indexed vectors
855 CoinIndexedVector
operator /(const CoinIndexedVector & op2)856 CoinIndexedVector::operator/(const CoinIndexedVector &op2)
857 {
858   assert(!packedMode_);
859   // I am treating 0.0/0.0 as 0.0
860   int i;
861   int nElements = nElements_;
862   int capacity = CoinMax(capacity_, op2.capacity_);
863   CoinIndexedVector newOne(*this);
864   newOne.reserve(capacity);
865   bool needClean = false;
866   // new one now can hold everything so just modify old and add new
867   for (i = 0; i < op2.nElements_; i++) {
868     int indexValue = op2.indices_[i];
869     double value = op2.elements_[indexValue];
870     double oldValue = elements_[indexValue];
871     if (oldValue) {
872       if (!value)
873         throw CoinError("zero divisor", "/", "CoinIndexedVector");
874       value = oldValue / value;
875       newOne.elements_[indexValue] = value;
876       if (fabs(value) < COIN_INDEXED_TINY_ELEMENT) {
877         needClean = true;
878       }
879     }
880   }
881 
882   newOne.nElements_ = nElements;
883 
884   if (needClean) {
885     // go through again
886     newOne.nElements_ = 0;
887     for (i = 0; i < nElements; i++) {
888       int indexValue = newOne.indices_[i];
889       double value = newOne.elements_[indexValue];
890       if (fabs(value) >= COIN_INDEXED_TINY_ELEMENT) {
891         newOne.indices_[newOne.nElements_++] = indexValue;
892       } else {
893         newOne.elements_[indexValue] = 0.0;
894       }
895     }
896   }
897   return newOne;
898 }
899 // The sum of two indexed vectors
operator +=(const CoinIndexedVector & op2)900 void CoinIndexedVector::operator+=(const CoinIndexedVector &op2)
901 {
902   // do slowly at first
903   *this = *this + op2;
904 }
905 
906 // The difference of two indexed vectors
operator -=(const CoinIndexedVector & op2)907 void CoinIndexedVector::operator-=(const CoinIndexedVector &op2)
908 {
909   // do slowly at first
910   *this = *this - op2;
911 }
912 
913 // The element-wise product of two indexed vectors
operator *=(const CoinIndexedVector & op2)914 void CoinIndexedVector::operator*=(const CoinIndexedVector &op2)
915 {
916   // do slowly at first
917   *this = *this * op2;
918 }
919 
920 // The element-wise ratio of two indexed vectors (0.0/0.0 => 0.0) (0 vanishes)
operator /=(const CoinIndexedVector & op2)921 void CoinIndexedVector::operator/=(const CoinIndexedVector &op2)
922 {
923   // do slowly at first
924   *this = *this / op2;
925 }
926 //#############################################################################
sortDecrIndex()927 void CoinIndexedVector::sortDecrIndex()
928 {
929   // Should replace with std sort
930   double *elements = new double[nElements_];
931   CoinZeroN(elements, nElements_);
932   CoinSort_2(indices_, indices_ + nElements_, elements,
933     CoinFirstGreater_2< int, double >());
934   delete[] elements;
935 }
936 
sortPacked()937 void CoinIndexedVector::sortPacked()
938 {
939   assert(packedMode_);
940   CoinSort_2(indices_, indices_ + nElements_, elements_);
941 }
942 
sortIncrElement()943 void CoinIndexedVector::sortIncrElement()
944 {
945   double *elements = new double[nElements_];
946   int i;
947   for (i = 0; i < nElements_; i++)
948     elements[i] = elements_[indices_[i]];
949   CoinSort_2(elements, elements + nElements_, indices_,
950     CoinFirstLess_2< double, int >());
951   delete[] elements;
952 }
953 
sortDecrElement()954 void CoinIndexedVector::sortDecrElement()
955 {
956   double *elements = new double[nElements_];
957   int i;
958   for (i = 0; i < nElements_; i++)
959     elements[i] = elements_[indices_[i]];
960   CoinSort_2(elements, elements + nElements_, indices_,
961     CoinFirstGreater_2< double, int >());
962   delete[] elements;
963 }
964 //#############################################################################
965 
gutsOfSetVector(int size,const int * inds,const double * elems)966 void CoinIndexedVector::gutsOfSetVector(int size,
967   const int *inds, const double *elems)
968 {
969 #ifndef COIN_FAST_CODE
970   if (size < 0)
971     throw CoinError("negative number of indices", "setVector", "CoinIndexedVector");
972 #endif
973   assert(!packedMode_);
974   // find largest
975   int i;
976   int maxIndex = -1;
977   for (i = 0; i < size; i++) {
978     int indexValue = inds[i];
979 #ifndef COIN_FAST_CODE
980     if (indexValue < 0)
981       throw CoinError("negative index", "setVector", "CoinIndexedVector");
982 #endif
983     if (maxIndex < indexValue)
984       maxIndex = indexValue;
985   }
986   reserve(maxIndex + 1);
987   nElements_ = 0;
988   // elements_ array is all zero
989   bool needClean = false;
990   int numberDuplicates = 0;
991   for (i = 0; i < size; i++) {
992     int indexValue = inds[i];
993     if (elements_[indexValue] == 0) {
994       if (fabs(elems[i]) >= COIN_INDEXED_TINY_ELEMENT) {
995         indices_[nElements_++] = indexValue;
996         elements_[indexValue] = elems[i];
997       }
998     } else {
999       numberDuplicates++;
1000       elements_[indexValue] += elems[i];
1001       if (fabs(elements_[indexValue]) < COIN_INDEXED_TINY_ELEMENT)
1002         needClean = true; // need to go through again
1003     }
1004   }
1005   if (needClean) {
1006     // go through again
1007     size = nElements_;
1008     nElements_ = 0;
1009     for (i = 0; i < size; i++) {
1010       int indexValue = indices_[i];
1011       double value = elements_[indexValue];
1012       if (fabs(value) >= COIN_INDEXED_TINY_ELEMENT) {
1013         indices_[nElements_++] = indexValue;
1014       } else {
1015         elements_[indexValue] = 0.0;
1016       }
1017     }
1018   }
1019   if (numberDuplicates)
1020     throw CoinError("duplicate index", "setVector", "CoinIndexedVector");
1021 }
1022 
1023 //#############################################################################
1024 
gutsOfSetVector(int size,int numberIndices,const int * inds,const double * elems)1025 void CoinIndexedVector::gutsOfSetVector(int size, int numberIndices,
1026   const int *inds, const double *elems)
1027 {
1028   assert(!packedMode_);
1029 
1030   int i;
1031   reserve(size);
1032 #ifndef COIN_FAST_CODE
1033   if (numberIndices < 0)
1034     throw CoinError("negative number of indices", "setVector", "CoinIndexedVector");
1035 #endif
1036   nElements_ = 0;
1037   // elements_ array is all zero
1038   bool needClean = false;
1039   int numberDuplicates = 0;
1040   for (i = 0; i < numberIndices; i++) {
1041     int indexValue = inds[i];
1042 #ifndef COIN_FAST_CODE
1043     if (indexValue < 0)
1044       throw CoinError("negative index", "setVector", "CoinIndexedVector");
1045     else if (indexValue >= size)
1046       throw CoinError("too large an index", "setVector", "CoinIndexedVector");
1047 #endif
1048     if (elements_[indexValue]) {
1049       numberDuplicates++;
1050       elements_[indexValue] += elems[indexValue];
1051       if (fabs(elements_[indexValue]) < COIN_INDEXED_TINY_ELEMENT)
1052         needClean = true; // need to go through again
1053     } else {
1054 #ifndef COIN_FAC_NEW
1055       if (fabs(elems[indexValue]) >= COIN_INDEXED_TINY_ELEMENT) {
1056 #endif
1057         elements_[indexValue] = elems[indexValue];
1058         indices_[nElements_++] = indexValue;
1059 #ifndef COIN_FAC_NEW
1060       }
1061 #endif
1062     }
1063   }
1064   if (needClean) {
1065     // go through again
1066     size = nElements_;
1067     nElements_ = 0;
1068     for (i = 0; i < size; i++) {
1069       int indexValue = indices_[i];
1070       double value = elements_[indexValue];
1071       if (fabs(value) >= COIN_INDEXED_TINY_ELEMENT) {
1072         indices_[nElements_++] = indexValue;
1073       } else {
1074         elements_[indexValue] = 0.0;
1075       }
1076     }
1077   }
1078   if (numberDuplicates)
1079     throw CoinError("duplicate index", "setVector", "CoinIndexedVector");
1080 }
1081 
1082 //-----------------------------------------------------------------------------
1083 
gutsOfSetPackedVector(int size,int numberIndices,const int * inds,const double * elems)1084 void CoinIndexedVector::gutsOfSetPackedVector(int size, int numberIndices,
1085   const int *inds, const double *elems)
1086 {
1087   packedMode_ = true;
1088 
1089   int i;
1090   reserve(size);
1091 #ifndef COIN_FAST_CODE
1092   if (numberIndices < 0)
1093     throw CoinError("negative number of indices", "setVector", "CoinIndexedVector");
1094 #endif
1095   nElements_ = 0;
1096   // elements_ array is all zero
1097   // Does not check for duplicates
1098   for (i = 0; i < numberIndices; i++) {
1099     int indexValue = inds[i];
1100 #ifndef COIN_FAST_CODE
1101     if (indexValue < 0)
1102       throw CoinError("negative index", "setVector", "CoinIndexedVector");
1103       //else if (indexValue>=size)
1104       //throw CoinError("too large an index", "setVector", "CoinIndexedVector");
1105 #endif
1106     if (fabs(elems[i]) >= COIN_INDEXED_TINY_ELEMENT) {
1107       elements_[nElements_] = elems[i];
1108       indices_[nElements_++] = indexValue;
1109     }
1110   }
1111 }
1112 
1113 //-----------------------------------------------------------------------------
1114 
gutsOfSetConstant(int size,const int * inds,double value)1115 void CoinIndexedVector::gutsOfSetConstant(int size,
1116   const int *inds, double value)
1117 {
1118 
1119   assert(!packedMode_);
1120 #ifndef COIN_FAST_CODE
1121   if (size < 0)
1122     throw CoinError("negative number of indices", "setConstant", "CoinIndexedVector");
1123 #endif
1124   // find largest
1125   int i;
1126   int maxIndex = -1;
1127   for (i = 0; i < size; i++) {
1128     int indexValue = inds[i];
1129 #ifndef COIN_FAST_CODE
1130     if (indexValue < 0)
1131       throw CoinError("negative index", "setConstant", "CoinIndexedVector");
1132 #endif
1133     if (maxIndex < indexValue)
1134       maxIndex = indexValue;
1135   }
1136 
1137   reserve(maxIndex + 1);
1138   nElements_ = 0;
1139   int numberDuplicates = 0;
1140   // elements_ array is all zero
1141   bool needClean = false;
1142   for (i = 0; i < size; i++) {
1143     int indexValue = inds[i];
1144     if (elements_[indexValue] == 0) {
1145       if (fabs(value) >= COIN_INDEXED_TINY_ELEMENT) {
1146         elements_[indexValue] += value;
1147         indices_[nElements_++] = indexValue;
1148       }
1149     } else {
1150       numberDuplicates++;
1151       elements_[indexValue] += value;
1152       if (fabs(elements_[indexValue]) < COIN_INDEXED_TINY_ELEMENT)
1153         needClean = true; // need to go through again
1154     }
1155   }
1156   if (needClean) {
1157     // go through again
1158     size = nElements_;
1159     nElements_ = 0;
1160     for (i = 0; i < size; i++) {
1161       int indexValue = indices_[i];
1162       double value = elements_[indexValue];
1163       if (fabs(value) >= COIN_INDEXED_TINY_ELEMENT) {
1164         indices_[nElements_++] = indexValue;
1165       } else {
1166         elements_[indexValue] = 0.0;
1167       }
1168     }
1169   }
1170   if (numberDuplicates)
1171     throw CoinError("duplicate index", "setConstant", "CoinIndexedVector");
1172 }
1173 
1174 //#############################################################################
1175 // Append a CoinIndexedVector to the end
append(const CoinIndexedVector & caboose)1176 void CoinIndexedVector::append(const CoinIndexedVector &caboose)
1177 {
1178   const int cs = caboose.getNumElements();
1179 
1180   const int *cind = caboose.getIndices();
1181   const double *celem = caboose.denseVector();
1182   int maxIndex = -1;
1183   int i;
1184   for (i = 0; i < cs; i++) {
1185     int indexValue = cind[i];
1186 #ifndef COIN_FAST_CODE
1187     if (indexValue < 0)
1188       throw CoinError("negative index", "append", "CoinIndexedVector");
1189 #endif
1190     if (maxIndex < indexValue)
1191       maxIndex = indexValue;
1192   }
1193   reserve(maxIndex + 1);
1194   bool needClean = false;
1195   int numberDuplicates = 0;
1196   for (i = 0; i < cs; i++) {
1197     int indexValue = cind[i];
1198     if (elements_[indexValue]) {
1199       numberDuplicates++;
1200       elements_[indexValue] += celem[indexValue];
1201       if (fabs(elements_[indexValue]) < COIN_INDEXED_TINY_ELEMENT)
1202         needClean = true; // need to go through again
1203     } else {
1204       if (fabs(celem[indexValue]) >= COIN_INDEXED_TINY_ELEMENT) {
1205         elements_[indexValue] = celem[indexValue];
1206         indices_[nElements_++] = indexValue;
1207       }
1208     }
1209   }
1210   if (needClean) {
1211     // go through again
1212     int size = nElements_;
1213     nElements_ = 0;
1214     for (i = 0; i < size; i++) {
1215       int indexValue = indices_[i];
1216       double value = elements_[indexValue];
1217       if (fabs(value) >= COIN_INDEXED_TINY_ELEMENT) {
1218         indices_[nElements_++] = indexValue;
1219       } else {
1220         elements_[indexValue] = 0.0;
1221       }
1222     }
1223   }
1224   if (numberDuplicates)
1225     throw CoinError("duplicate index", "append", "CoinIndexedVector");
1226 }
1227 // Append a CoinIndexedVector to the end and modify indices
append(CoinIndexedVector & other,int adjustIndex,bool zapElements)1228 void CoinIndexedVector::append(CoinIndexedVector &other, int adjustIndex, bool zapElements /*,double multiplier*/)
1229 {
1230   const int cs = other.nElements_;
1231   const int *cind = other.indices_;
1232   double *celem = other.elements_;
1233   int *newInd = indices_ + nElements_;
1234   if (packedMode_) {
1235     double *newEls = elements_ + nElements_;
1236     if (zapElements) {
1237       if (other.packedMode_) {
1238         for (int i = 0; i < cs; i++) {
1239           newInd[i] = cind[i] + adjustIndex;
1240           newEls[i] = celem[i] /* *multiplier*/;
1241           celem[i] = 0.0;
1242         }
1243       } else {
1244         for (int i = 0; i < cs; i++) {
1245           int k = cind[i];
1246           newInd[i] = k + adjustIndex;
1247           newEls[i] = celem[k] /* *multiplier*/;
1248           celem[k] = 0.0;
1249         }
1250       }
1251     } else {
1252       if (other.packedMode_) {
1253         for (int i = 0; i < cs; i++) {
1254           newEls[i] = celem[i] /* *multiplier*/;
1255           newInd[i] = cind[i] + adjustIndex;
1256         }
1257       } else {
1258         for (int i = 0; i < cs; i++) {
1259           int k = cind[i];
1260           newInd[i] = k + adjustIndex;
1261           newEls[i] = celem[k] /* *multiplier*/;
1262         }
1263       }
1264     }
1265   } else {
1266     double *newEls = elements_ + adjustIndex;
1267     if (zapElements) {
1268       if (other.packedMode_) {
1269         for (int i = 0; i < cs; i++) {
1270           int k = cind[i];
1271           newInd[i] = k + adjustIndex;
1272           newEls[k] = celem[i] /* *multiplier*/;
1273           celem[i] = 0.0;
1274         }
1275       } else {
1276         for (int i = 0; i < cs; i++) {
1277           int k = cind[i];
1278           newInd[i] = k + adjustIndex;
1279           newEls[k] = celem[k] /* *multiplier*/;
1280           celem[k] = 0.0;
1281         }
1282       }
1283     } else {
1284       if (other.packedMode_) {
1285         for (int i = 0; i < cs; i++) {
1286           int k = cind[i];
1287           newInd[i] = k + adjustIndex;
1288           newEls[k] = celem[i] /* *multiplier*/;
1289         }
1290       } else {
1291         for (int i = 0; i < cs; i++) {
1292           int k = cind[i];
1293           newInd[i] = k + adjustIndex;
1294           newEls[k] = celem[k] /* *multiplier*/;
1295         }
1296       }
1297     }
1298   }
1299   nElements_ += cs;
1300   if (zapElements)
1301     other.nElements_ = 0;
1302 }
1303 #ifndef CLP_NO_VECTOR
1304 /* Equal. Returns true if vectors have same length and corresponding
1305    element of each vector is equal. */
operator ==(const CoinPackedVectorBase & rhs) const1306 bool CoinIndexedVector::operator==(const CoinPackedVectorBase &rhs) const
1307 {
1308   const int cs = rhs.getNumElements();
1309 
1310   const int *cind = rhs.getIndices();
1311   const double *celem = rhs.getElements();
1312   if (nElements_ != cs)
1313     return false;
1314   int i;
1315   bool okay = true;
1316   for (i = 0; i < cs; i++) {
1317     int iRow = cind[i];
1318     if (celem[i] != elements_[iRow]) {
1319       okay = false;
1320       break;
1321     }
1322   }
1323   return okay;
1324 }
1325 // Not equal
operator !=(const CoinPackedVectorBase & rhs) const1326 bool CoinIndexedVector::operator!=(const CoinPackedVectorBase &rhs) const
1327 {
1328   const int cs = rhs.getNumElements();
1329 
1330   const int *cind = rhs.getIndices();
1331   const double *celem = rhs.getElements();
1332   if (nElements_ != cs)
1333     return true;
1334   int i;
1335   bool okay = false;
1336   for (i = 0; i < cs; i++) {
1337     int iRow = cind[i];
1338     if (celem[i] != elements_[iRow]) {
1339       okay = true;
1340       break;
1341     }
1342   }
1343   return okay;
1344 }
1345 #endif
1346 // Equal with a tolerance (returns -1 or position of inequality).
isApproximatelyEqual(const CoinIndexedVector & rhs,double tolerance) const1347 int CoinIndexedVector::isApproximatelyEqual(const CoinIndexedVector &rhs, double tolerance) const
1348 {
1349   CoinIndexedVector tempA(*this);
1350   CoinIndexedVector tempB(rhs);
1351   int *cind = tempB.indices_;
1352   double *celem = tempB.elements_;
1353   double *elem = tempA.elements_;
1354   int cs = tempB.nElements_;
1355   int bad = -1;
1356   CoinRelFltEq eq(tolerance);
1357   if (!packedMode_ && !tempB.packedMode_) {
1358     for (int i = 0; i < cs; i++) {
1359       int iRow = cind[i];
1360       if (!eq(celem[iRow], elem[iRow])) {
1361         bad = iRow;
1362         break;
1363       } else {
1364         celem[iRow] = elem[iRow] = 0.0;
1365       }
1366     }
1367     cs = tempA.nElements_;
1368     cind = tempA.indices_;
1369     for (int i = 0; i < cs; i++) {
1370       int iRow = cind[i];
1371       if (!eq(celem[iRow], elem[iRow])) {
1372         bad = iRow;
1373         break;
1374       } else {
1375         celem[iRow] = elem[iRow] = 0.0;
1376       }
1377     }
1378   } else if (packedMode_ && tempB.packedMode_) {
1379     double *celem2 = rhs.elements_;
1380     memset(celem, 0, CoinMin(capacity_, tempB.capacity_) * sizeof(double));
1381     for (int i = 0; i < cs; i++) {
1382       int iRow = cind[i];
1383       celem[iRow] = celem2[i];
1384     }
1385     for (int i = 0; i < cs; i++) {
1386       int iRow = cind[i];
1387       if (!eq(celem[iRow], elem[i])) {
1388         bad = iRow;
1389         break;
1390       } else {
1391         celem[iRow] = elem[i] = 0.0;
1392       }
1393     }
1394   } else {
1395     double *celem2 = elem;
1396     double *celem3 = celem;
1397     if (packedMode_) {
1398       celem2 = celem;
1399       celem3 = elem;
1400     }
1401     for (int i = 0; i < cs; i++) {
1402       int iRow = cind[i];
1403       if (!eq(celem2[iRow], celem3[i])) {
1404         bad = iRow;
1405         break;
1406       } else {
1407         celem2[iRow] = celem3[i] = 0.0;
1408       }
1409     }
1410   }
1411   if (bad < 0) {
1412     for (int i = 0; i < tempA.capacity_; i++) {
1413       if (elem[i]) {
1414         if (fabs(elem[i]) > tolerance) {
1415           bad = i;
1416           break;
1417         }
1418       }
1419     }
1420     for (int i = 0; i < tempB.capacity_; i++) {
1421       if (celem[i]) {
1422         if (fabs(celem[i]) > tolerance) {
1423           bad = i;
1424           break;
1425         }
1426       }
1427     }
1428   }
1429   return bad;
1430 }
1431 /* Equal. Returns true if vectors have same length and corresponding
1432    element of each vector is equal. */
operator ==(const CoinIndexedVector & rhs) const1433 bool CoinIndexedVector::operator==(const CoinIndexedVector &rhs) const
1434 {
1435   const int cs = rhs.nElements_;
1436 
1437   const int *cind = rhs.indices_;
1438   const double *celem = rhs.elements_;
1439   if (nElements_ != cs)
1440     return false;
1441   bool okay = true;
1442   CoinRelFltEq eq(1.0e-8);
1443   if (!packedMode_ && !rhs.packedMode_) {
1444     for (int i = 0; i < cs; i++) {
1445       int iRow = cind[i];
1446       if (!eq(celem[iRow], elements_[iRow])) {
1447         okay = false;
1448         break;
1449       }
1450     }
1451   } else if (packedMode_ && rhs.packedMode_) {
1452     double *temp = new double[CoinMax(capacity_, rhs.capacity_)];
1453     memset(temp, 0, CoinMax(capacity_, rhs.capacity_) * sizeof(double));
1454     for (int i = 0; i < cs; i++) {
1455       int iRow = cind[i];
1456       temp[iRow] = celem[i];
1457     }
1458     for (int i = 0; i < cs; i++) {
1459       int iRow = cind[i];
1460       if (!eq(temp[iRow], elements_[i])) {
1461         okay = false;
1462         break;
1463       }
1464     }
1465     delete[] temp;
1466   } else {
1467     const double *celem2 = elements_;
1468     if (packedMode_) {
1469       celem2 = celem;
1470       celem = elements_;
1471     }
1472     for (int i = 0; i < cs; i++) {
1473       int iRow = cind[i];
1474       if (!eq(celem2[iRow], celem[i])) {
1475         okay = false;
1476         break;
1477       }
1478     }
1479   }
1480   return okay;
1481 }
1482 /// Not equal
operator !=(const CoinIndexedVector & rhs) const1483 bool CoinIndexedVector::operator!=(const CoinIndexedVector &rhs) const
1484 {
1485   const int cs = rhs.nElements_;
1486 
1487   const int *cind = rhs.indices_;
1488   const double *celem = rhs.elements_;
1489   if (nElements_ != cs)
1490     return true;
1491   int i;
1492   bool okay = false;
1493   for (i = 0; i < cs; i++) {
1494     int iRow = cind[i];
1495     if (celem[iRow] != elements_[iRow]) {
1496       okay = true;
1497       break;
1498     }
1499   }
1500   return okay;
1501 }
1502 // Get value of maximum index
getMaxIndex() const1503 int CoinIndexedVector::getMaxIndex() const
1504 {
1505   int maxIndex = -COIN_INT_MAX;
1506   int i;
1507   for (i = 0; i < nElements_; i++)
1508     maxIndex = CoinMax(maxIndex, indices_[i]);
1509   return maxIndex;
1510 }
1511 // Get value of minimum index
getMinIndex() const1512 int CoinIndexedVector::getMinIndex() const
1513 {
1514   int minIndex = COIN_INT_MAX;
1515   int i;
1516   for (i = 0; i < nElements_; i++)
1517     minIndex = CoinMin(minIndex, indices_[i]);
1518   return minIndex;
1519 }
1520 // Scan dense region and set up indices
scan()1521 int CoinIndexedVector::scan()
1522 {
1523   nElements_ = 0;
1524   return scan(0, capacity_);
1525 }
1526 // Scan dense region from start to < end and set up indices
scan(int start,int end)1527 int CoinIndexedVector::scan(int start, int end)
1528 {
1529   assert(!packedMode_);
1530   end = CoinMin(end, capacity_);
1531   start = CoinMax(start, 0);
1532   int i;
1533   int number = 0;
1534   int *indices = indices_ + nElements_;
1535   for (i = start; i < end; i++)
1536     if (elements_[i])
1537       indices[number++] = i;
1538   nElements_ += number;
1539   return number;
1540 }
1541 // Scan dense region and set up indices with tolerance
scan(double tolerance)1542 int CoinIndexedVector::scan(double tolerance)
1543 {
1544   nElements_ = 0;
1545 #if ABOCA_LITE_FACTORIZATION == 0
1546   return scan(0, capacity_, tolerance);
1547 #else
1548   return scan(0, capacity_ & 0x7fffffff, tolerance);
1549 #endif
1550 }
1551 // Scan dense region from start to < end and set up indices with tolerance
scan(int start,int end,double tolerance)1552 int CoinIndexedVector::scan(int start, int end, double tolerance)
1553 {
1554   assert(!packedMode_);
1555 #if ABOCA_LITE_FACTORIZATION == 0
1556   end = CoinMin(end, capacity_);
1557 #else
1558   end = CoinMin(end, capacity_ & 0x7fffffff);
1559 #endif
1560   start = CoinMax(start, 0);
1561   int i;
1562   int number = 0;
1563   int *indices = indices_ + nElements_;
1564   for (i = start; i < end; i++) {
1565     double value = elements_[i];
1566     if (value) {
1567       if (fabs(value) >= tolerance)
1568         indices[number++] = i;
1569       else
1570         elements_[i] = 0.0;
1571     }
1572   }
1573   nElements_ += number;
1574   return number;
1575 }
1576 // These pack down
cleanAndPack(double tolerance)1577 int CoinIndexedVector::cleanAndPack(double tolerance)
1578 {
1579   if (!packedMode_) {
1580     int number = nElements_;
1581     int i;
1582     nElements_ = 0;
1583     for (i = 0; i < number; i++) {
1584       int indexValue = indices_[i];
1585       double value = elements_[indexValue];
1586       elements_[indexValue] = 0.0;
1587       if (fabs(value) >= tolerance) {
1588         elements_[nElements_] = value;
1589         indices_[nElements_++] = indexValue;
1590       }
1591     }
1592     packedMode_ = true;
1593   }
1594   return nElements_;
1595 }
1596 // These pack down
cleanAndPackSafe(double tolerance)1597 int CoinIndexedVector::cleanAndPackSafe(double tolerance)
1598 {
1599   int number = nElements_;
1600   if (number) {
1601     int i;
1602     nElements_ = 0;
1603     assert(!packedMode_);
1604     double *temp = NULL;
1605     bool gotMemory;
1606     if (number * 3 < capacity_ - 3 - 9999999) {
1607       // can find room without new
1608       gotMemory = false;
1609       // But may need to align on 8 byte boundary
1610       char *tempC = reinterpret_cast< char * >(indices_ + number);
1611       CoinInt64 xx = reinterpret_cast< CoinInt64 >(tempC);
1612       CoinInt64 iBottom = xx & 7;
1613       if (iBottom)
1614         tempC += 8 - iBottom;
1615       temp = reinterpret_cast< double * >(tempC);
1616       xx = reinterpret_cast< CoinInt64 >(temp);
1617       iBottom = xx & 7;
1618       assert(!iBottom);
1619     } else {
1620       // might be better to do complete scan
1621       gotMemory = true;
1622       temp = new double[number];
1623     }
1624     for (i = 0; i < number; i++) {
1625       int indexValue = indices_[i];
1626       double value = elements_[indexValue];
1627       elements_[indexValue] = 0.0;
1628       if (fabs(value) >= tolerance) {
1629         temp[nElements_] = value;
1630         indices_[nElements_++] = indexValue;
1631       }
1632     }
1633     CoinMemcpyN(temp, nElements_, elements_);
1634     if (gotMemory)
1635       delete[] temp;
1636     packedMode_ = true;
1637   }
1638   return nElements_;
1639 }
1640 // Scan dense region and set up indices
scanAndPack()1641 int CoinIndexedVector::scanAndPack()
1642 {
1643   nElements_ = 0;
1644   return scanAndPack(0, capacity_);
1645 }
1646 // Scan dense region from start to < end and set up indices
scanAndPack(int start,int end)1647 int CoinIndexedVector::scanAndPack(int start, int end)
1648 {
1649   assert(!packedMode_);
1650   end = CoinMin(end, capacity_);
1651   start = CoinMax(start, 0);
1652   int i;
1653   int number = 0;
1654   int *indices = indices_ + nElements_;
1655   for (i = start; i < end; i++) {
1656     double value = elements_[i];
1657     elements_[i] = 0.0;
1658     if (value) {
1659       elements_[number] = value;
1660       indices[number++] = i;
1661     }
1662   }
1663   nElements_ += number;
1664   packedMode_ = true;
1665   return number;
1666 }
1667 // Scan dense region and set up indices with tolerance
scanAndPack(double tolerance)1668 int CoinIndexedVector::scanAndPack(double tolerance)
1669 {
1670   nElements_ = 0;
1671   return scanAndPack(0, capacity_, tolerance);
1672 }
1673 // Scan dense region from start to < end and set up indices with tolerance
scanAndPack(int start,int end,double tolerance)1674 int CoinIndexedVector::scanAndPack(int start, int end, double tolerance)
1675 {
1676   assert(!packedMode_);
1677   end = CoinMin(end, capacity_);
1678   start = CoinMax(start, 0);
1679   int i;
1680   int number = 0;
1681   int *indices = indices_ + nElements_;
1682   for (i = start; i < end; i++) {
1683     double value = elements_[i];
1684     elements_[i] = 0.0;
1685     if (fabs(value) >= tolerance) {
1686       elements_[number] = value;
1687       indices[number++] = i;
1688     }
1689   }
1690   nElements_ += number;
1691   packedMode_ = true;
1692   return number;
1693 }
1694 // This is mainly for testing - goes from packed to indexed
expand()1695 void CoinIndexedVector::expand()
1696 {
1697   if (nElements_ && packedMode_) {
1698     double *temp = new double[capacity_];
1699     int i;
1700     for (i = 0; i < nElements_; i++)
1701       temp[indices_[i]] = elements_[i];
1702     CoinZeroN(elements_, nElements_);
1703     for (i = 0; i < nElements_; i++) {
1704       int iRow = indices_[i];
1705       elements_[iRow] = temp[iRow];
1706     }
1707     delete[] temp;
1708   }
1709   packedMode_ = false;
1710 }
1711 // Create packed array
createPacked(int number,const int * indices,const double * elements)1712 void CoinIndexedVector::createPacked(int number, const int *indices,
1713   const double *elements)
1714 {
1715   nElements_ = number;
1716   packedMode_ = true;
1717   CoinMemcpyN(indices, number, indices_);
1718   CoinMemcpyN(elements, number, elements_);
1719 }
1720 // Create packed array
createUnpacked(int number,const int * indices,const double * elements)1721 void CoinIndexedVector::createUnpacked(int number, const int *indices,
1722   const double *elements)
1723 {
1724   nElements_ = number;
1725   packedMode_ = false;
1726   for (int i = 0; i < nElements_; i++) {
1727     int iRow = indices[i];
1728     indices_[i] = iRow;
1729     elements_[iRow] = elements[i];
1730   }
1731 }
1732 // Create unpacked singleton
createOneUnpackedElement(int index,double element)1733 void CoinIndexedVector::createOneUnpackedElement(int index, double element)
1734 {
1735   nElements_ = 1;
1736   packedMode_ = false;
1737   indices_[0] = index;
1738   elements_[index] = element;
1739 }
1740 //  Print out
print() const1741 void CoinIndexedVector::print() const
1742 {
1743   printf("Vector has %d elements (%spacked mode)\n", nElements_, packedMode_ ? "" : "un");
1744   for (int i = 0; i < nElements_; i++) {
1745     if (i && (i % 5 == 0))
1746       printf("\n");
1747     int index = indices_[i];
1748     double value = packedMode_ ? elements_[i] : elements_[index];
1749     printf(" (%d,%g)", index, value);
1750   }
1751   printf("\n");
1752 }
1753 
1754 // Zero out array
clear()1755 void CoinArrayWithLength::clear()
1756 {
1757   assert((size_ > 0 && array_) || !array_);
1758   memset(array_, 0, size_);
1759 }
1760 // Get array with alignment
getArray(CoinBigIndex size)1761 void CoinArrayWithLength::getArray(CoinBigIndex size)
1762 {
1763   if (size > 0) {
1764     if (alignment_ > 2) {
1765       offset_ = 1 << alignment_;
1766     } else {
1767       offset_ = 0;
1768     }
1769     assert(size > 0);
1770     char *array = new char[size + offset_];
1771     if (offset_) {
1772       // need offset
1773       CoinInt64 xx = reinterpret_cast< CoinInt64 >(array);
1774       int iBottom = static_cast< int >(xx & ((offset_ - 1)));
1775       if (iBottom)
1776         offset_ = offset_ - iBottom;
1777       else
1778         offset_ = 0;
1779       array_ = array + offset_;
1780       ;
1781     } else {
1782       array_ = array;
1783     }
1784     if (size_ != -1)
1785       size_ = size;
1786   } else {
1787     array_ = NULL;
1788   }
1789 }
1790 // Get rid of array with alignment
conditionalDelete()1791 void CoinArrayWithLength::conditionalDelete()
1792 {
1793   if (size_ == -1) {
1794     char *charArray = reinterpret_cast< char * >(array_);
1795     if (charArray)
1796       delete[](charArray - offset_);
1797     array_ = NULL;
1798   } else if (size_ >= 0) {
1799     size_ = -size_ - 2;
1800   }
1801 }
1802 // Really get rid of array with alignment
reallyFreeArray()1803 void CoinArrayWithLength::reallyFreeArray()
1804 {
1805   char *charArray = reinterpret_cast< char * >(array_);
1806   if (charArray)
1807     delete[](charArray - offset_);
1808   array_ = NULL;
1809   size_ = -1;
1810 }
1811 // Get enough space
getCapacity(CoinBigIndex numberBytes,CoinBigIndex numberNeeded)1812 void CoinArrayWithLength::getCapacity(CoinBigIndex numberBytes, CoinBigIndex numberNeeded)
1813 {
1814   CoinBigIndex k = capacity();
1815   if (k < numberBytes) {
1816     CoinBigIndex saveSize = size_;
1817     reallyFreeArray();
1818     size_ = saveSize;
1819     getArray(CoinMax(numberBytes, numberNeeded));
1820   } else if (size_ < 0) {
1821     size_ = -size_ - 2;
1822   }
1823 }
1824 /* Alternate Constructor - length in bytes
1825    mode -  0 size_ set to size
1826    >0 size_ set to size and zeroed
1827    if size<=0 just does alignment
1828    If abs(mode) >2 then align on that as power of 2
1829 */
CoinArrayWithLength(CoinBigIndex size,int mode)1830 CoinArrayWithLength::CoinArrayWithLength(CoinBigIndex size, int mode)
1831 {
1832   alignment_ = abs(mode);
1833   size_ = size;
1834   getArray(size);
1835   if (mode > 0 && array_)
1836     memset(array_, 0, size);
1837 }
~CoinArrayWithLength()1838 CoinArrayWithLength::~CoinArrayWithLength()
1839 {
1840   if (array_)
1841     delete[](array_ - offset_);
1842 }
1843 // Conditionally gets new array
1844 char *
conditionalNew(CoinBigIndex sizeWanted)1845 CoinArrayWithLength::conditionalNew(CoinBigIndex sizeWanted)
1846 {
1847   if (size_ == -1) {
1848     getCapacity(static_cast< int >(sizeWanted));
1849   } else {
1850     int newSize = static_cast< int >(sizeWanted * 101 / 100) + 64;
1851     // round to multiple of 16
1852     newSize -= newSize & 15;
1853     getCapacity(static_cast< int >(sizeWanted), newSize);
1854   }
1855   return array_;
1856 }
1857 /* Copy constructor. */
CoinArrayWithLength(const CoinArrayWithLength & rhs)1858 CoinArrayWithLength::CoinArrayWithLength(const CoinArrayWithLength &rhs)
1859 {
1860   assert(capacity() >= 0);
1861   getArray(rhs.capacity());
1862   if (size_ > 0)
1863     CoinMemcpyN(rhs.array_, size_, array_);
1864 }
1865 
1866 /* Copy constructor.2 */
CoinArrayWithLength(const CoinArrayWithLength * rhs)1867 CoinArrayWithLength::CoinArrayWithLength(const CoinArrayWithLength *rhs)
1868 {
1869   assert(rhs->capacity() >= 0);
1870   size_ = rhs->size_;
1871   getArray(rhs->capacity());
1872   if (size_ > 0)
1873     CoinMemcpyN(rhs->array_, size_, array_);
1874 }
1875 /* Assignment operator. */
1876 CoinArrayWithLength &
operator =(const CoinArrayWithLength & rhs)1877 CoinArrayWithLength::operator=(const CoinArrayWithLength &rhs)
1878 {
1879   if (this != &rhs) {
1880     assert(rhs.size_ != -1 || !rhs.array_);
1881     if (rhs.size_ == -1) {
1882       reallyFreeArray();
1883     } else {
1884       getCapacity(rhs.size_);
1885       if (size_ > 0)
1886         CoinMemcpyN(rhs.array_, size_, array_);
1887     }
1888   }
1889   return *this;
1890 }
1891 /* Assignment with length (if -1 use internal length) */
copy(const CoinArrayWithLength & rhs,int numberBytes)1892 void CoinArrayWithLength::copy(const CoinArrayWithLength &rhs, int numberBytes)
1893 {
1894   if (numberBytes == -1 || numberBytes <= rhs.capacity()) {
1895     CoinArrayWithLength::operator=(rhs);
1896   } else {
1897     assert(numberBytes >= 0);
1898     getCapacity(numberBytes);
1899     if (rhs.array_)
1900       CoinMemcpyN(rhs.array_, numberBytes, array_);
1901   }
1902 }
1903 /* Assignment with length - does not copy */
allocate(const CoinArrayWithLength & rhs,CoinBigIndex numberBytes)1904 void CoinArrayWithLength::allocate(const CoinArrayWithLength &rhs, CoinBigIndex numberBytes)
1905 {
1906   if (numberBytes == -1 || numberBytes <= rhs.capacity()) {
1907     assert(rhs.size_ != -1 || !rhs.array_);
1908     if (rhs.size_ == -1) {
1909       reallyFreeArray();
1910     } else {
1911       getCapacity(rhs.size_);
1912     }
1913   } else {
1914     assert(numberBytes >= 0);
1915     if (size_ == -1) {
1916       delete[] array_;
1917       array_ = NULL;
1918     } else {
1919       size_ = -1;
1920     }
1921     if (rhs.size_ >= 0)
1922       size_ = numberBytes;
1923     assert(numberBytes >= 0);
1924     assert(!array_);
1925     if (numberBytes)
1926       array_ = new char[numberBytes];
1927   }
1928 }
1929 // Does what is needed to set persistence
setPersistence(int flag,int currentLength)1930 void CoinArrayWithLength::setPersistence(int flag, int currentLength)
1931 {
1932   if (flag) {
1933     if (size_ == -1) {
1934       if (currentLength && array_) {
1935         size_ = currentLength;
1936       } else {
1937         conditionalDelete();
1938         size_ = 0;
1939         array_ = NULL;
1940       }
1941     }
1942   } else {
1943     size_ = -1;
1944   }
1945 }
1946 // Swaps memory between two members
swap(CoinArrayWithLength & other)1947 void CoinArrayWithLength::swap(CoinArrayWithLength &other)
1948 {
1949 #ifdef COIN_DEVELOP
1950   if (!(size_ == other.size_ || size_ == -1 || other.size_ == -1))
1951     printf("Two arrays have sizes - %d and %d\n", size_, other.size_);
1952 #endif
1953   assert(alignment_ == other.alignment_);
1954   char *swapArray = other.array_;
1955   other.array_ = array_;
1956   array_ = swapArray;
1957   CoinBigIndex swapSize = other.size_;
1958   other.size_ = size_;
1959   size_ = swapSize;
1960   int swapOffset = other.offset_;
1961   other.offset_ = offset_;
1962   offset_ = swapOffset;
1963 }
1964 // Extend a persistent array keeping data (size in bytes)
extend(int newSize)1965 void CoinArrayWithLength::extend(int newSize)
1966 {
1967   //assert (newSize>=capacity()&&capacity()>=0);
1968   assert(size_ >= 0); // not much point otherwise
1969   if (newSize > size_) {
1970     char *temp = array_;
1971     getArray(newSize);
1972     if (temp) {
1973       CoinMemcpyN(array_, size_, temp);
1974       delete[](temp - offset_);
1975     }
1976     size_ = newSize;
1977   }
1978 }
1979 
1980 /* Default constructor */
CoinPartitionedVector()1981 CoinPartitionedVector::CoinPartitionedVector()
1982   : CoinIndexedVector()
1983 {
1984   memset(startPartition_, 0, ((&numberPartitions_ - startPartition_) + 1) * sizeof(int));
1985 }
1986 /* Copy constructor. */
CoinPartitionedVector(const CoinPartitionedVector & rhs)1987 CoinPartitionedVector::CoinPartitionedVector(const CoinPartitionedVector &rhs)
1988   : CoinIndexedVector(rhs)
1989 {
1990   memcpy(startPartition_, rhs.startPartition_, ((&numberPartitions_ - startPartition_) + 1) * sizeof(int));
1991 }
1992 /* Copy constructor.2 */
CoinPartitionedVector(const CoinPartitionedVector * rhs)1993 CoinPartitionedVector::CoinPartitionedVector(const CoinPartitionedVector *rhs)
1994   : CoinIndexedVector(rhs)
1995 {
1996   memcpy(startPartition_, rhs->startPartition_, ((&numberPartitions_ - startPartition_) + 1) * sizeof(int));
1997 }
1998 /* Assignment operator. */
operator =(const CoinPartitionedVector & rhs)1999 CoinPartitionedVector &CoinPartitionedVector::operator=(const CoinPartitionedVector &rhs)
2000 {
2001   if (this != &rhs) {
2002     CoinIndexedVector::operator=(rhs);
2003     memcpy(startPartition_, rhs.startPartition_, ((&numberPartitions_ - startPartition_) + 1) * sizeof(int));
2004   }
2005   return *this;
2006 }
2007 /* Destructor */
~CoinPartitionedVector()2008 CoinPartitionedVector::~CoinPartitionedVector()
2009 {
2010 }
2011 // Add up number of elements in partitions
computeNumberElements()2012 void CoinPartitionedVector::computeNumberElements()
2013 {
2014   if (numberPartitions_) {
2015     assert(packedMode_);
2016     int n = 0;
2017     for (int i = 0; i < numberPartitions_; i++)
2018       n += numberElementsPartition_[i];
2019     nElements_ = n;
2020   }
2021 }
2022 // Add up number of elements in partitions and pack and get rid of partitions
compact()2023 void CoinPartitionedVector::compact()
2024 {
2025   if (numberPartitions_) {
2026     int n = numberElementsPartition_[0];
2027     numberElementsPartition_[0] = 0;
2028     for (int i = 1; i < numberPartitions_; i++) {
2029       int nThis = numberElementsPartition_[i];
2030       int start = startPartition_[i];
2031       memmove(indices_ + n, indices_ + start, nThis * sizeof(int));
2032       memmove(elements_ + n, elements_ + start, nThis * sizeof(double));
2033       n += nThis;
2034     }
2035     nElements_ = n;
2036     // clean up
2037     for (int i = 1; i < numberPartitions_; i++) {
2038       int nThis = numberElementsPartition_[i];
2039       int start = startPartition_[i];
2040       numberElementsPartition_[i] = 0;
2041       int end = nThis + start;
2042       if (nElements_ < end) {
2043         int offset = CoinMax(nElements_ - start, 0);
2044         start += offset;
2045         nThis -= offset;
2046         memset(elements_ + start, 0, nThis * sizeof(double));
2047       }
2048     }
2049     packedMode_ = true;
2050     numberPartitions_ = 0;
2051   }
2052 }
2053 /* Reserve space.
2054  */
reserve(int n)2055 void CoinPartitionedVector::reserve(int n)
2056 {
2057   CoinIndexedVector::reserve(n);
2058   memset(startPartition_, 0, ((&numberPartitions_ - startPartition_) + 1) * sizeof(int));
2059   startPartition_[1] = capacity_; // for safety
2060 }
2061 // Setup partitions (needs end as well)
setPartitions(int number,const int * starts)2062 void CoinPartitionedVector::setPartitions(int number, const int *starts)
2063 {
2064   if (number) {
2065     packedMode_ = true;
2066     assert(number <= COIN_PARTITIONS);
2067     memcpy(startPartition_, starts, (number + 1) * sizeof(int));
2068     numberPartitions_ = number;
2069 #ifndef NDEBUG
2070     assert(startPartition_[0] == 0);
2071     int last = -1;
2072     for (int i = 0; i < numberPartitions_; i++) {
2073       assert(startPartition_[i] >= last);
2074       assert(numberElementsPartition_[i] == 0);
2075       last = startPartition_[i];
2076     }
2077     assert(startPartition_[numberPartitions_] >= last && startPartition_[numberPartitions_] <= capacity_);
2078 #endif
2079   } else {
2080     clearAndReset();
2081   }
2082 }
2083 // Reset the vector (as if were just created an empty vector). Gets rid of partitions
clearAndReset()2084 void CoinPartitionedVector::clearAndReset()
2085 {
2086   if (numberPartitions_) {
2087     assert(packedMode_ || !nElements_);
2088     for (int i = 0; i < numberPartitions_; i++) {
2089       int n = numberElementsPartition_[i];
2090       memset(elements_ + startPartition_[i], 0, n * sizeof(double));
2091       numberElementsPartition_[i] = 0;
2092     }
2093   } else {
2094     memset(elements_, 0, nElements_ * sizeof(double));
2095   }
2096   nElements_ = 0;
2097   numberPartitions_ = 0;
2098   startPartition_[1] = capacity_;
2099   packedMode_ = false;
2100 }
2101 // Reset the vector (as if were just created an empty vector). Keeps partitions
clearAndKeep()2102 void CoinPartitionedVector::clearAndKeep()
2103 {
2104   assert(packedMode_);
2105   for (int i = 0; i < numberPartitions_; i++) {
2106     int n = numberElementsPartition_[i];
2107     memset(elements_ + startPartition_[i], 0, n * sizeof(double));
2108     numberElementsPartition_[i] = 0;
2109   }
2110   nElements_ = 0;
2111 }
2112 // Clear a partition.
clearPartition(int partition)2113 void CoinPartitionedVector::clearPartition(int partition)
2114 {
2115   assert(packedMode_);
2116   assert(partition < COIN_PARTITIONS);
2117   int n = numberElementsPartition_[partition];
2118   memset(elements_ + startPartition_[partition], 0, n * sizeof(double));
2119   numberElementsPartition_[partition] = 0;
2120 }
2121 #ifndef NDEBUG
2122 // For debug check vector is clear i.e. no elements
checkClear()2123 void CoinPartitionedVector::checkClear()
2124 {
2125 #ifndef NDEBUG
2126   //printf("checkClear %p\n",this);
2127   assert(!nElements_);
2128   //assert(!packedMode_);
2129   int i;
2130   for (i = 0; i < capacity_; i++) {
2131     assert(!elements_[i]);
2132   }
2133 #else
2134   if (nElements_) {
2135     printf("%d nElements_ - checkClear\n", nElements_);
2136     abort();
2137   }
2138   int i;
2139   int n = 0;
2140   int k = -1;
2141   for (i = 0; i < capacity_; i++) {
2142     if (elements_[i]) {
2143       n++;
2144       if (k < 0)
2145         k = i;
2146     }
2147   }
2148   if (n) {
2149     printf("%d elements, first %d - checkClear\n", n, k);
2150     abort();
2151   }
2152 #endif
2153 }
2154 // For debug check vector is clean i.e. elements match indices
checkClean()2155 void CoinPartitionedVector::checkClean()
2156 {
2157   //printf("checkClean %p\n",this);
2158   int i;
2159   if (!nElements_) {
2160     checkClear();
2161   } else {
2162     if (packedMode_) {
2163       for (i = 0; i < nElements_; i++)
2164         assert(elements_[i]);
2165       for (; i < capacity_; i++)
2166         assert(!elements_[i]);
2167     } else {
2168       double *copy = new double[capacity_];
2169       CoinMemcpyN(elements_, capacity_, copy);
2170       for (i = 0; i < nElements_; i++) {
2171         int indexValue = indices_[i];
2172         assert(copy[indexValue]);
2173         copy[indexValue] = 0.0;
2174       }
2175       for (i = 0; i < capacity_; i++)
2176         assert(!copy[i]);
2177       delete[] copy;
2178     }
2179 #ifndef NDEBUG
2180     // check mark array zeroed
2181     char *mark = reinterpret_cast< char * >(indices_ + capacity_);
2182     for (i = 0; i < capacity_; i++) {
2183       assert(!mark[i]);
2184     }
2185 #endif
2186   }
2187 }
2188 #endif
2189 // Scan dense region and set up indices (returns number found)
scan(int partition,double tolerance)2190 int CoinPartitionedVector::scan(int partition, double tolerance)
2191 {
2192   assert(packedMode_);
2193   assert(partition < COIN_PARTITIONS);
2194   int n = 0;
2195   int start = startPartition_[partition];
2196   double *COIN_RESTRICT elements = elements_ + start;
2197   int *COIN_RESTRICT indices = indices_ + start;
2198   int sizePartition = startPartition_[partition + 1] - start;
2199   if (!tolerance) {
2200     for (int i = 0; i < sizePartition; i++) {
2201       double value = elements[i];
2202       if (value) {
2203         elements[i] = 0.0;
2204         elements[n] = value;
2205         indices[n++] = i + start;
2206       }
2207     }
2208   } else {
2209     for (int i = 0; i < sizePartition; i++) {
2210       double value = elements[i];
2211       if (value) {
2212         elements[i] = 0.0;
2213         if (fabs(value) > tolerance) {
2214           elements[n] = value;
2215           indices[n++] = i + start;
2216         }
2217       }
2218     }
2219   }
2220   numberElementsPartition_[partition] = n;
2221   return n;
2222 }
2223 //  Print out
print() const2224 void CoinPartitionedVector::print() const
2225 {
2226   printf("Vector has %d elements (%d partitions)\n", nElements_, numberPartitions_);
2227   if (!numberPartitions_) {
2228     CoinIndexedVector::print();
2229     return;
2230   }
2231   double *tempElements = CoinCopyOfArray(elements_, capacity_);
2232   int *tempIndices = CoinCopyOfArray(indices_, capacity_);
2233   for (int partition = 0; partition < numberPartitions_; partition++) {
2234     printf("Partition %d has %d elements\n", partition, numberElementsPartition_[partition]);
2235     int start = startPartition_[partition];
2236     double *COIN_RESTRICT elements = tempElements + start;
2237     int *COIN_RESTRICT indices = tempIndices + start;
2238     CoinSort_2(indices, indices + numberElementsPartition_[partition], elements);
2239     for (int i = 0; i < numberElementsPartition_[partition]; i++) {
2240       if (i && (i % 5 == 0))
2241         printf("\n");
2242       int index = indices[i];
2243       double value = elements[i];
2244       printf(" (%d,%g)", index, value);
2245     }
2246     printf("\n");
2247   }
2248   delete[] tempElements;
2249   delete[] tempIndices;
2250 }
2251 /* Sort the indexed storage vector (increasing indices). */
sort()2252 void CoinPartitionedVector::sort()
2253 {
2254   assert(packedMode_);
2255   for (int partition = 0; partition < numberPartitions_; partition++) {
2256     int start = startPartition_[partition];
2257     double *COIN_RESTRICT elements = elements_ + start;
2258     int *COIN_RESTRICT indices = indices_ + start;
2259     CoinSort_2(indices, indices + numberElementsPartition_[partition], elements);
2260   }
2261 }
2262 
2263 /* vi: softtabstop=2 shiftwidth=2 expandtab tabstop=2
2264 */
2265