dune-istl
2.2.0
|
00001 #ifndef DUNE_SUPERLU_HH 00002 #define DUNE_SUPERLU_HH 00003 00004 #if HAVE_SUPERLU 00005 #ifdef TRUE 00006 #undef TRUE 00007 #endif 00008 #ifdef FALSE 00009 #undef FALSE 00010 #endif 00011 #ifndef SUPERLU_NTYPE 00012 #define SUPERLU_NTYPE 1 00013 #endif 00014 #ifdef SUPERLU_POST_2005_VERSION 00015 00016 #if SUPERLU_NTYPE==0 00017 #include "slu_sdefs.h" 00018 #endif 00019 00020 #if SUPERLU_NTYPE==1 00021 #include "slu_ddefs.h" 00022 #endif 00023 00024 #if SUPERLU_NTYPE==2 00025 #include "slu_cdefs.h" 00026 #endif 00027 00028 #if SUPERLU_NTYPE>=3 00029 #include "slu_zdefs.h" 00030 #endif 00031 00032 #else 00033 00034 #if SUPERLU_NTYPE==0 00035 #include "ssp_defs.h" 00036 #endif 00037 00038 #if SUPERLU_NTYPE==1 00039 #include "dsp_defs.h" 00040 #warning Support for SuperLU older than SuperLU 3.0 from August 2005 is deprecated. 00041 #endif 00042 00043 #if SUPERLU_NTYPE==2 00044 #include "csp_defs.h" 00045 #endif 00046 00047 #if SUPERLU_NTYPE>=3 00048 #include "zsp_defs.h" 00049 #endif 00050 00051 #endif 00052 #include "solvers.hh" 00053 #include "supermatrix.hh" 00054 #include<algorithm> 00055 #include<functional> 00056 #include"bcrsmatrix.hh" 00057 #include"bvector.hh" 00058 #include"istlexception.hh" 00059 #include<dune/common/fmatrix.hh> 00060 #include<dune/common/fvector.hh> 00061 #include<dune/common/stdstreams.hh> 00062 #include <dune/istl/solvertype.hh> 00063 00064 namespace Dune 00065 { 00066 00077 template<class Matrix> 00078 class SuperLU 00079 { 00080 }; 00081 00082 template<class M, class T, class TM, class TD, class TA> 00083 class SeqOverlappingSchwarz; 00084 00085 template<class T> 00086 class SeqOverlappingSchwarzAssembler; 00087 00088 template<class T> 00089 struct SuperLUSolveChooser 00090 { 00091 }; 00092 00093 template<class T> 00094 struct SuperLUDenseMatChooser 00095 { 00096 }; 00097 00098 template<class T> 00099 struct SuperLUQueryChooser 00100 { 00101 }; 00102 00103 template<class T> 00104 struct QuerySpaceChooser 00105 { 00106 }; 00107 00108 #if SUPERLU_NTYPE==0 00109 template<> 00110 struct SuperLUDenseMatChooser<float> 00111 { 00112 static void create(SuperMatrix *mat, int n, int m, float *dat, int n1, 00113 Stype_t stype, Dtype_t dtype, Mtype_t mtype) 00114 { 00115 sCreate_Dense_Matrix(mat, n, m, dat, n1, stype, dtype, mtype); 00116 00117 } 00118 00119 static void destroy(SuperMatrix *m) 00120 { 00121 } 00122 00123 }; 00124 template<> 00125 struct SuperLUSolveChooser<float> 00126 { 00127 static void solve(superlu_options_t *options, SuperMatrix *mat, int *perm_c, int *perm_r, int *etree, 00128 char *equed, float *R, float *C, SuperMatrix *L, SuperMatrix *U, 00129 void *work, int lwork, SuperMatrix *B, SuperMatrix *X, 00130 float *rpg, float *rcond, float *ferr, float *berr, 00131 mem_usage_t *memusage, SuperLUStat_t *stat, int *info) 00132 { 00133 sgssvx(options, mat, perm_c, perm_r, etree, equed, R, C, 00134 L, U, work, lwork, B, X, rpg, rcond, ferr, berr, 00135 memusage, stat, info); 00136 } 00137 }; 00138 00139 template<> 00140 struct QuerySpaceChooser<float> 00141 { 00142 static void querySpace(SuperMatrix* L, SuperMatrix* U, mem_usage_t* memusage) 00143 { 00144 sQuerySpace(L,U,memusage); 00145 } 00146 }; 00147 00148 #endif 00149 00150 #if SUPERLU_NTYPE==1 00151 00152 template<> 00153 struct SuperLUDenseMatChooser<double> 00154 { 00155 static void create(SuperMatrix *mat, int n, int m, double *dat, int n1, 00156 Stype_t stype, Dtype_t dtype, Mtype_t mtype) 00157 { 00158 dCreate_Dense_Matrix(mat, n, m, dat, n1, stype, dtype, mtype); 00159 00160 } 00161 00162 static void destroy(SuperMatrix *mat) 00163 { 00164 } 00165 }; 00166 template<> 00167 struct SuperLUSolveChooser<double> 00168 { 00169 static void solve(superlu_options_t *options, SuperMatrix *mat, int *perm_c, int *perm_r, int *etree, 00170 char *equed, double *R, double *C, SuperMatrix *L, SuperMatrix *U, 00171 void *work, int lwork, SuperMatrix *B, SuperMatrix *X, 00172 double *rpg, double *rcond, double *ferr, double *berr, 00173 mem_usage_t *memusage, SuperLUStat_t *stat, int *info) 00174 { 00175 dgssvx(options, mat, perm_c, perm_r, etree, equed, R, C, 00176 L, U, work, lwork, B, X, rpg, rcond, ferr, berr, 00177 memusage, stat, info); 00178 } 00179 }; 00180 00181 template<> 00182 struct QuerySpaceChooser<double> 00183 { 00184 static void querySpace(SuperMatrix* L, SuperMatrix* U, mem_usage_t* memusage) 00185 { 00186 dQuerySpace(L,U,memusage); 00187 } 00188 }; 00189 #endif 00190 00191 #if SUPERLU_NTYPE>=3 00192 template<> 00193 struct SuperLUDenseMatChooser<std::complex<double> > 00194 { 00195 static void create(SuperMatrix *mat, int n, int m, std::complex<double> *dat, int n1, 00196 Stype_t stype, Dtype_t dtype, Mtype_t mtype) 00197 { 00198 zCreate_Dense_Matrix(mat, n, m, reinterpret_cast<doublecomplex*>(dat), n1, stype, dtype, mtype); 00199 00200 } 00201 00202 static void destroy(SuperMatrix *mat) 00203 { 00204 } 00205 }; 00206 00207 template<> 00208 struct SuperLUSolveChooser<std::complex<double> > 00209 { 00210 static void solve(superlu_options_t *options, SuperMatrix *mat, int *perm_c, int *perm_r, int *etree, 00211 char *equed, double *R, double *C, SuperMatrix *L, SuperMatrix *U, 00212 void *work, int lwork, SuperMatrix *B, SuperMatrix *X, 00213 double *rpg, double *rcond, double *ferr, double *berr, 00214 mem_usage_t *memusage, SuperLUStat_t *stat, int *info) 00215 { 00216 zgssvx(options, mat, perm_c, perm_r, etree, equed, R, C, 00217 L, U, work, lwork, B, X, rpg, rcond, ferr, berr, 00218 memusage, stat, info); 00219 } 00220 }; 00221 00222 template<> 00223 struct QuerySpaceChooser<std::complex<double> > 00224 { 00225 static void querySpace(SuperMatrix* L, SuperMatrix* U, mem_usage_t* memusage) 00226 { 00227 zQuerySpace(L,U,memusage); 00228 } 00229 }; 00230 #endif 00231 00232 #if SUPERLU_NTYPE==2 00233 template<> 00234 struct SuperLUDenseMatChooser<std::complex<float> > 00235 { 00236 static void create(SuperMatrix *mat, int n, int m, std::complex<float> *dat, int n1, 00237 Stype_t stype, Dtype_t dtype, Mtype_t mtype) 00238 { 00239 cCreate_Dense_Matrix(mat, n, m, reinterpret_cast< ::complex*>(dat), n1, stype, dtype, mtype); 00240 00241 } 00242 00243 static void destroy(SuperMatrix *mat) 00244 { 00245 } 00246 }; 00247 00248 template<> 00249 struct SuperLUSolveChooser<std::complex<float> > 00250 { 00251 static void solve(superlu_options_t *options, SuperMatrix *mat, int *perm_c, int *perm_r, int *etree, 00252 char *equed, float *R, float *C, SuperMatrix *L, SuperMatrix *U, 00253 void *work, int lwork, SuperMatrix *B, SuperMatrix *X, 00254 float *rpg, float *rcond, float *ferr, float *berr, 00255 mem_usage_t *memusage, SuperLUStat_t *stat, int *info) 00256 { 00257 cgssvx(options, mat, perm_c, perm_r, etree, equed, R, C, 00258 L, U, work, lwork, B, X, rpg, rcond, ferr, berr, 00259 memusage, stat, info); 00260 } 00261 }; 00262 00263 template<> 00264 struct QuerySpaceChooser<std::complex<float> > 00265 { 00266 static void querySpace(SuperMatrix* L, SuperMatrix* U, mem_usage_t* memusage) 00267 { 00268 cQuerySpace(L,U,memusage); 00269 } 00270 }; 00271 #endif 00272 00286 template<typename T, typename A, int n, int m> 00287 class SuperLU<BCRSMatrix<FieldMatrix<T,n,m>,A > > 00288 : public InverseOperator< 00289 BlockVector<FieldVector<T,m>, 00290 typename A::template rebind<FieldVector<T,m> >::other>, 00291 BlockVector<FieldVector<T,n>, 00292 typename A::template rebind<FieldVector<T,n> >::other> > 00293 { 00294 public: 00295 /* @brief The matrix type. */ 00296 typedef Dune::BCRSMatrix<FieldMatrix<T,n,m>,A> Matrix; 00297 /* @brief The corresponding SuperLU Matrix type.*/ 00298 typedef Dune::SuperLUMatrix<Matrix> SuperLUMatrix; 00300 typedef Dune::BlockVector< 00301 FieldVector<T,m>, 00302 typename A::template rebind<FieldVector<T,m> >::other> domain_type; 00304 typedef Dune::BlockVector< 00305 FieldVector<T,n>, 00306 typename A::template rebind<FieldVector<T,n> >::other> range_type; 00316 explicit SuperLU(const Matrix& mat, bool verbose=false); 00322 SuperLU(); 00323 00324 ~SuperLU(); 00325 00329 void apply(domain_type& x, range_type& b, InverseOperatorResult& res); 00330 00334 void apply (domain_type& x, range_type& b, double reduction, InverseOperatorResult& res) 00335 { 00336 apply(x,b,res); 00337 } 00338 00342 void apply(T* x, T* b); 00343 00345 void setMatrix(const Matrix& mat); 00346 00347 typename SuperLUMatrix::size_type nnz() const 00348 { 00349 return mat.nnz(); 00350 } 00351 00352 template<class S> 00353 void setSubMatrix(const Matrix& mat, const S& rowIndexSet); 00354 00355 void setVerbosity(bool v); 00356 00361 void free(); 00362 private: 00363 friend class std::mem_fun_ref_t<void,SuperLU>; 00364 template<class M,class X, class TM, class TD, class T1> 00365 friend class SeqOverlappingSchwarz; 00366 friend class SeqOverlappingSchwarzAssembler<SuperLU<Matrix> >; 00367 00369 void decompose(); 00370 00371 SuperLUMatrix mat; 00372 SuperMatrix L, U, B, X; 00373 int *perm_c, *perm_r, *etree; 00374 typename GetSuperLUType<T>::float_type *R, *C; 00375 T *bstore; 00376 superlu_options_t options; 00377 char equed; 00378 void *work; 00379 int lwork; 00380 bool first, verbose; 00381 }; 00382 00383 template<typename T, typename A, int n, int m> 00384 SuperLU<BCRSMatrix<FieldMatrix<T,n,m>,A> > 00385 ::~SuperLU() 00386 { 00387 if(mat.N()+mat.M()>0) 00388 free(); 00389 } 00390 00391 template<typename T, typename A, int n, int m> 00392 void SuperLU<BCRSMatrix<FieldMatrix<T,n,m>,A> >::free() 00393 { 00394 delete[] perm_c; 00395 delete[] perm_r; 00396 delete[] etree; 00397 delete[] R; 00398 delete[] C; 00399 if(lwork>=0){ 00400 Destroy_SuperNode_Matrix(&L); 00401 Destroy_CompCol_Matrix(&U); 00402 } 00403 lwork=0; 00404 if(!first){ 00405 SUPERLU_FREE(B.Store); 00406 SUPERLU_FREE(X.Store); 00407 } 00408 mat.free(); 00409 } 00410 00411 template<typename T, typename A, int n, int m> 00412 SuperLU<BCRSMatrix<FieldMatrix<T,n,m>,A> > 00413 ::SuperLU(const Matrix& mat_, bool verbose_) 00414 : work(0), lwork(0), first(true), verbose(verbose_) 00415 { 00416 setMatrix(mat_); 00417 00418 } 00419 template<typename T, typename A, int n, int m> 00420 SuperLU<BCRSMatrix<FieldMatrix<T,n,m>,A> >::SuperLU() 00421 : work(0), lwork(0),verbose(false) 00422 {} 00423 template<typename T, typename A, int n, int m> 00424 void SuperLU<BCRSMatrix<FieldMatrix<T,n,m>,A> >::setVerbosity(bool v) 00425 { 00426 verbose=v; 00427 } 00428 00429 template<typename T, typename A, int n, int m> 00430 void SuperLU<BCRSMatrix<FieldMatrix<T,n,m>,A> >::setMatrix(const Matrix& mat_) 00431 { 00432 if(mat.N()+mat.M()>0){ 00433 free(); 00434 } 00435 lwork=0; 00436 work=0; 00437 //a=&mat_; 00438 mat=mat_; 00439 decompose(); 00440 } 00441 00442 template<typename T, typename A, int n, int m> 00443 template<class S> 00444 void SuperLU<BCRSMatrix<FieldMatrix<T,n,m>,A> >::setSubMatrix(const Matrix& mat_, 00445 const S& mrs) 00446 { 00447 if(mat.N()+mat.M()>0){ 00448 free(); 00449 } 00450 lwork=0; 00451 work=0; 00452 //a=&mat_; 00453 mat.setMatrix(mat_,mrs); 00454 decompose(); 00455 } 00456 00457 template<typename T, typename A, int n, int m> 00458 void SuperLU<BCRSMatrix<FieldMatrix<T,n,m>,A> >::decompose() 00459 { 00460 00461 first = true; 00462 perm_c = new int[mat.M()]; 00463 perm_r = new int[mat.N()]; 00464 etree = new int[mat.M()]; 00465 R = new typename GetSuperLUType<T>::float_type[mat.N()]; 00466 C = new typename GetSuperLUType<T>::float_type[mat.M()]; 00467 00468 set_default_options(&options); 00469 // Do the factorization 00470 B.ncol=0; 00471 B.Stype=SLU_DN; 00472 B.Dtype=GetSuperLUType<T>::type; 00473 B.Mtype= SLU_GE; 00474 DNformat fakeFormat; 00475 fakeFormat.lda=mat.N(); 00476 B.Store=&fakeFormat; 00477 X.Stype=SLU_DN; 00478 X.Dtype=GetSuperLUType<T>::type; 00479 X.Mtype= SLU_GE; 00480 X.ncol=0; 00481 X.Store=&fakeFormat; 00482 00483 typename GetSuperLUType<T>::float_type rpg, rcond, ferr, berr; 00484 int info; 00485 mem_usage_t memusage; 00486 SuperLUStat_t stat; 00487 00488 StatInit(&stat); 00489 SuperLUSolveChooser<T>::solve(&options, &static_cast<SuperMatrix&>(mat), perm_c, perm_r, etree, &equed, R, C, 00490 &L, &U, work, lwork, &B, &X, &rpg, &rcond, &ferr, 00491 &berr, &memusage, &stat, &info); 00492 00493 if(verbose){ 00494 dinfo<<"LU factorization: dgssvx() returns info "<< info<<std::endl; 00495 00496 if ( info == 0 || info == n+1 ) { 00497 00498 if ( options.PivotGrowth ) 00499 dinfo<<"Recip. pivot growth = "<<rpg<<std::endl; 00500 if ( options.ConditionNumber ) 00501 dinfo<<"Recip. condition number = %e\n"<< rcond<<std::endl; 00502 SCformat* Lstore = (SCformat *) L.Store; 00503 NCformat* Ustore = (NCformat *) U.Store; 00504 dinfo<<"No of nonzeros in factor L = "<< Lstore->nnz<<std::endl; 00505 dinfo<<"No of nonzeros in factor U = "<< Ustore->nnz<<std::endl; 00506 dinfo<<"No of nonzeros in L+U = "<< Lstore->nnz + Ustore->nnz - n<<std::endl; 00507 QuerySpaceChooser<T>::querySpace(&L, &U, &memusage); 00508 dinfo<<"L\\U MB "<<memusage.for_lu/1e6<<" \ttotal MB needed "<<memusage.total_needed/1e6 00509 <<" \texpansions "; 00510 00511 #ifdef HAVE_MEM_USAGE_T_EXPANSIONS 00512 std::cout<<memusage.expansions<<std::endl; 00513 #else 00514 std::cout<<stat.expansions<<std::endl; 00515 #endif 00516 } else if ( info > 0 && lwork == -1 ) { 00517 dinfo<<"** Estimated memory: "<< info - n<<std::endl; 00518 } 00519 if ( options.PrintStat ) StatPrint(&stat); 00520 } 00521 StatFree(&stat); 00522 /* 00523 NCformat* Ustore = (NCformat *) U.Store; 00524 int k=0; 00525 dPrint_CompCol_Matrix("U", &U); 00526 for(int i=0; i < U.ncol; ++i, ++k){ 00527 std::cout<<i<<": "; 00528 for(int c=Ustore->colptr[i]; c < Ustore->colptr[i+1]; ++c) 00529 //if(Ustore->rowind[c]==i) 00530 std::cout<<Ustore->rowind[c]<<"->"<<((double*)Ustore->nzval)[c]<<" "; 00531 if(k==0){ 00532 // 00533 k=-1; 00534 }std::cout<<std::endl; 00535 } 00536 dPrint_SuperNode_Matrix("L", &L); 00537 for(int i=0; i < U.ncol; ++i, ++k){ 00538 std::cout<<i<<": "; 00539 for(int c=Ustore->colptr[i]; c < Ustore->colptr[i+1]; ++c) 00540 //if(Ustore->rowind[c]==i) 00541 std::cout<<Ustore->rowind[c]<<"->"<<((double*)Ustore->nzval)[c]<<" "; 00542 if(k==0){ 00543 // 00544 k=-1; 00545 }std::cout<<std::endl; 00546 } */ 00547 options.Fact = FACTORED; 00548 } 00549 00550 template<typename T, typename A, int n, int m> 00551 void SuperLU<BCRSMatrix<FieldMatrix<T,n,m>,A> > 00552 ::apply(domain_type& x, range_type& b, InverseOperatorResult& res) 00553 { 00554 if(mat.M()+mat.N()==0) 00555 DUNE_THROW(ISTLError, "Matrix of SuperLU is null!"); 00556 00557 if(first){ 00558 SuperLUDenseMatChooser<T>::create(&B, (int)mat.N(), 1, reinterpret_cast<T*>(&b[0]), (int)mat.N(), SLU_DN, GetSuperLUType<T>::type, SLU_GE); 00559 SuperLUDenseMatChooser<T>::create(&X, (int)mat.N(), 1, reinterpret_cast<T*>(&x[0]), (int)mat.N(), SLU_DN, GetSuperLUType<T>::type, SLU_GE); 00560 first=false; 00561 }else{ 00562 ((DNformat*) B.Store)->nzval=&b[0]; 00563 ((DNformat*)X.Store)->nzval=&x[0]; 00564 } 00565 typename GetSuperLUType<T>::float_type rpg, rcond, ferr, berr; 00566 int info; 00567 mem_usage_t memusage; 00568 SuperLUStat_t stat; 00569 /* Initialize the statistics variables. */ 00570 StatInit(&stat); 00571 /* 00572 range_type d=b; 00573 a->usmv(-1, x, d); 00574 00575 double def0=d.two_norm(); 00576 */ 00577 #ifdef SUPERLU_MIN_VERSION_4_3 00578 options.IterRefine=SLU_DOUBLE; 00579 #else 00580 options.IterRefine=DOUBLE; 00581 #endif 00582 00583 SuperLUSolveChooser<T>::solve(&options, &static_cast<SuperMatrix&>(mat), perm_c, perm_r, etree, &equed, R, C, 00584 &L, &U, work, lwork, &B, &X, &rpg, &rcond, &ferr, &berr, 00585 &memusage, &stat, &info); 00586 00587 res.iterations=1; 00588 00589 /* 00590 if(options.Equil==YES) 00591 // undo scaling of right hand side 00592 std::transform(reinterpret_cast<T*>(&b[0]),reinterpret_cast<T*>(&b[0])+mat.M(), 00593 C, reinterpret_cast<T*>(&d[0]), std::divides<T>()); 00594 else 00595 d=b; 00596 a->usmv(-1, x, d); 00597 res.reduction=d.two_norm()/def0; 00598 res.conv_rate = res.reduction; 00599 res.converged=(res.reduction<1e-10||d.two_norm()<1e-18); 00600 */ 00601 res.converged=true; 00602 00603 if(verbose){ 00604 00605 dinfo<<"Triangular solve: dgssvx() returns info "<< info<<std::endl; 00606 00607 if ( info == 0 || info == n+1 ) { 00608 00609 if ( options.IterRefine ) { 00610 std::cout<<"Iterative Refinement: steps=" 00611 <<stat.RefineSteps<<" FERR="<<ferr<<" BERR="<<berr<<std::endl; 00612 }else 00613 std::cout<<" FERR="<<ferr<<" BERR="<<berr<<std::endl; 00614 } else if ( info > 0 && lwork == -1 ) { 00615 std::cout<<"** Estimated memory: "<< info - n<<" bytes"<<std::endl; 00616 } 00617 00618 if ( options.PrintStat ) StatPrint(&stat); 00619 } 00620 StatFree(&stat); 00621 } 00622 00623 template<typename T, typename A, int n, int m> 00624 void SuperLU<BCRSMatrix<FieldMatrix<T,n,m>,A> > 00625 ::apply(T* x, T* b) 00626 { 00627 if(mat.N()+mat.M()==0) 00628 DUNE_THROW(ISTLError, "Matrix of SuperLU is null!"); 00629 00630 if(first){ 00631 SuperLUDenseMatChooser<T>::create(&B, mat.N(), 1, b, mat.N(), SLU_DN, GetSuperLUType<T>::type, SLU_GE); 00632 SuperLUDenseMatChooser<T>::create(&X, mat.N(), 1, x, mat.N(), SLU_DN, GetSuperLUType<T>::type, SLU_GE); 00633 first=false; 00634 }else{ 00635 ((DNformat*) B.Store)->nzval=b; 00636 ((DNformat*)X.Store)->nzval=x; 00637 } 00638 00639 typename GetSuperLUType<T>::float_type rpg, rcond, ferr, berr; 00640 int info; 00641 mem_usage_t memusage; 00642 SuperLUStat_t stat; 00643 /* Initialize the statistics variables. */ 00644 StatInit(&stat); 00645 00646 #ifdef SUPERLU_MIN_VERSION_4_3 00647 options.IterRefine=SLU_DOUBLE; 00648 #else 00649 options.IterRefine=DOUBLE; 00650 #endif 00651 00652 SuperLUSolveChooser<T>::solve(&options, &static_cast<SuperMatrix&>(mat), perm_c, perm_r, etree, &equed, R, C, 00653 &L, &U, work, lwork, &B, &X, &rpg, &rcond, &ferr, &berr, 00654 &memusage, &stat, &info); 00655 00656 if(verbose){ 00657 dinfo<<"Triangular solve: dgssvx() returns info "<< info<<std::endl; 00658 00659 if ( info == 0 || info == n+1 ) { 00660 00661 if ( options.IterRefine ) { 00662 dinfo<<"Iterative Refinement: steps=" 00663 <<stat.RefineSteps<<" FERR="<<ferr<<" BERR="<<berr<<std::endl; 00664 }else 00665 dinfo<<" FERR="<<ferr<<" BERR="<<berr<<std::endl; 00666 } else if ( info > 0 && lwork == -1 ) { 00667 dinfo<<"** Estimated memory: "<< info - n<<" bytes"<<std::endl; 00668 } 00669 if ( options.PrintStat ) StatPrint(&stat); 00670 } 00671 00672 StatFree(&stat); 00673 } 00676 template<typename T, typename A, int n, int m> 00677 struct IsDirectSolver<SuperLU<BCRSMatrix<FieldMatrix<T,n,m>,A> > > 00678 { 00679 enum{ value=true}; 00680 }; 00681 } 00682 00683 #endif // HAVE_SUPERLU 00684 #endif // DUNE_SUPERLU_HH