Add a reasonable (unfinished) petsc test.

git-svn-id: https://svn.dealii.org/trunk@31651 0785d39b-7218-0410-832d-ea1e28bc413d

diff --git a/deal.II/tests/quick_tests/CMakeLists.txt b/deal.II/tests/quick_tests/CMakeLists.txt
index 747e562420db6dd9a9b37490ef19001dc6bb1e56..727911cb15a917e748abd5569c70edeaa073bcf7 100644 (file)
--- a/deal.II/tests/quick_tests/CMakeLists.txt
+++ b/deal.II/tests/quick_tests/CMakeLists.txt
@@ -85,8 +85,13 @@ IF (DEAL_II_WITH_P4EST)
   make_quicktest("p4est" ${_mybuild} 10)
 ENDIF()
 
+# Test petsc
+IF (DEAL_II_WITH_PETSC)
+  make_quicktest("step-petsc" ${_mybuild} "")
+ENDIF()
+
 # Test slepc
-IF (DEAL_II_WITH_SLEPC)
+IF (DEAL_II_WITH_PETSC AND DEAL_II_WITH_SLEPC)
   make_quicktest("step-slepc" ${_mybuild} "")
 ENDIF()
 

diff --git a/deal.II/tests/quick_tests/step-petsc.cc b/deal.II/tests/quick_tests/step-petsc.cc
new file mode 100644 (file)
index 0000000..9290f55
--- /dev/null
+++ b/deal.II/tests/quick_tests/step-petsc.cc
@@ -0,0 +1,220 @@
+/* ---------------------------------------------------------------------\r
+ * $Id: step-petsc.cc $\r
+ *\r
+ * Copyright (C) 2013 by the deal.II authors\r
+ *\r
+ * This file is part of the deal.II library.\r
+ *\r
+ * The deal.II library is free software; you can use it, redistribute\r
+ * it, and/or modify it under the terms of the GNU Lesser General\r
+ * Public License as published by the Free Software Foundation; either\r
+ * version 2.1 of the License, or (at your option) any later version.\r
+ * The full text of the license can be found in the file LICENSE at\r
+ * the top level of the deal.II distribution.\r
+ *\r
+ * ---------------------------------------------------------------------\r
+ */\r
+\r
+#include <deal.II/base/logstream.h>\r
+#include <deal.II/base/table_handler.h>\r
+#include <deal.II/base/quadrature_lib.h>\r
+#include <deal.II/base/function.h>\r
+#include <deal.II/base/utilities.h>\r
+#include <deal.II/grid/tria.h>\r
+#include <deal.II/grid/grid_generator.h>\r
+#include <deal.II/grid/tria_accessor.h>\r
+#include <deal.II/grid/tria_iterator.h>\r
+#include <deal.II/dofs/dof_handler.h>\r
+#include <deal.II/dofs/dof_accessor.h>\r
+#include <deal.II/dofs/dof_tools.h>\r
+#include <deal.II/fe/fe_q.h>\r
+#include <deal.II/fe/fe_values.h>\r
+#include <deal.II/numerics/vector_tools.h>\r
+#include <deal.II/numerics/matrix_tools.h>\r
+#include <deal.II/numerics/data_out.h>\r
+#include <deal.II/lac/full_matrix.h>\r
+#include <deal.II/lac/petsc_sparse_matrix.h>\r
+#include <deal.II/lac/petsc_vector.h>\r
+#include <deal.II/lac/petsc_solver.h>\r
+#include <deal.II/lac/petsc_precondition.h>\r
+\r
+#include <fstream>\r
+#include <iostream>\r
+\r
+using namespace dealii;\r
+\r
+// Test that deal.II is working with PETSc by solving the Laplace's\r
+// problem in 2d.\r
+class LaplaceProblem\r
+{\r
+public:\r
+  LaplaceProblem ();\r
+  void run ();\r
+  \r
+private:\r
+  void setup_system ();\r
+  void assemble_system ();\r
+  void solve ();\r
+  \r
+  Triangulation<2> triangulation;\r
+  FE_Q<2>          fe;\r
+  DoFHandler<2>    dof_handler;\r
+  \r
+  PETScWrappers::SparseMatrix A;\r
+  PETScWrappers::Vector       b, x;  \r
+  ConstraintMatrix            constraints;\r
+\r
+  TableHandler output_table;\r
+};\r
+\r
+LaplaceProblem::LaplaceProblem ()\r
+  :\r
+  fe (1),\r
+  dof_handler (triangulation)\r
+{}\r
+\r
+void LaplaceProblem::setup_system ()\r
+{\r
+  dof_handler.distribute_dofs (fe);\r
+\r
+  constraints.clear ();\r
+  DoFTools::make_zero_boundary_constraints (dof_handler, constraints);\r
+  constraints.close ();\r
+  \r
+  A.reinit (dof_handler.n_dofs(), dof_handler.n_dofs(),\r
+           dof_handler.max_couplings_between_dofs());\r
+  b.reinit (dof_handler.n_dofs());\r
+  x.reinit (dof_handler.n_dofs());\r
+\r
+  // some output\r
+  output_table.add_value ("cells", triangulation.n_active_cells());\r
+  output_table.add_value ("dofs",  dof_handler.n_dofs());\r
+}\r
+\r
+void LaplaceProblem::assemble_system ()\r
+{\r
+  QGauss<2> quadrature_formula(2);\r
+  \r
+  FEValues<2> fe_values (fe, quadrature_formula,\r
+                        update_values            | \r
+                        update_gradients         |\r
+                        update_quadrature_points | \r
+                        update_JxW_values);\r
+  \r
+  const unsigned int dofs_per_cell = fe.dofs_per_cell;\r
+  const unsigned int n_q_points    = quadrature_formula.size();\r
+  \r
+  FullMatrix<double> cell_A (dofs_per_cell, dofs_per_cell);\r
+  Vector<double>     cell_b (dofs_per_cell);\r
+  \r
+  std::vector<types::global_dof_index> local_dof_indices (dofs_per_cell);\r
+  \r
+  typename DoFHandler<2>::active_cell_iterator\r
+    cell = dof_handler.begin_active (),\r
+    endc = dof_handler.end ();\r
+  \r
+  for (; cell!=endc; ++cell)\r
+    {\r
+      fe_values.reinit (cell);\r
+      cell_A = 0;\r
+      cell_b = 0;\r
+      \r
+      for (unsigned int q_point=0; q_point<n_q_points; ++q_point)\r
+       for (unsigned int i=0; i<dofs_per_cell; ++i)\r
+         {\r
+           for (unsigned int j=0; j<dofs_per_cell; ++j)\r
+             {\r
+               cell_A (i, j)\r
+                 += \r
+                 fe_values.shape_grad (i, q_point) *\r
+                 fe_values.shape_grad (j, q_point)\r
+                 * \r
+                 fe_values.JxW (q_point);\r
+             }\r
+\r
+           cell_b (i)\r
+             += \r
+             fe_values.shape_value (i, q_point)\r
+             * \r
+             fe_values.JxW (q_point);\r
+         }\r
+      \r
+      cell->get_dof_indices (local_dof_indices);\r
+      \r
+      constraints.distribute_local_to_global (cell_A, local_dof_indices, A);\r
+      constraints.distribute_local_to_global (cell_b, local_dof_indices, b);\r
+    }\r
+  \r
+  A.compress (VectorOperation::add);\r
+  b.compress (VectorOperation::add);\r
+}\r
+\r
+void LaplaceProblem::solve ()\r
+{\r
+  SolverControl solver_control (1e03, 1e-03);\r
+  PETScWrappers::SolverCG cg_solver (solver_control);\r
+  PETScWrappers::PreconditionBlockJacobi preconditioner (A);\r
+  cg_solver.solve (A, x, b, preconditioner);\r
+\r
+  // some output\r
+  // ?\r
+}\r
+\r
+void LaplaceProblem::run ()\r
+{\r
+  GridGenerator::hyper_cube (triangulation, -1, 1);\r
+\r
+  for (unsigned int c=0; c<5; ++c)\r
+    {\r
+      triangulation.refine_global (1);\r
+      setup_system ();\r
+      assemble_system ();\r
+      solve ();\r
+    }\r
+\r
+  // finialise output\r
+  output_table.write_text (std::cout);\r
+  deallog << std::endl;\r
+}\r
+\r
+\r
+int main (int argc, char **argv)\r
+{\r
+  try\r
+    {\r
+      dealii::Utilities::MPI::MPI_InitFinalize mpi_initialization (argc, argv, 1);\r
+      {\r
+        deallog.depth_console (0);\r
+        LaplaceProblem problem;\r
+        problem.run ();\r
+       deallog << "OK" << std::endl;\r
+      }\r
+    }\r
+\r
+  catch (std::exception &exc)\r
+    {\r
+      std::cerr << std::endl << std::endl\r
+                << "----------------------------------------------------"\r
+                << std::endl;\r
+      std::cerr << "Exception on processing: " << std::endl\r
+                << exc.what() << std::endl\r
+                << "Aborting!" << std::endl\r
+                << "----------------------------------------------------"\r
+                << std::endl;\r
+\r
+      return 1;\r
+    }\r
+  catch (...)\r
+    {\r
+      std::cerr << std::endl << std::endl\r
+                << "----------------------------------------------------"\r
+                << std::endl;\r
+      std::cerr << "Unknown exception!" << std::endl\r
+                << "Aborting!" << std::endl\r
+                << "----------------------------------------------------"\r
+                << std::endl;\r
+      return 1;\r
+    }\r
+\r
+  return 0;\r
+}\r