From: heister Date: Wed, 6 Nov 2013 14:47:38 +0000 (+0000) Subject: add slepc quicktest, thanks Toby X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=90591c71341ecc00d684731d2de27fe81b013850;p=dealii-svn.git add slepc quicktest, thanks Toby git-svn-id: https://svn.dealii.org/trunk@31563 0785d39b-7218-0410-832d-ea1e28bc413d --- diff --git a/deal.II/tests/quick_tests/CMakeLists.txt b/deal.II/tests/quick_tests/CMakeLists.txt index 57ef944d87..2a4031cfab 100644 --- a/deal.II/tests/quick_tests/CMakeLists.txt +++ b/deal.II/tests/quick_tests/CMakeLists.txt @@ -85,6 +85,10 @@ IF (DEAL_II_WITH_P4EST) make_quicktest("p4est" ${_mybuild} 10) ENDIF() +# Test slepc +IF (DEAL_II_WITH_SLEPC) + make_quicktest("step-slepc" ${_mybuild} "") +ENDIF() # A custom test target: ADD_CUSTOM_TARGET(test diff --git a/deal.II/tests/quick_tests/step-slepc.cc b/deal.II/tests/quick_tests/step-slepc.cc new file mode 100644 index 0000000000..51bc2a664f --- /dev/null +++ b/deal.II/tests/quick_tests/step-slepc.cc @@ -0,0 +1,205 @@ +/* --------------------------------------------------------------------- + * $Id$ + * + * Copyright (C) 2013 by the deal.II authors + * + * This file is part of the deal.II library. + * + * The deal.II library is free software; you can use it, redistribute + * it, and/or modify it under the terms of the GNU Lesser General + * Public License as published by the Free Software Foundation; either + * version 2.1 of the License, or (at your option) any later version. + * The full text of the license can be found in the file LICENSE at + * the top level of the deal.II distribution. + * + * --------------------------------------------------------------------- + */ + +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include + +#include +#include + +using namespace dealii; + +// Test that deal.II is working with SLEPc by solving the Laplace's +// eigenspectrum problem in 2d. +class LaplaceEigenspectrumProblem +{ +public: + LaplaceEigenspectrumProblem (); + void run (); + +private: + void setup_system (); + void assemble_system (); + void solve (); + + Triangulation<2> triangulation; + FE_Q<2> fe; + DoFHandler<2> dof_handler; + + PETScWrappers::SparseMatrix A, B; + std::vector x; + std::vector lambda; + + ConstraintMatrix constraints; +}; + +LaplaceEigenspectrumProblem::LaplaceEigenspectrumProblem () + : + fe (1), + dof_handler (triangulation) +{} + +void LaplaceEigenspectrumProblem::setup_system () +{ + GridGenerator::hyper_cube (triangulation, -1, 1); + triangulation.refine_global (3); + dof_handler.distribute_dofs (fe); + + DoFTools::make_zero_boundary_constraints (dof_handler, constraints); + constraints.close (); + + A.reinit (dof_handler.n_dofs(), dof_handler.n_dofs(), + dof_handler.max_couplings_between_dofs()); + B.reinit (dof_handler.n_dofs(), dof_handler.n_dofs(), + dof_handler.max_couplings_between_dofs()); + + x.resize (1); + x[0].reinit (dof_handler.n_dofs ()); + lambda.resize (1); + lambda[0] = 0.; +} + +void LaplaceEigenspectrumProblem::assemble_system () +{ + QGauss<2> quadrature_formula(2); + + FEValues<2> fe_values (fe, quadrature_formula, + update_values | + update_gradients | + update_quadrature_points | + update_JxW_values); + + const unsigned int dofs_per_cell = fe.dofs_per_cell; + const unsigned int n_q_points = quadrature_formula.size(); + + FullMatrix cell_A (dofs_per_cell, dofs_per_cell); + FullMatrix cell_B (dofs_per_cell, dofs_per_cell); + + std::vector local_dof_indices (dofs_per_cell); + + typename DoFHandler<2>::active_cell_iterator + cell = dof_handler.begin_active (), + endc = dof_handler.end (); + + for (; cell!=endc; ++cell) + { + fe_values.reinit (cell); + cell_A = 0; + cell_B = 0; + + for (unsigned int q_point=0; q_pointget_dof_indices (local_dof_indices); + + constraints.distribute_local_to_global (cell_A, local_dof_indices, A); + constraints.distribute_local_to_global (cell_B, local_dof_indices, B); + } + + A.compress (VectorOperation::add); + B.compress (VectorOperation::add); +} + +void LaplaceEigenspectrumProblem::solve () +{ + SolverControl solver_control (1., 1e-03); + SLEPcWrappers::SolverLAPACK eigensolver (solver_control); + eigensolver.solve (A, B, lambda, x, x.size()); +} + +void LaplaceEigenspectrumProblem::run () +{ + setup_system (); + assemble_system (); + solve (); +} + + +int main (int argc, char **argv) +{ + try + { + dealii::Utilities::MPI::MPI_InitFinalize mpi_initialization (argc, argv, 1); + { + deallog.depth_console (0); + LaplaceEigenspectrumProblem problem; + problem.run (); + deallog << "OK" << std::endl; + } + } + + catch (std::exception &exc) + { + std::cerr << std::endl << std::endl + << "----------------------------------------------------" + << std::endl; + std::cerr << "Exception on processing: " << std::endl + << exc.what() << std::endl + << "Aborting!" << std::endl + << "----------------------------------------------------" + << std::endl; + + return 1; + } + catch (...) + { + std::cerr << std::endl << std::endl + << "----------------------------------------------------" + << std::endl; + std::cerr << "Unknown exception!" << std::endl + << "Aborting!" << std::endl + << "----------------------------------------------------" + << std::endl; + return 1; + } + + return 0; +}