sparse matrix must be row major like Scilab 5.
[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.getArgs().begin())->getLocation());
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     GenericType* 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
369     /**
370        "in-place" cast into a sparse matrix of comlpex values
371      */
372     void toComplex();
373
374     /* coefficient wise relational operator < between *this sparse matrix and an other.
375        Matrices must have the same dimensions except if one of them is of size (1,1)
376        (i.e. a scalar) : it is then treated as a constant matrix of thre required dimensions.
377
378        @param o other sparse matrix
379
380        @return ptr to a new Sparse matrix where each element is the result of the logical operator
381         '<' between the elements of *this and those of o.
382      */
383     SparseBool* newLessThan(Sparse const&o) const;
384
385     /* coefficient wise relational operator > between *this sparse matrix and an other.
386        Matrices must have the same dimensions except if one of them is of size (1,1)
387        (i.e. a scalar) : it is then treated as a constant matrix of thre required dimensions.
388
389        @param o other sparse matrix
390
391        @return ptr to a new Sparse matrix where each element is the result of the logical operator
392         '>' between the elements of *this and those of o.
393      */
394     SparseBool* newGreaterThan(Sparse const&o) const;
395
396     /* coefficient wise relational operator != between *this sparse matrix and an other.
397        Matrices must have the same dimensions except if one of them is of size (1,1)
398        (i.e. a scalar) : it is then treated as a constant matrix of thre required dimensions.
399
400        @param o other sparse matrix
401
402        @return ptr to a new Sparse matrix where each element is the result of the logical operator
403         '!=' between the elements of *this and those of o.
404      */
405     SparseBool* newNotEqualTo(Sparse const&o) const;
406
407     /* coefficient wise relational operator <= between *this sparse matrix and an other.
408        Matrices must have the same dimensions except if one of them is of size (1,1)
409        (i.e. a scalar) : it is then treated as a constant matrix of thre required dimensions.
410
411        Do not use this function is possible as the result will be dense because
412        0. <= 0. is true, hence the result matrix will hold a non default value (i.e. true)
413        for each pair of default values (0.) of the sparse arguments !
414
415        @param o other sparse matrix
416
417        @return ptr to a new Sparse matrix where each element is the result of the logical operator
418         '<=' between the elements of *this and those of o.
419      */
420     SparseBool* newLessOrEqual(Sparse const&o) const;
421
422     /* coefficient wise relational operator >= between *this sparse matrix and an other.
423        Matrices must have the same dimensions except if one of them is of size (1,1)
424        (i.e. a scalar) : it is then treated as a constant matrix of thre required dimensions.
425
426        Do not use this function is possible as the result will be dense because
427        0. >= 0. is true, hence the result matrix will hold a non default value (i.e. true)
428        for each pair of default values (0.) of the sparse arguments !
429
430        @param o other sparse matrix
431
432        @return ptr to a new Sparse matrix where each element is the result of the logical operator
433         '>=' between the elements of *this and those of o.
434      */
435     SparseBool* newGreaterOrEqual(Sparse const&o) const;
436
437     /* coefficient wise relational operator == between *this sparse matrix and an other.
438        Matrices must have the same dimensions except if one of them is of size (1,1)
439        (i.e. a scalar) : it is then treated as a constant matrix of thre required dimensions.
440
441        Do not use this function is possible as the result will be dense because
442        0. == 0. is true, hence the result matrix will hold a non default value (i.e. true)
443        for each pair of default values (0.) of the sparse arguments !
444
445        @param o other sparse matrix
446
447        @return ptr to a new Sparse matrix where each element is the result of the logical operator
448         '==' between the elements of *this and those of o.
449      */
450     SparseBool* newEqualTo(Sparse const&o) const;
451
452     /**
453        output 1-base column numbers of the non zero elements
454        @param out : ptr used as an output iterator over double values
455        @return past-the-end output iterator after ouput is done
456      */
457     double* outputCols(double* out) const;
458
459     /**
460        output real and imaginary values of the non zero elements
461        @param outReal : ptr used as an output iterator over double values for real values
462        @param outImag : ptr used as an output iterator over double values for imaginary values if any
463        @return pair of past-the-end output iterators after ouput is done
464      */
465     std::pair<double*, double*> outputValues(double* outReal, double* outImag)const;
466
467     /**
468        ouput rows and afterwards columns of the non zero elements
469        @param out : ptr used as an output iterator over double values
470        @return past-the-end output iterators after ouput is done
471      */
472     int* outputRowCol(int* out)const;
473
474     /**
475        @param dest Double to be filled with values from the current sparse matrix.
476      */
477     void fill(Double& dest, int r = 0, int c = 0) SPARSE_CONST;
478
479
480     inline ScilabType  getType(void) SPARSE_CONST
481     {
482         return ScilabSparse;
483     }
484
485     inline ScilabId    getId(void) SPARSE_CONST
486     {
487         if (isComplex())
488         {
489             return IdSparseComplex;
490         }
491         return IdSparse;
492     }
493
494
495
496     SparseBool* newLesserThan(Sparse const&o);
497
498     typedef Eigen::SparseMatrix<double, Eigen::RowMajor>                 RealSparse_t;
499     typedef Eigen::SparseMatrix<std::complex<double>, Eigen::RowMajor>   CplxSparse_t;
500     /**
501        One and only one of the args should be 0.
502        @param realSp ptr to an Eigen sparse matrix of double values
503        @param cplxSp ptr to an Eigen sparse matrix of std::complex<double> elements
504      */
505     Sparse(RealSparse_t* realSp, CplxSparse_t* cplxSp);
506
507     RealSparse_t* matrixReal;
508     CplxSparse_t* matrixCplx;
509
510 protected :
511 private :
512
513     /** utility function used by constructors
514         @param rows : nb of rows
515         @param cols : nb of columns
516         @param src : Double matrix data source
517         @param : iterator (cf MatrixIterator.hxx) with indices
518         @param n : nb of elements to copy from data source.
519      */
520     template<typename DestIter>
521     void create(int rows, int cols, Double SPARSE_CONST& src, DestIter o, std::size_t n);
522
523     /** utility function used by insert functions conceptually : sp[destTrav]= src[srcTrav]
524         @param src : data source
525         @param SrcTrav : iterator over the data source
526         @param n : nb of elements to copy
527         @param sp : sparse destination matrix
528         @param destTrav : iterator over the data sink (i.e. sp)
529      */
530     template<typename Src, typename SrcTraversal, typename Sz, typename DestTraversal>
531     static bool copyToSparse(Src SPARSE_CONST& src, SrcTraversal srcTrav, Sz n, Sparse& sp, DestTraversal destTrav);
532 };
533 /*
534   Implement sparse boolean matrix
535  */
536 struct EXTERN_AST SparseBool : GenericType
537 {
538
539     /* @param src: Bool matrix to copy into a new sparse matrix
540     **/
541     SparseBool(Bool SPARSE_CONST& src);
542     /* @param src : Bool matrix to copy into a new sparse matrix
543        @param idx : Double matrix to use as indexes to get values from the src
544     **/
545     SparseBool(Bool SPARSE_CONST& src, Double SPARSE_CONST& idx);
546     /* @param src : Bool matrix to copy into a new sparse matrix
547        @param idx : Double matrix to use as indexes to get values from the src
548        @param dims : Double matrix containing the dimensions of the new matrix
549     **/
550     SparseBool(Bool SPARSE_CONST& src, Double SPARSE_CONST& idx, Double SPARSE_CONST& dims);
551     /*
552        @param rows : nb of rows of the new matrix
553        @param rows : nb of columns of the new matrix
554     */
555     SparseBool(int rows, int cols);
556
557     SparseBool(SparseBool const& o);
558
559     bool isSparseBool()
560     {
561         return true;
562     }
563     void finalize();
564
565     bool toString(std::wostringstream& ostr) const;
566     bool toString(std::wostringstream& ostr)
567     {
568         return const_cast<SparseBool const*>(this)->toString(ostr);
569     }
570
571     /* Config management and GenericType methods overrides */
572     SparseBool* clone(void) const;
573     SparseBool* clone(void)
574     {
575         return const_cast<SparseBool const*>(this)->clone();
576     }
577     bool resize(int _iNewRows, int _iNewCols);
578
579     bool reshape(int* _piNewDims, int _iNewDims);
580     bool reshape(int _iNewRows, int _iNewCols);
581
582     SparseBool* insert(typed_list* _pArgs, InternalType* _pSource);
583     SparseBool* insert(typed_list* _pArgs, SparseBool* _pSource);
584     SparseBool* remove(typed_list* _pArgs);
585
586     bool append(int _iRows, int _iCols, SparseBool SPARSE_CONST* _poSource);
587
588     static InternalType* insertNew(typed_list* _pArgs, InternalType* _pSource);
589     SparseBool* extract(int _iSeqCount, int* _piSeqCoord, int* _piMaxDim, int* _piDimSize, bool _bAsVector) SPARSE_CONST;
590     InternalType* extract(typed_list* _pArgs);
591
592     virtual bool invoke(typed_list & in, optional_list &/*opt*/, int /*_iRetCount*/, typed_list & out, ast::ConstVisitor & /*execFunc*/, const ast::CallExp & e)
593     {
594         if (in.size() == 0)
595         {
596             out.push_back(this);
597         }
598         else
599         {
600             InternalType * _out = extract(&in);
601             if (!_out)
602             {
603                 std::wostringstream os;
604                 os << _W("Invalid index.\n");
605                 throw ast::ScilabError(os.str(), 999, (*e.getArgs().begin())->getLocation());
606             }
607             out.push_back(_out);
608         }
609
610         return true;
611     }
612
613     virtual bool isInvokable() const
614     {
615         return true;
616     }
617
618     virtual bool hasInvokeOption() const
619     {
620         return false;
621     }
622
623     virtual int getInvokeNbIn()
624     {
625         return -1;
626     }
627
628     virtual int getInvokeNbOut()
629     {
630         return 1;
631     }
632
633     bool transpose(InternalType *& out)
634     {
635         out = new SparseBool(new BoolSparse_t(matrixBool->transpose()));
636         return true;
637     }
638
639     /** @return the nb of non zero values.
640      */
641     std::size_t nbTrue() const;
642     /* @param i row of the current sparse matrix
643        @return the nb of non zero values in row i
644      */
645     std::size_t nbTrue(std::size_t i) const;
646
647     int* getNbItemByRow(int* _piNbItemByRows);
648     int* getColPos(int* _piColPos);
649     /**
650        output 1-base column numbers of the non zero elements
651        @param out : ptr used as an output iterator over double values
652        @return past-the-end output iterator after ouput is done
653      */
654
655     int* outputRowCol(int* out)const;
656
657     bool operator==(const InternalType& it) SPARSE_CONST;
658     bool operator!=(const InternalType& it) SPARSE_CONST;
659
660     /* return type as string ( double, int, cell, list, ... )*/
661     virtual std::wstring getTypeStr() SPARSE_CONST {return std::wstring(L"boolean sparse");}
662     /* return type as short string ( s, i, ce, l, ... )*/
663     virtual std::wstring getShortTypeStr() SPARSE_CONST {return std::wstring(L"spb");}
664
665     inline ScilabType getType(void) SPARSE_CONST
666     {
667         return ScilabSparseBool;
668     }
669
670     inline ScilabId getId(void) SPARSE_CONST
671     {
672         return IdSparseBool;
673     }
674
675     inline bool isScalar()
676     {
677         return (getRows() == 1 && getCols() == 1);
678     }
679
680     bool neg(InternalType *& out)
681     {
682         SparseBool * _out = new SparseBool(getRows(), getCols());
683         type_traits::neg(getRows(), getCols(), matrixBool, _out->matrixBool);
684         _out->finalize();
685         out = _out;
686         return true;
687     }
688
689     void whoAmI() SPARSE_CONST;
690
691     bool get(int r, int c) SPARSE_CONST;
692     bool get(int _iIndex) SPARSE_CONST
693     {
694         return get(_iIndex % m_iRows, _iIndex / m_iRows);
695     }
696
697     bool set(int r, int c, bool b, bool _bFinalize = true) SPARSE_CONST;
698     bool set(int _iIndex, bool b, bool _bFinalize = true) SPARSE_CONST
699     {
700         return set(_iIndex % m_iRows, _iIndex / m_iRows, b, _bFinalize);
701     }
702
703     void fill(Bool& dest, int r = 0, int c = 0) SPARSE_CONST;
704
705     Sparse* newOnes() const;
706     SparseBool* newNotEqualTo(SparseBool const&o) const;
707     SparseBool* newEqualTo(SparseBool const&o) const;
708
709     SparseBool* newLogicalOr(SparseBool const&o) const;
710     SparseBool* newLogicalAnd(SparseBool const&o) const;
711
712     typedef Eigen::SparseMatrix<bool, Eigen::RowMajor> BoolSparse_t;
713     SparseBool(BoolSparse_t* o);
714     BoolSparse_t* matrixBool;
715
716 private:
717     template<typename DestIter>
718     void create(int rows, int cols, Bool SPARSE_CONST& src, DestIter o, std::size_t n);
719
720 };
721 template<typename T>
722 struct SparseTraits
723 {
724     typedef types::Sparse type;
725 };
726 template<>
727 struct SparseTraits<types::Bool>
728 {
729     typedef types::SparseBool type;
730 };
731 template<>
732 struct SparseTraits<types::SparseBool>
733 {
734     typedef types::SparseBool type;
735 };
736
737 }
738
739 #endif /* !__SPARSE_HH__ */