From: Niklas Fehn Date: Thu, 9 Mar 2017 10:21:34 +0000 (+0100) Subject: test has been added to verify the state of LAPACKFullMatrix after having computed... X-Git-Tag: v8.5.0-rc1~45^2~1 X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=3760edde5f023aece81e87e7943cd999fc765800;p=dealii.git test has been added to verify the state of LAPACKFullMatrix after having computed the LU factorization --- diff --git a/source/lac/lapack_full_matrix.cc b/source/lac/lapack_full_matrix.cc index d1e8fc612e..ae80e56929 100644 --- a/source/lac/lapack_full_matrix.cc +++ b/source/lac/lapack_full_matrix.cc @@ -504,10 +504,10 @@ LAPACKFullMatrix::compute_lu_factorization() getrf(&mm, &nn, values, &mm, &ipiv[0], &info); Assert(info >= 0, ExcInternalError()); - + // if info >= 0, the factorization has been completed state = lu; - + AssertThrow(info == 0, LACExceptions::ExcSingular()); } diff --git a/tests/lapack/full_matrix_14.cc b/tests/lapack/full_matrix_14.cc new file mode 100644 index 0000000000..77e3f1cec9 --- /dev/null +++ b/tests/lapack/full_matrix_14.cc @@ -0,0 +1,68 @@ +// --------------------------------------------------------------------- +// +// Copyright (C) 2014 - 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. +// +// --------------------------------------------------------------------- + + +// Tests LAPACKFullMatrix::operator*= and operator/= + +#include "../tests.h" +#include +#include + +#include +#include + + + +void test(const bool is_singular) +{ + const unsigned int n=10; + LAPACKFullMatrix A(n,n); + A = 0; + + // generate a singular matrix if is_singular == true + for (unsigned int i=0; i v(n); + v = 1.0; + v(n-1) = 0.0; + // LU factorization can only be applied if state == lu! + A.apply_lu_factorization(v,false); + + deallog << "apply lu factorization succeeded with norm " << v.l2_norm() << std::endl; +} + +int main() +{ + const std::string logname = "output"; + std::ofstream logfile(logname.c_str()); + logfile.precision(3); + deallog.attach(logfile); + deallog.threshold_double(1.e-10); + + test(false); + test(true); +} diff --git a/tests/lapack/full_matrix_14.output b/tests/lapack/full_matrix_14.output new file mode 100644 index 0000000000..566cdbebc4 --- /dev/null +++ b/tests/lapack/full_matrix_14.output @@ -0,0 +1,4 @@ + +DEAL::apply lu factorization succeeded with norm 3.00000 +DEAL::matrix is singular +DEAL::apply lu factorization succeeded with norm 3.00000