00001
00002
00003
00004
00005
00006
00007
00008
00009
00010
00011
00012
00013
00014
00015
00016
00017
00018 #include <vector>
00019
00020
00021
00022
00023
00024
00025
00026
00027
00028
00029
00030
00031
00032
00041 class SparseView {
00042 public:
00043
00044 Index nrows() const;
00045 Index ncols() const;
00046 Index nnz() const;
00047
00048
00049 Numeric& rw(Index r, Index c);
00050 Numeric ro(Index r, Index c) const;
00051
00052
00053 friend std::ostream& operator<<(std::ostream& os, const SparseView& v);
00054
00055 protected:
00056
00057 SparseView();
00058 SparseView(std::vector<Numeric> *data,
00059 std::vector<Index> *rowind,
00060 std::vector<Index> *colptr,
00061 const Range& rr,
00062 const Range& cr );
00063
00064
00066 std::vector<Numeric> *mdata;
00068 std::vector<Index> *mrowind;
00070 std::vector<Index> *mcolptr;
00072 Index mnr;
00073 };
00074
00080 class Sparse : public SparseView {
00081 public:
00082
00083 Sparse();
00084 Sparse(Index r, Index c);
00085 Sparse(const Sparse& m);
00086
00087
00088 ~Sparse();
00089 };
00090
00091
00092
00093
00094
00096 inline Index SparseView::nrows() const
00097 {
00098 return mrr.mextent;
00099 }
00100
00102 inline Index SparseView::ncols() const
00103 {
00104 return mcr.mextent;
00105 }
00106
00108 inline Index SparseView::nnz() const
00109 {
00110 return *(mcolptr->end()-1);
00111 }
00112
00120 inline Numeric& SparseView::rw(Index r, Index c)
00121 {
00122
00123 assert( 0<=r );
00124 assert( 0<=c );
00125 assert( r<mrr.mextent );
00126 assert( c<mcr.mextent );
00127
00128
00129
00130 r = mrr.mstart + r*mrr.mstride;
00131 c = mcr.mstart + c*mcr.mstride;
00132
00133
00134 Index i = (*mcolptr)[c];
00135
00136
00137
00138
00139 const Index end = (*mcolptr)[c+1];
00140
00141
00142
00143
00144
00145 for ( ; i<end; ++i )
00146 {
00147 Index rowi = (*mrowind)[i];
00148 if ( r < rowi )
00149 break;
00150 if ( r == rowi )
00151 return (*mdata)[i];
00152 }
00153
00154
00155
00156
00157
00158
00159 for ( std::vector<Index>::iterator j = mcolptr->begin() + c + 1;
00160 j < mcolptr->end();
00161 ++j )
00162 ++(*j);
00163
00164
00165
00166
00167
00168 mrowind->insert( mrowind->begin()+i, r );
00169 return *( mdata->insert( mdata->begin()+i, 0 ) );
00170 }
00171
00173 inline Numeric SparseView::ro(Index r, Index c) const
00174 {
00175
00176 assert( 0<=r );
00177 assert( 0<=c );
00178 assert( r<mrr.mextent );
00179 assert( c<mcr.mextent );
00180
00181
00182 r = mrr.mstart + r*mrr.mstride;
00183 c = mcr.mstart + c*mcr.mstride;
00184
00185
00186 Index begin = (*mcolptr)[c];
00187
00188
00189
00190
00191 const Index end = (*mcolptr)[c+1];
00192
00193
00194
00195
00196 for ( Index i=begin; i<end; ++i )
00197 {
00198 Index rowi = (*mrowind)[i];
00199 if ( r < rowi )
00200 return 0;
00201 if ( r == rowi )
00202 return (*mdata)[i];
00203 }
00204 return 0;
00205 }
00206
00209 inline SparseView::SparseView() :
00210 mdata(NULL),
00211 mrowind(NULL),
00212 mcolptr(NULL),
00213 mrr(0,0),
00214 mcr(0,0)
00215 {
00216
00217 }
00218
00222 inline SparseView::SparseView(std::vector<Numeric> *data,
00223 std::vector<Index> *rowind,
00224 std::vector<Index> *colptr,
00225 const Range& rr,
00226 const Range& cr ) :
00227 mdata(data),
00228 mrowind(rowind),
00229 mcolptr(colptr),
00230 mrr(rr),
00231 mcr(cr)
00232 {
00233
00234 }
00235
00236
00237
00238
00239
00240
00242 inline Sparse::Sparse()
00243 {
00244
00245 }
00246
00260 inline Sparse::Sparse(Index r, Index c) :
00261 SparseView( new std::vector<Numeric>,
00262 new std::vector<Index>,
00263 new std::vector<Index>(c+1,0),
00264 Range(0,r),
00265 Range(0,c) )
00266 {
00267
00268 }
00269
00272 inline Sparse::Sparse(const Sparse& m) :
00273 SparseView( new std::vector<Numeric>(*m.mdata),
00274 new std::vector<Index>(*m.mrowind),
00275 new std::vector<Index>(*m.mcolptr),
00276 m.mrr,
00277 m.mcr )
00278 {
00279
00280 }
00281
00282
00285 inline Sparse::~Sparse()
00286 {
00287 delete mdata;
00288 delete mrowind;
00289 delete mcolptr;
00290 }
00291
00292
00293
00294 inline std::ostream& operator<<(std::ostream& os, const SparseView& v)
00295 {
00296 for (size_t c=0; c<v.mcolptr->size()-1; ++c)
00297 {
00298
00299 Index begin = (*v.mcolptr)[c];
00300
00301
00302
00303
00304 const Index end = (*v.mcolptr)[c+1];
00305
00306
00307 for ( Index i=begin; i<end; ++i )
00308 {
00309 Index r = (*v.mrowind)[i];
00310
00311
00312
00313
00314 Index ra = (r-v.mrr.mstart)/v.mrr.mstride;
00315 Index ca = (c-v.mcr.mstart)/v.mcr.mstride;
00316
00317
00318 if ( 0 <= ra &&
00319 ra < v.mrr.mextent &&
00320 0 <= ca &&
00321 ca < v.mcr.mextent )
00322 {
00323
00324 os << setw(3) << ra << " "
00325 << setw(3) << ca << " "
00326 << setw(3) << (*v.mdata)[i] << "\n";
00327 }
00328 }
00329 }
00330
00331 return os;
00332 }