From 86325b9ac02095cb96718651885ab02f7e9a4150 Mon Sep 17 00:00:00 2001 From: Denis Davydov Date: Thu, 13 Oct 2016 11:29:19 +0200 Subject: [PATCH] add a quick test for Arpack --- tests/quick_tests/CMakeLists.txt | 5 ++ tests/quick_tests/arpack.cc | 106 +++++++++++++++++++++++++++++++ 2 files changed, 111 insertions(+) create mode 100644 tests/quick_tests/arpack.cc diff --git a/tests/quick_tests/CMakeLists.txt b/tests/quick_tests/CMakeLists.txt index 97bae28fca..84fcbfb053 100644 --- a/tests/quick_tests/CMakeLists.txt +++ b/tests/quick_tests/CMakeLists.txt @@ -136,6 +136,11 @@ IF (DEAL_II_WITH_GSL) make_quicktest("gsl" ${_mybuild} "") ENDIF() +# Test Arpack +IF (DEAL_II_WITH_ARPACK AND DEAL_II_WITH_UMFPACK) + make_quicktest("arpack" ${_mybuild} "") +ENDIF() + # A custom test target: diff --git a/tests/quick_tests/arpack.cc b/tests/quick_tests/arpack.cc new file mode 100644 index 0000000000..087d061d27 --- /dev/null +++ b/tests/quick_tests/arpack.cc @@ -0,0 +1,106 @@ +// --------------------------------------------------------------------- +// +// Copyright (C) 2016 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. +// +// --------------------------------------------------------------------- + + +// test Arpack by calculating eigenvalues of Laplace matrix. + +#include "../tests.h" +#include +#include +#include + +#include +#include + +#include +#include + +#include + +#include +#include +#include +#include +#include +#include +#include + +#include +#include + +#include +#include + + + +template +void test () +{ + const unsigned int n_eigenvalues=3; + + Triangulation tria; + GridGenerator::hyper_cube (tria,0,1); + tria.refine_global (3); + + FE_Q fe (1); + DoFHandler dof_handler (tria); + dof_handler.distribute_dofs (fe); + + deallog << "Number of dofs = " << dof_handler.n_dofs() << std::endl; + + SparsityPattern sparsity_pattern; + sparsity_pattern.reinit (dof_handler.n_dofs(), + dof_handler.n_dofs(), + dof_handler.max_couplings_between_dofs()); + DoFTools::make_sparsity_pattern (dof_handler, sparsity_pattern); + sparsity_pattern.compress (); + + SparseMatrix A,B; + A.reinit (sparsity_pattern); + B.reinit (sparsity_pattern); + + std::vector > eigenvectors(n_eigenvalues, + Vector(dof_handler.n_dofs())); + std::vector > eigenvalues(n_eigenvalues); + + QGauss qr (2); + MatrixTools::create_laplace_matrix (dof_handler, qr, A); + MatrixTools::create_mass_matrix (dof_handler, qr, B); + + SolverControl solver_control (dof_handler.n_dofs(), 1e-10); + SparseDirectUMFPACK inverse; + inverse.initialize (A); + const unsigned int num_arnoldi_vectors = 2*eigenvalues.size() + 2; + ArpackSolver::AdditionalData additional_data(num_arnoldi_vectors, + ArpackSolver::largest_magnitude, + true); + ArpackSolver eigensolver (solver_control, additional_data); + eigensolver.solve (A, + B, + inverse, + eigenvalues, + eigenvectors, + eigenvalues.size()); +} + + +int main () +{ + std::ofstream logfile("output"); + deallog.attach(logfile); + deallog.threshold_double(1.e-9); + + test<2> (); +} -- 2.39.5