<ol>
+ <li> <p>New: There is now a function <code
+ class="member">FullMatrix::trace</code> that does what its name suggests
+ it does.
+ <br>
+ (WB 2006/06/16)
+ </p>
+
<li> <p>Improved: <code class="class">PointerMatrixAux</code> now has a
- default constructor and a function for setting the <code
- class="class">VectorMemory</code> object.
- <br>
- (GK 2006/06/14)
- </p>
+ default constructor and a function for setting the <code
+ class="class">VectorMemory</code> object.
+ <br>
+ (GK 2006/06/14)
+ </p>
<li> <p>Fixed: <code
class="member">FullMatrix::print</code> would yield a link error
*/
double determinant () const;
+ /**
+ * Return the trace of the matrix,
+ * i.e. the sum of the diagonal values
+ * (which happens to also equal the sum
+ * of the eigenvalues of a matrix).
+ * Obviously, the matrix needs to
+ * be quadratic for this function.
+ */
+ double trace () const;
+
/**
* Output of the matrix in
* user-defined format.
}
+
+template <typename number>
+double
+FullMatrix<number>::trace () const
+{
+ Assert (!this->empty(), ExcEmptyMatrix());
+
+ Assert (this->n_cols() == this->n_rows(),
+ ExcDimensionMismatch(this->n_cols(), this->n_rows()));
+
+ double tr = 0;
+ for (unsigned int i=0; i<this->n_rows(); ++i)
+ tr += this->el(i,i);
+
+ return tr;
+}
+
+
+
template <typename number>
number
FullMatrix<number>::frobenius_norm () const
matrix_lib matrix_out \
solver eigen \
sparse_ilu sparse_ilu_t lapack lapack_fill block_minres \
- identity_matrix*
+ identity_matrix* trace
# from above list of regular expressions, generate the real set of
--- /dev/null
+//---------------------------------------------------------------------------
+// $Id$
+// Version: $Name$
+//
+// Copyright (C) 2005, 2006 by the deal.II authors
+//
+// This file is subject to QPL and may not be distributed
+// without copyright and license information. Please refer
+// to the file deal.II/doc/license.html for the text and
+// further information on this license.
+//
+//---------------------------------------------------------------------------
+
+// check FullMatrix::trace
+
+
+#include "../tests.h"
+#include <base/logstream.h>
+#include <lac/full_matrix.h>
+
+#include <fstream>
+#include <iostream>
+#include <cmath>
+
+
+int main()
+{
+ std::ofstream logfile("trace/output");
+ logfile.setf(std::ios::fixed);
+ logfile.precision(0);
+ deallog.attach(logfile);
+ deallog.depth_console(0);
+ deallog.threshold_double(1.e-10);
+
+ const unsigned int N=20;
+ FullMatrix<double> m (N,N);
+
+ double tr = 0;
+ for (unsigned int i=0; i<N; ++i)
+ for (unsigned int j=0; j<N; ++j)
+ {
+ m(i,j) = i+j;
+ if (i==j)
+ tr += i+j;
+ }
+
+ deallog << "Trace=" << m.trace() << std::endl;
+ Assert (m.trace() == tr, ExcInternalError());
+}
--- /dev/null
+obj.* exe OK output
--- /dev/null
+
+DEAL::Trace=380.