3#ifndef DUNE_ISTL_ILU_HH
4#define DUNE_ISTL_ILU_HH
14#include <dune/common/fmatrix.hh>
40 typedef typename M::block_type block;
51 for (
ij=(*i).begin();
ij.index()<i.index(); ++
ij)
57 (*ij).rightmultiply(*
jj);
64 if (
ik.index()==
jk.index())
73 if (
ik.index()<
jk.index())
81 if (
ij.index()!=i.index())
86 catch (Dune::FMatrixError &
e) {
88 << i.index() <<
"][" <<
ij.index() <<
"]" <<
e.what();
95 template<
class M,
class X,
class Y>
101 typedef typename Y::block_type
dblock;
102 typedef typename X::block_type
vblock;
109 for (
coliterator j=(*i).begin(); j.index()<i.index(); ++j)
110 (*j).mmv(v[j.index()],rhs);
120 for (j=(*i).beforeEnd(); j.index()>i.index(); --j)
121 (*j).mmv(v[j.index()],rhs);
123 (*j).umv(rhs,v[i.index()]);
136 template<
class K,
int n,
int m>
163 typedef typename M::field_type K;
164 typedef std::map<size_t, int> map;
211 ci.insert((*ik).first);
231 ILUij = ILU[i.index()].begin();
void bilu_backsolve(const M &A, X &v, const Y &d)
LU backsolve with stored inverse.
Definition ilu.hh:96
M::field_type & firstmatrixelement(M &A)
Definition ilu.hh:131
void bilu_decomposition(const M &A, int n, M &ILU)
Definition ilu.hh:156
void bilu0_decomposition(M &A)
compute ILU decomposition of A. A is overwritten by its decomposition
Definition ilu.hh:35
Definition basearray.hh:19
iterator begin()
begin iterator
Definition basearray.hh:170
Statistics about compression achieved in implicit mode.
Definition bcrsmatrix.hh:81
int c
Definition ilu.hh:29
int r
Definition ilu.hh:29
derive error class from the base class in common
Definition istlexception.hh:16
RowIterator beforeBegin()
Definition matrix.hh:629
RowIterator beforeEnd()
Definition matrix.hh:622
RowIterator end()
Get iterator to one beyond last row.
Definition matrix.hh:615
RowIterator begin()
Get iterator to first row.
Definition matrix.hh:609