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