dune-istl
2.1.1
|
00001 // -*- tab-width: 8; indent-tabs-mode: nil; c-basic-offset: 2 -*- 00002 // vi: set et ts=8 sw=2 sts=2: 00003 #ifndef DUNE_SUPERLUMATRIX_HH 00004 #define DUNE_SUPERLUMATRIX_HH 00005 00006 #if HAVE_SUPERLU 00007 #ifdef SUPERLU_POST_2005_VERSION 00008 #include "slu_ddefs.h" 00009 #else 00010 #include "dsp_defs.h" 00011 #endif 00012 #include"bcrsmatrix.hh" 00013 #include"bvector.hh" 00014 #include<dune/common/fmatrix.hh> 00015 #include<dune/common/fvector.hh> 00016 #include<limits> 00017 00018 namespace Dune 00019 { 00020 00026 template<class M> 00027 class MatrixRowSet 00028 { 00029 public: 00030 // @brief The type of the matrix. 00031 typedef M Matrix; 00032 // @brief The matrix row iterator type. 00033 typedef typename Matrix::ConstRowIterator const_iterator; 00034 00039 MatrixRowSet(const Matrix& m) 00040 : m_(m) 00041 {} 00042 00043 // @brief Get the row iterator at the first row. 00044 const_iterator begin() const 00045 { 00046 return m_.begin(); 00047 } 00048 //@brief Get the row iterator at the end of all rows. 00049 const_iterator end() const 00050 { 00051 return m_.end(); 00052 } 00053 private: 00054 const Matrix& m_; 00055 }; 00056 00064 template<class M, class S> 00065 class MatrixRowSubset 00066 { 00067 public: 00068 /* @brief the type of the matrix class. */ 00069 typedef M Matrix; 00070 /* @brief the type of the set of valid row indices. */ 00071 typedef S RowIndexSet; 00072 00078 MatrixRowSubset(const Matrix& m, const RowIndexSet& s) 00079 : m_(m), s_(s) 00080 {} 00081 00082 const Matrix& matrix() const 00083 { 00084 return m_; 00085 } 00086 00087 const RowIndexSet& rowIndexSet() const 00088 { 00089 return s_; 00090 } 00091 00092 // @brief The matrix row iterator type. 00093 class const_iterator 00094 : public ForwardIteratorFacade<const_iterator, const typename Matrix::row_type> 00095 { 00096 public: 00097 const_iterator(typename Matrix::const_iterator firstRow, 00098 typename RowIndexSet::const_iterator pos) 00099 : firstRow_(firstRow), pos_(pos) 00100 {} 00101 00102 00103 const typename Matrix::row_type& dereference() const 00104 { 00105 return *(firstRow_+ *pos_); 00106 } 00107 bool equals(const const_iterator& o) const 00108 { 00109 return pos_==o.pos_; 00110 } 00111 void increment() 00112 { 00113 ++pos_; 00114 } 00115 typename RowIndexSet::value_type index() const 00116 { 00117 return *pos_; 00118 } 00119 00120 private: 00121 // @brief Iterator pointing to the first row of the matrix. 00122 typename Matrix::const_iterator firstRow_; 00123 // @brief Iterator pointing to the current row index. 00124 typename RowIndexSet::const_iterator pos_; 00125 }; 00126 00127 // @brief Get the row iterator at the first row. 00128 const_iterator begin() const 00129 { 00130 return const_iterator(m_.begin(), s_.begin()); 00131 } 00132 //@brief Get the row iterator at the end of all rows. 00133 const_iterator end() const 00134 { 00135 return const_iterator(m_.begin(), s_.end()); 00136 } 00137 00138 private: 00139 // @brief The matrix for which we manage the row subset. 00140 const Matrix& m_; 00141 // @brief The set of row indices we manage. 00142 const RowIndexSet& s_; 00143 00144 }; 00145 00146 template<class T> 00147 struct GetSuperLUType 00148 { 00149 }; 00150 00151 template<> 00152 struct GetSuperLUType<double> 00153 { 00154 enum{ type = SLU_D}; 00155 00156 }; 00157 00158 template<> 00159 struct GetSuperLUType<float> 00160 { 00161 enum{ type = SLU_S}; 00162 00163 }; 00164 00165 template<> 00166 struct GetSuperLUType<std::complex<double> > 00167 { 00168 enum{ type = SLU_Z}; 00169 00170 }; 00171 00172 template<> 00173 struct GetSuperLUType<std::complex<float> > 00174 { 00175 enum{ type = SLU_C}; 00176 00177 }; 00178 00183 template<class M> 00184 struct SuperLUMatrix 00185 { 00186 00187 }; 00188 00189 template<class M> 00190 struct SuperMatrixInitializer 00191 {}; 00192 00193 template<class M, class X, class TM, class TD, class T1> 00194 class SeqOverlappingSchwarz; 00195 00196 00197 template<class T> 00198 class SeqOverlappingSchwarzAssembler; 00199 00200 template<class T> 00201 class SuperLU; 00202 00206 template<class B, class TA, int n, int m> 00207 class SuperLUMatrix<BCRSMatrix<FieldMatrix<B,n,m>,TA> > 00208 { 00209 template<class M, class X, class TM, class TD, class T1> 00210 friend class SeqOverlappingSchwarz; 00211 friend class SuperMatrixInitializer<BCRSMatrix<FieldMatrix<B,n,m>,TA> >; 00212 00213 public: 00215 typedef BCRSMatrix<FieldMatrix<B,n,m>,TA> Matrix; 00216 00217 friend class SeqOverlappingSchwarzAssembler<SuperLU<Matrix> >; 00218 00219 typedef typename Matrix::size_type size_type; 00220 00225 SuperLUMatrix(const Matrix& mat); 00226 00227 SuperLUMatrix(); 00228 00230 ~SuperLUMatrix(); 00231 00233 operator SuperMatrix&() 00234 { 00235 return A; 00236 } 00237 00239 operator const SuperMatrix&()const 00240 { 00241 return A; 00242 } 00243 bool operator==(const Matrix& mat) const; 00244 00249 size_type N() const 00250 { 00251 return N_; 00252 } 00253 00254 size_type nnz() const 00255 { 00256 return Nnz_/n/m; 00257 } 00258 00263 size_type M() const 00264 { 00265 return M_; 00266 } 00267 00268 SuperLUMatrix& operator=(const Matrix& mat); 00269 00270 SuperLUMatrix& operator=(const SuperLUMatrix& mat); 00271 00278 template<class S> 00279 void setMatrix(const Matrix& mat, const S& mrs); 00281 void free(); 00282 private: 00284 void setMatrix(const Matrix& mat); 00285 00286 int N_, M_, Nnz_; 00287 B* values; 00288 int* rowindex; 00289 int* colstart; 00290 SuperMatrix A; 00291 }; 00292 00293 template<class T, class A, int n, int m> 00294 void writeCompColMatrixToMatlab(const SuperLUMatrix<BCRSMatrix<FieldMatrix<T,n,m>,A> >& mat, 00295 std::ostream& os) 00296 { 00297 const SuperMatrix a=static_cast<const SuperMatrix&>(mat); 00298 NCformat *astore = (NCformat *) a.Store; 00299 double* dp = (double*)astore->nzval; 00300 00301 // remember old flags 00302 std::ios_base::fmtflags oldflags = os.flags(); 00303 // set the output format 00304 //os.setf(std::ios_base::scientific, std::ios_base::floatfield); 00305 int oldprec = os.precision(); 00306 //os.precision(10); 00307 //dPrint_CompCol_Matrix("A", const_cast<SuperMatrix*>(&a)); 00308 00309 os <<"["; 00310 for(int row=0; row<a.nrow; ++row){ 00311 for(int col= 0; col < a.ncol; ++col){ 00312 // linear search for col 00313 int i; 00314 for(i=astore->colptr[col]; i < astore->colptr[col+1]; ++i) 00315 if(astore->rowind[i]==row){ 00316 os<<dp[i]<<" "; 00317 break; 00318 } 00319 if(i==astore->colptr[col+1]) 00320 // entry not present 00321 os<<0<<" "; 00322 } 00323 if(row==a.nrow-1) 00324 os<<"]"; 00325 os<<std::endl; 00326 } 00327 // reset the output format 00328 os.flags(oldflags); 00329 os.precision(oldprec); 00330 } 00331 00332 00333 template<class T, class A, int n, int m> 00334 class SuperMatrixInitializer<BCRSMatrix<FieldMatrix<T,n,m>,A> > 00335 { 00336 template<class I, class S, class D> 00337 friend class OverlappingSchwarzInitializer; 00338 public: 00339 typedef Dune::BCRSMatrix<FieldMatrix<T,n,m>,A> Matrix; 00340 typedef Dune::SuperLUMatrix<Matrix> SuperLUMatrix; 00341 typedef typename Matrix::row_type::const_iterator CIter; 00342 typedef typename Matrix::size_type size_type; 00343 00344 SuperMatrixInitializer(SuperLUMatrix& lum); 00345 00346 SuperMatrixInitializer(); 00347 00348 ~SuperMatrixInitializer(); 00349 00350 template<typename Iter> 00351 void addRowNnz(const Iter& row)const; 00352 00353 template<typename Iter, typename Set> 00354 void addRowNnz(const Iter& row, const Set& s)const; 00355 00356 void allocate(); 00357 00358 template<typename Iter> 00359 void countEntries(const Iter& row, const CIter& col) const; 00360 00361 void countEntries(size_type colidx)const; 00362 00363 void calcColstart()const; 00364 00365 template<typename Iter> 00366 void copyValue(const Iter& row, const CIter& col) const; 00367 00368 void copyValue(const CIter& col, size_type rowindex, size_type colidx)const; 00369 00370 void createMatrix() const; 00371 00372 private: 00373 00374 void allocateMatrixStorage()const; 00375 00376 void allocateMarker(); 00377 00378 SuperLUMatrix* mat; 00379 size_type cols; 00380 mutable size_type *marker; 00381 }; 00382 00383 template<class T, class A, int n, int m> 00384 SuperMatrixInitializer<BCRSMatrix<FieldMatrix<T,n,m>,A> >::SuperMatrixInitializer(SuperLUMatrix& mat_) 00385 : mat(&mat_), cols(mat_.N()), marker(0) 00386 { 00387 mat->Nnz_=0; 00388 } 00389 00390 template<class T, class A, int n, int m> 00391 SuperMatrixInitializer<BCRSMatrix<FieldMatrix<T,n,m>,A> >::SuperMatrixInitializer() 00392 : mat(0), cols(0), marker(0) 00393 {} 00394 00395 template<class T, class A, int n, int m> 00396 SuperMatrixInitializer<BCRSMatrix<FieldMatrix<T,n,m>,A> >::~SuperMatrixInitializer() 00397 { 00398 if(marker) 00399 delete[] marker; 00400 } 00401 00402 template<class T, class A, int n, int m> 00403 template<typename Iter> 00404 void SuperMatrixInitializer<BCRSMatrix<FieldMatrix<T,n,m>,A> >::addRowNnz(const Iter& row)const 00405 { 00406 mat->Nnz_+=row->getsize(); 00407 } 00408 00409 template<class T, class A, int n, int m> 00410 template<typename Iter, typename Map> 00411 void SuperMatrixInitializer<BCRSMatrix<FieldMatrix<T,n,m>,A> >::addRowNnz(const Iter& row, 00412 const Map& indices)const 00413 { 00414 typedef typename Iter::value_type::const_iterator RIter; 00415 typedef typename Map::const_iterator MIter; 00416 MIter siter =indices.begin(); 00417 for(RIter entry=row->begin();entry!=row->end();++entry) 00418 { 00419 for(;siter!=indices.end() && *siter<entry.index(); ++siter); 00420 if(siter==indices.end()) 00421 break; 00422 if(*siter==entry.index()) 00423 // index is in subdomain 00424 ++mat->Nnz_; 00425 } 00426 } 00427 00428 template<class T, class A, int n, int m> 00429 void SuperMatrixInitializer<BCRSMatrix<FieldMatrix<T,n,m>,A> >::allocate() 00430 { 00431 allocateMatrixStorage(); 00432 allocateMarker(); 00433 } 00434 00435 template<class T, class A, int n, int m> 00436 void SuperMatrixInitializer<BCRSMatrix<FieldMatrix<T,n,m>,A> >::allocateMatrixStorage()const 00437 { 00438 mat->Nnz_*=n*m; 00439 // initialize data 00440 mat->values=new T[mat->Nnz_]; 00441 mat->rowindex=new int[mat->Nnz_]; 00442 mat->colstart=new int[cols+1]; 00443 } 00444 00445 template<class T, class A, int n, int m> 00446 void SuperMatrixInitializer<BCRSMatrix<FieldMatrix<T,n,m>,A> >::allocateMarker() 00447 { 00448 marker = new typename Matrix::size_type[cols]; 00449 00450 for(size_type i=0; i < cols; ++i) 00451 marker[i]=0; 00452 } 00453 00454 template<class T, class A, int n, int m> 00455 template<typename Iter> 00456 void SuperMatrixInitializer<BCRSMatrix<FieldMatrix<T,n,m>,A> >::countEntries(const Iter& row, const CIter& col)const 00457 { 00458 countEntries(col.index()); 00459 00460 } 00461 00462 template<class T, class A, int n, int m> 00463 void SuperMatrixInitializer<BCRSMatrix<FieldMatrix<T,n,m>,A> >::countEntries(size_type colindex)const 00464 { 00465 for(size_type i=0; i < m; ++i){ 00466 assert(colindex*m+i<cols); 00467 marker[colindex*m+i]+=n; 00468 } 00469 00470 } 00471 00472 template<class T, class A, int n, int m> 00473 void SuperMatrixInitializer<BCRSMatrix<FieldMatrix<T,n,m>,A> >::calcColstart()const 00474 { 00475 mat->colstart[0]=0; 00476 for(size_type i=0; i < cols; ++i){ 00477 assert(i<cols); 00478 mat->colstart[i+1]=mat->colstart[i]+marker[i]; 00479 marker[i]=mat->colstart[i]; 00480 } 00481 } 00482 00483 template<class T, class A, int n, int m> 00484 template<typename Iter> 00485 void SuperMatrixInitializer<BCRSMatrix<FieldMatrix<T,n,m>,A> >::copyValue(const Iter& row, const CIter& col)const 00486 { 00487 copyValue(col, row.index(), col.index()); 00488 } 00489 00490 template<class T, class A, int n, int m> 00491 void SuperMatrixInitializer<BCRSMatrix<FieldMatrix<T,n,m>,A> >::copyValue(const CIter& col, size_type rowindex, size_type colindex)const 00492 { 00493 for(size_type i=0; i<n;i++){ 00494 for(size_type j=0; j<m; j++){ 00495 assert(colindex*m+j<cols-1 || (int)marker[colindex*m+j]<mat->colstart[colindex*m+j+1]); 00496 assert((int)marker[colindex*m+j]<mat->Nnz_); 00497 mat->rowindex[marker[colindex*m+j]]=rowindex*n+i; 00498 mat->values[marker[colindex*m+j]]=(*col)[i][j]; 00499 ++marker[colindex*m+j]; // index for next entry in column 00500 } 00501 } 00502 } 00503 00504 template<class T, class A, int n, int m> 00505 void SuperMatrixInitializer<BCRSMatrix<FieldMatrix<T,n,m>,A> >::createMatrix() const 00506 { 00507 delete[] marker; 00508 marker=0; 00509 dCreate_CompCol_Matrix(&mat->A, mat->N_, mat->M_, mat->colstart[cols], 00510 mat->values, mat->rowindex, mat->colstart, SLU_NC, static_cast<Dtype_t>(GetSuperLUType<T>::type), SLU_GE); 00511 } 00512 00513 template<class F, class MRS> 00514 void copyToSuperMatrix(F& initializer, const MRS& mrs) 00515 { 00516 typedef typename MRS::const_iterator Iter; 00517 typedef typename std::iterator_traits<Iter>::value_type::const_iterator CIter; 00518 for(Iter row=mrs.begin(); row!= mrs.end(); ++row) 00519 initializer.addRowNnz(row); 00520 00521 initializer.allocate(); 00522 00523 for(Iter row=mrs.begin(); row!= mrs.end(); ++row){ 00524 00525 for(CIter col=row->begin(); col != row->end(); ++col) 00526 initializer.countEntries(row, col); 00527 } 00528 00529 initializer.calcColstart(); 00530 00531 for(Iter row=mrs.begin(); row!= mrs.end(); ++row){ 00532 for(CIter col=row->begin(); col != row->end(); ++col){ 00533 initializer.copyValue(row, col); 00534 } 00535 00536 } 00537 initializer.createMatrix(); 00538 } 00539 00540 template<class F, class M,class S> 00541 void copyToSuperMatrix(F& initializer, const MatrixRowSubset<M,S>& mrs) 00542 { 00543 typedef MatrixRowSubset<M,S> MRS; 00544 typedef typename MRS::RowIndexSet SIS; 00545 typedef typename SIS::const_iterator SIter; 00546 typedef typename MRS::const_iterator Iter; 00547 typedef typename std::iterator_traits<Iter>::value_type row_type; 00548 typedef typename row_type::const_iterator CIter; 00549 00550 // Calculate upper Bound for nonzeros 00551 for(Iter row=mrs.begin(); row!= mrs.end(); ++row) 00552 initializer.addRowNnz(row, mrs.rowIndexSet()); 00553 00554 initializer.allocate(); 00555 00556 typedef typename MRS::Matrix::size_type size_type; 00557 00558 // A vector containing the corresponding indices in 00559 // the to create submatrix. 00560 // If an entry is the maximum of size_type then this index will not appear in 00561 // the submatrix. 00562 std::vector<size_type> subMatrixIndex(mrs.matrix().N(), 00563 std::numeric_limits<size_type>::max()); 00564 size_type s=0; 00565 for(SIter index = mrs.rowIndexSet().begin(); index!=mrs.rowIndexSet().end(); ++index) 00566 subMatrixIndex[*index]=s++; 00567 00568 for(Iter row=mrs.begin(); row!= mrs.end(); ++row) 00569 for(CIter col=row->begin(); col != row->end(); ++col){ 00570 if(subMatrixIndex[col.index()]!=std::numeric_limits<size_type>::max()) 00571 // This column is in our subset (use submatrix column index) 00572 initializer.countEntries(subMatrixIndex[col.index()]); 00573 } 00574 00575 initializer.calcColstart(); 00576 00577 for(Iter row=mrs.begin(); row!= mrs.end(); ++row) 00578 for(CIter col=row->begin(); col != row->end(); ++col){ 00579 if(subMatrixIndex[col.index()]!=std::numeric_limits<size_type>::max()) 00580 // This value is in our submatrix -> copy (use submatrix indices 00581 initializer.copyValue(col, subMatrixIndex[row.index()], subMatrixIndex[col.index()]); 00582 } 00583 initializer.createMatrix(); 00584 } 00585 00586 #ifndef DOXYGEN 00587 00588 template<class B, class TA, int n, int m> 00589 bool SuperLUMatrix<BCRSMatrix<FieldMatrix<B,n,m>,TA> >::operator==(const BCRSMatrix<FieldMatrix<B,n,m>,TA>& mat) const 00590 { 00591 const NCformat* S=static_cast<const NCformat *>(A.Store); 00592 for(size_type col=0; col < M(); ++col){ 00593 for(int j=S->colptr[col]; j < S->colptr[col+1]; ++j){ 00594 int row=S->rowind[j]; 00595 if((mat[row/n][col/m])[row%n][col%m]!=reinterpret_cast<B*>(S->nzval)[j]){ 00596 std::cerr<<" bcrs["<<row/n<<"]["<<col/m<<"]["<<row%n<<"]["<<row%m 00597 <<"]="<<(mat[row/n][col/m])[row%n][col%m]<<" super["<<row<<"]["<<col<<"]="<<reinterpret_cast<B*>(S->nzval)[j]; 00598 00599 return false; 00600 } 00601 } 00602 } 00603 00604 return true; 00605 } 00606 00607 #endif // DOYXGEN 00608 00609 template<class B, class TA, int n, int m> 00610 bool operator==(SuperLUMatrix<BCRSMatrix<FieldMatrix<B,n,m>,TA> >& sla, BCRSMatrix<FieldMatrix<B,n,m>,TA>& a) 00611 { 00612 return a==sla; 00613 } 00614 00615 #ifndef DOXYGEN 00616 00617 template<class B, class TA, int n, int m> 00618 SuperLUMatrix<BCRSMatrix<FieldMatrix<B,n,m>,TA> >::SuperLUMatrix() 00619 : N_(0), M_(0), Nnz_(0), values(0), rowindex(0), colstart(0) 00620 {} 00621 00622 template<class B, class TA, int n, int m> 00623 SuperLUMatrix<BCRSMatrix<FieldMatrix<B,n,m>,TA> > 00624 ::SuperLUMatrix(const Matrix& mat) 00625 : N_(n*mat.N()), M_(m*mat.M()), Nnz_(n*m*mat.Nnz()) 00626 { 00627 00628 } 00629 00630 template<class B, class TA, int n, int m> 00631 SuperLUMatrix<BCRSMatrix<FieldMatrix<B,n,m>,TA> >& 00632 SuperLUMatrix<BCRSMatrix<FieldMatrix<B,n,m>,TA> >::operator=(const Matrix& mat) 00633 { 00634 if(N_+M_+Nnz_!=0) 00635 free(); 00636 setMatrix(mat); 00637 return *this; 00638 } 00639 00640 template<class B, class TA, int n, int m> 00641 SuperLUMatrix<BCRSMatrix<FieldMatrix<B,n,m>,TA> >& 00642 SuperLUMatrix<BCRSMatrix<FieldMatrix<B,n,m>,TA> >::operator=(const SuperLUMatrix& mat) 00643 { 00644 if(N_+M_+Nnz_!=0) 00645 free(); 00646 N_=mat.N_; 00647 M_=mat.M_; 00648 Nnz_= mat.Nnz_; 00649 if(M_>0){ 00650 colstart=new int[M_+1]; 00651 for(int i=0; i<=M_; ++i) 00652 colstart[i]=mat.colstart[i]; 00653 } 00654 00655 if(Nnz_>0){ 00656 values = new B[Nnz_]; 00657 rowindex = new int[Nnz_]; 00658 00659 for(int i=0; i<Nnz_; ++i) 00660 values[i]=mat.values[i]; 00661 00662 for(int i=0; i<Nnz_; ++i) 00663 rowindex[i]=mat.rowindex[i]; 00664 } 00665 if(M_+Nnz_>0) 00666 dCreate_CompCol_Matrix(&A, N_, M_, Nnz_, 00667 values, rowindex, colstart, SLU_NC, static_cast<Dtype_t>(GetSuperLUType<B>::type), SLU_GE); 00668 return *this; 00669 } 00670 00671 template<class B, class TA, int n, int m> 00672 void SuperLUMatrix<BCRSMatrix<FieldMatrix<B,n,m>,TA> > 00673 ::setMatrix(const Matrix& mat) 00674 { 00675 N_=n*mat.N(); 00676 M_=m*mat.M(); 00677 SuperMatrixInitializer<Matrix> initializer(*this); 00678 00679 copyToSuperMatrix(initializer, MatrixRowSet<Matrix>(mat)); 00680 00681 #ifdef DUNE_ISTL_WITH_CHECKING 00682 char name[] = {'A',0}; 00683 if(N_<0) 00684 dPrint_CompCol_Matrix(name,&A); 00685 assert(*this==mat); 00686 #endif 00687 } 00688 00689 template<class B, class TA, int n, int m> 00690 template<class S> 00691 void SuperLUMatrix<BCRSMatrix<FieldMatrix<B,n,m>,TA> > 00692 ::setMatrix(const Matrix& mat, const S& mrs) 00693 { 00694 if(N_+M_+Nnz_!=0) 00695 free(); 00696 N_=mrs.size()*n; 00697 M_=mrs.size()*m; 00698 SuperMatrixInitializer<Matrix> initializer(*this); 00699 00700 copyToSuperMatrix(initializer, MatrixRowSubset<Matrix,S>(mat,mrs)); 00701 00702 #ifdef DUNE_ISTL_WITH_CHECKING 00703 char name[] = {'A',0}; 00704 if(N_<0) 00705 dPrint_CompCol_Matrix(name,&A); 00706 #endif 00707 } 00708 00709 template<class B, class TA, int n, int m> 00710 SuperLUMatrix<BCRSMatrix<FieldMatrix<B,n,m>,TA> >::~SuperLUMatrix() 00711 { 00712 if(N_+M_+Nnz_!=0) 00713 free(); 00714 } 00715 00716 template<class B, class TA, int n, int m> 00717 void SuperLUMatrix<BCRSMatrix<FieldMatrix<B,n,m>,TA> >::free() 00718 { 00719 delete[] values; 00720 delete[] rowindex; 00721 delete[] colstart; 00722 SUPERLU_FREE(A.Store); 00723 N_=M_=Nnz_=0; 00724 } 00725 00726 #endif // DOXYGEN 00727 00728 } 00729 #endif 00730 #endif