1    	/* * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * */
2    	/*                                                                           */
3    	/*                  This file is part of the class library                   */
4    	/*       SoPlex --- the Sequential object-oriented simPlex.                  */
5    	/*                                                                           */
6    	/*  Copyright (c) 1996-2023 Zuse Institute Berlin (ZIB)                      */
7    	/*                                                                           */
8    	/*  Licensed under the Apache License, Version 2.0 (the "License");          */
9    	/*  you may not use this file except in compliance with the License.         */
10   	/*  You may obtain a copy of the License at                                  */
11   	/*                                                                           */
12   	/*      http://www.apache.org/licenses/LICENSE-2.0                           */
13   	/*                                                                           */
14   	/*  Unless required by applicable law or agreed to in writing, software      */
15   	/*  distributed under the License is distributed on an "AS IS" BASIS,        */
16   	/*  WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied. */
17   	/*  See the License for the specific language governing permissions and      */
18   	/*  limitations under the License.                                           */
19   	/*                                                                           */
20   	/*  You should have received a copy of the Apache-2.0 license                */
21   	/*  along with SoPlex; see the file LICENSE. If not email to soplex@zib.de.  */
22   	/*                                                                           */
23   	/* * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * */
24   	
25   	/**@file  ssvectorbase.h
26   	 * @brief Semi sparse vector.
27   	 */
28   	#ifndef _SSVECTORBASE_H_
29   	#define _SSVECTORBASE_H_
30   	
31   	#include <assert.h>
32   	
33   	#include "soplex/spxdefines.h"
34   	#include "soplex/vectorbase.h"
35   	#include "soplex/idxset.h"
36   	#include "soplex/spxalloc.h"
37   	#include "soplex/timer.h"
38   	#include "soplex/stablesum.h"
39   	
40   	namespace soplex
41   	{
42   	template < class R > class SVectorBase;
43   	template < class R > class SVSetBase;
44   	
45   	/**@brief   Semi sparse vector.
46   	 * @ingroup Algebra
47   	 *
48   	 *  This class implements semi-sparse vectors.  Such are #VectorBase%s where the indices of its nonzero elements can be
49   	 *  stored in an extra IdxSet.  Only elements with absolute value > #epsilon are considered to be nonzero.  Since really
50   	 *  storing the nonzeros is not always convenient, an SSVectorBase provides two different stati: setup and not setup.
51   	 *  An SSVectorBase being setup means that the nonzero indices are available, otherwise an SSVectorBase is just an
52   	 *  ordinary VectorBase with an empty IdxSet.  Note that due to arithmetic operation, zeros can slip in, i.e., it is
53   	 *  only guaranteed that at least every non-zero is in the IdxSet.
54   	 */
55   	template < class R >
56   	class SSVectorBase : public VectorBase<R>, protected IdxSet
57   	{
58   	private:
59   	
60   	   friend class VectorBase<R>;
61   	   template < class S > friend class DSVectorBase;
62   	
63   	   // ------------------------------------------------------------------------------------------------------------------
64   	   /**@name Data */
65   	   ///@{
66   	
67   	   /// Is the SSVectorBase set up?
68   	   bool setupStatus;
69   	
70   	   /// Allocates enough space to accommodate \p newmax values.
71   	   void setMax(int newmax)
72   	   {
73   	      assert(idx != 0);
74   	      assert(newmax != 0);
75   	      assert(newmax >= IdxSet::size());
76   	
77   	      len = newmax;
78   	      spx_realloc(idx, len);
79   	   }
80   	
81   	   ///@}
82   	
83   	protected:
84   	   std::shared_ptr<Tolerances> _tolerances;
85   	
86   	public:
87   	
88   	   // ------------------------------------------------------------------------------------------------------------------
89   	   /**@name Status of an SSVectorBase
90   	    *
91   	    *  An SSVectorBase can be set up or not. In case it is set up, its IdxSet correctly contains all indices of nonzero
92   	    *  elements of the SSVectorBase.  Otherwise, it does not contain any useful data. Whether or not an SSVectorBase is
93   	    *  setup can be determined with the method \ref soplex::SSVectorBase::isSetup() "isSetup()".
94   	    *
95   	    *  There are three methods for directly affecting the setup status of an SSVectorBase:
96   	    *
97   	    *  - unSetup():    This method sets the status to ``not setup''.
98   	    *
99   	    *  - setup():      This method initializes the IdxSet to the SSVectorBase's nonzero indices and sets the status to
100  	    *                  ``setup''.
101  	    *
102  	    *  - forceSetup(): This method sets the status to ``setup'' without verifying that the IdxSet correctly contains all
103  	    *                  nonzero indices. It may be used when the nonzero indices have been computed externally.
104  	    */
105  	   ///@{
106  	
107  	   /// Only used in slufactor.hpp.
108  	   R* get_ptr()
109  	   {
110  	      return VectorBase<R>::get_ptr();
111  	   }
112  	
113  	   /// set the _tolerances member variable
114  	   virtual void setTolerances(std::shared_ptr<Tolerances> newTolerances)
115  	   {
116  	      this->_tolerances = newTolerances;
117  	   }
118  	
119  	   /// returns current tolerances
120  	   const std::shared_ptr<Tolerances>& tolerances() const
121  	   {
122  	      assert(this->_tolerances != nullptr);
123  	      return this->_tolerances;
124  	   }
125  	
126  	   /// Returns setup status.
127  	   bool isSetup() const
128  	   {
129  	      return setupStatus;
130  	   }
131  	
132  	   R getEpsilon() const
133  	   {
134  	      assert(this->_tolerances != nullptr);
135  	      return this->_tolerances == nullptr ? R(0) : R(this->tolerances()->epsilon());
136  	   }
137  	
138  	   /// Makes SSVectorBase not setup.
139  	   void unSetup()
140  	   {
141  	      setupStatus = false;
142  	   }
143  	
144  	   /// Initializes nonzero indices for elements with absolute values above #epsilon and sets all other elements to 0.
145  	   void setup()
146  	   {
147  	      if(!isSetup())
148  	      {
149  	         IdxSet::clear();
150  	
151  	         int d = dim();
152  	         num = 0;
153  	
154  	         for(int i = 0; i < d; ++i)
155  	         {
156  	            if(VectorBase<R>::val[i] != R(0))
157  	            {
158  	               if(spxAbs(VectorBase<R>::val[i]) <= this->getEpsilon())
159  	                  VectorBase<R>::val[i] = R(0);
160  	               else
161  	               {
162  	                  idx[num] = i;
163  	                  num++;
164  	               }
165  	            }
166  	         }
167  	
168  	         setupStatus = true;
169  	
170  	         assert(isConsistent());
171  	      }
172  	   }
173  	
174  	   /// Forces setup status.
175  	   void forceSetup()
176  	   {
177  	      setupStatus = true;
178  	   }
179  	
180  	   ///@}
181  	
182  	   // ------------------------------------------------------------------------------------------------------------------
183  	   /**@name Methods for setup SSVectorBases */
184  	   ///@{
185  	
186  	   /// Returns index of the \p n 'th nonzero element.
187  	   int index(int n) const
188  	   {
189  	      assert(isSetup());
190  	
(1) Event deref_parm_in_call: Function "index" dereferences "this->idx". [details]
191  	      return IdxSet::index(n);
192  	   }
193  	
194  	   /// Returns value of the \p n 'th nonzero element.
195  	   R value(int n) const
196  	   {
197  	      assert(isSetup());
198  	      assert(n >= 0 && n < size());
199  	
200  	      return VectorBase<R>::val[idx[n]];
201  	   }
202  	
203  	   /// Finds the position of index \p i in the #IdxSet, or -1 if \p i doesn't exist.
204  	   int pos(int i) const
205  	   {
206  	      assert(isSetup());
207  	
208  	      return IdxSet::pos(i);
209  	   }
210  	
211  	   /// Returns the number of nonzeros.
212  	   int size() const
213  	   {
214  	      assert(isSetup());
215  	
216  	      return IdxSet::size();
217  	   }
218  	
219  	   /// Adds nonzero (\p i, \p x) to SSVectorBase.
220  	   /** No nonzero with index \p i must exist in the SSVectorBase. */
221  	   void add(int i, R x)
222  	   {
223  	      assert(VectorBase<R>::val[i] == R(0));
224  	      assert(pos(i) < 0);
225  	
226  	      addIdx(i);
227  	      VectorBase<R>::val[i] = x;
228  	   }
229  	
230  	   /// Sets \p i 'th element to \p x.
231  	   void setValue(int i, R x)
232  	   {
233  	      assert(i >= 0);
234  	      assert(i < VectorBase<R>::dim());
235  	
236  	      if(isSetup())
237  	      {
238  	         int n = pos(i);
239  	
240  	         if(n < 0)
241  	         {
242  	            if(spxAbs(x) > this->getEpsilon())
243  	               IdxSet::add(1, &i);
244  	         }
245  	         else if(x == R(0))
246  	            clearNum(n);
247  	      }
248  	
249  	      VectorBase<R>::val[i] = x;
250  	
251  	      assert(isConsistent());
252  	   }
253  	
254  	   /// Scale \p i 'th element by a
255  	   void scaleValue(int i, int scaleExp)
256  	   {
257  	      assert(i >= 0);
258  	      assert(i < VectorBase<R>::dim());
259  	
260  	      VectorBase<R>::val[i] = spxLdexp(VectorBase<R>::val[i], scaleExp);
261  	
262  	      assert(isConsistent());
263  	   }
264  	
265  	   /// Clears element \p i.
266  	   void clearIdx(int i)
267  	   {
268  	      if(isSetup())
269  	      {
270  	         int n = pos(i);
271  	
272  	         if(n >= 0)
273  	            remove(n);
274  	      }
275  	
276  	      VectorBase<R>::val[i] = 0;
277  	
278  	      assert(isConsistent());
279  	   }
280  	
281  	   /// Sets \p n 'th nonzero element to 0 (index \p n must exist).
282  	   void clearNum(int n)
283  	   {
284  	      assert(isSetup());
285  	      assert(index(n) >= 0);
286  	
287  	      VectorBase<R>::val[index(n)] = 0;
288  	      remove(n);
289  	
290  	      assert(isConsistent());
291  	   }
292  	
293  	   ///@}
294  	
295  	   // ------------------------------------------------------------------------------------------------------------------
296  	   /**@name Methods independent of the Status */
297  	   ///@{
298  	
299  	   /// Returns \p i 'th value.
300  	   R operator[](int i) const
301  	   {
302  	      return VectorBase<R>::val[i];
303  	   }
304  	
305  	   /// Returns array indices.
306  	   const int* indexMem() const
307  	   {
308  	      return idx;
309  	   }
310  	
311  	   /// Returns array values.
312  	   const R* values() const
313  	   {
314  	      return VectorBase<R>::val.data();
315  	   }
316  	
317  	   /// Returns indices.
318  	   const IdxSet& indices() const
319  	   {
320  	      return *this;
321  	   }
322  	
323  	   /// Returns array indices.
324  	   int* altIndexMem()
325  	   {
326  	      unSetup();
327  	      return idx;
328  	   }
329  	
330  	   /// Returns array values.
331  	   R* altValues()
332  	   {
333  	      unSetup();
334  	      return VectorBase<R>::val.data();
335  	   }
336  	
337  	   /// Returns indices.
338  	   IdxSet& altIndices()
339  	   {
340  	      unSetup();
341  	      return *this;
342  	   }
343  	
344  	   ///@}
345  	
346  	   // ------------------------------------------------------------------------------------------------------------------
347  	   /**@name Arithmetic operations */
348  	   ///@{
349  	
350  	   /// Addition.
351  	   template < class S >
352  	   SSVectorBase<R>& operator+=(const VectorBase<S>& vec)
353  	   {
354  	      VectorBase<S>::operator+=(vec);
355  	
356  	      if(isSetup())
357  	      {
358  	         setupStatus = false;
359  	         setup();
360  	      }
361  	
362  	      return *this;
363  	   }
364  	
365  	   /// Addition.
366  	   template < class S >
367  	   SSVectorBase<R>& operator+=(const SVectorBase<S>& vec);
368  	
369  	   /// Addition.
370  	   template < class S >
371  	   SSVectorBase<R>& operator+=(const SSVectorBase<S>& vec)
372  	   {
373  	      assert(vec.isSetup());
374  	
375  	      for(int i = vec.size() - 1; i >= 0; --i)
376  	         VectorBase<R>::val[vec.index(i)] += vec.value(i);
377  	
378  	      if(isSetup())
379  	      {
380  	         setupStatus = false;
381  	         setup();
382  	      }
383  	
384  	      return *this;
385  	   }
386  	
387  	   /// Subtraction.
388  	   template < class S >
389  	   SSVectorBase<R>& operator-=(const VectorBase<S>& vec)
390  	   {
391  	      VectorBase<R>::operator-=(vec);
392  	
393  	      if(isSetup())
394  	      {
395  	         setupStatus = false;
396  	         setup();
397  	      }
398  	
399  	      return *this;
400  	   }
401  	
402  	   /// Subtraction.
403  	   template < class S >
404  	   SSVectorBase<R>& operator-=(const SVectorBase<S>& vec);
405  	
406  	   /// Subtraction.
407  	   template < class S >
408  	   SSVectorBase<R>& operator-=(const SSVectorBase<S>& vec)
409  	   {
410  	      if(vec.isSetup())
411  	      {
412  	         for(int i = vec.size() - 1; i >= 0; --i)
413  	            VectorBase<R>::val[vec.index(i)] -= vec.value(i);
414  	      }
415  	      else
416  	         VectorBase<R>::operator-=(VectorBase<S>(vec));
417  	
418  	      if(isSetup())
419  	      {
420  	         setupStatus = false;
421  	         setup();
422  	      }
423  	
424  	      return *this;
425  	   }
426  	
427  	   /// Scaling.
428  	   template < class S >
429  	   SSVectorBase<R>& operator*=(S x)
430  	   {
431  	      assert(isSetup());
432  	      assert(x != S(0));
433  	
434  	      for(int i = size() - 1; i >= 0; --i)
435  	         VectorBase<R>::val[index(i)] *= x;
436  	
437  	      assert(isConsistent());
438  	
439  	      return *this;
440  	   }
441  	
442  	   // Inner product.
443  	   template < class S >
444  	   R operator*(const SSVectorBase<S>& w)
445  	   {
446  	      setup();
447  	
448  	      StableSum<R> x;
449  	      int i = size() - 1;
450  	      int j = w.size() - 1;
451  	
452  	      // both *this and w non-zero vectors?
453  	      if(i >= 0 && j >= 0)
454  	      {
455  	         int vi = index(i);
456  	         int wj = w.index(j);
457  	
458  	         while(i != 0 && j != 0)
459  	         {
460  	            if(vi == wj)
461  	            {
462  	               x += VectorBase<R>::val[vi] * R(w.val[wj]);
463  	               vi = index(--i);
464  	               wj = w.index(--j);
465  	            }
466  	            else if(vi > wj)
467  	               vi = index(--i);
468  	            else
469  	               wj = w.index(--j);
470  	         }
471  	
472  	         /* check remaining indices */
473  	
474  	         while(i != 0 && vi != wj)
475  	            vi = index(--i);
476  	
477  	         while(j != 0 && vi != wj)
478  	            wj = w.index(--j);
479  	
480  	         if(vi == wj)
481  	            x += VectorBase<R>::val[vi] * R(w.val[wj]);
482  	      }
483  	
484  	      return x;
485  	   }
486  	
487  	   /// Addition of a scaled vector.
488  	   ///@todo SSVectorBase::multAdd() should be rewritten without pointer arithmetic.
489  	   template < class S, class T >
490  	   SSVectorBase<R>& multAdd(S xx, const SVectorBase<T>& vec);
491  	
492  	   /// Addition of a scaled vector.
493  	   template < class S, class T >
494  	   SSVectorBase<R>& multAdd(S x, const VectorBase<T>& vec)
495  	   {
496  	      VectorBase<R>::multAdd(x, vec);
497  	
498  	      if(isSetup())
499  	      {
500  	         setupStatus = false;
501  	         setup();
502  	      }
503  	
504  	      return *this;
505  	   }
506  	
507  	   /// Assigns pair wise vector product to SSVectorBase.
508  	   template < class S, class T >
509  	   SSVectorBase<R>& assignPWproduct4setup(const SSVectorBase<S>& x, const SSVectorBase<T>& y);
510  	
511  	   /// Assigns \f$x^T \cdot A\f$ to SSVectorBase.
512  	   template < class S, class T >
513  	   SSVectorBase<R>& assign2product(const SSVectorBase<S>& x, const SVSetBase<T>& A);
514  	
515  	   /// Assigns SSVectorBase to \f$A \cdot x\f$ for a setup \p x.
516  	   template < class S, class T >
517  	   SSVectorBase<R>& assign2product4setup(const SVSetBase<S>& A, const SSVectorBase<T>& x,
518  	                                         Timer* timeSparse, Timer* timeFull, int& nCallsSparse, int& nCallsFull);
519  	
520  	public:
521  	
522  	   /// Assigns SSVectorBase to \f$A \cdot x\f$ thereby setting up \p x.
523  	   template < class S, class T >
524  	   SSVectorBase<R>& assign2productAndSetup(const SVSetBase<S>& A, SSVectorBase<T>& x);
525  	
526  	   /// Maximum absolute value, i.e., infinity norm.
527  	   R maxAbs() const
528  	   {
529  	      if(isSetup())
530  	      {
531  	         R maxabs = 0;
532  	
533  	         for(int i = 0; i < num; ++i)
534  	         {
535  	            R x = spxAbs(VectorBase<R>::val[idx[i]]);
536  	
537  	            if(x > maxabs)
538  	               maxabs = x;
539  	         }
540  	
541  	         return maxabs;
542  	      }
543  	      else
544  	         return VectorBase<R>::maxAbs();
545  	   }
546  	
547  	   /// Squared euclidian norm.
548  	   R length2() const
549  	   {
550  	      R x = 0;
551  	
552  	      if(isSetup())
553  	      {
554  	         for(int i = 0; i < num; ++i)
555  	            x += VectorBase<R>::val[idx[i]] * VectorBase<R>::val[idx[i]];
556  	      }
557  	      else
558  	         x = VectorBase<R>::length2();
559  	
560  	      return x;
561  	   }
562  	
563  	   /// Floating point approximation of euclidian norm (without any approximation guarantee).
564  	   R length() const
565  	   {
566  	      return spxSqrt(R(length2()));
567  	   }
568  	
569  	   ///@}
570  	
571  	   // ------------------------------------------------------------------------------------------------------------------
572  	   /**@name Miscellaneous */
573  	   ///@{
574  	
575  	   /// Dimension of VectorBase.
576  	   int dim() const
577  	   {
578  	      return VectorBase<R>::dim();
579  	   }
580  	
581  	   /// Resets dimension to \p newdim.
582  	   void reDim(int newdim)
583  	   {
584  	      for(int i = IdxSet::size() - 1; i >= 0; --i)
585  	      {
586  	         if(index(i) >= newdim)
587  	            remove(i);
588  	      }
589  	
590  	      VectorBase<R>::reDim(newdim);
591  	      setMax(VectorBase<R>::memSize() + 1);
592  	
593  	      assert(isConsistent());
594  	   }
595  	
596  	   /// Sets number of nonzeros (thereby unSetup SSVectorBase).
597  	   void setSize(int n)
598  	   {
599  	      assert(n >= 0);
600  	      assert(n <= IdxSet::max());
601  	
602  	      unSetup();
603  	      num = n;
604  	   }
605  	
606  	   /// Resets memory consumption to \p newsize.
607  	   void reMem(int newsize)
608  	   {
609  	      VectorBase<R>::reSize(newsize);
610  	      assert(isConsistent());
611  	
612  	      setMax(VectorBase<R>::memSize() + 1);
613  	   }
614  	
615  	   /// Clears vector.
616  	   void clear()
617  	   {
618  	      if(isSetup())
619  	      {
620  	         for(int i = 0; i < num; ++i)
621  	            VectorBase<R>::val[idx[i]] = 0;
622  	      }
623  	      else
624  	         VectorBase<R>::clear();
625  	
626  	      IdxSet::clear();
627  	      setupStatus = true;
628  	
629  	      assert(isConsistent());
630  	   }
631  	
632  	   /// consistency check.
633  	   bool isConsistent() const
634  	   {
635  	#ifdef ENABLE_CONSISTENCY_CHECKS
636  	
637  	      if(VectorBase<R>::dim() > IdxSet::max())
638  	         return SPX_MSG_INCONSISTENT("SSVectorBase");
639  	
640  	      if(VectorBase<R>::dim() < IdxSet::dim())
641  	         return SPX_MSG_INCONSISTENT("SSVectorBase");
642  	
643  	      if(isSetup())
644  	      {
645  	         for(int i = 0; i < VectorBase<R>::dim(); ++i)
646  	         {
647  	            int j = pos(i);
648  	
649  	            if(j < 0 && spxAbs(VectorBase<R>::val[i]) > 0)
650  	            {
651  	               SPX_MSG_ERROR(std::cerr << "ESSVEC01 i = " << i
652  	                             << "\tidx = " << j
653  	                             << "\tval = " << std::setprecision(16) << VectorBase<R>::val[i]
654  	                             << std::endl;)
655  	
656  	               return SPX_MSG_INCONSISTENT("SSVectorBase");
657  	            }
658  	         }
659  	      }
660  	
661  	      return VectorBase<R>::isConsistent() && IdxSet::isConsistent();
662  	#else
663  	      return true;
664  	#endif
665  	   }
666  	
667  	   ///@}
668  	
669  	   // ------------------------------------------------------------------------------------------------------------------
670  	   /**@name Constructors / Destructors */
671  	   ///@{
672  	
673  	   /// Default constructor.
674  	   explicit SSVectorBase<R>(int p_dim, std::shared_ptr<Tolerances> tol = nullptr)
675  	      : VectorBase<R>(p_dim)
(1) Event write_constant_to_parm_in_call: Called function writes 0 to a dereference of parameter "this". [details]
676  	      , IdxSet()
677  	      , setupStatus(true)
678  	   {
(2) Event cond_true: Condition "p_dim < 1", taking true branch.
679  	      len = (p_dim < 1) ? 1 : p_dim;
680  	      spx_alloc(idx, len);
681  	      VectorBase<R>::clear();
682  	      _tolerances = tol;
683  	
684  	      assert(isConsistent());
685  	   }
686  	
687  	   /// Copy constructor.
688  	   template < class S >
689  	   SSVectorBase(const SSVectorBase<S>& vec)
690  	      : VectorBase<R>(vec)
691  	      , IdxSet()
692  	      , setupStatus(vec.setupStatus)
693  	   {
694  	      len = (vec.dim() < 1) ? 1 : vec.dim();
695  	      spx_alloc(idx, len);
696  	      IdxSet::operator=(vec);
697  	      _tolerances = vec._tolerances;
698  	
699  	      assert(isConsistent());
700  	   }
701  	
702  	   /// Copy constructor.
703  	   /** The redundancy with the copy constructor below is necessary since otherwise the compiler doesn't realize that it
704  	    *  could use the more general one with S = R and generates a shallow copy constructor.
705  	    */
706  	   SSVectorBase(const SSVectorBase<R>& vec)
707  	      : VectorBase<R>(vec)
708  	      , IdxSet()
709  	      , setupStatus(vec.setupStatus)
710  	   {
711  	      len = (vec.dim() < 1) ? 1 : vec.dim();
712  	      spx_alloc(idx, len);
713  	      IdxSet::operator=(vec);
714  	      _tolerances = vec._tolerances;
715  	
716  	      assert(isConsistent());
717  	   }
718  	
719  	   /// Constructs nonsetup copy of \p vec.
720  	   template < class S >
721  	   explicit SSVectorBase<R>(const VectorBase<S>& vec)
722  	      : VectorBase<R>(vec)
723  	      , IdxSet()
724  	      , setupStatus(false)
725  	   {
726  	      len = (vec.dim() < 1) ? 1 : vec.dim();
727  	      spx_alloc(idx, len);
728  	
729  	      assert(isConsistent());
730  	   }
731  	
732  	   /// Sets up \p rhs vector, and assigns it.
733  	   template < class S >
734  	   void setup_and_assign(SSVectorBase<S>& rhs)
735  	   {
736  	      clear();
737  	      setMax(rhs.max());
738  	      VectorBase<R>::reDim(rhs.dim());
739  	      _tolerances = rhs.tolerances();
740  	
741  	      if(rhs.isSetup())
742  	      {
743  	         IdxSet::operator=(rhs);
744  	
745  	         for(int i = size() - 1; i >= 0; --i)
746  	         {
747  	            int j  = index(i);
748  	            VectorBase<R>::val[j] = rhs.val[j];
749  	         }
750  	      }
751  	      else
752  	      {
753  	         int d = rhs.dim();
754  	         num = 0;
755  	
756  	         for(int i = 0; i < d; ++i)
757  	         {
758  	            if(rhs.val[i] != 0)
759  	            {
760  	               if(spxAbs(rhs.val[i]) > this->getEpsilon())
761  	               {
762  	                  rhs.idx[num] = i;
763  	                  idx[num] = i;
764  	                  VectorBase<R>::val[i] = rhs.val[i];
765  	                  num++;
766  	               }
767  	               else
768  	                  rhs.val[i] = 0;
769  	            }
770  	         }
771  	
772  	         rhs.num = num;
773  	         rhs.setupStatus = true;
774  	      }
775  	
776  	      setupStatus = true;
777  	
778  	      assert(rhs.isConsistent());
779  	      assert(isConsistent());
780  	   }
781  	
782  	   /// Assigns only the elements of \p rhs.
783  	   template < class S >
784  	   SSVectorBase<R>& assign(const SVectorBase<S>& rhs);
785  	
786  	   /// Assignment operator.
787  	   template < class S >
788  	   SSVectorBase<R>& operator=(const SSVectorBase<S>& rhs)
789  	   {
790  	      assert(rhs.isConsistent());
791  	
792  	      if(this != &rhs)
793  	      {
794  	         clear();
795  	         _tolerances = rhs._tolerances;
796  	         setMax(rhs.max());
797  	         VectorBase<R>::reDim(rhs.dim());
798  	
799  	         if(rhs.isSetup())
800  	         {
801  	            IdxSet::operator=(rhs);
802  	
803  	            for(int i = size() - 1; i >= 0; --i)
804  	            {
805  	               int j = index(i);
806  	               VectorBase<R>::val[j] = rhs.val[j];
807  	            }
808  	         }
809  	         else
810  	         {
811  	            int d = rhs.dim();
812  	            num = 0;
813  	
814  	            for(int i = 0; i < d; ++i)
815  	            {
816  	               if(spxAbs(rhs.val[i]) > this->getEpsilon())
817  	               {
818  	                  VectorBase<R>::val[i] = rhs.val[i];
819  	                  idx[num] = i;
820  	                  num++;
821  	               }
822  	            }
823  	         }
824  	
825  	         setupStatus = true;
826  	      }
827  	
828  	      assert(isConsistent());
829  	
830  	      return *this;
831  	   }
832  	
833  	   /// Assignment operator.
834  	   SSVectorBase<R>& operator=(const SSVectorBase<R>& rhs)
835  	   {
836  	      assert(rhs.isConsistent());
837  	
838  	      if(this != &rhs)
839  	      {
840  	         clear();
841  	         _tolerances = rhs._tolerances;
842  	         setMax(rhs.max());
843  	         VectorBase<R>::reDim(rhs.dim());
844  	
845  	         if(rhs.isSetup())
846  	         {
847  	            IdxSet::operator=(rhs);
848  	
849  	            for(int i = size() - 1; i >= 0; --i)
850  	            {
851  	               int j = index(i);
852  	               VectorBase<R>::val[j] = rhs.val[j];
853  	            }
854  	         }
855  	         else
856  	         {
857  	            num = 0;
858  	
859  	            for(int i = 0; i < rhs.dim(); ++i)
860  	            {
861  	               if(spxAbs(rhs.val[i]) > this->getEpsilon())
862  	               {
863  	                  VectorBase<R>::val[i] = rhs.val[i];
864  	                  idx[num] = i;
865  	                  num++;
866  	               }
867  	            }
868  	         }
869  	
870  	         setupStatus = true;
871  	      }
872  	
873  	      assert(isConsistent());
874  	
875  	      return *this;
876  	   }
877  	
878  	   /// Assignment operator.
879  	   template < class S >
880  	   SSVectorBase<R>& operator=(const SVectorBase<S>& rhs);
881  	
882  	   /// Assignment operator.
883  	   template < class S >
884  	   SSVectorBase<R>& operator=(const VectorBase<S>& rhs)
885  	   {
886  	      unSetup();
887  	      VectorBase<R>::operator=(rhs);
888  	
889  	      assert(isConsistent());
890  	
891  	      return *this;
892  	   }
893  	
894  	   /// destructor
895  	   ~SSVectorBase()
896  	   {
897  	      if(idx)
898  	         spx_free(idx);
899  	   }
900  	
901  	   ///@}
902  	
903  	private:
904  	
905  	   // ------------------------------------------------------------------------------------------------------------------
906  	   /**@name Private helpers */
907  	   ///@{
908  	
909  	   /// Assignment helper.
910  	   template < class S, class T >
911  	   SSVectorBase<R>& assign2product1(const SVSetBase<S>& A, const SSVectorBase<T>& x);
912  	
913  	   /// Assignment helper.
914  	   template < class S, class T >
915  	   SSVectorBase<R>& assign2productShort(const SVSetBase<S>& A, const SSVectorBase<T>& x);
916  	
917  	   /// Assignment helper.
918  	   template < class S, class T >
919  	   SSVectorBase<R>& assign2productFull(const SVSetBase<S>& A, const SSVectorBase<T>& x);
920  	
921  	   ///@}
922  	};
923  	
924  	} // namespace soplex
925  	#endif // _SSVECTORBASE_H_
926