]> https://gitweb.dealii.org/ - dealii.git/commitdiff
provide a doxygen module, finish documentation 750/head
authorMatthias Maier <tamiko@43-1.org>
Sat, 18 Apr 2015 15:07:36 +0000 (17:07 +0200)
committerMatthias Maier <tamiko@43-1.org>
Sun, 19 Apr 2015 20:55:07 +0000 (22:55 +0200)
doc/doxygen/headers/laoperators.h [new file with mode: 0644]
include/deal.II/lac/block_matrix.h
include/deal.II/lac/iterative_inverse.h
include/deal.II/lac/linear_operator.h
include/deal.II/lac/matrix_lib.h
include/deal.II/lac/shifted_matrix.h
include/deal.II/lac/transpose_matrix.h

diff --git a/doc/doxygen/headers/laoperators.h b/doc/doxygen/headers/laoperators.h
new file mode 100644 (file)
index 0000000..85b8176
--- /dev/null
@@ -0,0 +1,82 @@
+// ---------------------------------------------------------------------
+//
+// Copyright (C) 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.
+//
+// ---------------------------------------------------------------------
+
+
+
+/**
+ * @defgroup LAOperators Linear Operators
+ *
+ * <h3>Linear Operator</h3>
+ *
+ * If deal.II is configured with C++11 support (i.e.,
+ * <code>DEAL_II_WITH_CXX11=on</code> during configuration) a versatile
+ * mechanism for storing the concept of a linear operator is available.
+ *
+ * This is done with a LinearOperator class that, similarly to the abstract
+ * MATRIX interface, defines a minimal interface for <i>applying</i> a
+ * linear operation on a vector.
+ * @code
+ *   std::function<void(Range &, const Domain &)> vmult;
+ *   std::function<void(Range &, const Domain &)> vmult_add;
+ *   std::function<void(Domain &, const Range &)> Tvmult;
+ *   std::function<void(Domain &, const Range &)> Tvmult_add;
+ * @endcode
+ *
+ * Thus, such an object can be used as a matrix object in all
+ * @ref Solvers "iterative solver" classes, either as a matrix object, or as
+ * @ref Preconditioners "preconditioner".
+ *
+ * The big advantage of the LinearOperator class is that it provides
+ * syntactic sugar for complex matrix-vector operations. As an example
+ * consider the operation $(A+k\,B)\,C$, where $A$, $B$ and $C$ denote
+ * (possibly different) SparseMatrix objects. In order to construct a
+ * LinearOperator <code>op</code> that performs above computation when
+ * applied on a vector, one can write:
+ * @code
+ * dealii::SparseMatrix<double> A, B, C;
+ * // Setup and assembly...
+ * const double k = ...;
+ *
+ * const auto op_a = linear_operator(A);
+ * const auto op_b = linear_operator(B);
+ * const auto op_c = linear_operator(C);
+ *
+ * const auto op = (op_a + k * op_b) * op_c;
+ * @endcode
+ * Now, <code>op</code> can be used as a matrix object for further
+ * computation.
+ *
+ * The linear_operator() function can be used to wrap an ordinary matrix or
+ * preconditioner object into a LinearOperator. A linear operator can be
+ * transposed with transpose_operator(), or inverted by using the
+ * inverse_operator() together with an iterative solver.
+ *
+ * For objects of type LinearOperator, all vector space operations, i.e.,
+ * addition and subtraction, scalar multiplication and composition (of
+ * compatible linear operators) are implemented.
+ *
+ * block_operator() and block_diagonal_operator() provide further
+ * encapsulation of individual linear operators into blocked linear
+ * operator variants.
+ *
+ * @note The LinearOperator facility obsoletes some of the @ref Matrix2
+ * "derived matrix" classes, such as BlockDiagonalMatrix, IterativeInverse,
+ * ProductMatrix, ScaledMatrix, ProductSparseMatrix,
+ * InverseMatrixRichardson, SchurMatrix, ShiftedMatrix,
+ * ShiftedMatrixGeneralized, TransposeMatrix
+ *
+ * @ingroup LAC
+ * @ingroup MATRICES
+ */
index 05ff3570a13153a5b4cf6d930508e952b68031b0..60788236b4ca358e6e8aa3a96d331cb54f507c02 100644 (file)
@@ -1,6 +1,6 @@
 // ---------------------------------------------------------------------
 //
-// Copyright (C) 2000 - 2014 by the deal.II authors
+// Copyright (C) 2000 - 2015 by the deal.II authors
 //
 // This file is part of the deal.II library.
 //
@@ -38,6 +38,10 @@ DEAL_II_NAMESPACE_OPEN
  * One special application is a one by one block matrix, allowing to apply the
  * @p vmult of the original matrix (or preconditioner) to a block vector.
  *
+ * @deprecated If deal.II was configured with C++11 support, use the
+ * LinearOperator class instead, see the module on
+ * @ref LAOperators "linear operators" for further details.
+ *
  * @see
  * @ref GlossBlockLA "Block (linear algebra)"
  * @author Guido Kanschat, 2000
index b1574ed10e8c55755f890a0af110090588ebdeb5..f47ca91a11fcffaae864c71fcbd411d7ffa9ed49 100644 (file)
@@ -1,6 +1,6 @@
 // ---------------------------------------------------------------------
 //
-// Copyright (C) 1999 - 2014 by the deal.II authors
+// Copyright (C) 1999 - 2015 by the deal.II authors
 //
 // This file is part of the deal.II library.
 //
@@ -66,6 +66,10 @@ DEAL_II_NAMESPACE_OPEN
  * from machine accuracy, even if the errors of the outer loop are in the
  * range of machine accuracy.
  *
+ * @deprecated If deal.II was configured with C++11 support, use the
+ * LinearOperator class instead, see the module on
+ * @ref LAOperators "linear operators" for further details.
+ *
  * @ingroup Matrix2
  * @author Guido Kanschat
  * @date 2010
index 26d3ba2133a976ff8800372becf4e31256c68821..7ee7bcda82a6ac1c964b2a5f13c482bd155197bd 100644 (file)
@@ -53,10 +53,10 @@ LinearOperator<Range, Domain> linear_operator (const Matrix &);
  * that store the knowledge of how to apply the linear operator by
  * implementing the abstract @ref Matrix interface:
  * @code
- * std::function<const void(Range &, const Domain &)> vmult;
- * std::function<const void(Range &, const Domain &)> vmult_add;
- * std::function<const void(Domain &, const Range &)> Tvmult;
- * std::function<const void(Domain &, const Range &)> Tvmult_add;
+ *   std::function<void(Range &, const Domain &)> vmult;
+ *   std::function<void(Range &, const Domain &)> vmult_add;
+ *   std::function<void(Domain &, const Range &)> Tvmult;
+ *   std::function<void(Domain &, const Range &)> Tvmult_add;
  * @endcode
  *
  * But, in contrast to a usual matrix object, the domain and range of the
@@ -92,11 +92,13 @@ LinearOperator<Range, Domain> linear_operator (const Matrix &);
  * const auto op = (op_a + k * op_b) * op_c;
  * @endcode
  *
- * @note: This class is only available if deal.II was configured with C++11
+ * @note This class is only available if deal.II was configured with C++11
  * support, i.e., if <code>DEAL_II_WITH_CXX11</code> is enabled during
  * cmake configure.
  *
  * @author Luca Heltai, Matthias Maier, 2015
+ *
+ * @ingroup LAOperators
  */
 template <typename Range, typename Domain> class LinearOperator
 {
@@ -264,6 +266,8 @@ public:
  *
  * Addition of two linear operators @p first_op and @p second_op given
  * by $(\text{first_op}+\text{second_op})x:=\text{first_op}(x)+\text{second_op}(x)$
+ *
+ * @ingroup LAOperators
  */
 template <typename Range, typename Domain>
 LinearOperator<Range, Domain>
@@ -311,6 +315,8 @@ operator+(const LinearOperator<Range, Domain> &first_op,
  *
  * Subtraction of two linear operators @p first_op and @p second_op given
  * by $(\text{first_op}-\text{second_op})x:=\text{first_op}(x)-\text{second_op}(x)$
+ *
+ * @ingroup LAOperators
  */
 template <typename Range, typename Domain>
 LinearOperator<Range, Domain>
@@ -327,6 +333,8 @@ operator-(const LinearOperator<Range, Domain> &first_op,
  *
  * Concatenation of two linear operators @p first_op and @p second_op given
  * by $(\text{first_op}*\text{second_op})x:=\text{first_op}(\text{second_op}(x))$
+ *
+ * @ingroup LAOperators
  */
 template <typename Range, typename Intermediate, typename Domain>
 LinearOperator<Range, Domain>
@@ -403,6 +411,8 @@ operator*(const LinearOperator<Range, Intermediate> &first_op,
  * Domain & operator *=(Domain::value_type);
  * Range & operator *=(Range::value_type);
  * @endcode
+ *
+ * @ingroup LAOperators
  */
 template <typename Range, typename Domain>
 LinearOperator<Range, Domain>
@@ -460,6 +470,8 @@ operator*(typename Range::value_type  number,
  * Domain & operator *=(Domain::value_type);
  * Range & operator *=(Range::value_type);
  * @endcode
+ *
+ * @ingroup LAOperators
  */
 template <typename Range, typename Domain>
 LinearOperator<Range, Domain>
@@ -478,6 +490,8 @@ operator*(const LinearOperator<Range, Domain> &op,
  * \relates LinearOperator
  *
  * Returns the transpose linear operations of @ref op.
+ *
+ * @ingroup LAOperators
  */
 template <typename Range, typename Domain>
 LinearOperator<Domain, Range>
@@ -507,6 +521,8 @@ transpose_operator(const LinearOperator<Range, Domain> &op)
  * reinit_vector as an argument to initialize the
  * <code>reinit_range_vector</code> and <code>reinit_domain_vector</code>
  * objects of the LinearOperator object.
+ *
+ * @ingroup LAOperators
  */
 template <typename Range>
 LinearOperator<Range, Range>
@@ -556,6 +572,8 @@ identity_operator(const std::function<void(Range &, bool)> &reinit_vector)
  * reference for the whole lifetime of the LinearOperator object. Internal
  * data structures of the @p solver object will be modified upon
  * invocation of <code>vmult</code> or <code>Tvmult</code>.
+ *
+ * @ingroup LAOperators
  */
 template <typename Solver, typename Preconditioner>
 LinearOperator<typename Solver::vector_type, typename Solver::vector_type>
@@ -632,6 +650,8 @@ inverse_operator(const LinearOperator<typename Solver::vector_type, typename Sol
  * @code
  * block_operator<2, 2, BlockVector<double>>({op_a00, op_a01, op_a10, op_a11});
  * @endcode
+ *
+ * @ingroup LAOperators
  */
 template <unsigned int m, unsigned int n,
           typename Range = BlockVector<double>,
@@ -732,6 +752,8 @@ block_operator(const std::array<std::array<LinearOperator<typename Range::BlockT
  * @code
  * block_diagonal_operator<m, BlockVector<double>>({op_00, op_a1, ..., op_am});
  * @endcode
+ *
+ * @ingroup LAOperators
  */
 template <unsigned int m,
           typename Range = BlockVector<double>,
@@ -814,6 +836,8 @@ block_diagonal_operator(const std::array<LinearOperator<typename Range::BlockTyp
  * A variant of above function that only takes a single LinearOperator
  * argument @p op and creates a blockdiagonal linear operator with @p m
  * copies of it.
+ *
+ * @ingroup LAOperators
  */
 template <unsigned int m,
           typename Range = BlockVector<double>,
@@ -1117,6 +1141,8 @@ namespace
  * storage).
  *
  * @author Matthias Maier, 2015
+ *
+ * @ingroup LAOperators
  */
 template <typename Range = Vector<double>,
           typename Domain = Range,
@@ -1142,6 +1168,8 @@ LinearOperator<Range, Domain> linear_operator(const Matrix &matrix)
  * matrix).
  *
  * @author Matthias Maier, 2015
+ *
+ * @ingroup LAOperators
  */
 template <typename Range = Vector<double>,
           typename Domain = Range,
index 405902c27bc3a517bbce265e9cbd1f08e464631e..a776be95d71a1b8de889ce08c8430140cdb18c4f 100644 (file)
@@ -1,6 +1,6 @@
 // ---------------------------------------------------------------------
 //
-// Copyright (C) 2002 - 2014 by the deal.II authors
+// Copyright (C) 2002 - 2015 by the deal.II authors
 //
 // This file is part of the deal.II library.
 //
@@ -44,6 +44,10 @@ template<typename number> class SparseMatrix;
  * Here is an example multiplying two different FullMatrix objects:
  * @include product_matrix.cc
  *
+ * @deprecated If deal.II was configured with C++11 support, use the
+ * LinearOperator class instead, see the module on
+ * @ref LAOperators "linear operators" for further details.
+ *
  * @author Guido Kanschat, 2000, 2001, 2002, 2005
  */
 template<class VECTOR>
@@ -142,6 +146,10 @@ private:
  * Matrix-vector products of this matrix are composed of those of the original
  * matrix with the vector and then scaling of the result by a constant factor.
  *
+ * @deprecated If deal.II was configured with C++11 support, use the
+ * LinearOperator class instead, see the module on
+ * @ref LAOperators "linear operators" for further details.
+ *
  * @author Guido Kanschat, 2007
  */
 template<class VECTOR>
@@ -205,6 +213,10 @@ private:
  * The documentation of ProductMatrix applies with exception that these
  * matrices here may be rectangular.
  *
+ * @deprecated If deal.II was configured with C++11 support, use the
+ * LinearOperator class instead, see the module on
+ * @ref LAOperators "linear operators" for further details.
+ *
  * @author Guido Kanschat, 2000, 2001, 2002, 2005
  */
 template<typename number, typename vector_number>
@@ -407,6 +419,10 @@ private:
  * @ref Instantiations
  * in the manual).
  *
+ * @deprecated If deal.II was configured with C++11 support, use the
+ * LinearOperator class instead, see the module on
+ * @ref LAOperators "linear operators" for further details.
+ *
  * @author Guido Kanschat, 2005
  */
 template<class VECTOR>
index b46edd0bfb3bd7dc81104269568b2f191833cfe2..0eb4259dd360db78c95527c7705c974ee5874ae0 100644 (file)
@@ -1,6 +1,6 @@
 // ---------------------------------------------------------------------
 //
-// Copyright (C) 2001 - 2014 by the deal.II authors
+// Copyright (C) 2001 - 2015 by the deal.II authors
 //
 // This file is part of the deal.II library.
 //
@@ -32,6 +32,10 @@ DEAL_II_NAMESPACE_OPEN
  * Given a matrix <tt>A</tt>, this class implements a matrix-vector product
  * with <i>A+s I</i>, where <i>s</i> is a provided shift parameter.
  *
+ * @deprecated If deal.II was configured with C++11 support, use the
+ * LinearOperator class instead, see the module on @ref LAOperators "linear
+ * operators" for further details.
+ *
  * @author Guido Kanschat, 2000, 2001
  */
 template<class MATRIX>
@@ -91,6 +95,10 @@ private:
  * with <i>A+s M</i>, where <i>s</i> is a provided shift parameter and
  * <tt>M</tt> is the matrix representing the identity
  *
+ * @deprecated If deal.II was configured with C++11 support, use the
+ * LinearOperator class instead, see the module on
+ * @ref LAOperators "linear operators" for further details.
+ *
  * @author Guido Kanschat, 2001
  */
 template<class MATRIX, class MASSMATRIX, class VECTOR>
index b941bc9c4cbb2bb63a96578451468997a44ff015..f018f9e028f062d0a46d6bc8b0947493c3b7b3cf 100644 (file)
@@ -1,6 +1,6 @@
 // ---------------------------------------------------------------------
 //
-// Copyright (C) 2002 - 2014 by the deal.II authors
+// Copyright (C) 2002 - 2015 by the deal.II authors
 //
 // This file is part of the deal.II library.
 //
@@ -32,6 +32,10 @@ DEAL_II_NAMESPACE_OPEN
  * @note The transposed matrix is never actually assembled. Instead, only the
  * matrix vector multiplication is performed in a transposed way.
  *
+ * @deprecated If deal.II was configured with C++11 support, use the
+ * LinearOperator class instead, see the module on
+ * @ref LAOperators "linear operators" for further details.
+ *
  * @ingroup Matrix2
  * @author Guido Kanschat, 2006
  */

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.