]> https://gitweb.dealii.org/ - dealii.git/commitdiff
add LAPACKFullMatrix::trace() 5211/head
authorDenis Davydov <davydden@gmail.com>
Mon, 9 Oct 2017 09:06:34 +0000 (11:06 +0200)
committerDenis Davydov <davydden@gmail.com>
Mon, 9 Oct 2017 09:06:34 +0000 (11:06 +0200)
doc/news/changes/minor/20171009DenisDavydov [new file with mode: 0644]
include/deal.II/lac/lapack_full_matrix.h
source/lac/lapack_full_matrix.cc
tests/lapack/full_matrix_20.cc [new file with mode: 0644]
tests/lapack/full_matrix_20.output [new file with mode: 0644]

diff --git a/doc/news/changes/minor/20171009DenisDavydov b/doc/news/changes/minor/20171009DenisDavydov
new file mode 100644 (file)
index 0000000..10c2aa6
--- /dev/null
@@ -0,0 +1,3 @@
+New: add LAPACKFullMatrix::trace() to calculate trace of a matrix.
+<br>
+(Denis Davydov, 2017/10/09)
index a1a2c413a58694bf9abeae396472aa041e4dd5ec..ede5e3e48c08b601986fb85460bee4824f638e53 100644 (file)
@@ -459,6 +459,12 @@ public:
    */
   number frobenius_norm() const;
 
+  /**
+   * Compute trace of the matrix, i.e. the sum of the diagonal values.
+   * Obviously, the matrix needs to be quadratic for this function.
+   */
+  number trace() const;
+
   /**
    * Invert the matrix by first computing an LU/Cholesky factorization with the LAPACK
    * function Xgetrf/Xpotrf and then building the actual inverse using Xgetri/Xpotri.
index dc781474164b8d2680adacd218b55e98265f224a..e298960064eba90ad7b0fec20150c3acc4214e9d 100644 (file)
@@ -594,6 +594,25 @@ number LAPACKFullMatrix<number>::norm(const char type) const
 }
 
 
+
+template <typename number>
+number LAPACKFullMatrix<number>::trace() const
+{
+  Assert (state == LAPACKSupport::matrix ||
+          state == LAPACKSupport::inverse_matrix,
+          ExcMessage("Trace can be called in matrix state only."));
+  Assert (this->n_cols() == this->n_rows(),
+          ExcDimensionMismatch(this->n_cols(), this->n_rows()));
+
+  number tr = 0;
+  for (size_type i=0; i<this->n_rows(); ++i)
+    tr += (*this)(i,i);
+
+  return tr;
+}
+
+
+
 template <typename number>
 void
 LAPACKFullMatrix<number>::compute_cholesky_factorization()
diff --git a/tests/lapack/full_matrix_20.cc b/tests/lapack/full_matrix_20.cc
new file mode 100644 (file)
index 0000000..363acef
--- /dev/null
@@ -0,0 +1,57 @@
+// ---------------------------------------------------------------------
+//
+// 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::trace() by comparing to FullMatrix
+
+#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)
+{
+  // Full matrix:
+  FullMatrix<NumberType> F(size);
+  create_random(F);
+
+  // Lapack:
+  LAPACKFullMatrix<NumberType> M(size);
+  M = F;
+
+  deallog << "trace difference: " << (F.trace() - M.trace())<< std::endl;
+}
+
+
+int main()
+{
+  const std::string logname = "output";
+  std::ofstream logfile(logname.c_str());
+  logfile.precision(3);
+  deallog.attach(logfile);
+
+  const std::vector<unsigned int> sizes = {{17,391}};
+  for (const auto &s : sizes)
+    {
+      deallog << "size=" << s << std::endl;
+      test<double>(s);
+    }
+}
diff --git a/tests/lapack/full_matrix_20.output b/tests/lapack/full_matrix_20.output
new file mode 100644 (file)
index 0000000..2ea8f8c
--- /dev/null
@@ -0,0 +1,5 @@
+
+DEAL::size=17
+DEAL::trace difference: 0.00000
+DEAL::size=391
+DEAL::trace difference: 0.00000

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.