From: Wolfgang Bangerth Date: Tue, 9 Jan 2024 22:26:22 +0000 (-0700) Subject: Also test MUMPS symmetric mode. X-Git-Tag: relicensing~152^2 X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=refs%2Fpull%2F16444%2Fhead;p=dealii.git Also test MUMPS symmetric mode. --- diff --git a/tests/petsc/sparse_direct_mumps.cc b/tests/petsc/sparse_direct_mumps.cc index 58910e3e1a..a07ad2d6c9 100644 --- a/tests/petsc/sparse_direct_mumps.cc +++ b/tests/petsc/sparse_direct_mumps.cc @@ -40,30 +40,52 @@ main(int argc, char **argv) Utilities::MPI::MPI_InitFinalize mpi_initialization(argc, argv, 1); { - const unsigned int size = 32; - unsigned int dim = (size - 1) * (size - 1); + const unsigned int n_points = 32; + unsigned int size = (n_points - 1) * (n_points - 1); - deallog << "Size " << size << " Unknowns " << dim << std::endl; + deallog << "Size " << n_points << " Unknowns " << size << std::endl; // Make matrix - FDMatrix testproblem(size, size); - PETScWrappers::SparseMatrix A(dim, dim, 5); + FDMatrix testproblem(n_points, n_points); + PETScWrappers::SparseMatrix A(size, size, 5); testproblem.five_point(A); - IndexSet indices(dim); - indices.add_range(0, dim); + IndexSet indices(size); + indices.add_range(0, size); PETScWrappers::MPI::Vector f(indices, MPI_COMM_WORLD); PETScWrappers::MPI::Vector u(indices, MPI_COMM_WORLD); - u = 0.; - f = 1.; A.compress(VectorOperation::insert); - SolverControl cn; - PETScWrappers::SparseDirectMUMPS solver(cn); - // solver.set_symmetric_mode(true); - solver.solve(A, u, f); + // First test non-symmetric mode: + { + deallog << "Testing non-symmetric mode" << std::endl; - PETScWrappers::MPI::Vector tmp(indices, MPI_COMM_WORLD); - deallog << "residual = " << A.residual(tmp, u, f) << std::endl; + u = 0.; + f = 1.; + + SolverControl cn; + PETScWrappers::SparseDirectMUMPS solver(cn); + solver.solve(A, u, f); + + PETScWrappers::MPI::Vector tmp(indices, MPI_COMM_WORLD); + deallog << "residual = " << A.residual(tmp, u, f) << std::endl; + } + + // Now also test the case where the matrix can be assumed to be + // symmetric (which of course it is here): + { + deallog << "Testing symmetric mode" << std::endl; + + u = 0.; + f = 1.; + + SolverControl cn; + PETScWrappers::SparseDirectMUMPS solver(cn); + solver.set_symmetric_mode(true); + solver.solve(A, u, f); + + PETScWrappers::MPI::Vector tmp(indices, MPI_COMM_WORLD); + deallog << "residual = " << A.residual(tmp, u, f) << std::endl; + } } } diff --git a/tests/petsc/sparse_direct_mumps.with_petsc_with_mumps=on.with_64bit_indices=off.output b/tests/petsc/sparse_direct_mumps.with_petsc_with_mumps=on.with_64bit_indices=off.output index 932860edcc..f3462ba980 100644 --- a/tests/petsc/sparse_direct_mumps.with_petsc_with_mumps=on.with_64bit_indices=off.output +++ b/tests/petsc/sparse_direct_mumps.with_petsc_with_mumps=on.with_64bit_indices=off.output @@ -1,4 +1,8 @@ DEAL::Size 32 Unknowns 961 -DEAL::Convergence step 1 value 0 -DEAL::residual = 0 +DEAL::Testing non-symmetric mode +DEAL::Convergence step 1 value 0.000 +DEAL::residual = 7.366e-13 +DEAL::Testing symmetric mode +DEAL::Convergence step 1 value 0.000 +DEAL::residual = 5.412e-13