From 62d0ec667eff3054270e39c5fadfcc14b5e535c7 Mon Sep 17 00:00:00 2001 From: Denis Davydov Date: Thu, 13 Oct 2016 00:20:22 +0200 Subject: [PATCH] add a quick test for trilinos --- tests/quick_tests/CMakeLists.txt | 6 + tests/quick_tests/step-trilinos.cc | 222 +++++++++++++++++++++++++++++ 2 files changed, 228 insertions(+) create mode 100644 tests/quick_tests/step-trilinos.cc diff --git a/tests/quick_tests/CMakeLists.txt b/tests/quick_tests/CMakeLists.txt index 8ed1d64706..07b30b05cf 100644 --- a/tests/quick_tests/CMakeLists.txt +++ b/tests/quick_tests/CMakeLists.txt @@ -111,6 +111,12 @@ IF (DEAL_II_WITH_PETSC AND DEAL_II_WITH_SLEPC) make_quicktest("step-slepc" ${_mybuild} "") ENDIF() +# Test trilinos +IF (DEAL_II_WITH_TRILINOS) + make_quicktest("step-trilinos" ${_mybuild} "") +ENDIF() + + # A custom test target: ADD_CUSTOM_TARGET(test COMMAND ${CMAKE_COMMAND} -D ALL_TESTS="${ALL_TESTS}" -P ${CMAKE_CURRENT_SOURCE_DIR}/run.cmake diff --git a/tests/quick_tests/step-trilinos.cc b/tests/quick_tests/step-trilinos.cc new file mode 100644 index 0000000000..df87c68ddb --- /dev/null +++ b/tests/quick_tests/step-trilinos.cc @@ -0,0 +1,222 @@ +/* --------------------------------------------------------------------- + * + * Copyright (C) 2013 - 2015 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 +#include + +#include +#include + +using namespace dealii; + +// Test that deal.II is working with Trilinos by solving the Laplace's +// problem in 2d. +class LaplaceProblem +{ +public: + LaplaceProblem (); + void run (); + +private: + void setup_system (); + void assemble_system (); + void solve (); + + Triangulation<2> triangulation; + FE_Q<2> fe; + DoFHandler<2> dof_handler; + + TrilinosWrappers::SparseMatrix A; + TrilinosWrappers::Vector b, x; + ConstraintMatrix constraints; + + TableHandler output_table; +}; + +LaplaceProblem::LaplaceProblem () + : + fe (1), + dof_handler (triangulation) +{} + +void LaplaceProblem::setup_system () +{ + dof_handler.distribute_dofs (fe); + + constraints.clear (); + DoFTools::make_zero_boundary_constraints (dof_handler, constraints); + constraints.close (); + + DynamicSparsityPattern dsp(dof_handler.n_dofs()); + DoFTools::make_sparsity_pattern (dof_handler, dsp, constraints, false); + + A.reinit (dsp); + b.reinit (dof_handler.n_dofs()); + x.reinit (dof_handler.n_dofs()); + + // some output + output_table.add_value ("cells", triangulation.n_active_cells()); + output_table.add_value ("dofs", dof_handler.n_dofs()); +} + +void LaplaceProblem::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); + Vector cell_b (dofs_per_cell); + + std::vector local_dof_indices (dofs_per_cell); + + 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 LaplaceProblem::solve () +{ + SolverControl solver_control (1e03, 1e-03); + TrilinosWrappers::SolverCG cg_solver (solver_control); + TrilinosWrappers::PreconditionBlockJacobi preconditioner; + preconditioner.initialize (A); + cg_solver.solve (A, x, b, preconditioner); + + // some output + // ? +} + +void LaplaceProblem::run () +{ + GridGenerator::hyper_cube (triangulation, -1, 1); + + for (unsigned int c=0; c<5; ++c) + { + triangulation.refine_global (1); + setup_system (); + assemble_system (); + solve (); + } + + // finialise output + output_table.write_text (std::cout); + deallog << std::endl; +} + + +int main (int argc, char **argv) +{ + try + { + Utilities::MPI::MPI_InitFinalize mpi_initialization (argc, argv, 1); + { + LaplaceProblem 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; +} -- 2.39.5