From: Matthias Maier Date: Sat, 18 Apr 2015 15:07:36 +0000 (+0200) Subject: provide a doxygen module, finish documentation X-Git-Tag: v8.3.0-rc1~242^2 X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=60c6f4276c6331ea2b60b51415060871ceb0a7ce;p=dealii.git provide a doxygen module, finish documentation --- diff --git a/doc/doxygen/headers/laoperators.h b/doc/doxygen/headers/laoperators.h new file mode 100644 index 0000000000..85b8176087 --- /dev/null +++ b/doc/doxygen/headers/laoperators.h @@ -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 + * + *

Linear Operator

+ * + * If deal.II is configured with C++11 support (i.e., + * DEAL_II_WITH_CXX11=on 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 applying a + * linear operation on a vector. + * @code + * std::function vmult; + * std::function vmult_add; + * std::function Tvmult; + * std::function 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 op that performs above computation when + * applied on a vector, one can write: + * @code + * dealii::SparseMatrix 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, op 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 + */ diff --git a/include/deal.II/lac/block_matrix.h b/include/deal.II/lac/block_matrix.h index 05ff3570a1..60788236b4 100644 --- a/include/deal.II/lac/block_matrix.h +++ b/include/deal.II/lac/block_matrix.h @@ -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 diff --git a/include/deal.II/lac/iterative_inverse.h b/include/deal.II/lac/iterative_inverse.h index b1574ed10e..f47ca91a11 100644 --- a/include/deal.II/lac/iterative_inverse.h +++ b/include/deal.II/lac/iterative_inverse.h @@ -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 diff --git a/include/deal.II/lac/linear_operator.h b/include/deal.II/lac/linear_operator.h index 26d3ba2133..7ee7bcda82 100644 --- a/include/deal.II/lac/linear_operator.h +++ b/include/deal.II/lac/linear_operator.h @@ -53,10 +53,10 @@ LinearOperator 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 vmult; - * std::function vmult_add; - * std::function Tvmult; - * std::function Tvmult_add; + * std::function vmult; + * std::function vmult_add; + * std::function Tvmult; + * std::function Tvmult_add; * @endcode * * But, in contrast to a usual matrix object, the domain and range of the @@ -92,11 +92,13 @@ LinearOperator 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 DEAL_II_WITH_CXX11 is enabled during * cmake configure. * * @author Luca Heltai, Matthias Maier, 2015 + * + * @ingroup LAOperators */ template 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 LinearOperator @@ -311,6 +315,8 @@ operator+(const LinearOperator &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 LinearOperator @@ -327,6 +333,8 @@ operator-(const LinearOperator &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 LinearOperator @@ -403,6 +411,8 @@ operator*(const LinearOperator &first_op, * Domain & operator *=(Domain::value_type); * Range & operator *=(Range::value_type); * @endcode + * + * @ingroup LAOperators */ template LinearOperator @@ -460,6 +470,8 @@ operator*(typename Range::value_type number, * Domain & operator *=(Domain::value_type); * Range & operator *=(Range::value_type); * @endcode + * + * @ingroup LAOperators */ template LinearOperator @@ -478,6 +490,8 @@ operator*(const LinearOperator &op, * \relates LinearOperator * * Returns the transpose linear operations of @ref op. + * + * @ingroup LAOperators */ template LinearOperator @@ -507,6 +521,8 @@ transpose_operator(const LinearOperator &op) * reinit_vector as an argument to initialize the * reinit_range_vector and reinit_domain_vector * objects of the LinearOperator object. + * + * @ingroup LAOperators */ template LinearOperator @@ -556,6 +572,8 @@ identity_operator(const std::function &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 vmult or Tvmult. + * + * @ingroup LAOperators */ template LinearOperator @@ -632,6 +650,8 @@ inverse_operator(const LinearOperator>({op_a00, op_a01, op_a10, op_a11}); * @endcode + * + * @ingroup LAOperators */ template , @@ -732,6 +752,8 @@ block_operator(const std::array>({op_00, op_a1, ..., op_am}); * @endcode + * + * @ingroup LAOperators */ template , @@ -814,6 +836,8 @@ block_diagonal_operator(const std::array, @@ -1117,6 +1141,8 @@ namespace * storage). * * @author Matthias Maier, 2015 + * + * @ingroup LAOperators */ template , typename Domain = Range, @@ -1142,6 +1168,8 @@ LinearOperator linear_operator(const Matrix &matrix) * matrix). * * @author Matthias Maier, 2015 + * + * @ingroup LAOperators */ template , typename Domain = Range, diff --git a/include/deal.II/lac/matrix_lib.h b/include/deal.II/lac/matrix_lib.h index 405902c27b..a776be95d7 100644 --- a/include/deal.II/lac/matrix_lib.h +++ b/include/deal.II/lac/matrix_lib.h @@ -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 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 @@ -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 @@ -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 @@ -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 diff --git a/include/deal.II/lac/shifted_matrix.h b/include/deal.II/lac/shifted_matrix.h index b46edd0bfb..0eb4259dd3 100644 --- a/include/deal.II/lac/shifted_matrix.h +++ b/include/deal.II/lac/shifted_matrix.h @@ -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 A, this class implements a matrix-vector product * with A+s I, where s 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 @@ -91,6 +95,10 @@ private: * with A+s M, where s is a provided shift parameter and * M 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 diff --git a/include/deal.II/lac/transpose_matrix.h b/include/deal.II/lac/transpose_matrix.h index b941bc9c4c..f018f9e028 100644 --- a/include/deal.II/lac/transpose_matrix.h +++ b/include/deal.II/lac/transpose_matrix.h @@ -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 */