1#ifndef DUNE_PDELAB_BACKEND_ISTL_OVLP_AMG_DG_BACKEND_HH
2#define DUNE_PDELAB_BACKEND_ISTL_OVLP_AMG_DG_BACKEND_HH
4#include <dune/common/parametertree.hh>
5#include <dune/common/power.hh>
7#include <dune/istl/matrixmatrix.hh>
9#include <dune/grid/common/datahandleif.hh>
32 :
public Dune::CommDataHandleIF<LocalGlobalMapDataHandle<GFS>,int>
36 typedef typename GFS::Traits::GridView
GV;
39 typedef typename GV::Grid
Grid;
41 typedef typename GlobalIdSet::IdType
IdType;
66 template<
class EntityType>
67 size_t size (EntityType& e)
const
73 template<
class MessageBuffer,
class EntityType>
74 void gather (MessageBuffer& buff,
const EntityType& e)
const
78 IdType myid = globalIdSet.id(e);
90 template<
class MessageBuffer,
class EntityType>
91 void scatter (MessageBuffer& buff,
const EntityType& e,
size_t n)
97 IdType myid = globalIdSet.id(e);
105 : gfs(gfs_), gv(gfs.gridView()), indexSet(gv.indexSet()),
106 grid(gv.grid()), globalIdSet(grid.globalIdSet()),
123 template<
class GFS,
class M>
125 :
public Dune::CommDataHandleIF<MatrixExchangeDataHandle<GFS,M>,
126 std::pair<typename GFS::Traits::GridView::Grid::Traits::GlobalIdSet::IdType,
127 typename M::block_type> >
131 typedef typename GFS::Traits::GridView
GV;
134 typedef typename GV::Grid
Grid;
136 typedef typename GlobalIdSet::IdType
IdType;
139 typedef typename M::block_type
B;
164 template<
class EntityType>
165 size_t size (EntityType& e)
const
168 typename M::row_type myrow = m[myindex];
169 typename M::row_type::iterator endit=myrow.end();
171 for (
typename M::row_type::iterator it=myrow.begin(); it!=endit; ++it)
173 typename LocalToGlobalMap::const_iterator find=l2g.find(it.index());
184 template<
class MessageBuffer,
class EntityType>
185 void gather (MessageBuffer& buff,
const EntityType& e)
const
189 typename M::row_type myrow = m[myindex];
190 typename M::row_type::iterator endit=myrow.end();
192 for (
typename M::row_type::iterator it=myrow.begin(); it!=endit; ++it)
194 typename LocalToGlobalMap::const_iterator find=l2g.find(it.index());
197 buff.write(std::make_pair(find->second,*it));
209 template<
class MessageBuffer,
class EntityType>
210 void scatter (MessageBuffer& buff,
const EntityType& e,
size_t n)
213 std::cout << gv.comm().rank() <<
": begin scatter local=" << myindex <<
" size=" << n << std::endl;
216 for (
size_t i=0; i<n; ++i)
219 std::cout << gv.comm().rank() <<
": --> received global=" << x.first << std::endl;
220 typename GlobalToLocalMap::const_iterator find=g2l.find(x.first);
224 if (m.exists(myindex,nbindex))
226 m[myindex][nbindex] = x.second;
228 block -= m2[myindex][nbindex];
229 std::cout << gv.comm().rank() <<
": compare i=" << myindex <<
" j=" << nbindex
230 <<
" norm=" << block.infinity_norm() << std::endl;
242 : gfs(gfs_), m(m_), gv(gfs.gridView()), indexSet(gv.indexSet()),
243 grid(gv.grid()), globalIdSet(grid.globalIdSet()),
244 l2g(l2g_), g2l(g2l_), m2(m2_)
263 template<
class GFS,
class T,
class A,
int n,
int m>
265 Dune::BCRSMatrix<Dune::FieldMatrix<T,n,m>,A>& matrix2)
267 typedef Dune::FieldMatrix<T,n,m> B;
268 typedef Dune::BCRSMatrix<B,A> M;
273 LocalToGlobalMap l2g;
274 GlobalToLocalMap g2l;
276 if (gfs.gridView().comm().size()>1)
277 gfs.gridView().communicate(lgdh,Dune::InteriorBorder_All_Interface,Dune::ForwardCommunication);
281 if (gfs.gridView().comm().size()>1)
282 gfs.gridView().communicate(mexdh,Dune::InteriorBorder_All_Interface,Dune::ForwardCommunication);
298 template<
class DGGFS,
class DGMatrix,
class DGPrec,
class DGCC,
299 class CGGFS,
class CGPrec,
class CGCC,
300 class P,
class DGHelper,
class Comm>
302 :
public Dune::Preconditioner<Dune::PDELab::Backend::Vector<DGGFS,typename DGPrec::domain_type::field_type>,
303 Dune::PDELab::Backend::Vector<DGGFS,typename DGPrec::range_type::field_type>>
313 const DGHelper& dghelper;
328 return SolverCategory::overlapping;
338 OvlpDGAMGPrec (
const DGGFS& dggfs_, DGMatrix& dgmatrix_, DGPrec& dgprec_,
const DGCC& dgcc_,
339 const CGGFS& cggfs_, CGPrec& cgprec_,
const CGCC& cgcc_, P& p_,
340 const DGHelper& dghelper_,
const Comm& comm_,
int n1_,
int n2_)
341 : dggfs(dggfs_), dgmatrix(dgmatrix_), dgprec(dgprec_), dgcc(dgcc_),
342 cggfs(cggfs_), cgprec(cgprec_), cgcc(cgcc_),
p(p_), dghelper(dghelper_),
343 comm(comm_), n1(n1_), n2(n2_)
352 virtual void pre (
V& x,
W& b)
override
355 dgprec.pre(native(x),native(b));
358 cgprec.pre(native(cgv),native(cgd));
366 virtual void apply (
V& x,
const W& b)
override
375 for (
int i=0; i<n1; i++)
379 dgprec.apply(native(v),native(d));
381 if (dggfs.gridView().comm().size()>1)
382 dggfs.gridView().communicate(adddh,Dune::All_All_Interface,Dune::ForwardCommunication);
383 dgmatrix.mmv(native(v),native(d));
389 dghelper.maskForeignDOFs(d);
391 p.mtv(native(d),native(cgd));
393 if (cggfs.gridView().comm().size()>1)
394 cggfs.gridView().communicate(adddh,Dune::All_All_Interface,Dune::ForwardCommunication);
396 comm.project(native(cgd));
401 cgprec.apply(native(cgv),native(cgd));
404 p.mv(native(cgv),native(v));
405 dgmatrix.mmv(native(v),native(d));
410 for (
int i=0; i<n2; i++)
413 dgprec.apply(native(v),native(d));
415 if (dggfs.gridView().comm().size()>1)
416 dggfs.gridView().communicate(adddh,Dune::All_All_Interface,Dune::ForwardCommunication);
417 dgmatrix.mmv(native(v),native(d));
463template<
class DGGO,
class DGCC,
class CGGFS,
class CGCC,
class TransferLOP,
464 template<
class,
class,
class,
int>
class DGPrec,
template<
class>
class Solver,
int s=96>
471 typedef typename DGGO::Traits::TrialGridFunctionSpace
GFS;
474 typedef typename DGGO::Traits::Jacobian
M;
475 typedef typename DGGO::Traits::Domain
V;
493 typedef typename Dune::TransposedMatMultMatResult<P,Matrix>::type
PTADG;
494 typedef typename Dune::MatMultMatResult<PTADG,P>::type
CGMatrix;
498 typedef Dune::OverlappingSchwarzOperator<CGMatrix,CGVector,CGVector,Comm>
ParCGOperator;
499 typedef Dune::SeqSSOR<CGMatrix,CGVector,CGVector,1>
Smoother;
500 typedef Dune::BlockPreconditioner<CGVector,CGVector,Comm,Smoother>
ParSmoother;
501 typedef Dune::Amg::AMG<ParCGOperator,CGVector,ParSmoother,Comm>
AMG;
511 std::shared_ptr<AMG> amg;
518 std::size_t low_order_space_entries_per_row;
521 double smoother_relaxation = 0.92;
558 amg_parameters = amg_parameters_;
570 return amg_parameters;
588 unsigned maxiter_=5000,
int verbose_=1,
bool reuse_=
false,
589 bool usesuperlu_=
true)
591 , gfs(dggo_.trialGridFunctionSpace())
596 , amg_parameters(15,2000)
601 , usesuperlu(usesuperlu_)
602 , low_order_space_entries_per_row(StaticPower<3,
GFS::Traits::GridView::dimension>::power)
604 , pgo(cggfs,dggo.trialGridFunctionSpace(),cgtodglop,
MBE(low_order_space_entries_per_row))
606 , acg(Backend::attached_container())
608 amg_parameters.setDefaultValuesIsotropic(GFS::Traits::GridViewType::Traits::Grid::dimension);
609 amg_parameters.setDebugLevel(verbose_);
611 if (usesuperlu ==
true)
613 if (gfs.gridView().comm().rank()==0)
614 std::cout <<
"WARNING: You are using AMG without SuperLU!"
615 <<
" Please consider installing SuperLU,"
616 <<
" or set the usesuperlu flag to false"
617 <<
" to suppress this warning." << std::endl;
623 if (verbose>0 && gfs.gridView().comm().rank()==0) std::cout <<
"allocated prolongation matrix of size " << pmatrix.N() <<
" x " << pmatrix.M() << std::endl;
632 const ParameterTree& params)
634 , gfs(dggo_.trialGridFunctionSpace())
639 , maxiter(params.get<int>(
"max_iterations",5000))
640 , amg_parameters(15,2000)
641 , verbose(params.get<int>(
"verbose",1))
642 , reuse(params.get<bool>(
"reuse",false))
644 , usesuperlu(params.get<bool>(
"use_superlu",true))
645 , low_order_space_entries_per_row(params.get<
std::size_t>(
"low_order_space.entries_per_row",StaticPower<3,
GFS::Traits::GridView::dimension>::power))
647 , pgo(cggfs,dggo.trialGridFunctionSpace(),cgtodglop,
MBE(low_order_space_entries_per_row))
649 , acg(Backend::attached_container())
651 amg_parameters.setDefaultValuesIsotropic(GFS::Traits::GridViewType::Traits::Grid::dimension);
652 amg_parameters.setDebugLevel(params.get<
int>(
"verbose",1));
654 if (usesuperlu ==
true)
656 if (gfs.gridView().comm().rank()==0)
657 std::cout <<
"WARNING: You are using AMG without SuperLU!"
658 <<
" Please consider installing SuperLU,"
659 <<
" or set the usesuperlu flag to false"
660 <<
" to suppress this warning." << std::endl;
666 if (verbose>0 && gfs.gridView().comm().rank()==0) std::cout <<
"allocated prolongation matrix of size " << pmatrix.N() <<
" x " << pmatrix.M() << std::endl;
675 smoother_relaxation = relaxation_;
698 void apply (
M& A,
V& z,
V& r,
typename Dune::template FieldTraits<typename V::ElementType >::real_type reduction)
710 CGGO cggo(cggfs,cgcc,cggfs,cgcc,emptylop,
MBE(low_order_space_entries_per_row));
716 double triple_product_time = 0.0;
718 if(reuse ==
false || firstapply ==
true) {
720 Dune::transposeMatMultMat(ptadg,native(pmatrix),native(A));
722 Dune::matMultMat(native(acg),ptadg,native(pmatrix));
723 triple_product_time = watch.elapsed();
724 if (verbose>0 && gfs.gridView().comm().rank()==0)
725 std::cout <<
"=== triple matrix product " << triple_product_time <<
" s" << std::endl;
731 else if(verbose>0 && gfs.gridView().comm().rank()==0)
732 std::cout <<
"=== reuse CG matrix, SKIPPING triple matrix product " << std::endl;
736 DGGOEmpty dggoempty(gfs,dgcc,gfs,dgcc,emptylop,
MBE(1 << GFS::Traits::GridView::dimension));
737 dggoempty.jacobian(z,A);
743 Comm oocc(gfs.gridView().comm());
745 CGHELPER cghelper(cggfs,2);
746 cghelper.createIndexSetAndProjectForAMG(acg,oocc);
748 Dune::OverlappingSchwarzScalarProduct<CGVector,Comm> sp(oocc);
750 typedef typename Dune::Amg::SmootherTraits<ParSmoother>::Arguments SmootherArgs;
751 SmootherArgs smootherArgs;
752 smootherArgs.iterations = 1;
753 smootherArgs.relaxationFactor = 1.0;
754 typedef Dune::Amg::CoarsenCriterion<Dune::Amg::SymmetricCriterion<CGMatrix,Dune::Amg::FirstDiagonal> > Criterion;
755 Criterion criterion(amg_parameters);
759 double amg_setup_time = 0.0;
760 if(reuse ==
false || firstapply ==
true) {
761 amg.reset(
new AMG(paroop,criterion,smootherArgs,oocc));
763 amg_setup_time = watch.elapsed();
764 if (verbose>0 && gfs.gridView().comm().rank()==0)
765 std::cout <<
"=== AMG setup " <<amg_setup_time <<
" s" << std::endl;
767 else if (verbose>0 && gfs.gridView().comm().rank()==0)
768 std::cout <<
"=== reuse CG matrix, SKIPPING AMG setup " << std::endl;
771 typedef DGPrec<Matrix,Vector,Vector,1> DGPrecType;
772 DGPrecType dgprec(native(A),1,smoother_relaxation);
775 HybridPrec hybridprec(gfs,native(A),dgprec,dgcc,cggfs,*amg,cgcc,native(pmatrix),
790 if (gfs.gridView().comm().rank()>0) verb=0;
791 Solver<V> solver(pop,psp,hybridprec,reduction,maxiter,verb);
794 Dune::InverseOperatorResult stat;
796 solver.apply(z,r,stat);
797 double amg_solve_time = watch.elapsed();
798 if (verbose>0 && gfs.gridView().comm().rank()==0) std::cout <<
"=== Hybrid total solve time " << amg_solve_time+amg_setup_time+triple_product_time <<
" s" << std::endl;
801 res.
elapsed = amg_solve_time+amg_setup_time+triple_product_time;
const std::string s
Definition: function.hh:843
static const int dim
Definition: adaptivity.hh:84
const P & p
Definition: constraints.hh:148
void set_constrained_dofs(const CG &cg, typename XG::ElementType x, XG &xg)
construct constraints from given boundary condition function
Definition: constraints.hh:796
For backward compatibility – Do not use this!
Definition: adaptivity.hh:28
void restore_overlap_entries(const GFS &gfs, Dune::BCRSMatrix< Dune::FieldMatrix< T, n, m >, A > &matrix, Dune::BCRSMatrix< Dune::FieldMatrix< T, n, m >, A > &matrix2)
Definition: ovlp_amg_dg_backend.hh:264
typename impl::BackendVectorSelector< GridFunctionSpace, FieldType >::Type Vector
alias of the return type of BackendVectorSelector
Definition: backend/interface.hh:106
std::enable_if< std::is_base_of< impl::WrapperBase, T >::value, Native< T > & >::type native(T &t)
Definition: backend/interface.hh:192
typename native_type< T >::type Native
Alias of the native container type associated with T or T itself if it is not a backend wrapper.
Definition: backend/interface.hh:176
Backend using (possibly nested) ISTL BCRSMatrices.
Definition: bcrsmatrixbackend.hh:188
Definition: ovlp_amg_dg_backend.hh:33
GV::IndexSet IndexSet
Definition: ovlp_amg_dg_backend.hh:37
std::map< IndexType, IdType > LocalToGlobalMap
Definition: ovlp_amg_dg_backend.hh:47
LocalGlobalMapDataHandle(const GFS &gfs_, LocalToGlobalMap &l2g_, GlobalToLocalMap &g2l_)
constructor
Definition: ovlp_amg_dg_backend.hh:104
bool fixedSize(int dim, int codim) const
returns true if size per entity of given dim and codim is a constant
Definition: ovlp_amg_dg_backend.hh:57
GlobalIdSet::IdType IdType
Definition: ovlp_amg_dg_backend.hh:41
GFS::Traits::GridView GV
Definition: ovlp_amg_dg_backend.hh:36
int DataType
export type of data for message buffer
Definition: ovlp_amg_dg_backend.hh:44
std::map< IdType, IndexType > GlobalToLocalMap
Definition: ovlp_amg_dg_backend.hh:48
void gather(MessageBuffer &buff, const EntityType &e) const
pack data from user to message buffer
Definition: ovlp_amg_dg_backend.hh:74
size_t size(EntityType &e) const
Definition: ovlp_amg_dg_backend.hh:67
bool contains(int dim, int codim) const
returns true if data for this codim should be communicated
Definition: ovlp_amg_dg_backend.hh:51
Grid::Traits::GlobalIdSet GlobalIdSet
Definition: ovlp_amg_dg_backend.hh:40
GV::Grid Grid
Definition: ovlp_amg_dg_backend.hh:39
IndexSet::IndexType IndexType
Definition: ovlp_amg_dg_backend.hh:38
void scatter(MessageBuffer &buff, const EntityType &e, size_t n)
Definition: ovlp_amg_dg_backend.hh:91
Definition: ovlp_amg_dg_backend.hh:128
bool contains(int dim, int codim) const
returns true if data for this codim should be communicated
Definition: ovlp_amg_dg_backend.hh:149
MatrixExchangeDataHandle(const GFS &gfs_, M &m_, const LocalToGlobalMap &l2g_, const GlobalToLocalMap &g2l_, M &m2_)
constructor
Definition: ovlp_amg_dg_backend.hh:241
std::map< IndexType, IdType > LocalToGlobalMap
Definition: ovlp_amg_dg_backend.hh:145
M::block_type B
Definition: ovlp_amg_dg_backend.hh:139
GV::Grid Grid
Definition: ovlp_amg_dg_backend.hh:134
void scatter(MessageBuffer &buff, const EntityType &e, size_t n)
Definition: ovlp_amg_dg_backend.hh:210
void gather(MessageBuffer &buff, const EntityType &e) const
pack data from user to message buffer
Definition: ovlp_amg_dg_backend.hh:185
Grid::Traits::GlobalIdSet GlobalIdSet
Definition: ovlp_amg_dg_backend.hh:135
std::map< IdType, IndexType > GlobalToLocalMap
Definition: ovlp_amg_dg_backend.hh:146
std::pair< IdType, B > DataType
export type of data for message buffer
Definition: ovlp_amg_dg_backend.hh:142
GV::IndexSet IndexSet
Definition: ovlp_amg_dg_backend.hh:132
IndexSet::IndexType IndexType
Definition: ovlp_amg_dg_backend.hh:133
GlobalIdSet::IdType IdType
Definition: ovlp_amg_dg_backend.hh:136
bool fixedSize(int dim, int codim) const
returns true if size per entity of given dim and codim is a constant
Definition: ovlp_amg_dg_backend.hh:155
GFS::Traits::GridView GV
Definition: ovlp_amg_dg_backend.hh:131
size_t size(EntityType &e) const
Definition: ovlp_amg_dg_backend.hh:165
Definition: ovlp_amg_dg_backend.hh:304
virtual void pre(V &x, W &b) override
Prepare the preconditioner.
Definition: ovlp_amg_dg_backend.hh:352
Dune::PDELab::Backend::Vector< CGGFS, typename CGPrec::domain_type::field_type > CGV
Definition: ovlp_amg_dg_backend.hh:322
virtual void apply(V &x, const W &b) override
Apply the precondioner.
Definition: ovlp_amg_dg_backend.hh:366
SolverCategory::Category category() const override
Definition: ovlp_amg_dg_backend.hh:326
Dune::PDELab::Backend::Vector< DGGFS, typename DGPrec::range_type::field_type > W
Definition: ovlp_amg_dg_backend.hh:319
Backend::Native< V > X
Definition: ovlp_amg_dg_backend.hh:320
Backend::Native< W > Y
Definition: ovlp_amg_dg_backend.hh:321
Dune::PDELab::Backend::Vector< CGGFS, typename CGPrec::range_type::field_type > CGW
Definition: ovlp_amg_dg_backend.hh:323
virtual void post(V &x) override
Clean up.
Definition: ovlp_amg_dg_backend.hh:428
Dune::PDELab::Backend::Vector< DGGFS, typename DGPrec::domain_type::field_type > V
Definition: ovlp_amg_dg_backend.hh:318
OvlpDGAMGPrec(const DGGFS &dggfs_, DGMatrix &dgmatrix_, DGPrec &dgprec_, const DGCC &dgcc_, const CGGFS &cggfs_, CGPrec &cgprec_, const CGCC &cgcc_, P &p_, const DGHelper &dghelper_, const Comm &comm_, int n1_, int n2_)
Constructor.
Definition: ovlp_amg_dg_backend.hh:338
Definition: ovlp_amg_dg_backend.hh:468
Dune::SeqSSOR< CGMatrix, CGVector, CGVector, 1 > Smoother
Definition: ovlp_amg_dg_backend.hh:499
Dune::TransposedMatMultMatResult< P, Matrix >::type PTADG
Definition: ovlp_amg_dg_backend.hh:493
Dune::BlockPreconditioner< CGVector, CGVector, Comm, Smoother > ParSmoother
Definition: ovlp_amg_dg_backend.hh:500
Dune::Amg::AMG< ParCGOperator, CGVector, ParSmoother, Comm > AMG
Definition: ovlp_amg_dg_backend.hh:501
DGGO::Traits::TrialGridFunctionSpace GFS
Definition: ovlp_amg_dg_backend.hh:471
DGGO::Traits::Jacobian M
Definition: ovlp_amg_dg_backend.hh:474
Backend::Native< CGV > CGVector
Definition: ovlp_amg_dg_backend.hh:482
void setDGSmootherRelaxation(double relaxation_)
set number of presmoothing steps on the DG level
Definition: ovlp_amg_dg_backend.hh:673
Dune::PDELab::ISTL::BCRSMatrixBackend MBE
Definition: ovlp_amg_dg_backend.hh:485
Backend::Native< V > Vector
Definition: ovlp_amg_dg_backend.hh:477
Dune::PDELab::Backend::Vector< CGGFS, field_type > CGV
Definition: ovlp_amg_dg_backend.hh:481
Backend::Native< M > Matrix
Definition: ovlp_amg_dg_backend.hh:476
void setNoDGPostSmoothSteps(int n2_)
set number of postsmoothing steps on the DG level
Definition: ovlp_amg_dg_backend.hh:685
bool getReuse() const
Return whether the AMG is reused during call to apply()
Definition: ovlp_amg_dg_backend.hh:580
Vector::field_type field_type
Definition: ovlp_amg_dg_backend.hh:478
Backend::Native< PMatrix > P
Definition: ovlp_amg_dg_backend.hh:490
DGGO::Traits::Domain V
Definition: ovlp_amg_dg_backend.hh:475
Dune::OverlappingSchwarzOperator< CGMatrix, CGVector, CGVector, Comm > ParCGOperator
Definition: ovlp_amg_dg_backend.hh:498
ISTLBackend_OVLP_AMG_4_DG(DGGO &dggo_, const DGCC &dgcc_, CGGFS &cggfs_, const CGCC &cgcc_, unsigned maxiter_=5000, int verbose_=1, bool reuse_=false, bool usesuperlu_=true)
Definition: ovlp_amg_dg_backend.hh:587
Dune::Amg::Parameters Parameters
Definition: ovlp_amg_dg_backend.hh:502
void setNoDGPreSmoothSteps(int n1_)
set number of presmoothing steps on the DG level
Definition: ovlp_amg_dg_backend.hh:679
Dune::PDELab::GridOperator< CGGFS, GFS, CGTODGLOP, MBE, field_type, field_type, field_type, CC, CC > PGO
Definition: ovlp_amg_dg_backend.hh:488
void apply(M &A, V &z, V &r, typename Dune::template FieldTraits< typename V::ElementType >::real_type reduction)
solve the given linear system
Definition: ovlp_amg_dg_backend.hh:698
Dune::MatMultMatResult< PTADG, P >::type CGMatrix
Definition: ovlp_amg_dg_backend.hh:494
PMatrix & prolongation_matrix()
Definition: ovlp_amg_dg_backend.hh:547
void setReuse(bool reuse_)
Set whether the AMG should be reused again during call to apply().
Definition: ovlp_amg_dg_backend.hh:574
TransferLOP CGTODGLOP
Definition: ovlp_amg_dg_backend.hh:487
void setParameters(const Parameters &amg_parameters_)
set AMG parameters
Definition: ovlp_amg_dg_backend.hh:556
Dune::PDELab::ISTL::CommSelector< s, Dune::MPIHelper::isFake >::type Comm
Definition: ovlp_amg_dg_backend.hh:497
Dune::PDELab::EmptyTransformation CC
Definition: ovlp_amg_dg_backend.hh:486
const Parameters & parameters() const
Get the parameters describing the behaviuour of AMG.
Definition: ovlp_amg_dg_backend.hh:568
PGO::Jacobian PMatrix
Definition: ovlp_amg_dg_backend.hh:489
ISTLBackend_OVLP_AMG_4_DG(DGGO &dggo_, const DGCC &dgcc_, CGGFS &cggfs_, const CGCC &cgcc_, const ParameterTree ¶ms)
Definition: ovlp_amg_dg_backend.hh:631
Definition: ovlpistlsolverbackend.hh:43
Definition: ovlpistlsolverbackend.hh:384
const ISTL::ParallelHelper< DGGO::Traits::TrialGridFunctionSpace > & parallelHelper() const
Definition: ovlpistlsolverbackend.hh:414
Definition: ovlpistlsolverbackend.hh:434
Definition: parallelhelper.hh:51
Dune::Amg::SequentialInformation type
Definition: parallelhelper.hh:429
RFType conv_rate
Definition: solver.hh:36
bool converged
Definition: solver.hh:32
RFType reduction
Definition: solver.hh:35
unsigned int iterations
Definition: solver.hh:33
double elapsed
Definition: solver.hh:34
Dune::PDELab::LinearSolverResult< double > res
Definition: solver.hh:63
Definition: constraintstransformation.hh:112
Definition: genericdatahandle.hh:667
Standard grid operator implementation.
Definition: gridoperator.hh:36
void jacobian(const Domain &x, Jacobian &a) const
Assembler jacobian.
Definition: gridoperator.hh:180
Dune::PDELab::Backend::Matrix< MBE, Domain, Range, field_type > Jacobian
The type of the jacobian.
Definition: gridoperator.hh:47
Default flags for all local operators.
Definition: flags.hh:19
Default class for additional methods in instationary local operators.
Definition: idefault.hh:90
Implements linear and nonlinear versions of jacobian_apply_volume() based on alpha_volume()
Definition: numericaljacobianapply.hh:34
sparsity pattern generator
Definition: pattern.hh:14
Definition: recipe-operator-splitting.cc:108