]> https://gitweb.dealii.org/ - dealii.git/commitdiff
lapack: test solve() with Cholesky factorization, improve documentation 5585/head
authorDenis Davydov <davydden@gmail.com>
Wed, 6 Dec 2017 17:23:18 +0000 (18:23 +0100)
committerDenis Davydov <davydden@gmail.com>
Wed, 6 Dec 2017 17:23:18 +0000 (18:23 +0100)
include/deal.II/lac/lapack_full_matrix.h
tests/lapack/full_matrix_25.cc [new file with mode: 0644]
tests/lapack/full_matrix_25.output [new file with mode: 0644]

index daf294a1f1406ed8d264b2c0aa33f31ebd85e591..a1be34d185b241aeaeef99256f2e14047f846518 100644 (file)
@@ -445,8 +445,10 @@ public:
   number reciprocal_condition_number(const number l1_norm) const;
 
   /**
-   * Same as above but for triangular matrices. The matrix has to have
-   * the LAPACKSupport::Property set, see set_property().
+   * Estimate the reciprocal of the condition number $1/k(A)$ in $L_1$ norm
+   * for triangular matrices. The matrix has to have
+   * the LAPACKSupport::Property set to either LAPACKSupport::Property::upper_triangular
+   * or LAPACKSupport::Property::lower_triangular, see set_property().
    */
   number reciprocal_condition_number() const;
 
@@ -485,21 +487,22 @@ public:
   void invert ();
 
   /**
-   * Solve the linear system with right hand side. The matrix should
-   * be either triangular or LU factorization should be previously computed.
+   * Solve the linear system with right hand side @p v and put the solution
+   * back to @p v. The matrix should be either triangular or LU/Cholesky
+   * factorization should be previously computed.
    *
    * The flag transposed indicates whether the solution of the transposed
    * system is to be performed.
    */
   void solve(Vector<number> &v,
-             const bool transposed) const;
+             const bool transposed = false) const;
 
   /**
    * Same as above but for multiple right hand sides (as many as there
    * are columns in the matrix @p B).
    */
   void solve(LAPACKFullMatrix<number> &B,
-             const bool transposed) const;
+             const bool transposed = false) const;
 
   /**
    * Solve the linear system with right hand side given by applying
diff --git a/tests/lapack/full_matrix_25.cc b/tests/lapack/full_matrix_25.cc
new file mode 100644 (file)
index 0000000..0b5fc92
--- /dev/null
@@ -0,0 +1,79 @@
+// ---------------------------------------------------------------------
+//
+// Copyright (C) 2005 - 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.
+//
+// ---------------------------------------------------------------------
+
+
+// test LAPACKFullMatrix::solve() for Cholesky factorization
+
+/* MWE for size=3 in Octave:
+R = [10,2,3; 2, 20, 6; 3, 6, 90]
+x = [2; -7; 1]
+
+> R\x
+ans =
+
+   0.268693
+  -0.385220
+   0.027836
+*/
+
+#include "../tests.h"
+#include "create_matrix.h"
+#include <deal.II/lac/lapack_full_matrix.h>
+#include <deal.II/lac/full_matrix.h>
+#include <deal.II/lac/vector.h>
+
+#include <iostream>
+
+
+template <typename NumberType>
+void
+test()
+{
+  const unsigned int size = 3;
+  LAPACKFullMatrix<NumberType> M(size);
+  M.set_property(LAPACKSupport::symmetric);
+
+  M(0,0) = 10;
+  M(0,1) = 2;
+  M(1,0) = 2;
+  M(0,2) = 3;
+  M(2,0) = 3;
+  M(1,1) = 20;
+  M(1,2) = 6;
+  M(2,1) = 6;
+  M(2,2) = 90;
+
+  Vector<NumberType> x(size), y(size);
+  x[0] = 2;
+  x[1] = -7;
+  x[2] = 1;
+
+  y = x;
+  M.compute_cholesky_factorization();
+  M.solve(y);
+  y.print(deallog.get_file_stream(), 6, false);
+}
+
+
+int main()
+{
+  const std::string logname = "output";
+  std::ofstream logfile(logname.c_str());
+  logfile.precision(3);
+  deallog.attach(logfile);
+
+  test<double>();
+
+}
diff --git a/tests/lapack/full_matrix_25.output b/tests/lapack/full_matrix_25.output
new file mode 100644 (file)
index 0000000..4f8485e
--- /dev/null
@@ -0,0 +1,2 @@
+
+0.268693 -0.385220 0.027836 

In the beginning the Universe was created. This has made a lot of people very angry and has been widely regarded as a bad move.

Douglas Adams


Typeset in Trocchi and Trocchi Bold Sans Serif.