From d3c0fecbbc7ccde552fc80d831f3edfe4ba5a38a Mon Sep 17 00:00:00 2001 From: Guido Kanschat Date: Wed, 21 Dec 2005 13:30:36 +0000 Subject: [PATCH] replace petsc preconditioner by SparseILU git-svn-id: https://svn.dealii.org/trunk@11909 0785d39b-7218-0410-832d-ea1e28bc413d --- deal.II/examples/step-21/Makefile | 17 ---- deal.II/examples/step-21/step-21.cc | 141 ++++++++++++++-------------- 2 files changed, 68 insertions(+), 90 deletions(-) diff --git a/deal.II/examples/step-21/Makefile b/deal.II/examples/step-21/Makefile index f807c0589d..ad64734a77 100644 --- a/deal.II/examples/step-21/Makefile +++ b/deal.II/examples/step-21/Makefile @@ -46,21 +46,6 @@ clean-up-files = *gmv *gnuplot *gpl *eps *pov include $D/common/Make.global_options -################################################################ -# This example program will only work if PETSc is installed. If this -# is not the case, then simply redefine the main targets to to nothing -ifneq ($(USE_CONTRIB_PETSC),yes) -default run clean: - @echo - @echo "===========================================================" - @echo "= This program cannot be compiled without PETSc. Make =" - @echo "= sure you have PETSc installed and detected during =" - @echo "= configuration of deal.II =" - @echo "===========================================================" - @echo - -else -# ################################################################ # Since the whole project consists of only one file, we need not # consider difficult dependencies. We only have to declare the @@ -167,5 +152,3 @@ Makefile.dep: $(target).cc Makefile \ # them: include Makefile.dep - -endif # CONTRIB_USE_PETSC diff --git a/deal.II/examples/step-21/step-21.cc b/deal.II/examples/step-21/step-21.cc index 3223c43b58..3cc8396412 100644 --- a/deal.II/examples/step-21/step-21.cc +++ b/deal.II/examples/step-21/step-21.cc @@ -21,6 +21,9 @@ #include #include #include +#include +#include +#include #include #include #include @@ -42,21 +45,12 @@ // include files. #include #include +#include #include #include #include - // The PETSc vector and matrix - // classes are used to have access to - // the PETSc preconditioners. -#include -#include - // Then we also need interfaces for solvers - // and preconditioners that PETSc provides: -#include -#include - // We are going to use gradients as // refinement indicator. #include @@ -599,7 +593,7 @@ class DGMethod void setup_system (); void assemble_system1 (); void assemble_system2 (); - void solve (PETScWrappers::Vector &solution); + void solve (Vector &solution); void refine_grid (); void output_results (const unsigned int cycle) const; @@ -628,7 +622,8 @@ class DGMethod // replaced by a hpDoFHandler. hpDoFHandler dof_handler; - PETScWrappers::SparseMatrix system_matrix; + SparsityPattern sparsity; + SparseMatrix system_matrix; // We define the quadrature // formulae for the cell and the @@ -651,9 +646,9 @@ class DGMethod // different assembling routines // ``assemble_system1'' and // ``assemble_system2''; - PETScWrappers::Vector solution1; - PETScWrappers::Vector solution2; - PETScWrappers::Vector right_hand_side; + Vector solution1; + Vector solution2; + Vector right_hand_side; // Finally this class includes an // object of the @@ -732,7 +727,22 @@ void DGMethod::setup_system () // First we need to distribute the // DoFs. dof_handler.distribute_dofs (fe_collection); - + // In order to get a good + // preconditioner, the degrees of + // freedom should be ordered in + // downstream direction. First, we + // initalize a vector fairly close + // to the real vector field; since + // this is for preconditioning + // only, a rough approximation is + // sufficient. + Point sorting_direction; + for (unsigned int d=0;d::setup_system () DoFTools::make_sparsity_pattern (dof_handler, compressed_pattern); DoFTools::make_flux_sparsity_pattern (dof_handler, compressed_pattern); - system_matrix.reinit (compressed_pattern); + sparsity.copy_from(compressed_pattern); + system_matrix.reinit (sparsity); solution1.reinit (dof_handler.n_dofs()); solution2.reinit (dof_handler.n_dofs()); @@ -1259,9 +1270,6 @@ void DGMethod::assemble_system1 () for (unsigned int i=0; i::assemble_system2 () for (unsigned int i=0; i -void DGMethod::solve (PETScWrappers::Vector &solution) +void DGMethod::solve (Vector &solution) { - SolverControl solver_control (solution.size(), - 1e-13, false, false); - PETScWrappers::SolverBiCG bicg (solver_control); - PETScWrappers::PreconditionILU preconditioner(system_matrix); - - // Then solve the system: - bicg.solve (system_matrix, solution, right_hand_side, + SolverControl solver_control (10000, 1e-12, false, true); + SolverGMRES > solver (solver_control); + // Initialize the ILU + // preconditioner. We decide for + // two additional off diagonals in + // order to enhance its + // performance. + SparseILU::AdditionalData data(0., 2); + SparseILU preconditioner; + preconditioner.initialize (system_matrix, data); + // Then solve the system: + solver.solve (system_matrix, solution, right_hand_side, preconditioner); - - solution.compress (); - - std::cout << "Iterations : " << solver_control.last_step() << std::endl; } @@ -1655,7 +1656,7 @@ void DGMethod::output_results (const unsigned int cycle) const Assert (cycle < 10, ExcInternalError()); filename += ".eps"; - std::cout << "Writing grid to <" << filename << ">..." << std::endl; + deallog << "Writing grid to <" << filename << ">..." << std::endl; std::ofstream eps_output (filename.c_str()); GridOut grid_out; @@ -1669,7 +1670,7 @@ void DGMethod::output_results (const unsigned int cycle) const // filename += ".gnuplot"; filename += ".vtk"; - std::cout << "Writing solution to <" << filename << ">..." + deallog << "Writing solution to <" << filename << ">..." << std::endl << std::endl; std::ofstream gnuplot_output (filename.c_str()); @@ -1702,7 +1703,7 @@ void DGMethod::run () { for (unsigned int cycle=0; cycle<5; ++cycle) { - std::cout << "Cycle " << cycle << ':' << std::endl; + deallog << "Cycle " << cycle << ':' << std::endl; if (cycle == 0) { @@ -1714,13 +1715,13 @@ void DGMethod::run () refine_grid (); - std::cout << " Number of active cells: " + deallog << " Number of active cells: " << triangulation.n_active_cells() << std::endl; setup_system (); - std::cout << " Number of degrees of freedom: " + deallog << " Number of degrees of freedom: " << dof_handler.n_dofs() << std::endl; @@ -1734,7 +1735,7 @@ void DGMethod::run () // current time without // disturbing the time // measurement. - std::cout << "Time of assemble_system1: " + deallog << "Time of assemble_system1: " << assemble_timer() << std::endl; solve (solution1); @@ -1754,7 +1755,7 @@ void DGMethod::run () // call the second assembling routine assemble_system2 (); // and access the current time. - std::cout << "Time of assemble_system2: " + deallog << "Time of assemble_system2: " << assemble_timer() << std::endl; solve (solution2); @@ -1768,10 +1769,10 @@ void DGMethod::run () solution1-=solution2; const double difference=solution1.linfty_norm(); - if (difference>1e-13) - std::cout << "solution1 and solution2 differ!!" << std::endl; + if (difference>1e-12) + deallog << "solution1 and solution2 differ!!" << std::endl; else - std::cout << "solution1 and solution2 coincide." << std::endl; + deallog << "solution1 and solution2 coincide." << std::endl; // Finally we perform the // output. @@ -1782,18 +1783,12 @@ void DGMethod::run () // The following ``main'' function is // similar to previous examples and // need not to be commented on. -int main (int argc, char **argv) +int main () { try { - PetscInitialize(&argc,&argv,0,0); - - { DGMethod<2> dgmethod; dgmethod.run (); - } - - PetscFinalize (); } catch (std::exception &exc) { -- 2.39.5