Dune Core Modules (2.6.0)

superlu.hh
Go to the documentation of this file.
1// -*- tab-width: 4; indent-tabs-mode: nil; c-basic-offset: 2 -*-
2// vi: set et ts=4 sw=2 sts=2:
3#ifndef DUNE_ISTL_SUPERLU_HH
4#define DUNE_ISTL_SUPERLU_HH
5
6#if HAVE_SUPERLU
7
8#include "superlufunctions.hh"
9#include "solvers.hh"
10#include "supermatrix.hh"
11#include <algorithm>
12#include <functional>
13#include "bcrsmatrix.hh"
14#include "bvector.hh"
15#include "istlexception.hh"
20
21namespace Dune
22{
23
34 template<class Matrix>
35 class SuperLU
36 {};
37
38 template<class M, class T, class TM, class TD, class TA>
39 class SeqOverlappingSchwarz;
40
41 template<class T, bool tag>
42 struct SeqOverlappingSchwarzAssemblerHelper;
43
44 template<class T>
45 struct SuperLUSolveChooser
46 {};
47
48 template<class T>
49 struct SuperLUDenseMatChooser
50 {};
51
52 template<class T>
53 struct SuperLUQueryChooser
54 {};
55
56 template<class T>
57 struct QuerySpaceChooser
58 {};
59
60#if HAVE_SLU_SDEFS_H
61 template<>
62 struct SuperLUDenseMatChooser<float>
63 {
64 static void create(SuperMatrix *mat, int n, int m, float *dat, int n1,
65 Stype_t stype, Dtype_t dtype, Mtype_t mtype)
66 {
67 sCreate_Dense_Matrix(mat, n, m, dat, n1, stype, dtype, mtype);
68
69 }
70
71 static void destroy(SuperMatrix*)
72 {}
73
74 };
75 template<>
76 struct SuperLUSolveChooser<float>
77 {
78 static void solve(superlu_options_t *options, SuperMatrix *mat, int *perm_c, int *perm_r, int *etree,
79 char *equed, float *R, float *C, SuperMatrix *L, SuperMatrix *U,
80 void *work, int lwork, SuperMatrix *B, SuperMatrix *X,
81 float *rpg, float *rcond, float *ferr, float *berr,
82 mem_usage_t *memusage, SuperLUStat_t *stat, int *info)
83 {
84#if SUPERLU_MIN_VERSION_5
85 GlobalLU_t gLU;
86 sgssvx(options, mat, perm_c, perm_r, etree, equed, R, C,
87 L, U, work, lwork, B, X, rpg, rcond, ferr, berr,
88 &gLU, memusage, stat, info);
89#else
90 sgssvx(options, mat, perm_c, perm_r, etree, equed, R, C,
91 L, U, work, lwork, B, X, rpg, rcond, ferr, berr,
92 memusage, stat, info);
93#endif
94 }
95 };
96
97 template<>
98 struct QuerySpaceChooser<float>
99 {
100 static void querySpace(SuperMatrix* L, SuperMatrix* U, mem_usage_t* memusage)
101 {
102 sQuerySpace(L,U,memusage);
103 }
104 };
105
106#endif
107
108#if HAVE_SLU_DDEFS_H
109
110 template<>
111 struct SuperLUDenseMatChooser<double>
112 {
113 static void create(SuperMatrix *mat, int n, int m, double *dat, int n1,
114 Stype_t stype, Dtype_t dtype, Mtype_t mtype)
115 {
116 dCreate_Dense_Matrix(mat, n, m, dat, n1, stype, dtype, mtype);
117
118 }
119
120 static void destroy(SuperMatrix * /* mat */)
121 {}
122 };
123 template<>
124 struct SuperLUSolveChooser<double>
125 {
126 static void solve(superlu_options_t *options, SuperMatrix *mat, int *perm_c, int *perm_r, int *etree,
127 char *equed, double *R, double *C, SuperMatrix *L, SuperMatrix *U,
128 void *work, int lwork, SuperMatrix *B, SuperMatrix *X,
129 double *rpg, double *rcond, double *ferr, double *berr,
130 mem_usage_t *memusage, SuperLUStat_t *stat, int *info)
131 {
132#if SUPERLU_MIN_VERSION_5
133 GlobalLU_t gLU;
134 dgssvx(options, mat, perm_c, perm_r, etree, equed, R, C,
135 L, U, work, lwork, B, X, rpg, rcond, ferr, berr,
136 &gLU, memusage, stat, info);
137#else
138 dgssvx(options, mat, perm_c, perm_r, etree, equed, R, C,
139 L, U, work, lwork, B, X, rpg, rcond, ferr, berr,
140 memusage, stat, info);
141#endif
142 }
143 };
144
145 template<>
146 struct QuerySpaceChooser<double>
147 {
148 static void querySpace(SuperMatrix* L, SuperMatrix* U, mem_usage_t* memusage)
149 {
150 dQuerySpace(L,U,memusage);
151 }
152 };
153#endif
154
155#if HAVE_SLU_ZDEFS_H
156 template<>
157 struct SuperLUDenseMatChooser<std::complex<double> >
158 {
159 static void create(SuperMatrix *mat, int n, int m, std::complex<double> *dat, int n1,
160 Stype_t stype, Dtype_t dtype, Mtype_t mtype)
161 {
162 zCreate_Dense_Matrix(mat, n, m, reinterpret_cast<doublecomplex*>(dat), n1, stype, dtype, mtype);
163
164 }
165
166 static void destroy(SuperMatrix*)
167 {}
168 };
169
170 template<>
171 struct SuperLUSolveChooser<std::complex<double> >
172 {
173 static void solve(superlu_options_t *options, SuperMatrix *mat, int *perm_c, int *perm_r, int *etree,
174 char *equed, double *R, double *C, SuperMatrix *L, SuperMatrix *U,
175 void *work, int lwork, SuperMatrix *B, SuperMatrix *X,
176 double *rpg, double *rcond, double *ferr, double *berr,
177 mem_usage_t *memusage, SuperLUStat_t *stat, int *info)
178 {
179#if SUPERLU_MIN_VERSION_5
180 GlobalLU_t gLU;
181 zgssvx(options, mat, perm_c, perm_r, etree, equed, R, C,
182 L, U, work, lwork, B, X, rpg, rcond, ferr, berr,
183 &gLU, memusage, stat, info);
184#else
185 zgssvx(options, mat, perm_c, perm_r, etree, equed, R, C,
186 L, U, work, lwork, B, X, rpg, rcond, ferr, berr,
187 memusage, stat, info);
188#endif
189 }
190 };
191
192 template<>
193 struct QuerySpaceChooser<std::complex<double> >
194 {
195 static void querySpace(SuperMatrix* L, SuperMatrix* U, mem_usage_t* memusage)
196 {
197 zQuerySpace(L,U,memusage);
198 }
199 };
200#endif
201
202#if HAVE_SLU_CDEFS_H
203 template<>
204 struct SuperLUDenseMatChooser<std::complex<float> >
205 {
206 static void create(SuperMatrix *mat, int n, int m, std::complex<float> *dat, int n1,
207 Stype_t stype, Dtype_t dtype, Mtype_t mtype)
208 {
209 cCreate_Dense_Matrix(mat, n, m, reinterpret_cast< ::complex*>(dat), n1, stype, dtype, mtype);
210
211 }
212
213 static void destroy(SuperMatrix* /* mat */)
214 {}
215 };
216
217 template<>
218 struct SuperLUSolveChooser<std::complex<float> >
219 {
220 static void solve(superlu_options_t *options, SuperMatrix *mat, int *perm_c, int *perm_r, int *etree,
221 char *equed, float *R, float *C, SuperMatrix *L, SuperMatrix *U,
222 void *work, int lwork, SuperMatrix *B, SuperMatrix *X,
223 float *rpg, float *rcond, float *ferr, float *berr,
224 mem_usage_t *memusage, SuperLUStat_t *stat, int *info)
225 {
226#if SUPERLU_MIN_VERSION_5
227 GlobalLU_t gLU;
228 cgssvx(options, mat, perm_c, perm_r, etree, equed, R, C,
229 L, U, work, lwork, B, X, rpg, rcond, ferr, berr,
230 &gLU, memusage, stat, info);
231#else
232 cgssvx(options, mat, perm_c, perm_r, etree, equed, R, C,
233 L, U, work, lwork, B, X, rpg, rcond, ferr, berr,
234 memusage, stat, info);
235#endif
236 }
237 };
238
239 template<>
240 struct QuerySpaceChooser<std::complex<float> >
241 {
242 static void querySpace(SuperMatrix* L, SuperMatrix* U, mem_usage_t* memusage)
243 {
244 cQuerySpace(L,U,memusage);
245 }
246 };
247#endif
248
262 template<typename T, typename A, int n, int m>
263 class SuperLU<BCRSMatrix<FieldMatrix<T,n,m>,A > >
264 : public InverseOperator<
265 BlockVector<FieldVector<T,m>,
266 typename A::template rebind<FieldVector<T,m> >::other>,
267 BlockVector<FieldVector<T,n>,
268 typename A::template rebind<FieldVector<T,n> >::other> >
269 {
270 public:
277 typedef SuperMatrixInitializer<BCRSMatrix<FieldMatrix<T,n,m>,A> > MatrixInitializer;
279 typedef Dune::BlockVector<
281 typename A::template rebind<FieldVector<T,m> >::other> domain_type;
283 typedef Dune::BlockVector<
285 typename A::template rebind<FieldVector<T,n> >::other> range_type;
286
289 {
290 return SolverCategory::Category::sequential;
291 }
292
307 explicit SuperLU(const Matrix& mat, bool verbose=false,
308 bool reusevector=true);
315 SuperLU();
316
317 ~SuperLU();
318
322 void apply(domain_type& x, range_type& b, InverseOperatorResult& res);
323
327 void apply (domain_type& x, range_type& b, double reduction, InverseOperatorResult& res)
328 {
329 DUNE_UNUSED_PARAMETER(reduction);
330 apply(x,b,res);
331 }
332
336 void apply(T* x, T* b);
337
339 void setMatrix(const Matrix& mat);
340
341 typename SuperLUMatrix::size_type nnz() const
342 {
343 return mat.nnz();
344 }
345
346 template<class S>
347 void setSubMatrix(const Matrix& mat, const S& rowIndexSet);
348
349 void setVerbosity(bool v);
350
355 void free();
356
357 const char* name() { return "SuperLU"; }
358 private:
359 template<class M,class X, class TM, class TD, class T1>
360 friend class SeqOverlappingSchwarz;
361 friend struct SeqOverlappingSchwarzAssemblerHelper<SuperLU<Matrix>,true>;
362
363 SuperLUMatrix& getInternalMatrix() { return mat; }
364
366 void decompose();
367
368 SuperLUMatrix mat;
369 SuperMatrix L, U, B, X;
370 int *perm_c, *perm_r, *etree;
371 typename GetSuperLUType<T>::float_type *R, *C;
372 T *bstore;
373 superlu_options_t options;
374 char equed;
375 void *work;
376 int lwork;
377 bool first, verbose, reusevector;
378 };
379
380 template<typename T, typename A, int n, int m>
381 SuperLU<BCRSMatrix<FieldMatrix<T,n,m>,A> >
382 ::~SuperLU()
383 {
384 if(mat.N()+mat.M()>0)
385 free();
386 }
387
388 template<typename T, typename A, int n, int m>
389 void SuperLU<BCRSMatrix<FieldMatrix<T,n,m>,A> >::free()
390 {
391 delete[] perm_c;
392 delete[] perm_r;
393 delete[] etree;
394 delete[] R;
395 delete[] C;
396 if(lwork>=0) {
397 Destroy_SuperNode_Matrix(&L);
398 Destroy_CompCol_Matrix(&U);
399 }
400 lwork=0;
401 if(!first && reusevector) {
402 SUPERLU_FREE(B.Store);
403 SUPERLU_FREE(X.Store);
404 }
405 mat.free();
406 }
407
408 template<typename T, typename A, int n, int m>
409 SuperLU<BCRSMatrix<FieldMatrix<T,n,m>,A> >
410 ::SuperLU(const Matrix& mat_, bool verbose_, bool reusevector_)
411 : work(0), lwork(0), first(true), verbose(verbose_),
412 reusevector(reusevector_)
413 {
414 setMatrix(mat_);
415
416 }
417 template<typename T, typename A, int n, int m>
418 SuperLU<BCRSMatrix<FieldMatrix<T,n,m>,A> >::SuperLU()
419 : work(0), lwork(0),verbose(false),
420 reusevector(false)
421 {}
422 template<typename T, typename A, int n, int m>
423 void SuperLU<BCRSMatrix<FieldMatrix<T,n,m>,A> >::setVerbosity(bool v)
424 {
425 verbose=v;
426 }
427
428 template<typename T, typename A, int n, int m>
429 void SuperLU<BCRSMatrix<FieldMatrix<T,n,m>,A> >::setMatrix(const Matrix& mat_)
430 {
431 if(mat.N()+mat.M()>0) {
432 free();
433 }
434 lwork=0;
435 work=0;
436 //a=&mat_;
437 mat=mat_;
438 decompose();
439 }
440
441 template<typename T, typename A, int n, int m>
442 template<class S>
443 void SuperLU<BCRSMatrix<FieldMatrix<T,n,m>,A> >::setSubMatrix(const Matrix& mat_,
444 const S& mrs)
445 {
446 if(mat.N()+mat.M()>0) {
447 free();
448 }
449 lwork=0;
450 work=0;
451 //a=&mat_;
452 mat.setMatrix(mat_,mrs);
453 decompose();
454 }
455
456 template<typename T, typename A, int n, int m>
457 void SuperLU<BCRSMatrix<FieldMatrix<T,n,m>,A> >::decompose()
458 {
459
460 first = true;
461 perm_c = new int[mat.M()];
462 perm_r = new int[mat.N()];
463 etree = new int[mat.M()];
464 R = new typename GetSuperLUType<T>::float_type[mat.N()];
465 C = new typename GetSuperLUType<T>::float_type[mat.M()];
466
467 set_default_options(&options);
468 // Do the factorization
469 B.ncol=0;
470 B.Stype=SLU_DN;
471 B.Dtype=GetSuperLUType<T>::type;
472 B.Mtype= SLU_GE;
473 DNformat fakeFormat;
474 fakeFormat.lda=mat.N();
475 B.Store=&fakeFormat;
476 X.Stype=SLU_DN;
477 X.Dtype=GetSuperLUType<T>::type;
478 X.Mtype= SLU_GE;
479 X.ncol=0;
480 X.Store=&fakeFormat;
481
482 typename GetSuperLUType<T>::float_type rpg, rcond, ferr=1e10, berr=1e10;
483 int info;
484 mem_usage_t memusage;
485 SuperLUStat_t stat;
486
487 StatInit(&stat);
488 SuperLUSolveChooser<T>::solve(&options, &static_cast<SuperMatrix&>(mat), perm_c, perm_r, etree, &equed, R, C,
489 &L, &U, work, lwork, &B, &X, &rpg, &rcond, &ferr,
490 &berr, &memusage, &stat, &info);
491
492 if(verbose) {
493 dinfo<<"LU factorization: dgssvx() returns info "<< info<<std::endl;
494
495 if ( info == 0 || info == n+1 ) {
496
497 if ( options.PivotGrowth )
498 dinfo<<"Recip. pivot growth = "<<rpg<<std::endl;
499 if ( options.ConditionNumber )
500 dinfo<<"Recip. condition number = %e\n"<< rcond<<std::endl;
501 SCformat* Lstore = (SCformat *) L.Store;
502 NCformat* Ustore = (NCformat *) U.Store;
503 dinfo<<"No of nonzeros in factor L = "<< Lstore->nnz<<std::endl;
504 dinfo<<"No of nonzeros in factor U = "<< Ustore->nnz<<std::endl;
505 dinfo<<"No of nonzeros in L+U = "<< Lstore->nnz + Ustore->nnz - n<<std::endl;
506 QuerySpaceChooser<T>::querySpace(&L, &U, &memusage);
507 dinfo<<"L\\U MB "<<memusage.for_lu/1e6<<" \ttotal MB needed "<<memusage.total_needed/1e6
508 <<" \texpansions ";
509 std::cout<<stat.expansions<<std::endl;
510
511 } else if ( info > 0 && lwork == -1 ) {
512 dinfo<<"** Estimated memory: "<< info - n<<std::endl;
513 }
514 if ( options.PrintStat ) StatPrint(&stat);
515 }
516 StatFree(&stat);
517 /*
518 NCformat* Ustore = (NCformat *) U.Store;
519 int k=0;
520 dPrint_CompCol_Matrix("U", &U);
521 for(int i=0; i < U.ncol; ++i, ++k){
522 std::cout<<i<<": ";
523 for(int c=Ustore->colptr[i]; c < Ustore->colptr[i+1]; ++c)
524 //if(Ustore->rowind[c]==i)
525 std::cout<<Ustore->rowind[c]<<"->"<<((double*)Ustore->nzval)[c]<<" ";
526 if(k==0){
527 //
528 k=-1;
529 }std::cout<<std::endl;
530 }
531 dPrint_SuperNode_Matrix("L", &L);
532 for(int i=0; i < U.ncol; ++i, ++k){
533 std::cout<<i<<": ";
534 for(int c=Ustore->colptr[i]; c < Ustore->colptr[i+1]; ++c)
535 //if(Ustore->rowind[c]==i)
536 std::cout<<Ustore->rowind[c]<<"->"<<((double*)Ustore->nzval)[c]<<" ";
537 if(k==0){
538 //
539 k=-1;
540 }std::cout<<std::endl;
541 } */
542 options.Fact = FACTORED;
543 }
544
545 template<typename T, typename A, int n, int m>
546 void SuperLU<BCRSMatrix<FieldMatrix<T,n,m>,A> >
548 {
549 if (mat.N() != b.dim())
550 DUNE_THROW(ISTLError, "Size of right-hand-side vector b does not match the number of matrix rows!");
551 if (mat.M() != x.dim())
552 DUNE_THROW(ISTLError, "Size of solution vector x does not match the number of matrix columns!");
553 if (mat.M()+mat.N()==0)
554 DUNE_THROW(ISTLError, "Matrix of SuperLU is null!");
555
556 SuperMatrix* mB = &B;
557 SuperMatrix* mX = &X;
558 SuperMatrix rB, rX;
559 if (reusevector) {
560 if(first) {
561 SuperLUDenseMatChooser<T>::create(&B, (int)mat.N(), 1, reinterpret_cast<T*>(&b[0]), (int)mat.N(), SLU_DN, GetSuperLUType<T>::type, SLU_GE);
562 SuperLUDenseMatChooser<T>::create(&X, (int)mat.N(), 1, reinterpret_cast<T*>(&x[0]), (int)mat.N(), SLU_DN, GetSuperLUType<T>::type, SLU_GE);
563 first=false;
564 }else{
565 ((DNformat*)B.Store)->nzval=&b[0];
566 ((DNformat*)X.Store)->nzval=&x[0];
567 }
568 } else {
569 SuperLUDenseMatChooser<T>::create(&rB, (int)mat.N(), 1, reinterpret_cast<T*>(&b[0]), (int)mat.N(), SLU_DN, GetSuperLUType<T>::type, SLU_GE);
570 SuperLUDenseMatChooser<T>::create(&rX, (int)mat.N(), 1, reinterpret_cast<T*>(&x[0]), (int)mat.N(), SLU_DN, GetSuperLUType<T>::type, SLU_GE);
571 mB = &rB;
572 mX = &rX;
573 }
574 typename GetSuperLUType<T>::float_type rpg, rcond, ferr=1e10, berr;
575 int info;
576 mem_usage_t memusage;
577 SuperLUStat_t stat;
578 /* Initialize the statistics variables. */
579 StatInit(&stat);
580 /*
581 range_type d=b;
582 a->usmv(-1, x, d);
583
584 double def0=d.two_norm();
585 */
586#ifdef SUPERLU_MIN_VERSION_4_3
587 options.IterRefine=SLU_DOUBLE;
588#else
589 options.IterRefine=DOUBLE;
590#endif
591
592 SuperLUSolveChooser<T>::solve(&options, &static_cast<SuperMatrix&>(mat), perm_c, perm_r, etree, &equed, R, C,
593 &L, &U, work, lwork, mB, mX, &rpg, &rcond, &ferr, &berr,
594 &memusage, &stat, &info);
595
596 res.iterations=1;
597
598 /*
599 if(options.Equil==YES)
600 // undo scaling of right hand side
601 std::transform(reinterpret_cast<T*>(&b[0]),reinterpret_cast<T*>(&b[0])+mat.M(),
602 C, reinterpret_cast<T*>(&d[0]), std::divides<T>());
603 else
604 d=b;
605 a->usmv(-1, x, d);
606 res.reduction=d.two_norm()/def0;
607 res.conv_rate = res.reduction;
608 res.converged=(res.reduction<1e-10||d.two_norm()<1e-18);
609 */
610 res.converged=true;
611
612 if(verbose) {
613
614 dinfo<<"Triangular solve: dgssvx() returns info "<< info<<std::endl;
615
616 if ( info == 0 || info == n+1 ) {
617
618 if ( options.IterRefine ) {
619 std::cout<<"Iterative Refinement: steps="
620 <<stat.RefineSteps<<" FERR="<<ferr<<" BERR="<<berr<<std::endl;
621 }else
622 std::cout<<" FERR="<<ferr<<" BERR="<<berr<<std::endl;
623 } else if ( info > 0 && lwork == -1 ) {
624 std::cout<<"** Estimated memory: "<< info - n<<" bytes"<<std::endl;
625 }
626
627 if ( options.PrintStat ) StatPrint(&stat);
628 }
629 StatFree(&stat);
630 if (!reusevector) {
631 SUPERLU_FREE(rB.Store);
632 SUPERLU_FREE(rX.Store);
633 }
634 }
635
636 template<typename T, typename A, int n, int m>
637 void SuperLU<BCRSMatrix<FieldMatrix<T,n,m>,A> >
638 ::apply(T* x, T* b)
639 {
640 if(mat.N()+mat.M()==0)
641 DUNE_THROW(ISTLError, "Matrix of SuperLU is null!");
642
643 SuperMatrix* mB = &B;
644 SuperMatrix* mX = &X;
645 SuperMatrix rB, rX;
646 if (reusevector) {
647 if(first) {
648 SuperLUDenseMatChooser<T>::create(&B, mat.N(), 1, b, mat.N(), SLU_DN, GetSuperLUType<T>::type, SLU_GE);
649 SuperLUDenseMatChooser<T>::create(&X, mat.N(), 1, x, mat.N(), SLU_DN, GetSuperLUType<T>::type, SLU_GE);
650 first=false;
651 }else{
652 ((DNformat*) B.Store)->nzval=b;
653 ((DNformat*)X.Store)->nzval=x;
654 }
655 } else {
656 SuperLUDenseMatChooser<T>::create(&rB, mat.N(), 1, b, mat.N(), SLU_DN, GetSuperLUType<T>::type, SLU_GE);
657 SuperLUDenseMatChooser<T>::create(&rX, mat.N(), 1, x, mat.N(), SLU_DN, GetSuperLUType<T>::type, SLU_GE);
658 mB = &rB;
659 mX = &rX;
660 }
661
662 typename GetSuperLUType<T>::float_type rpg, rcond, ferr=1e10, berr;
663 int info;
664 mem_usage_t memusage;
665 SuperLUStat_t stat;
666 /* Initialize the statistics variables. */
667 StatInit(&stat);
668
669#ifdef SUPERLU_MIN_VERSION_4_3
670 options.IterRefine=SLU_DOUBLE;
671#else
672 options.IterRefine=DOUBLE;
673#endif
674
675 SuperLUSolveChooser<T>::solve(&options, &static_cast<SuperMatrix&>(mat), perm_c, perm_r, etree, &equed, R, C,
676 &L, &U, work, lwork, mB, mX, &rpg, &rcond, &ferr, &berr,
677 &memusage, &stat, &info);
678
679 if(verbose) {
680 dinfo<<"Triangular solve: dgssvx() returns info "<< info<<std::endl;
681
682 if ( info == 0 || info == n+1 ) {
683
684 if ( options.IterRefine ) {
685 dinfo<<"Iterative Refinement: steps="
686 <<stat.RefineSteps<<" FERR="<<ferr<<" BERR="<<berr<<std::endl;
687 }else
688 dinfo<<" FERR="<<ferr<<" BERR="<<berr<<std::endl;
689 } else if ( info > 0 && lwork == -1 ) {
690 dinfo<<"** Estimated memory: "<< info - n<<" bytes"<<std::endl;
691 }
692 if ( options.PrintStat ) StatPrint(&stat);
693 }
694
695 StatFree(&stat);
696 if (!reusevector) {
697 SUPERLU_FREE(rB.Store);
698 SUPERLU_FREE(rX.Store);
699 }
700 }
703 template<typename T, typename A, int n, int m>
704 struct IsDirectSolver<SuperLU<BCRSMatrix<FieldMatrix<T,n,m>,A> > >
705 {
706 enum { value=true};
707 };
708
709 template<typename T, typename A, int n, int m>
710 struct StoresColumnCompressed<SuperLU<BCRSMatrix<FieldMatrix<T,n,m>,A> > >
711 {
712 enum { value = true };
713 };
714}
715
716// undefine macros from SuperLU's slu_util.h
717#undef FIRSTCOL_OF_SNODE
718#undef NO_MARKER
719#undef NUM_TEMPV
720#undef USER_ABORT
721#undef USER_MALLOC
722#undef SUPERLU_MALLOC
723#undef USER_FREE
724#undef SUPERLU_FREE
725#undef CHECK_MALLOC
726#undef SUPERLU_MAX
727#undef SUPERLU_MIN
728#undef L_SUB_START
729#undef L_SUB
730#undef L_NZ_START
731#undef L_FST_SUPC
732#undef U_NZ_START
733#undef U_SUB
734#undef TRUE
735#undef FALSE
736#undef EMPTY
737#undef NODROP
738#undef DROP_BASIC
739#undef DROP_PROWS
740#undef DROP_COLUMN
741#undef DROP_AREA
742#undef DROP_SECONDARY
743#undef DROP_DYNAMIC
744#undef DROP_INTERP
745#undef MILU_ALPHA
746
747#endif // HAVE_SUPERLU
748#endif // DUNE_SUPERLU_HH
Implementation of the BCRSMatrix class.
This file implements a vector space as a tensor product of a given vector space. The number of compon...
A sparse block matrix with compressed row storage.
Definition: bcrsmatrix.hh:423
A vector of blocks with memory management.
Definition: bvector.hh:317
A dense n x m matrix.
Definition: fmatrix.hh:68
vector space out of a tensor product of fields.
Definition: fvector.hh:93
derive error class from the base class in common
Definition: istlexception.hh:16
Abstract base class for all solvers.
Definition: solver.hh:91
A generic dynamic dense matrix.
Definition: matrix.hh:555
size_type M() const
Return the number of columns.
Definition: matrix.hh:695
size_type N() const
Return the number of rows.
Definition: matrix.hh:690
Dune::BlockVector< FieldVector< T, n >, typename A::template rebind< FieldVector< T, n > >::other > range_type
The type of the range of the solver.
Definition: superlu.hh:285
virtual SolverCategory::Category category() const
Category of the solver (see SolverCategory::Category)
Definition: superlu.hh:288
Dune::SuperLUMatrix< Matrix > SuperLUMatrix
The corresponding SuperLU Matrix type.
Definition: superlu.hh:275
Dune::BlockVector< FieldVector< T, m >, typename A::template rebind< FieldVector< T, m > >::other > domain_type
The type of the domain of the solver.
Definition: superlu.hh:281
void apply(domain_type &x, range_type &b, double reduction, InverseOperatorResult &res)
apply inverse operator, with given convergence criteria.
Definition: superlu.hh:327
SuperMatrixInitializer< BCRSMatrix< FieldMatrix< T, n, m >, A > > MatrixInitializer
Type of an associated initializer class.
Definition: superlu.hh:277
Implements a matrix constructed from a given type representing a field and compile-time given number ...
Implements a vector constructed from a given type representing a field and a compile-time given size.
decltype(auto) apply(F &&f, ArgTuple &&args)
Apply function with arguments given as tuple.
Definition: apply.hh:58
#define DUNE_UNUSED_PARAMETER(parm)
A macro to mark intentionally unused function parameters with.
Definition: unused.hh:25
#define DUNE_THROW(E, m)
Definition: exceptions.hh:216
DInfoType dinfo(std::cout)
Stream for informative output.
Definition: stdstreams.hh:138
Dune namespace.
Definition: alignedallocator.hh:10
STL namespace.
Implementations of the inverse operator interface.
Templates characterizing the type of a solver.
Standard Dune debug streams.
Statistics about the application of an inverse operator.
Definition: solver.hh:41
int iterations
Number of iterations.
Definition: solver.hh:59
bool converged
True if convergence criterion has been met.
Definition: solver.hh:65
Category
Definition: solvercategory.hh:21
Creative Commons License   |  Legal Statements / Impressum  |  Hosted by TU Dresden  |  generated with Hugo v0.111.3 (Nov 24, 23:30, 2024)