From: Luca Heltai Date: Tue, 13 Mar 2018 09:03:45 +0000 (+0100) Subject: Added distributed test. X-Git-Tag: v9.0.0-rc1~325^2~1 X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=1145fa70025c51582b42b432121f5dfaa7d7b886;p=dealii.git Added distributed test. --- diff --git a/tests/non_matching/coupling_05.cc b/tests/non_matching/coupling_05.cc new file mode 100644 index 0000000000..4b43d03781 --- /dev/null +++ b/tests/non_matching/coupling_05.cc @@ -0,0 +1,175 @@ +// --------------------------------------------------------------------- +// +// Copyright (C) 2018 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 "../tests.h" + +using namespace dealii; + +// Test that a coupling matrix can be constructed for each pair of dimension +// and immersed dimension, and check that quadratic functions are correctly +// projected. + +template +void test() +{ + deallog << "dim: " << dim << ", spacedim: " + << spacedim << std::endl; + + const auto &comm = MPI_COMM_WORLD; + + parallel::shared::Triangulation tria(comm); + parallel::distributed::Triangulation space_tria(comm); + + GridGenerator::hyper_cube(tria,-.4,.3); + GridGenerator::hyper_cube(space_tria,-1,1); + + tria.refine_global(3); + space_tria.refine_global(3); + + FE_Q fe(2); + FE_Q space_fe(2); + + deallog << "FE : " << fe.get_name() << std::endl + << "Space FE: " << space_fe.get_name() << std::endl; + + DoFHandler dh(tria); + DoFHandler space_dh(space_tria); + + dh.distribute_dofs(fe); + space_dh.distribute_dofs(space_fe); + + auto space_locally_owned_dofs = space_dh.locally_owned_dofs (); + auto locally_owned_dofs = dh.locally_owned_dofs (); + + + deallog << "Dofs : " << dh.n_dofs() << std::endl + << "Space dofs: " << space_dh.n_dofs() << std::endl; + + deallog << "Local dofs : " << locally_owned_dofs.n_elements() + << std::endl; + deallog << "Local space dofs: " << space_locally_owned_dofs.n_elements() + << std::endl; + QGauss quad(3); // Quadrature for coupling + + + TrilinosWrappers::SparsityPattern sparsity(space_locally_owned_dofs, + locally_owned_dofs, comm); + NonMatching::create_coupling_sparsity_pattern(space_dh, dh, + quad, sparsity); + sparsity.compress(); + + TrilinosWrappers::SparseMatrix coupling(sparsity); + NonMatching::create_coupling_mass_matrix(space_dh, dh, quad, coupling); + coupling.compress(VectorOperation::add); + + TrilinosWrappers::SparsityPattern mass_sparsity(locally_owned_dofs, comm); + DoFTools::make_sparsity_pattern(dh,mass_sparsity); + mass_sparsity.compress(); + + TrilinosWrappers::SparseMatrix mass_matrix(mass_sparsity); + + { + QGauss quad(4); + FEValues fev(fe, quad, update_values|update_JxW_values); + std::vector dofs(fe.dofs_per_cell); + FullMatrix cell_matrix(fe.dofs_per_cell, fe.dofs_per_cell); + ConstraintMatrix constraints; + + for (auto cell : dh.active_cell_iterators()) + if (cell->is_locally_owned()) + { + cell_matrix = 0; + fev.reinit(cell); + cell->get_dof_indices(dofs); + for (unsigned int i=0; i(), + space_square); + VectorTools::interpolate(dh, Functions::SquareFunction(), + squares); + + coupling.Tvmult(Mprojected_squares, space_square); + + SolverControl cn(100, 1e-12); + TrilinosWrappers::SolverCG solver(cn); + TrilinosWrappers::PreconditionILU prec; + prec.initialize(mass_matrix); + + solver.solve(mass_matrix, projected_squares, Mprojected_squares, prec); + + deallog << "Squares norm : " << projected_squares.l2_norm() << std::endl; + + projected_squares -= squares; + + deallog << "Error on squares: " << projected_squares.l2_norm() << std::endl; +} + + + +int main(int argc, char **argv) +{ + auto init = Utilities::MPI::MPI_InitFinalize(argc, argv, 1); + MPILogInitAll log(true); + test<1,2>(); + test<2,2>(); + test<2,3>(); + test<3,3>(); +} diff --git a/tests/non_matching/coupling_05.with_trilinos=true,mpirun=2.output b/tests/non_matching/coupling_05.with_trilinos=true,mpirun=2.output new file mode 100644 index 0000000000..0f3760084c --- /dev/null +++ b/tests/non_matching/coupling_05.with_trilinos=true,mpirun=2.output @@ -0,0 +1,83 @@ + +DEAL:0::dim: 1, spacedim: 2 +DEAL:0::FE : FE_Q<1,2>(2) +DEAL:0::Space FE: FE_Q<2>(2) +DEAL:0::Dofs : 17 +DEAL:0::Space dofs: 289 +DEAL:0::Local dofs : 9 +DEAL:0::Local space dofs: 153 +DEAL:0::Convergence step 9 value 4.51115e-14 +DEAL:0::Squares norm : 0.275853 +DEAL:0::Error on squares: 1.07673e-12 +DEAL:0::dim: 2, spacedim: 2 +DEAL:0::FE : FE_Q<2>(2) +DEAL:0::Space FE: FE_Q<2>(2) +DEAL:0::Dofs : 289 +DEAL:0::Space dofs: 289 +DEAL:0::Local dofs : 153 +DEAL:0::Local space dofs: 153 +DEAL:0::Convergence step 12 value 2.66570e-13 +DEAL:0::Squares norm : 1.98578 +DEAL:0::Error on squares: 4.16199e-10 +DEAL:0::dim: 2, spacedim: 3 +DEAL:0::FE : FE_Q<2,3>(2) +DEAL:0::Space FE: FE_Q<3>(2) +DEAL:0::Dofs : 289 +DEAL:0::Space dofs: 4913 +DEAL:0::Local dofs : 153 +DEAL:0::Local space dofs: 2601 +DEAL:0::Convergence step 12 value 2.66570e-13 +DEAL:0::Squares norm : 1.98578 +DEAL:0::Error on squares: 4.16199e-10 +DEAL:0::dim: 3, spacedim: 3 +DEAL:0::FE : FE_Q<3>(2) +DEAL:0::Space FE: FE_Q<3>(2) +DEAL:0::Dofs : 4913 +DEAL:0::Space dofs: 4913 +DEAL:0::Local dofs : 2601 +DEAL:0::Local space dofs: 2601 +DEAL:0::Convergence step 12 value 6.65914e-13 +DEAL:0::Squares norm : 11.6248 +DEAL:0::Error on squares: 4.61425e-08 + +DEAL:1::dim: 1, spacedim: 2 +DEAL:1::FE : FE_Q<1,2>(2) +DEAL:1::Space FE: FE_Q<2>(2) +DEAL:1::Dofs : 17 +DEAL:1::Space dofs: 289 +DEAL:1::Local dofs : 8 +DEAL:1::Local space dofs: 136 +DEAL:1::Convergence step 9 value 4.51115e-14 +DEAL:1::Squares norm : 0.275853 +DEAL:1::Error on squares: 1.07673e-12 +DEAL:1::dim: 2, spacedim: 2 +DEAL:1::FE : FE_Q<2>(2) +DEAL:1::Space FE: FE_Q<2>(2) +DEAL:1::Dofs : 289 +DEAL:1::Space dofs: 289 +DEAL:1::Local dofs : 136 +DEAL:1::Local space dofs: 136 +DEAL:1::Convergence step 12 value 2.66570e-13 +DEAL:1::Squares norm : 1.98578 +DEAL:1::Error on squares: 4.16199e-10 +DEAL:1::dim: 2, spacedim: 3 +DEAL:1::FE : FE_Q<2,3>(2) +DEAL:1::Space FE: FE_Q<3>(2) +DEAL:1::Dofs : 289 +DEAL:1::Space dofs: 4913 +DEAL:1::Local dofs : 136 +DEAL:1::Local space dofs: 2312 +DEAL:1::Convergence step 12 value 2.66570e-13 +DEAL:1::Squares norm : 1.98578 +DEAL:1::Error on squares: 4.16199e-10 +DEAL:1::dim: 3, spacedim: 3 +DEAL:1::FE : FE_Q<3>(2) +DEAL:1::Space FE: FE_Q<3>(2) +DEAL:1::Dofs : 4913 +DEAL:1::Space dofs: 4913 +DEAL:1::Local dofs : 2312 +DEAL:1::Local space dofs: 2312 +DEAL:1::Convergence step 12 value 6.65914e-13 +DEAL:1::Squares norm : 11.6248 +DEAL:1::Error on squares: 4.61425e-08 +