390d38bbc25f884c2ba8775140eed0495109f2aa
[scilab.git] / scilab / modules / ast / includes / types / sparse.hxx
1 /*
2  *  Scilab ( http://www.scilab.org/ ) - This file is part of Scilab
3  *  Copyright (C) 2010-2010 - DIGITEO - Bernard Hugueney
4  *
5  *  This file must be used under the terms of the CeCILL.
6  *  This source file is licensed as described in the file COPYING, which
7  *  you should have received as part of this distribution.  The terms
8  *  are also available at
9  *  http://www.cecill.info/licences/Licence_CeCILL_V2-en.txt
10  *
11  */
12
13 #ifndef __SPARSE_HH__
14 #define __SPARSE_HH__
15
16 #include <Eigen/Sparse>
17 #include <complex>
18 #include "double.hxx"
19 #include "bool.hxx"
20
21 #define SPARSE_CONST
22
23 namespace types
24 {
25 /* Utility function to create a new var on the heap from another type
26  */
27 template<typename Dest, typename Arg>
28 Dest* create_new(Arg const&);
29
30 struct SparseBool;
31
32 /**
33    Sparse is a wrapper over Eigen sparse matrices templates for either double or std::complex<double> values.
34  */
35 struct EXTERN_AST Sparse : GenericType
36 {
37     virtual ~Sparse();
38     /* @param src: Double matrix to copy into a new sparse matrix
39     **/
40     Sparse(Double SPARSE_CONST& src);
41     /* @param src : Double matrix to copy into a new sparse matrix
42        @param idx : Double matrix to use as indexes to get values from the src
43     **/
44     Sparse(Double SPARSE_CONST& src, Double SPARSE_CONST& idx);
45     /* @param src : Double matrix to copy into a new sparse matrix
46        @param idx : Double matrix to use as indexes to get values from the src
47        @param dims : Double matrix containing the dimensions of the new matrix
48     **/
49     Sparse(Double SPARSE_CONST& src, Double SPARSE_CONST& idx, Double SPARSE_CONST& dims);
50     /*
51       @param rows : nb of rows of the new matrix
52       @param rows : nb of columns of the new matrix
53       @param cplx : if the new matrix contains complex numbers
54     **/
55     Sparse(int rows, int cols, bool cplx = false);
56     Sparse(Sparse const& o);
57     /* cf. adj2sp()
58       @param xadj : adjacency matrix for the new matrix
59       @param adjncy : adjacency matrix (row indexes) for the new matrix
60       @param src : data for the new matrix
61       @param r : nb of rows for the new matrix
62       @param c : nb of columns for the new matrix
63     **/
64     Sparse(Double SPARSE_CONST& xadj, Double SPARSE_CONST& adjncy, Double SPARSE_CONST& src, std::size_t r, std::size_t c);
65
66
67     bool isSparse()
68     {
69         return true;
70     }
71     void finalize();
72
73     /*data management member function defined for compatibility with the Double API*/
74     bool set(int _iRows, int _iCols, double _dblReal, bool _bFinalize = true);
75     bool set(int _iIndex, double _dblReal, bool _bFinalize = true)
76     {
77         return set(_iIndex % m_iRows, _iIndex / m_iRows, _dblReal, _bFinalize);
78     }
79
80     bool set(int _iRows, int _iCols, std::complex<double> v, bool _bFinalize = true);
81     bool set(int _iIndex, std::complex<double> v, bool _bFinalize = true)
82     {
83         return set(_iIndex % m_iRows, _iIndex / m_iRows, v, _bFinalize);
84     }
85     /*
86       set non zero values to 1.
87     **/
88     bool one_set();
89     /* get real value at coords (r,c)
90     **/
91     double getReal(int r, int c) const;
92     double getReal(int _iIndex) const
93     {
94         return getReal(_iIndex % m_iRows, _iIndex / m_iRows);
95     }
96
97     double get(int r, int c) const;
98     double get(int _iIndex) const
99     {
100         return get(_iIndex % m_iRows, _iIndex / m_iRows);
101     }
102
103     std::complex<double> getImg(int r, int c) const;
104     std::complex<double> getImg(int _iIndex) const
105     {
106         return getImg(_iIndex % m_iRows, _iIndex / m_iRows);
107     }
108
109     /* return true if matrix contains complex numbers, false otherwise.
110     **/
111     bool isComplex() const;
112     // overload of GenericType methode.
113     bool isComplex()
114     {
115         // force const to call isComplex const method.
116         const Sparse* sp = this;
117         return sp->isComplex();
118     }
119
120     inline bool isScalar()
121     {
122         return (getRows() == 1 && getCols() == 1);
123     }
124     /* clear all the values of the matrix to 0. (or 0.+0.i if complex)
125     **/
126     bool zero_set();
127
128     /*
129       Config management and GenericType methods overrides
130     */
131     void whoAmI() SPARSE_CONST;
132     bool isExtract() const;
133     Sparse* clone(void) const;
134     Sparse* clone(void)
135     {
136         return const_cast<Sparse const*>(this)->clone();
137     }
138     bool toString(std::wostringstream& ostr) const;
139     bool toString(std::wostringstream& ostr)
140     {
141         return const_cast<Sparse const*>(this)->toString(ostr);
142     }
143
144     /* post condition: dimensions are at least _iNewRows, _iNewCols
145        preserving existing data.
146        If dimensions where already >=, this is a no-op.
147
148        @param _iNewRows new minimum nb of rows
149        @param _iNewCols new minimum nb of cols
150        @return true upon succes, false otherwise.
151      */
152     bool resize(int _iNewRows, int _iNewCols);
153     /* post condition: new total size must be equal to the old size.
154                        Two dimensions maximum.
155
156        @param _iNewRows new nb of rows
157        @param _iNewCols new nb of cols
158        @param _piNewDims new nb of dimension
159        @param _iNewDims new size for each dimension
160        @return true upon succes, false otherwise.
161     */
162     bool reshape(int* _piNewDims, int _iNewDims);
163     bool reshape(int _iNewRows, int _iNewCols);
164     /*
165       insert _iSeqCount elements from _poSource at coords given by _piSeqCoord (max in _piMaxDim).
166       coords are considered 1D if _bAsVector, 2D otherwise.
167       @param _iSeqCount nb of elts to insert
168       @param _piSeqCoord dest coords
169       @param _poSource src
170       @param  _bAsVector if _piSeqCoord contains 1D coords.
171      */
172     Sparse* insert(typed_list* _pArgs, InternalType* _pSource);
173     Sparse* insert(typed_list* _pArgs, Sparse* _pSource);
174
175     Sparse* remove(typed_list* _pArgs);
176
177     static InternalType* insertNew(typed_list* _pArgs, InternalType* _pSource);
178
179     /* append _poSource from coords _iRows, _iCols
180        @param _iRows row to append from
181        @param _iCols col to append from
182        @param _poSource src data to append
183      */
184     bool append(int r, int c, types::Sparse SPARSE_CONST* src);
185
186     /*
187       extract a submatrix
188       @param _iSeqCount nb of elts to extract
189       @param _piSeqCoord src coords
190       @param _piMaxDim max coords
191       @param _piDimSize size of the extracted matrix
192       @param  _bAsVector if _piSeqCoord contains 1D coords.
193
194      */
195     InternalType* extract(typed_list* _pArgs);
196
197     virtual bool invoke(typed_list & in, optional_list & /*opt*/, int /*_iRetCount*/, typed_list & out, ast::ConstVisitor & /*execFunc*/, const ast::CallExp & e)
198     {
199         if (in.size() == 0)
200         {
201             out.push_back(this);
202         }
203         else
204         {
205             InternalType * _out = extract(&in);
206             if (!_out)
207             {
208                 std::wostringstream os;
209                 os << _W("Invalid index.\n");
210                 throw ast::ScilabError(os.str(), 999, e.getFirstLocation());
211             }
212             out.push_back(_out);
213         }
214
215         return true;
216     }
217
218     virtual bool isInvokable() const
219     {
220         return true;
221     }
222
223     virtual bool hasInvokeOption() const
224     {
225         return false;
226     }
227
228     virtual int getInvokeNbIn()
229     {
230         return -1;
231     }
232
233     virtual int getInvokeNbOut()
234     {
235         return 1;
236     }
237
238     Sparse* extract(int _iSeqCount, int* _piSeqCoord, int* _piMaxDim, int* _piDimSize, bool _bAsVector) SPARSE_CONST;
239
240     /*
241        change the sign (inplace).
242      */
243     void opposite();
244
245     /*
246       compares with an other value for equality (same nb of elts, with same values)
247       TODO: should it handle other types ?
248      */
249     bool operator==(const InternalType& it) SPARSE_CONST;
250     /*
251       compares with an other value for inequality (same nb of elts, with same values)
252       TODO: should it handle other types ?
253      */
254     bool operator!=(const InternalType& it) SPARSE_CONST
255     {
256         return !(*this == it);
257     }
258
259
260     /* return type as string ( double, int, cell, list, ... )*/
261     virtual std::wstring         getTypeStr() SPARSE_CONST {return std::wstring(L"sparse");}
262     /* return type as short string ( s, i, ce, l, ... ), as in overloading macros*/
263     virtual std::wstring         getShortTypeStr() SPARSE_CONST {return std::wstring(L"sp");}
264
265     /* create a new sparse matrix containing the result of an addition
266        @param o other matrix to add
267        @return ptr to the new matrix, 0 in case of failure
268      */
269     Sparse* add(Sparse const& o) const;
270
271     /* create a new sparse matrix containing the result of an addition
272        @param d scalar to add
273        @return ptr to the new matrix, 0 in case of failure
274      */
275     Double* add(double d) const;
276
277     /* create a new sparse matrix containing the result of an addition
278        @param c complex  to add
279        @return ptr to the new matrix, 0 in case of failure
280      */
281     Double* add(std::complex<double> c) const;
282
283
284
285
286     /* create a new sparse matrix containing the result of a substraction
287        @param o other matrix to substract
288        @return ptr to the new matrix, 0 in case of failure
289      */
290     Sparse* substract(Sparse const& o) const;
291
292     /* create a new sparse matrix containing the result of an subtraction
293        @param d scalar to subtract
294        @return ptr to the new matrix, 0 in case of failure
295      */
296     Double* substract(double d) const;
297
298     /* create a new sparse matrix containing the result of an subtraction
299        @param c scalar to subtract
300        @return ptr to the new matrix, 0 in case of failure
301      */
302     Double* substract(std::complex<double> c) const;
303
304
305     /* create a new sparse matrix containing the result of a multiplication
306        @param o other matrix to substract
307        @return ptr to the new matrix, 0 in case of failure
308      */
309     Sparse* multiply(Sparse const& o) const;
310
311     /* create a new sparse matrix containing the result of an multiplication
312        @param s scalar to multiply by
313        @return ptr to the new matrix, 0 in case of failure
314      */
315     Sparse* multiply(double s) const;
316
317     /* create a new sparse matrix containing the result of an multiplication
318        @param c scalar to subtract
319        @return ptr to the new matrix, 0 in case of failure
320      */
321     Sparse* multiply(std::complex<double> c) const;
322
323     /* create a new matrix containing the result of an .*
324        @param o sparse matrix to .*
325        @return ptr to the new matrix, 0 in case of failure
326      */
327     Sparse* dotMultiply(Sparse SPARSE_CONST& o) const;
328
329     /* create a new matrix containing the result of an ./
330       @param o sparse matrix to ./
331       @return ptr to the new matrix, 0 in case of failure
332     */
333     Sparse* dotDivide(Sparse SPARSE_CONST& o) const;
334
335     bool neg(InternalType *& out);
336
337     bool transpose(InternalType *& out)
338     {
339         out = new Sparse(matrixReal ? new RealSparse_t(matrixReal->transpose()) : 0, matrixCplx ? new CplxSparse_t(matrixCplx->transpose()) : 0);
340         return true;
341     }
342
343     bool adjoint(InternalType *& out)
344     {
345         out = new Sparse(matrixReal ? new RealSparse_t(matrixReal->adjoint()) : 0, matrixCplx ? new CplxSparse_t(matrixCplx->adjoint()) : 0);
346         return true;
347     }
348
349     /** create a new sparse matrix containing the non zero values set to 1.
350        equivalent but faster than calling one_set() on a new copy of the
351        current matrix.
352        @return ptr to the new matrix, 0 in case of failure
353      */
354     Sparse* newOnes() const;
355
356     Sparse* newReal() const;
357     /** @return the nb of non zero values.
358      */
359     std::size_t nonZeros() const;
360
361     /* @param i row of the current sparse matrix
362        @return the nb of non zero values in row i
363      */
364     std::size_t nonZeros(std::size_t i) const;
365
366     int* getNbItemByRow(int* _piNbItemByRows);
367     int* getColPos(int* _piColPos);
368     int  getNbItemByCol(int* _piNbItemByCols, int* _piRowPos);
369
370     /**
371        "in-place" cast into a sparse matrix of comlpex values
372      */
373     void toComplex();
374
375     /* coefficient wise relational operator < between *this sparse matrix and an other.
376        Matrices must have the same dimensions except if one of them is of size (1,1)
377        (i.e. a scalar) : it is then treated as a constant matrix of thre required dimensions.
378
379        @param o other sparse matrix
380
381        @return ptr to a new Sparse matrix where each element is the result of the logical operator
382         '<' between the elements of *this and those of o.
383      */
384     SparseBool* newLessThan(Sparse const&o) const;
385
386     /* coefficient wise relational operator > between *this sparse matrix and an other.
387        Matrices must have the same dimensions except if one of them is of size (1,1)
388        (i.e. a scalar) : it is then treated as a constant matrix of thre required dimensions.
389
390        @param o other sparse matrix
391
392        @return ptr to a new Sparse matrix where each element is the result of the logical operator
393         '>' between the elements of *this and those of o.
394      */
395     SparseBool* newGreaterThan(Sparse const&o) const;
396
397     /* coefficient wise relational operator != between *this sparse matrix and an other.
398        Matrices must have the same dimensions except if one of them is of size (1,1)
399        (i.e. a scalar) : it is then treated as a constant matrix of thre required dimensions.
400
401        @param o other sparse matrix
402
403        @return ptr to a new Sparse matrix where each element is the result of the logical operator
404         '!=' between the elements of *this and those of o.
405      */
406     SparseBool* newNotEqualTo(Sparse const&o) const;
407
408     /* coefficient wise relational operator <= between *this sparse matrix and an other.
409        Matrices must have the same dimensions except if one of them is of size (1,1)
410        (i.e. a scalar) : it is then treated as a constant matrix of thre required dimensions.
411
412        Do not use this function is possible as the result will be dense because
413        0. <= 0. is true, hence the result matrix will hold a non default value (i.e. true)
414        for each pair of default values (0.) of the sparse arguments !
415
416        @param o other sparse matrix
417
418        @return ptr to a new Sparse matrix where each element is the result of the logical operator
419         '<=' between the elements of *this and those of o.
420      */
421     SparseBool* newLessOrEqual(Sparse const&o) const;
422
423     /* coefficient wise relational operator >= between *this sparse matrix and an other.
424        Matrices must have the same dimensions except if one of them is of size (1,1)
425        (i.e. a scalar) : it is then treated as a constant matrix of thre required dimensions.
426
427        Do not use this function is possible as the result will be dense because
428        0. >= 0. is true, hence the result matrix will hold a non default value (i.e. true)
429        for each pair of default values (0.) of the sparse arguments !
430
431        @param o other sparse matrix
432
433        @return ptr to a new Sparse matrix where each element is the result of the logical operator
434         '>=' between the elements of *this and those of o.
435      */
436     SparseBool* newGreaterOrEqual(Sparse const&o) const;
437
438     /* coefficient wise relational operator == between *this sparse matrix and an other.
439        Matrices must have the same dimensions except if one of them is of size (1,1)
440        (i.e. a scalar) : it is then treated as a constant matrix of thre required dimensions.
441
442        Do not use this function is possible as the result will be dense because
443        0. == 0. is true, hence the result matrix will hold a non default value (i.e. true)
444        for each pair of default values (0.) of the sparse arguments !
445
446        @param o other sparse matrix
447
448        @return ptr to a new Sparse matrix where each element is the result of the logical operator
449         '==' between the elements of *this and those of o.
450      */
451     SparseBool* newEqualTo(Sparse const&o) const;
452
453     /**
454        output 1-base column numbers of the non zero elements
455        @param out : ptr used as an output iterator over double values
456        @return past-the-end output iterator after ouput is done
457      */
458     double* outputCols(double* out) const;
459
460     /**
461        output real and imaginary values of the non zero elements
462        @param outReal : ptr used as an output iterator over double values for real values
463        @param outImag : ptr used as an output iterator over double values for imaginary values if any
464        @return pair of past-the-end output iterators after ouput is done
465      */
466     std::pair<double*, double*> outputValues(double* outReal, double* outImag)const;
467
468     /**
469        ouput rows and afterwards columns of the non zero elements
470        @param out : ptr used as an output iterator over double values
471        @return past-the-end output iterators after ouput is done
472      */
473     int* outputRowCol(int* out)const;
474
475     /**
476        @param dest Double to be filled with values from the current sparse matrix.
477      */
478     void fill(Double& dest, int r = 0, int c = 0) SPARSE_CONST;
479
480
481     inline ScilabType  getType(void) SPARSE_CONST
482     {
483         return ScilabSparse;
484     }
485
486     inline ScilabId    getId(void) SPARSE_CONST
487     {
488         if (isComplex())
489         {
490             return IdSparseComplex;
491         }
492         return IdSparse;
493     }
494
495
496
497     SparseBool* newLesserThan(Sparse const&o);
498
499     typedef Eigen::SparseMatrix<double, Eigen::RowMajor>                 RealSparse_t;
500     typedef Eigen::SparseMatrix<std::complex<double>, Eigen::RowMajor>   CplxSparse_t;
501     /**
502        One and only one of the args should be 0.
503        @param realSp ptr to an Eigen sparse matrix of double values
504        @param cplxSp ptr to an Eigen sparse matrix of std::complex<double> elements
505      */
506     Sparse(RealSparse_t* realSp, CplxSparse_t* cplxSp);
507
508     RealSparse_t* matrixReal;
509     CplxSparse_t* matrixCplx;
510
511 protected :
512 private :
513
514     /** utility function used by constructors
515         @param rows : nb of rows
516         @param cols : nb of columns
517         @param src : Double matrix data source
518         @param : iterator (cf MatrixIterator.hxx) with indices
519         @param n : nb of elements to copy from data source.
520      */
521     template<typename DestIter>
522     void create(int rows, int cols, Double SPARSE_CONST& src, DestIter o, std::size_t n);
523     void create2(int rows, int cols, Double SPARSE_CONST& src, Double SPARSE_CONST& idx);
524
525     /** utility function used by insert functions conceptually : sp[destTrav]= src[srcTrav]
526         @param src : data source
527         @param SrcTrav : iterator over the data source
528         @param n : nb of elements to copy
529         @param sp : sparse destination matrix
530         @param destTrav : iterator over the data sink (i.e. sp)
531      */
532     template<typename Src, typename SrcTraversal, typename Sz, typename DestTraversal>
533     static bool copyToSparse(Src SPARSE_CONST& src, SrcTraversal srcTrav, Sz n, Sparse& sp, DestTraversal destTrav);
534 };
535 /*
536   Implement sparse boolean matrix
537  */
538 struct EXTERN_AST SparseBool : GenericType
539 {
540
541     /* @param src: Bool matrix to copy into a new sparse matrix
542     **/
543     SparseBool(Bool SPARSE_CONST& src);
544     /* @param src : Bool matrix to copy into a new sparse matrix
545        @param idx : Double matrix to use as indexes to get values from the src
546     **/
547     SparseBool(Bool SPARSE_CONST& src, Double SPARSE_CONST& idx);
548     /* @param src : Bool matrix to copy into a new sparse matrix
549        @param idx : Double matrix to use as indexes to get values from the src
550        @param dims : Double matrix containing the dimensions of the new matrix
551     **/
552     SparseBool(Bool SPARSE_CONST& src, Double SPARSE_CONST& idx, Double SPARSE_CONST& dims);
553     /*
554        @param rows : nb of rows of the new matrix
555        @param rows : nb of columns of the new matrix
556     */
557     SparseBool(int rows, int cols);
558
559     SparseBool(SparseBool const& o);
560
561     bool isSparseBool()
562     {
563         return true;
564     }
565     void finalize();
566
567     bool toString(std::wostringstream& ostr) const;
568     bool toString(std::wostringstream& ostr)
569     {
570         return const_cast<SparseBool const*>(this)->toString(ostr);
571     }
572
573     /* Config management and GenericType methods overrides */
574     SparseBool* clone(void) const;
575     SparseBool* clone(void)
576     {
577         return const_cast<SparseBool const*>(this)->clone();
578     }
579     bool resize(int _iNewRows, int _iNewCols);
580
581     bool reshape(int* _piNewDims, int _iNewDims);
582     bool reshape(int _iNewRows, int _iNewCols);
583
584     SparseBool* insert(typed_list* _pArgs, InternalType* _pSource);
585     SparseBool* insert(typed_list* _pArgs, SparseBool* _pSource);
586     SparseBool* remove(typed_list* _pArgs);
587
588     bool append(int _iRows, int _iCols, SparseBool SPARSE_CONST* _poSource);
589
590     static InternalType* insertNew(typed_list* _pArgs, InternalType* _pSource);
591     SparseBool* extract(int _iSeqCount, int* _piSeqCoord, int* _piMaxDim, int* _piDimSize, bool _bAsVector) SPARSE_CONST;
592     InternalType* extract(typed_list* _pArgs);
593
594     virtual bool invoke(typed_list & in, optional_list &/*opt*/, int /*_iRetCount*/, typed_list & out, ast::ConstVisitor & /*execFunc*/, const ast::CallExp & e)
595     {
596         if (in.size() == 0)
597         {
598             out.push_back(this);
599         }
600         else
601         {
602             InternalType * _out = extract(&in);
603             if (!_out)
604             {
605                 std::wostringstream os;
606                 os << _W("Invalid index.\n");
607                 throw ast::ScilabError(os.str(), 999, e.getFirstLocation());
608             }
609             out.push_back(_out);
610         }
611
612         return true;
613     }
614
615     virtual bool isInvokable() const
616     {
617         return true;
618     }
619
620     virtual bool hasInvokeOption() const
621     {
622         return false;
623     }
624
625     virtual int getInvokeNbIn()
626     {
627         return -1;
628     }
629
630     virtual int getInvokeNbOut()
631     {
632         return 1;
633     }
634
635     bool transpose(InternalType *& out)
636     {
637         out = new SparseBool(new BoolSparse_t(matrixBool->transpose()));
638         return true;
639     }
640
641     /** @return the nb of non zero values.
642      */
643     std::size_t nbTrue() const;
644     /* @param i row of the current sparse matrix
645        @return the nb of non zero values in row i
646      */
647     std::size_t nbTrue(std::size_t i) const;
648
649     int* getNbItemByRow(int* _piNbItemByRows);
650     int* getColPos(int* _piColPos);
651     /**
652        output 1-base column numbers of the non zero elements
653        @param out : ptr used as an output iterator over double values
654        @return past-the-end output iterator after ouput is done
655      */
656
657     int* outputRowCol(int* out)const;
658
659     bool operator==(const InternalType& it) SPARSE_CONST;
660     bool operator!=(const InternalType& it) SPARSE_CONST;
661
662     /* return type as string ( double, int, cell, list, ... )*/
663     virtual std::wstring getTypeStr() SPARSE_CONST {return std::wstring(L"boolean sparse");}
664     /* return type as short string ( s, i, ce, l, ... )*/
665     virtual std::wstring getShortTypeStr() SPARSE_CONST {return std::wstring(L"spb");}
666
667     inline ScilabType getType(void) SPARSE_CONST
668     {
669         return ScilabSparseBool;
670     }
671
672     inline ScilabId getId(void) SPARSE_CONST
673     {
674         return IdSparseBool;
675     }
676
677     inline bool isScalar()
678     {
679         return (getRows() == 1 && getCols() == 1);
680     }
681
682     bool isTrue()
683     {
684         if (nbTrue() == m_iSize)
685         {
686             return true;
687         }
688         return false;
689     }
690
691     bool neg(InternalType *& out)
692     {
693         SparseBool * _out = new SparseBool(getRows(), getCols());
694         type_traits::neg(getRows(), getCols(), matrixBool, _out->matrixBool);
695         _out->finalize();
696         out = _out;
697         return true;
698     }
699
700     void whoAmI() SPARSE_CONST;
701
702     bool get(int r, int c) SPARSE_CONST;
703     bool get(int _iIndex) SPARSE_CONST
704     {
705         return get(_iIndex % m_iRows, _iIndex / m_iRows);
706     }
707
708     bool set(int r, int c, bool b, bool _bFinalize = true) SPARSE_CONST;
709     bool set(int _iIndex, bool b, bool _bFinalize = true) SPARSE_CONST
710     {
711         return set(_iIndex % m_iRows, _iIndex / m_iRows, b, _bFinalize);
712     }
713
714     void fill(Bool& dest, int r = 0, int c = 0) SPARSE_CONST;
715
716     Sparse* newOnes() const;
717     SparseBool* newNotEqualTo(SparseBool const&o) const;
718     SparseBool* newEqualTo(SparseBool const&o) const;
719
720     SparseBool* newLogicalOr(SparseBool const&o) const;
721     SparseBool* newLogicalAnd(SparseBool const&o) const;
722
723     typedef Eigen::SparseMatrix<bool, Eigen::RowMajor> BoolSparse_t;
724     SparseBool(BoolSparse_t* o);
725     BoolSparse_t* matrixBool;
726
727 private:
728     void create2(int rows, int cols, Bool SPARSE_CONST& src, Double SPARSE_CONST& idx);
729 };
730 template<typename T>
731 struct SparseTraits
732 {
733     typedef types::Sparse type;
734 };
735 template<>
736 struct SparseTraits<types::Bool>
737 {
738     typedef types::SparseBool type;
739 };
740 template<>
741 struct SparseTraits<types::SparseBool>
742 {
743     typedef types::SparseBool type;
744 };
745
746 }
747
748 #endif /* !__SPARSE_HH__ */