From c0049fe39eae975f988e4c1c0dce59131094d9c6 Mon Sep 17 00:00:00 2001 From: guido Date: Mon, 10 May 1999 12:57:34 +0000 Subject: [PATCH] Multigrid on meshes without hanging nodes git-svn-id: https://svn.dealii.org/trunk@1307 0785d39b-7218-0410-832d-ea1e28bc413d --- deal.II/deal.II/include/multigrid/mg_base.h | 8 +- deal.II/deal.II/include/multigrid/multigrid.h | 49 ++---- .../include/multigrid/multigrid.templates.h | 66 +++++--- deal.II/deal.II/include/numerics/multigrid.h | 49 ++---- .../include/numerics/multigrid.templates.h | 66 +++++--- deal.II/deal.II/source/multigrid/mg_base.cc | 3 +- .../deal.II/source/multigrid/mg_smoother.cc | 5 - deal.II/deal.II/source/multigrid/multigrid.cc | 8 +- .../deal.II/source/numerics/mg_smoother.cc | 5 - deal.II/deal.II/source/numerics/multigrid.cc | 8 +- deal.II/lac/include/lac/mgbase.h | 8 +- .../lac/include/lac/sparse_matrix.templates.h | 4 +- deal.II/lac/source/mgbase.cc | 3 +- tests/deal.II/gradients.cc | 2 +- tests/deal.II/helmholtz.h | 26 +++ tests/deal.II/helmholtz1.th | 70 ++++++++ tests/deal.II/helmholtz1mg.th | 66 ++++++++ tests/deal.II/mg.cc | 160 ++++++++++++++++++ tests/deal.II/mg2.cc | 23 --- tests/deal.II/second_derivatives.cc | 2 +- tests/lac/mg.cc | 2 +- 21 files changed, 472 insertions(+), 161 deletions(-) create mode 100644 tests/deal.II/helmholtz.h create mode 100644 tests/deal.II/helmholtz1.th create mode 100644 tests/deal.II/helmholtz1mg.th create mode 100644 tests/deal.II/mg.cc delete mode 100644 tests/deal.II/mg2.cc diff --git a/deal.II/deal.II/include/multigrid/mg_base.h b/deal.II/deal.II/include/multigrid/mg_base.h index 5f16cbf0e9..41c746197e 100644 --- a/deal.II/deal.II/include/multigrid/mg_base.h +++ b/deal.II/deal.II/include/multigrid/mg_base.h @@ -318,6 +318,12 @@ class MGBase * Prolongation and restriction object. */ SmartPointer transfer; + /** + * Exception. + */ + DeclException2(ExcSwitchedLevels, int, int, + << "minlevel and maxlevel switched, should be: " + << arg1 << "<=" << arg2); protected: @@ -352,7 +358,7 @@ class MGBase * Constructor, subject to change. */ MGBase(const MGTransferBase& transfer, - unsigned int maxlevel, unsigned int minlevel); + unsigned int minlevel, unsigned int maxlevel); /** * Virtual destructor. */ diff --git a/deal.II/deal.II/include/multigrid/multigrid.h b/deal.II/deal.II/include/multigrid/multigrid.h index 916e2e7c6b..d9587dadfe 100644 --- a/deal.II/deal.II/include/multigrid/multigrid.h +++ b/deal.II/deal.II/include/multigrid/multigrid.h @@ -28,7 +28,9 @@ class MG public MGBase { public: - MG(const MGDoFHandler&, const MGTransferBase& transfer); + MG(const MGDoFHandler&, + const MGMatrix >&, + const MGTransferBase& transfer); /** Transfer from dVector to * MGVector. @@ -40,8 +42,7 @@ class MG * MGVector. */ template - void copy_to_mg(vector >& dst, - const Vector& src) const; + void copy_to_mg(const Vector& src); /** * Transfer from multi-grid vector to @@ -54,23 +55,15 @@ class MG * #src# are zero. */ template - void copy_from_mg(Vector& dst, - const vector >& src) const; + void copy_from_mg(Vector& dst) const; /** - * Access to matrix structure. + * Residual on a level. */ - SparseMatrixStruct& get_sparsity_pattern(unsigned int level); - - /** - * Read-only access to matrix structure. - */ - const SparseMatrixStruct& get_sparsity_pattern(unsigned int level) const; - - /** - * Access to level matrices. - */ - SparseMatrix& get_matrix(unsigned int level); + virtual void level_residual(unsigned int, + Vector &, + const Vector &, + const Vector &); /** * Read-only access to level matrices. @@ -83,21 +76,13 @@ class MG */ SmartPointer > dofs; - /** - * Matrix structures for each level. - * These are generated by the constructor - * of #MG# and can be accessed - * later. - */ - vector matrix_structures; - /** * Matrices for each level. * The matrices are prepared by * the constructor of #MG# and can * be accessed for assembling. */ - vector > matrices; + SmartPointer > > matrices; }; @@ -175,16 +160,6 @@ class MGTransferPrebuilt : public MGTransferBase }; -template -inline -SparseMatrix& -MG::get_matrix(unsigned int level) -{ - Assert((level>=minlevel) && (level inline const SparseMatrix& @@ -192,7 +167,7 @@ MG::get_matrix(unsigned int level) const { Assert((level>=minlevel) && (level -/* template -MultiGrid:: -*/ - -template -MG::MG(const MGDoFHandler& dofs, const MGTransferBase& transfer) +MG::MG(const MGDoFHandler& dofs, + const MGMatrix >& matrices, + const MGTransferBase& transfer) : MGBase(transfer, 0, dofs.get_tria().n_levels()-1), dofs(&dofs), - matrix_structures(maxlevel-minlevel), - matrices(maxlevel-minlevel) + matrices(&matrices) { - for(unsigned l = minlevel; l < maxlevel; ++l) + for (unsigned int level = minlevel; level<=maxlevel ; ++level) { - dofs.make_sparsity_pattern(l, matrix_structures[l-minlevel]); - matrices[l-minlevel].reinit(matrix_structures[l-minlevel]); + s[level].reinit(matrices[level].m()); + d[level].reinit(matrices[level].m()); } } -template template +template +template void -MG::copy_to_mg(vector >& dst, - const Vector& src) const +MG::copy_to_mg(const Vector& src) { + const unsigned int fe_dofs = dofs->get_fe().total_dofs; + + vector index(fe_dofs); + vector mgindex(fe_dofs); + DoFHandler::active_cell_iterator dc = dofs->DoFHandler::begin_active(); + for (MGDoFHandler::active_cell_iterator c = dofs->begin_active() + ; c != dofs->end() ; ++c, ++dc) + { + dc->get_dof_indices(index); + c->get_mg_dof_indices(mgindex); + for (unsigned int i=0;i template void -MG::copy_from_mg(Vector& dst, - const vector >& src) const +MG::copy_from_mg(Vector& dst) const { -// active iterator + const unsigned int fe_dofs = dofs->get_fe().total_dofs; + + vector index(fe_dofs); + vector mgindex(fe_dofs); + + DoFHandler::active_cell_iterator dc = dofs->DoFHandler::begin_active(); + for (MGDoFHandler::active_cell_iterator c = dofs->begin_active() + ; c != dofs->end() ; ++c, ++dc) + { + dc->get_dof_indices(index); + c->get_mg_dof_indices(mgindex); + for (unsigned int i=0;i +void +MG::level_residual(unsigned int l, + Vector &r, + const Vector &u, + const Vector &b) +{ + (*matrices)[l].residual(r,u,b); +} /* template @@ -46,7 +76,5 @@ template MG:: template MG:: -template -MG:: */ diff --git a/deal.II/deal.II/include/numerics/multigrid.h b/deal.II/deal.II/include/numerics/multigrid.h index 916e2e7c6b..d9587dadfe 100644 --- a/deal.II/deal.II/include/numerics/multigrid.h +++ b/deal.II/deal.II/include/numerics/multigrid.h @@ -28,7 +28,9 @@ class MG public MGBase { public: - MG(const MGDoFHandler&, const MGTransferBase& transfer); + MG(const MGDoFHandler&, + const MGMatrix >&, + const MGTransferBase& transfer); /** Transfer from dVector to * MGVector. @@ -40,8 +42,7 @@ class MG * MGVector. */ template - void copy_to_mg(vector >& dst, - const Vector& src) const; + void copy_to_mg(const Vector& src); /** * Transfer from multi-grid vector to @@ -54,23 +55,15 @@ class MG * #src# are zero. */ template - void copy_from_mg(Vector& dst, - const vector >& src) const; + void copy_from_mg(Vector& dst) const; /** - * Access to matrix structure. + * Residual on a level. */ - SparseMatrixStruct& get_sparsity_pattern(unsigned int level); - - /** - * Read-only access to matrix structure. - */ - const SparseMatrixStruct& get_sparsity_pattern(unsigned int level) const; - - /** - * Access to level matrices. - */ - SparseMatrix& get_matrix(unsigned int level); + virtual void level_residual(unsigned int, + Vector &, + const Vector &, + const Vector &); /** * Read-only access to level matrices. @@ -83,21 +76,13 @@ class MG */ SmartPointer > dofs; - /** - * Matrix structures for each level. - * These are generated by the constructor - * of #MG# and can be accessed - * later. - */ - vector matrix_structures; - /** * Matrices for each level. * The matrices are prepared by * the constructor of #MG# and can * be accessed for assembling. */ - vector > matrices; + SmartPointer > > matrices; }; @@ -175,16 +160,6 @@ class MGTransferPrebuilt : public MGTransferBase }; -template -inline -SparseMatrix& -MG::get_matrix(unsigned int level) -{ - Assert((level>=minlevel) && (level inline const SparseMatrix& @@ -192,7 +167,7 @@ MG::get_matrix(unsigned int level) const { Assert((level>=minlevel) && (level -/* template -MultiGrid:: -*/ - -template -MG::MG(const MGDoFHandler& dofs, const MGTransferBase& transfer) +MG::MG(const MGDoFHandler& dofs, + const MGMatrix >& matrices, + const MGTransferBase& transfer) : MGBase(transfer, 0, dofs.get_tria().n_levels()-1), dofs(&dofs), - matrix_structures(maxlevel-minlevel), - matrices(maxlevel-minlevel) + matrices(&matrices) { - for(unsigned l = minlevel; l < maxlevel; ++l) + for (unsigned int level = minlevel; level<=maxlevel ; ++level) { - dofs.make_sparsity_pattern(l, matrix_structures[l-minlevel]); - matrices[l-minlevel].reinit(matrix_structures[l-minlevel]); + s[level].reinit(matrices[level].m()); + d[level].reinit(matrices[level].m()); } } -template template +template +template void -MG::copy_to_mg(vector >& dst, - const Vector& src) const +MG::copy_to_mg(const Vector& src) { + const unsigned int fe_dofs = dofs->get_fe().total_dofs; + + vector index(fe_dofs); + vector mgindex(fe_dofs); + DoFHandler::active_cell_iterator dc = dofs->DoFHandler::begin_active(); + for (MGDoFHandler::active_cell_iterator c = dofs->begin_active() + ; c != dofs->end() ; ++c, ++dc) + { + dc->get_dof_indices(index); + c->get_mg_dof_indices(mgindex); + for (unsigned int i=0;i template void -MG::copy_from_mg(Vector& dst, - const vector >& src) const +MG::copy_from_mg(Vector& dst) const { -// active iterator + const unsigned int fe_dofs = dofs->get_fe().total_dofs; + + vector index(fe_dofs); + vector mgindex(fe_dofs); + + DoFHandler::active_cell_iterator dc = dofs->DoFHandler::begin_active(); + for (MGDoFHandler::active_cell_iterator c = dofs->begin_active() + ; c != dofs->end() ; ++c, ++dc) + { + dc->get_dof_indices(index); + c->get_mg_dof_indices(mgindex); + for (unsigned int i=0;i +void +MG::level_residual(unsigned int l, + Vector &r, + const Vector &u, + const Vector &b) +{ + (*matrices)[l].residual(r,u,b); +} /* template @@ -46,7 +76,5 @@ template MG:: template MG:: -template -MG:: */ diff --git a/deal.II/deal.II/source/multigrid/mg_base.cc b/deal.II/deal.II/source/multigrid/mg_base.cc index b08359efb4..e6340762e8 100644 --- a/deal.II/deal.II/source/multigrid/mg_base.cc +++ b/deal.II/deal.II/source/multigrid/mg_base.cc @@ -26,13 +26,14 @@ MGBase::~MGBase () MGBase::MGBase(const MGTransferBase& transfer, - unsigned maxlevel, unsigned minlevel) + unsigned minlevel, unsigned maxlevel) : maxlevel(maxlevel), minlevel(minlevel), d(minlevel,maxlevel), s(minlevel,maxlevel), transfer(&transfer) { + Assert(minlevel <= maxlevel, ExcSwitchedLevels(minlevel, maxlevel)); } void diff --git a/deal.II/deal.II/source/multigrid/mg_smoother.cc b/deal.II/deal.II/source/multigrid/mg_smoother.cc index 71b80cb4b3..23bea7acf9 100644 --- a/deal.II/deal.II/source/multigrid/mg_smoother.cc +++ b/deal.II/deal.II/source/multigrid/mg_smoother.cc @@ -90,11 +90,6 @@ MGSmoother::MGSmoother (const MGDoFHandler &mg_dof, unsigned int steps) ////////////////////////////////////////////////////////////////////// -MGSmootherBase::~MGSmootherBase () -{}; - -////////////////////////////////////////////////////////////////////// - void MGSmoother::set_zero_interior_boundary (const unsigned int level, Vector &u) const diff --git a/deal.II/deal.II/source/multigrid/multigrid.cc b/deal.II/deal.II/source/multigrid/multigrid.cc index b5db616dde..429f929615 100644 --- a/deal.II/deal.II/source/multigrid/multigrid.cc +++ b/deal.II/deal.II/source/multigrid/multigrid.cc @@ -137,8 +137,8 @@ void MGTransferPrebuilt::prolongate (const unsigned int to_level, Vector &dst, const Vector &src) const { - Assert ((to_level >= 1) && (to_level= 1) && (to_level<=prolongation_matrices.size()), + ExcIndexRange (to_level, 1, prolongation_matrices.size()+1)); prolongation_matrices[to_level-1].vmult (dst, src); }; @@ -149,8 +149,8 @@ void MGTransferPrebuilt::restrict (const unsigned int from_level, Vector &dst, const Vector &src) const { - Assert ((from_level >= 1) && (from_level= 1) && (from_level<=prolongation_matrices.size()), + ExcIndexRange (from_level, 1, prolongation_matrices.size()+1)); prolongation_matrices[from_level-1].Tvmult (dst, src); }; diff --git a/deal.II/deal.II/source/numerics/mg_smoother.cc b/deal.II/deal.II/source/numerics/mg_smoother.cc index 71b80cb4b3..23bea7acf9 100644 --- a/deal.II/deal.II/source/numerics/mg_smoother.cc +++ b/deal.II/deal.II/source/numerics/mg_smoother.cc @@ -90,11 +90,6 @@ MGSmoother::MGSmoother (const MGDoFHandler &mg_dof, unsigned int steps) ////////////////////////////////////////////////////////////////////// -MGSmootherBase::~MGSmootherBase () -{}; - -////////////////////////////////////////////////////////////////////// - void MGSmoother::set_zero_interior_boundary (const unsigned int level, Vector &u) const diff --git a/deal.II/deal.II/source/numerics/multigrid.cc b/deal.II/deal.II/source/numerics/multigrid.cc index b5db616dde..429f929615 100644 --- a/deal.II/deal.II/source/numerics/multigrid.cc +++ b/deal.II/deal.II/source/numerics/multigrid.cc @@ -137,8 +137,8 @@ void MGTransferPrebuilt::prolongate (const unsigned int to_level, Vector &dst, const Vector &src) const { - Assert ((to_level >= 1) && (to_level= 1) && (to_level<=prolongation_matrices.size()), + ExcIndexRange (to_level, 1, prolongation_matrices.size()+1)); prolongation_matrices[to_level-1].vmult (dst, src); }; @@ -149,8 +149,8 @@ void MGTransferPrebuilt::restrict (const unsigned int from_level, Vector &dst, const Vector &src) const { - Assert ((from_level >= 1) && (from_level= 1) && (from_level<=prolongation_matrices.size()), + ExcIndexRange (from_level, 1, prolongation_matrices.size()+1)); prolongation_matrices[from_level-1].Tvmult (dst, src); }; diff --git a/deal.II/lac/include/lac/mgbase.h b/deal.II/lac/include/lac/mgbase.h index 3e9396f8f0..cc71bc98f4 100644 --- a/deal.II/lac/include/lac/mgbase.h +++ b/deal.II/lac/include/lac/mgbase.h @@ -318,6 +318,12 @@ class MGBase * Prolongation and restriction object. */ SmartPointer transfer; + /** + * Exception. + */ + DeclException2(ExcSwitchedLevels, int, int, + << "minlevel and maxlevel switched, should be: " + << arg1 << "<=" << arg2); protected: @@ -352,7 +358,7 @@ class MGBase * Constructor, subject to change. */ MGBase(const MGTransferBase& transfer, - unsigned int maxlevel, unsigned int minlevel); + unsigned int minlevel, unsigned int maxlevel); /** * Virtual destructor. */ diff --git a/deal.II/lac/include/lac/sparse_matrix.templates.h b/deal.II/lac/include/lac/sparse_matrix.templates.h index a2342846fd..5218fd84b0 100644 --- a/deal.II/lac/include/lac/sparse_matrix.templates.h +++ b/deal.II/lac/include/lac/sparse_matrix.templates.h @@ -489,7 +489,9 @@ SparseMatrix::SSOR (Vector& dst, const number om) const template -const SparseMatrixStruct & SparseMatrix::get_sparsity_pattern () const { +const SparseMatrixStruct & +SparseMatrix::get_sparsity_pattern () const +{ Assert (cols != 0, ExcMatrixNotInitialized()); return *cols; }; diff --git a/deal.II/lac/source/mgbase.cc b/deal.II/lac/source/mgbase.cc index b08359efb4..e6340762e8 100644 --- a/deal.II/lac/source/mgbase.cc +++ b/deal.II/lac/source/mgbase.cc @@ -26,13 +26,14 @@ MGBase::~MGBase () MGBase::MGBase(const MGTransferBase& transfer, - unsigned maxlevel, unsigned minlevel) + unsigned minlevel, unsigned maxlevel) : maxlevel(maxlevel), minlevel(minlevel), d(minlevel,maxlevel), s(minlevel,maxlevel), transfer(&transfer) { + Assert(minlevel <= maxlevel, ExcSwitchedLevels(minlevel, maxlevel)); } void diff --git a/tests/deal.II/gradients.cc b/tests/deal.II/gradients.cc index 3963071d73..1b000fc53e 100644 --- a/tests/deal.II/gradients.cc +++ b/tests/deal.II/gradients.cc @@ -36,7 +36,7 @@ int main () { dof.distribute_dofs(fe); QTrapez<2> q; - FEValues<2> fevalues(fe,q,update_gradients,tria.get_boundary()); + FEValues<2> fevalues(fe,q,update_gradients); fevalues.reinit (dof.begin_active()); diff --git a/tests/deal.II/helmholtz.h b/tests/deal.II/helmholtz.h new file mode 100644 index 0000000000..f51156022a --- /dev/null +++ b/tests/deal.II/helmholtz.h @@ -0,0 +1,26 @@ +// $Id$ + +#include + +/** + * Generator for system matrix and right hand side. + */ +class Helmholtz +{ + public: + template + void build_all(MATRIX& A, + VECTOR& f, + const DoFHandler& dof, + const Quadrature& quadrature, + const Function& rhs); + + template + void build_mgmatrix(MGMatrix& A, + const MGDoFHandler& dof, + const Quadrature& quadrature); +}; + +#include "helmholtz1.th" +#include "helmholtz1mg.th" + diff --git a/tests/deal.II/helmholtz1.th b/tests/deal.II/helmholtz1.th new file mode 100644 index 0000000000..5abed88b92 --- /dev/null +++ b/tests/deal.II/helmholtz1.th @@ -0,0 +1,70 @@ +// $Id$ + +template +void +Helmholtz::build_all(MATRIX& A, + VECTOR& f, + const DoFHandler& dof, + const Quadrature& quadrature, + const Function& rhs) +{ + const unsigned int fe_dofs = dof.get_fe().total_dofs; + FEValues fe(dof.get_fe(), quadrature, + UpdateFlags(update_gradients | update_JxW_values | + update_q_points)); + + Vector elvec(fe_dofs); + FullMatrix elmat(fe_dofs); + vector indices(fe_dofs); + +//////////////////////////////////////////////////// +// Loop for computing the element matrices +//////////////////////////////////////////////////// + + for (DoFHandler::active_cell_iterator c = dof.begin_active() + ; c != dof.end() ; ++c) + { + elmat.clear(); + elvec.clear(); + fe.reinit(c); + + // Quadrature loop + + for (unsigned k=0;k loc = fe.quadrature_point(k); + + // Test function loop + for (unsigned i=0;i dv = fe.shape_grad(i,k); + double v = fe.shape_value(i,k); + + elvec(i) += fe.JxW(k) * + rhs(loc) * v; + + //Trial function loop + for (unsigned j=0;j du = fe.shape_grad(j,k); + double u = fe.shape_value(j,k); + + elmat(i,j) += fe.JxW(k) * + (1.*u*v + du*dv) + ; + } + } + } + + c->get_dof_indices(indices); + for (unsigned i=0;i +void +Helmholtz::build_mgmatrix(MGMatrix& A, + const MGDoFHandler& dof, + const Quadrature& quadrature) +{ + const unsigned int fe_dofs = dof.get_fe().total_dofs; + FEValues fe(dof.get_fe(), quadrature, + UpdateFlags(update_gradients | update_JxW_values | + update_q_points)); + + Vector elvec(fe_dofs); + FullMatrix elmat(fe_dofs); + vector indices(fe_dofs); + +//////////////////////////////////////////////////// +// Loop for computing the element matrices +//////////////////////////////////////////////////// + + DoFHandler::cell_iterator dc = dof.DoFHandler::begin(); + for (MGDoFHandler::cell_iterator c = dof.begin() + ; c != dof.end() ; ++c, ++dc) + { + elmat.clear(); + elvec.clear(); + fe.reinit(dc); + + // Quadrature loop + + for (unsigned k=0;k loc = fe.quadrature_point(k); + + // Test function loop + for (unsigned i=0;i dv = fe.shape_grad(i,k); + double v = fe.shape_value(i,k); + + //Trial function loop + for (unsigned j=0;j du = fe.shape_grad(j,k); + double u = fe.shape_value(j,k); + + elmat(i,j) += fe.JxW(k) * + (1.*u*v + du*dv) + ; + } + } + } + + const unsigned int level = c->level(); + + c->get_mg_dof_indices(indices); + for (unsigned i=0;i +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include + +#include "helmholtz.h" + +template +class RHSFunction + : + public Function +{ + public: + virtual double operator() (const Point&) const; +// virtual Tensor<1,dim> gradient (const Point &p) const; +}; + +class MGSmootherLAC + : + public MGSmootherBase +{ + private: + SmartPointer > >matrices; + public: + MGSmootherLAC(MGMatrix >&); + + virtual void smooth (const unsigned int level, + Vector &u, + const Vector &rhs) const; + +}; + +main(int argc, const char** argv) +{ + unsigned int refine; + if (argc<=1) + { + cerr << "Usage:" << argv[0] << " " << endl; + refine = 0; + } + else + refine = atoi(argv[1]); + + + Helmholtz equation; + RHSFunction<2> rhs; + QGauss3<2> quadrature; + + FELinear<2> fe; + Triangulation<2> tr; + MGDoFHandler<2> mgdof(&tr); + DoFHandler<2>& dof(mgdof); + + GridGenerator::hyper_cube(tr,-1.,1.); + tr.refine_global(refine); + dof.distribute_dofs(fe); + + const unsigned int size = dof.n_dofs(); + + SparseMatrixStruct structure(size, dof.max_couplings_between_dofs()); + dof.make_sparsity_pattern(structure); + structure.compress(); + + SparseMatrix A(structure); + Vector f(size); + + equation.build_all(A,f,dof, quadrature, rhs); + + Vector u; + u.reinit(f); + PrimitiveVectorMemory > mem; + SolverControl control(100, 1.e-15); + SolverCG, Vector >solver(control, mem); + + PreconditionIdentity > precondition; + solver.solve(A,u,f,precondition); + + u = 0.; + vector mgstruct(tr.n_levels()); + MGMatrix > mgA(0,tr.n_levels()-1); + for (unsigned int i=0;i > cgmem; + SolverCG, Vector > cgsolver(cgcontrol, cgmem); + PreconditionIdentity > cgprec; + MGCoarseGridLACIteration, Vector >, + SparseMatrix, PreconditionIdentity > > + coarse(cgsolver, mgA[0], cgprec); + + MGSmootherLAC smoother(mgA); + MGTransferPrebuilt transfer; + transfer.build_matrices(mgdof); + + deallog << "Levels: " << tr.n_levels() << endl; + + MG<2> multigrid(mgdof, mgA, transfer); + PreconditionMG, Vector > + mgprecondition(multigrid, smoother, smoother, coarse); + + solver.solve(A, u, f, mgprecondition); + + cerr << "OK" << endl; +} + +template +double +RHSFunction::operator() (const Point&) const +{ + return 1.; +} + +MGSmootherLAC::MGSmootherLAC(MGMatrix >& matrix) + : + matrices(&matrix) +{} + +void +MGSmootherLAC::smooth (const unsigned int level, + Vector &u, + const Vector &rhs) const +{ + SolverControl control(1,1.e-300,false,false); + PrimitiveVectorMemory > mem; + SolverRichardson , Vector > rich(control, mem); + PreconditionRelaxation , Vector > + prec((*matrices)[level], &SparseMatrix ::template precondition_SSOR, 1.); + + rich.solve((*matrices)[level], u, rhs, prec); +} diff --git a/tests/deal.II/mg2.cc b/tests/deal.II/mg2.cc deleted file mode 100644 index 44dd281b4a..0000000000 --- a/tests/deal.II/mg2.cc +++ /dev/null @@ -1,23 +0,0 @@ -// $Id$ - -// deal_II_libraries.g=-ldeal_II_2d.g -// deal_II_libraries=-ldeal_II_2d - -#include -#include -#include -#include -#include -#include - -main() -{ - Triangulation<2> tr; - MGDoFHandler<2> dof(&tr); - FELinear<2> fe; - - GridGenerator::hyper_cube(tr,-1.,1.); - tr.refine_global(3); - dof.distribute_dofs(fe); -} - diff --git a/tests/deal.II/second_derivatives.cc b/tests/deal.II/second_derivatives.cc index 2ad62be46c..0c4d5064cb 100644 --- a/tests/deal.II/second_derivatives.cc +++ b/tests/deal.II/second_derivatives.cc @@ -28,7 +28,7 @@ int main () { StraightBoundary<2> b; QTrapez<2> q; - FEValues<2> fevalues(fe,q,update_second_derivatives,b); + FEValues<2> fevalues(fe,q,update_second_derivatives); Vector val(4); diff --git a/tests/lac/mg.cc b/tests/lac/mg.cc index a8df6709d6..c667839d0b 100644 --- a/tests/lac/mg.cc +++ b/tests/lac/mg.cc @@ -160,7 +160,7 @@ main() FDMG::FDMG(unsigned int maxlevel, MGMatrix >& matrices, FDMGTransfer& transfer) : - MGBase(transfer, maxlevel, 0), + MGBase(transfer, 0, maxlevel), matrices(&matrices) { for (unsigned int level = minlevel; level<=maxlevel ; ++level) -- 2.39.5