From e074163179e0bc02c1b5eac36157a9adff194de9 Mon Sep 17 00:00:00 2001 From: Daniel Arndt Date: Thu, 26 Mar 2020 12:08:55 -0400 Subject: [PATCH] Delete BlockMatrixArray and MeanValueFilter --- .../incompatibilities/20200328DanielArndt | 4 + examples/doxygen/block_matrix_array.cc | 111 --- include/deal.II/lac/block_matrix_array.h | 660 ------------------ include/deal.II/lac/matrix_lib.h | 166 ----- include/deal.II/lac/matrix_lib.templates.h | 121 ---- include/deal.II/multigrid/mg_coarse.h | 1 - include/deal.II/multigrid/mg_transfer_block.h | 1 - source/lac/CMakeLists.txt | 2 - source/lac/block_matrix_array.cc | 463 ------------ source/lac/matrix_lib.cc | 54 -- tests/bits/block_matrix_array_01.cc | 48 -- tests/bits/block_matrix_array_01.output | 5 - tests/bits/block_sparse_matrix_1.cc | 1 - tests/lac/filtered_matrix.cc | 1 - tests/multigrid/smoother_block.cc | 222 ------ tests/multigrid/smoother_block.output | 40 -- 16 files changed, 4 insertions(+), 1896 deletions(-) create mode 100644 doc/news/changes/incompatibilities/20200328DanielArndt delete mode 100644 examples/doxygen/block_matrix_array.cc delete mode 100644 include/deal.II/lac/block_matrix_array.h delete mode 100644 include/deal.II/lac/matrix_lib.h delete mode 100644 include/deal.II/lac/matrix_lib.templates.h delete mode 100644 source/lac/block_matrix_array.cc delete mode 100644 source/lac/matrix_lib.cc delete mode 100644 tests/bits/block_matrix_array_01.cc delete mode 100644 tests/bits/block_matrix_array_01.output delete mode 100644 tests/multigrid/smoother_block.cc delete mode 100644 tests/multigrid/smoother_block.output diff --git a/doc/news/changes/incompatibilities/20200328DanielArndt b/doc/news/changes/incompatibilities/20200328DanielArndt new file mode 100644 index 0000000000..176788ff39 --- /dev/null +++ b/doc/news/changes/incompatibilities/20200328DanielArndt @@ -0,0 +1,4 @@ +Removed: The deprecated classes BlockMatrixArray (and derived classes) and MeanValueFilter +have been removed. Use BlockLinearOperator instead. +
+(Daniel Arndt, 2020/03/28) diff --git a/examples/doxygen/block_matrix_array.cc b/examples/doxygen/block_matrix_array.cc deleted file mode 100644 index 8b8ffc34a1..0000000000 --- a/examples/doxygen/block_matrix_array.cc +++ /dev/null @@ -1,111 +0,0 @@ -// --------------------------------------------------------------------- -// -// Copyright (C) 2005 - 2018 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.md at -// the top level directory of deal.II. -// -// --------------------------------------------------------------------- - - -// See documentation of BlockMatrixArray for documentation of this example - -#include -#include -#include -#include -#include -#include -#include -#include - -#include -#include - -using namespace dealii; - -double Adata[4][4] = {{4., .5, .1, 0.}, - {.5, 4., .5, .1}, - {.1, .5, 4., .5}, - {0., .1, .5, 4.}}; - -double B1data[4][2] = {{.5, .1}, {.4, .2}, {.3, .3}, {.2, .4}}; - -double B2data[2][4] = {{.3, 0., -.3, 0.}, {-.3, 0., .3, 0.}}; - -double Cdata[2][2] = {{8., 1.}, {1., 8.}}; - -int main() -{ - FullMatrix A(4, 4); - FullMatrix B1(4, 2); - FullMatrix B2(2, 4); - FullMatrix C(2, 2); - - A.fill(&Adata[0][0]); - B1.fill(&B1data[0][0]); - B2.fill(&B2data[0][0]); - C.fill(&Cdata[0][0]); - - BlockMatrixArray matrix(2, 2); - - matrix.enter(A, 0, 0, 2.); - matrix.enter(B1, 0, 1, -1.); - matrix.enter(B2, 0, 1, 1., true); - matrix.enter(B2, 1, 0, 1.); - matrix.enter(B1, 1, 0, -1., true); - matrix.enter(C, 1, 1); - matrix.print_latex(deallog); - - std::vector block_sizes(2); - block_sizes[0] = 4; - block_sizes[1] = 2; - - BlockVector result(block_sizes); - BlockVector x(block_sizes); - BlockVector y(block_sizes); - for (unsigned int i = 0; i < result.size(); ++i) - result(i) = i; - - matrix.vmult(y, result); - - SolverControl control(100, 1.e-10); - PreconditionIdentity id; - - SolverCG> cg(control); - cg.solve(matrix, x, y, id); - x.add(-1., result); - deallog << "Error " << x.l2_norm() << std::endl; - - deallog << "Error A-norm " << std::sqrt(matrix.matrix_norm_square(x)) - << std::endl; - - FullMatrix Ainv(4, 4); - Ainv.invert(A); - FullMatrix Cinv(2, 2); - Cinv.invert(C); - - BlockTrianglePrecondition precondition(2); - precondition.enter(Ainv, 0, 0, .5); - precondition.enter(Cinv, 1, 1); - - cg.solve(matrix, x, y, precondition); - x.add(-1., result); - deallog << "Error " << x.l2_norm() << std::endl; - - precondition.enter(B1, 1, 0, -1., true); - precondition.enter(B2, 1, 0, 1.); - - SolverGMRES> gmres(control); - gmres.solve(matrix, x, y, precondition); - x.add(-1., result); - deallog << "Error " << x.l2_norm() << std::endl; - - return 0; -} diff --git a/include/deal.II/lac/block_matrix_array.h b/include/deal.II/lac/block_matrix_array.h deleted file mode 100644 index e6ca1f4807..0000000000 --- a/include/deal.II/lac/block_matrix_array.h +++ /dev/null @@ -1,660 +0,0 @@ -// --------------------------------------------------------------------- -// -// Copyright (C) 2001 - 2019 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.md at -// the top level directory of deal.II. -// -// --------------------------------------------------------------------- - -#ifndef dealii_block_matrix_array_h -#define dealii_block_matrix_array_h - -#include - -#include -#include - -#include -#include -#include - -#include -#include -#include -#include -#include - -DEAL_II_NAMESPACE_OPEN - -/*! @addtogroup Matrix2 - *@{ - */ - - -/** - * Block matrix composed of different single matrices; these matrices may even - * be of different types. - * - * Given a set of arbitrary matrices Ai, this class - * implements a block matrix with block entries of the form Mjk - * = sjkAi. Each Ai may be used - * several times with different prefix. The matrices are not copied into the - * BlockMatrixArray object, but rather a PointerMatrix referencing each of - * them will be stored along with factors and transposition flags. - * - * Non-zero entries are registered by the function enter(), zero entries are - * not stored at all. Using enter() with the same location (i,j) - * several times will add the corresponding matrices in matrix-vector - * multiplications. These matrices will not be actually added, but the - * multiplications with them will be summed up. - * - * @note This mechanism makes it impossible to access single entries of - * BlockMatrixArray. In particular, (block) relaxation preconditioners based - * on PreconditionRelaxation or PreconditionBlock cannot be used with - * this class. If you need a preconditioner for a BlockMatrixArray object, use - * BlockTrianglePrecondition. - * - *

Requirements on MatrixType

- * - * The template argument MatrixType is a class providing the matrix- - * vector multiplication functions vmult(), Tvmult(), vmult_add() and - * Tvmult_add() used in this class, but with arguments of type - * Vector<number> instead of BlockVector<number>. Every matrix - * which can be used by PointerMatrix is allowed, in particular SparseMatrix - * is a possible entry type. - * - *

Example program

We document the relevant parts of - * examples/doxygen/block_matrix_array.cc. - * - * @dontinclude block_matrix_array.cc - * - * Obviously, we have to include the header file containing the definition of - * BlockMatrixArray: - * @skipline block_matrix_array.h - * - * First, we set up some matrices to be entered into the blocks. - * @skip main - * @until C.fill - * - * Now, we are ready to build a 2x2 BlockMatrixArray. - * @line Block - * First, we enter the matrix A multiplied by 2 in the upper left - * block - * @line enter - * Now -1 times B1 in the upper right block. - * @line enter - * We add the transpose of B2 to the upper right block and continue - * in a similar fashion. In the end, the block matrix structure is printed - * into an LaTeX table. - * @until latex - * - * Now, we set up vectors to be multiplied with this matrix and do a - * multiplication. - * @until vmult - * - * Finally, we solve a linear system with BlockMatrixArray, using no - * preconditioning and the conjugate gradient method. - * @until Error - * - * The remaining code of this sample program concerns preconditioning and is - * described in the documentation of BlockTrianglePrecondition. - * - * @deprecated This class has been superseded by BlockLinearOperator. - * - * @see - * @ref GlossBlockLA "Block (linear algebra)" - * @author Guido Kanschat - * @date 2000-2005, 2010 - */ -template > -class DEAL_II_DEPRECATED BlockMatrixArray : public Subscriptor -{ -public: - /** - * Declare the type for container size. - */ - using size_type = types::global_dof_index; - - /** - * Default constructor creating a useless object. initialize() must be - * called before using it. - */ - BlockMatrixArray(); - - /** - * Constructor fixing the dimensions. - */ - BlockMatrixArray(const unsigned int n_block_rows, - const unsigned int n_block_cols); - - /** - * Initialize object completely. This is the function to call for an object - * created by the default constructor. - */ - void - initialize(const unsigned int n_block_rows, const unsigned int n_block_cols); - - /** - * Adjust the matrix to a new size and delete all blocks. - */ - void - reinit(const unsigned int n_block_rows, const unsigned int n_block_cols); - - /** - * Add a block matrix entry. The matrix is entered into a list of - * blocks for multiplication, together with its coordinates row and - * col as well as optional multiplication factor prefix - * and transpose flag transpose. - * - * @note No check for consistency of block sizes is made. Therefore, - * entering a block of wrong dimension here will only lead to a - * ExcDimensionMismatch in one of the multiplication functions. - */ - template - void - enter(const MatrixType & matrix, - const unsigned int row, - const unsigned int col, - const number prefix = 1., - const bool transpose = false); - - /** - * Delete all entries, i.e. reset the matrix to an empty state. - */ - void - clear(); - - /** - * Number of block-entries per column. - */ - unsigned int - n_block_rows() const; - - /** - * Number of block-entries per row. - */ - unsigned int - n_block_cols() const; - - /** - * Matrix-vector multiplication. - */ - void - vmult(BlockVectorType &dst, const BlockVectorType &src) const; - - /** - * Matrix-vector multiplication adding to dst. - */ - void - vmult_add(BlockVectorType &dst, const BlockVectorType &src) const; - - /** - * Transposed matrix-vector multiplication. - */ - void - Tvmult(BlockVectorType &dst, const BlockVectorType &src) const; - - /** - * Transposed matrix-vector multiplication adding to dst. - */ - void - Tvmult_add(BlockVectorType &dst, const BlockVectorType &src) const; - - /** - * Matrix scalar product between two vectors (at least for a symmetric - * matrix). - */ - number - matrix_scalar_product(const BlockVectorType &u, - const BlockVectorType &v) const; - - /** - * Compute $u^T M u$. This is the square of the norm induced by the matrix - * assuming the matrix is symmetric positive definitive. - */ - number - matrix_norm_square(const BlockVectorType &u) const; - - /** - * Print the block structure as a LaTeX-array. This output will not be very - * intuitive, since the current object lacks knowledge about what the - * individual blocks represent or how they should be named. Instead, what - * you will see is an entry for each block showing all the matrices with - * their multiplication factors and possibly transpose marks. The matrices - * itself are named successively as they are encountered. If the same matrix - * is entered several times, it will be listed with different names. - * - * As an example, consider the following code: - * @code - * FullMatrix A1(4,4); - * FullMatrix A2(4,4); - * FullMatrix B(4,3); - * FullMatrix C(3,3); - * - * BlockMatrixArray block(2,2); - * - * block.enter(A1,0,0); - * block.enter(A2,0,0,2,true); - * block.enter(B,0,1,-3.); - * block.enter(B,0,1,-3.,true); - * block.enter(C,1,1,1.,true); - * - * block.print_latex(std::cout); - * @endcode - * The current function will then produce output of the following kind: - * @code - * \begin{array}{cc} - * M0+2xM1^T & -3xM2-3xM3^T\\ - * & M4^T - * \end{array} - * @endcode - * Note how the individual blocks here are just numbered successively as - * M0 to M4 and that the output misses the fact - * that M2 and M3 are, in fact, the same matrix. - * Nevertheless, the output at least gives some kind of idea of the block - * structure of this matrix. - */ - template - void - print_latex(StreamType &out) const; - -protected: - /** - * Internal data structure. - * - * For each entry of a BlockMatrixArray, its position, matrix, prefix and - * optional transposition must be stored. This structure encapsulates all of - * them. - * - * @author Guido Kanschat, 2000, 2001 - */ - class Entry - { - public: - /** - * Constructor initializing all data fields. A PointerMatrix object is - * generated for matrix. - */ - template - Entry(const MatrixType &matrix, - size_type row, - size_type col, - number prefix, - bool transpose); - - /** - * Copy constructor invalidating the old object. Since it is only used for - * entering temporary objects into a vector, this is ok. - * - * For a deep copy, we would need a reproduction operator in - * PointerMatixBase. - */ - Entry(const Entry &); - - /** - * Destructor, where we delete the PointerMatrix created by the - * constructor. - */ - ~Entry(); - - /** - * Row number in the block matrix. - */ - size_type row; - - /** - * Column number in the block matrix. - */ - size_type col; - - /** - * Factor in front of the matrix block. - */ - number prefix; - - /** - * Indicates that matrix block must be transposed for multiplication. - */ - bool transpose; - - /** - * The matrix block itself. - */ - PointerMatrixBase *matrix; - - /** - * Assignment operator. - * - * @note Since the copy constructor is destructive (see its documentation) - * and only exists for convenience there is no reasonable way to implement - * this, so it is explicitly deleted. - */ - Entry & - operator=(const Entry &) = delete; - }; - - /** - * Array of block entries in the matrix. - */ - std::vector entries; - -private: - /** - * Number of blocks per column. - */ - unsigned int block_rows; - /** - * number of blocks per row. - */ - unsigned int block_cols; -}; - -/*@}*/ - - -/** - * Inversion of a block-triangular matrix. - * - * In this block matrix, the inverses of the diagonal blocks are stored - * together with the off-diagonal blocks of a block matrix. Then, forward or - * backward insertion is performed block-wise. The diagonal blocks are NOT - * inverted for this purpose! - * - * Like for all preconditioners, the preconditioning operation is performed by - * the vmult() member function. - * - * @note While block indices may be duplicated (see BlockMatrixArray) to add - * blocks, this has to be used with caution, since summing up the inverse of - * two blocks does not yield the inverse of the sum. While the latter would be - * desirable, we can only perform the first. - * - * The implementation may be a little clumsy, but it should be sufficient as - * long as the block sizes are much larger than the number of blocks. - * - *

Example

Here, we document the second part of - * examples/doxygen/block_matrix_array.cc. For the beginning of this - * file, see BlockMatrixArray. - * - * In order to set up the preconditioner, we have to compute the inverses of - * the diagonal blocks ourselves. Since we used FullMatrix objects, this is - * fairly easy. - * @dontinclude block_matrix_array.cc - * @skip Error - * @until Cinv.invert - * - * After creating a 2x2 BlockTrianglePrecondition object, we only fill - * its diagonals. The scaling factor 1/2 used for A is the - * reciprocal of the scaling factor used for the matrix itself. - * Remember, this preconditioner actually multiplies with the diagonal - * blocks. - * @until Cinv - * - * Now, we have a block Jacobi preconditioner, which is still symmetric, since - * the blocks are symmetric. Therefore, we can still use the preconditioned - * conjugate gradient method. - * @until Error - * - * Now, we enter the subdiagonal block. This is the same as in - * matrix. - * @until B2 - * - * Since the preconditioner is not symmetric anymore, we use the GMRES method - * for solving. - * @until Error - * - * @deprecated This class has been superseded by block_back_substitution and - * block_forward_substitution, which build on BlockLinearOperator. - * - * @ingroup Preconditioners - * @author Guido Kanschat, 2001, 2005 - */ -template > -class DEAL_II_DEPRECATED BlockTrianglePrecondition - : private BlockMatrixArray -{ -public: - /** - * Declare type for container size. - */ - using size_type = types::global_dof_index; - - /** - * Default constructor creating a useless object. initialize() must be - * called before using it. - */ - BlockTrianglePrecondition(); - - /** - * Constructor. This matrix must be block-quadratic, and n_blocks - * is the number of blocks in each direction. - */ - BlockTrianglePrecondition(const unsigned int n_blocks); - - /** - * Resize preconditioner to a new size and clear all blocks. - */ - void - reinit(const unsigned int n_block_rows); - - - /** - * Enter a block. This calls BlockMatrixArray::enter(). Remember that the - * diagonal blocks should actually be inverse matrices or preconditioners. - */ - template - void - enter(const MatrixType &matrix, - const size_type row, - const size_type col, - const number prefix = 1., - const bool transpose = false); - - /** - * Preconditioning. - */ - void - vmult(BlockVectorType &dst, const BlockVectorType &src) const; - - /** - * Preconditioning adding to dst. - */ - void - vmult_add(BlockVectorType &dst, const BlockVectorType &src) const; - - /** - * Transposed preconditioning - */ - void - Tvmult(BlockVectorType &dst, const BlockVectorType &src) const; - - /** - * Transposed preconditioning adding to dst. - */ - void - Tvmult_add(BlockVectorType &dst, const BlockVectorType &src) const; - - /** - * Make function of base class available. - */ - using BlockMatrixArray::print_latex; - - /** - * Make function of base class available. - */ - using BlockMatrixArray::n_block_rows; - - /** - * Make function of base class available. - */ - using BlockMatrixArray::n_block_cols; - using BlockMatrixArray::clear; - using BlockMatrixArray::Subscriptor::subscribe; - using BlockMatrixArray::Subscriptor::unsubscribe; - - /** - * @addtogroup Exceptions - * @{ - */ - - /** - * Each diagonal block must contain one and only one matrix. If this - * exception is thrown, you did not enter a matrix here. - */ - DeclException1(ExcNoDiagonal, - size_type, - << "No diagonal entry was added for block " << arg1); - - /** - * Each diagonal block must contain one and only one matrix. If this - * exception is thrown, you entered a second matrix here. - */ - DeclException1(ExcMultipleDiagonal, - size_type, - << "Inverse diagonal entries may not be added in block " - << arg1); - //@} -private: - /** - * Add all off-diagonal contributions and return the entry of the diagonal - * element for one row. - */ - void - do_row(BlockVectorType &dst, size_type row_num) const; - - /** - * Flag for backward insertion. - */ - bool backward; -}; - - -#ifndef DOXYGEN -//--------------------------------------------------------------------------- - -template -template -inline BlockMatrixArray::Entry::Entry( - const MatrixType &m, - size_type row, - size_type col, - number prefix, - bool transpose) - : row(row) - , col(col) - , prefix(prefix) - , transpose(transpose) - , matrix(new_pointer_matrix_base(m, - typename BlockVectorType::BlockType(), - typeid(*this).name())) -{} - - - -template -template -inline void -BlockMatrixArray::enter(const MatrixType &matrix, - unsigned int row, - unsigned int col, - number prefix, - bool transpose) -{ - AssertIndexRange(row, n_block_rows()); - AssertIndexRange(col, n_block_cols()); - entries.push_back(Entry(matrix, row, col, prefix, transpose)); -} - - -template -template -inline void -BlockMatrixArray::print_latex(StreamType &out) const -{ - out << "\\begin{array}{" << std::string(n_block_cols(), 'c') << "}" - << std::endl; - - Table<2, std::string> array(n_block_rows(), n_block_cols()); - - using NameMap = - std::map *, - std::string>; - NameMap matrix_names; - - typename std::vector::const_iterator m = entries.begin(); - typename std::vector::const_iterator end = entries.end(); - - size_type matrix_number = 0; - for (; m != end; ++m) - { - if (matrix_names.find(m->matrix) == matrix_names.end()) - { - std::pair x = matrix_names.insert( - std::pair< - const PointerMatrixBase *, - std::string>(m->matrix, std::string("M"))); - std::ostringstream stream; - stream << matrix_number++; - - x.first->second += stream.str(); - } - - std::ostringstream stream; - - if (!array(m->row, m->col).empty() && m->prefix >= 0) - stream << "+"; - if (m->prefix != 1.) - stream << m->prefix << 'x'; - stream << matrix_names.find(m->matrix)->second; - // stream << '(' << m->matrix << ')'; - if (m->transpose) - stream << "^T"; - - array(m->row, m->col) += stream.str(); - } - for (unsigned int i = 0; i < n_block_rows(); ++i) - for (unsigned int j = 0; j < n_block_cols(); ++j) - { - out << '\t' << array(i, j); - if (j == n_block_cols() - 1) - { - if (i != n_block_rows() - 1) - out << "\\\\" << std::endl; - else - out << std::endl; - } - else - out << " &"; - } - out << "\\end{array}" << std::endl; -} - -template -template -inline void -BlockTrianglePrecondition::enter( - const MatrixType &matrix, - size_type row, - size_type col, - number prefix, - bool transpose) -{ - BlockMatrixArray::enter( - matrix, row, col, prefix, transpose); -} - - - -#endif // DOXYGEN - -DEAL_II_NAMESPACE_CLOSE - -#endif diff --git a/include/deal.II/lac/matrix_lib.h b/include/deal.II/lac/matrix_lib.h deleted file mode 100644 index 730377acc5..0000000000 --- a/include/deal.II/lac/matrix_lib.h +++ /dev/null @@ -1,166 +0,0 @@ -// --------------------------------------------------------------------- -// -// Copyright (C) 2002 - 2018 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.md at -// the top level directory of deal.II. -// -// --------------------------------------------------------------------- - -#ifndef dealii_matrix_lib_h -#define dealii_matrix_lib_h - -#include - -#include - -#include -#include - -DEAL_II_NAMESPACE_OPEN - -// Forward declarations -#ifndef DOXYGEN -template -class Vector; -template -class BlockVector; -template -class SparseMatrix; -#endif - -/*! @addtogroup Matrix2 - *@{ - */ - - -/** - * Mean value filter. The vmult() functions of this matrix filter out mean - * values of the vector. If the vector is of type BlockVector, then an - * additional parameter selects a single component for this operation. - * - * In mathematical terms, this class acts as if it was the matrix $I-\frac - * 1n{\mathbf 1}_n{\mathbf 1}_n^T$ where ${\mathbf 1}_n$ is a vector of size - * $n$ that has only ones as its entries. Thus, taking the dot product between - * a vector $\mathbf v$ and $\frac 1n {\mathbf 1}_n$ yields the mean - * value of the entries of ${\mathbf v}$. Consequently, $ \left[I-\frac - * 1n{\mathbf 1}_n{\mathbf 1}_n^T\right] \mathbf v = \mathbf v - \left[\frac - * 1n {\mathbf v} \cdot {\mathbf 1}_n\right]{\mathbf 1}_n$ subtracts from every - * vector element the mean value of all elements. - * - * @deprecated Use a LinearOperator, or a BlockLinearOperator instead. you - * can construct such a filter by using mean_value_filter, or - * block_diagonal_operator (with a mean_value_filter block), respectively. - * - * @author Guido Kanschat, 2002, 2003 - */ -class DEAL_II_DEPRECATED MeanValueFilter : public Subscriptor -{ -public: - /** - * Declare type for container size. - */ - using size_type = types::global_dof_index; - - /** - * Constructor, optionally selecting a component. - */ - MeanValueFilter(const size_type component = numbers::invalid_size_type); - - /** - * Subtract mean value from @p v. - */ - template - void - filter(Vector &v) const; - - /** - * Subtract mean value from @p v. - */ - template - void - filter(BlockVector &v) const; - - /** - * Return the source vector with subtracted mean value. - */ - template - void - vmult(Vector &dst, const Vector &src) const; - - /** - * Add source vector with subtracted mean value to dest. - */ - template - void - vmult_add(Vector &dst, const Vector &src) const; - - /** - * Return the source vector with subtracted mean value in selected - * component. - */ - template - void - vmult(BlockVector &dst, const BlockVector &src) const; - - /** - * Add a source to dest, where the mean value in the selected component is - * subtracted. - */ - template - void - vmult_add(BlockVector &dst, const BlockVector &src) const; - - - /** - * Not implemented. - */ - template - void - Tvmult(VectorType &, const VectorType &) const; - - /** - * Not implemented. - */ - template - void - Tvmult_add(VectorType &, const VectorType &) const; - -private: - /** - * Component for filtering block vectors. - */ - const size_type component; -}; - - - -/*@}*/ -//--------------------------------------------------------------------------- - - -template -inline void -MeanValueFilter::Tvmult(VectorType &, const VectorType &) const -{ - Assert(false, ExcNotImplemented()); -} - - -template -inline void -MeanValueFilter::Tvmult_add(VectorType &, const VectorType &) const -{ - Assert(false, ExcNotImplemented()); -} - - -DEAL_II_NAMESPACE_CLOSE - -#endif diff --git a/include/deal.II/lac/matrix_lib.templates.h b/include/deal.II/lac/matrix_lib.templates.h deleted file mode 100644 index 38885c3f17..0000000000 --- a/include/deal.II/lac/matrix_lib.templates.h +++ /dev/null @@ -1,121 +0,0 @@ -// --------------------------------------------------------------------- -// -// Copyright (C) 2002 - 2018 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.md at -// the top level directory of deal.II. -// -// --------------------------------------------------------------------- - -#ifndef dealii_matrix_lib_templates_h -#define dealii_matrix_lib_templates_h - -#include - -#include -#include -#include - -DEAL_II_NAMESPACE_OPEN - -template -void -MeanValueFilter::filter(Vector &v) const -{ - number mean = v.mean_value(); - - for (size_type i = 0; i < v.size(); ++i) - v(i) -= mean; -} - - - -template -void -MeanValueFilter::vmult(Vector &dst, const Vector &src) const -{ - Assert(dst.size() == src.size(), - ExcDimensionMismatch(dst.size(), src.size())); - - number mean = src.mean_value(); - - for (size_type i = 0; i < dst.size(); ++i) - dst(i) = src(i) - mean; -} - - - -template -void -MeanValueFilter::vmult_add(Vector &dst, const Vector &src) const -{ - Assert(dst.size() == src.size(), - ExcDimensionMismatch(dst.size(), src.size())); - - number mean = src.mean_value(); - - for (size_type i = 0; i < dst.size(); ++i) - dst(i) += src(i) - mean; -} - - - -template -void -MeanValueFilter::filter(BlockVector &v) const -{ - Assert(component != numbers::invalid_unsigned_int, ExcNotInitialized()); - - for (unsigned int i = 0; i < v.n_blocks(); ++i) - if (i == component) - vmult(v.block(i), v.block(i)); -} - - - -template -void -MeanValueFilter::vmult(BlockVector & dst, - const BlockVector &src) const -{ - Assert(component != numbers::invalid_unsigned_int, ExcNotInitialized()); - - Assert(dst.n_blocks() == src.n_blocks(), - ExcDimensionMismatch(dst.n_blocks(), src.n_blocks())); - - for (unsigned int i = 0; i < dst.n_blocks(); ++i) - if (i == component) - vmult(dst.block(i), src.block(i)); - else - dst.block(i) = src.block(i); -} - - - -template -void -MeanValueFilter::vmult_add(BlockVector & dst, - const BlockVector &src) const -{ - Assert(component != numbers::invalid_unsigned_int, ExcNotInitialized()); - - Assert(dst.n_blocks() == src.n_blocks(), - ExcDimensionMismatch(dst.n_blocks(), src.n_blocks())); - - for (unsigned int i = 0; i < dst.n_blocks(); ++i) - if (i == component) - vmult_add(dst.block(i), src.block(i)); - else - dst.block(i) += src.block(i); -} - - -DEAL_II_NAMESPACE_CLOSE - -#endif diff --git a/include/deal.II/multigrid/mg_coarse.h b/include/deal.II/multigrid/mg_coarse.h index 309f3c937d..6b1a8afdc2 100644 --- a/include/deal.II/multigrid/mg_coarse.h +++ b/include/deal.II/multigrid/mg_coarse.h @@ -22,7 +22,6 @@ #include #include #include -#include #include diff --git a/include/deal.II/multigrid/mg_transfer_block.h b/include/deal.II/multigrid/mg_transfer_block.h index 03a28cde9f..f9792fb310 100644 --- a/include/deal.II/multigrid/mg_transfer_block.h +++ b/include/deal.II/multigrid/mg_transfer_block.h @@ -22,7 +22,6 @@ #include -#include #include #include #include diff --git a/source/lac/CMakeLists.txt b/source/lac/CMakeLists.txt index 7d0862977a..7afbc26a98 100644 --- a/source/lac/CMakeLists.txt +++ b/source/lac/CMakeLists.txt @@ -17,7 +17,6 @@ INCLUDE_DIRECTORIES(BEFORE ${CMAKE_CURRENT_BINARY_DIR}) SET(_unity_include_src affine_constraints.cc - block_matrix_array.cc block_sparse_matrix.cc block_sparse_matrix_ez.cc block_sparsity_pattern.cc @@ -30,7 +29,6 @@ SET(_unity_include_src la_vector.cc la_parallel_vector.cc la_parallel_block_vector.cc - matrix_lib.cc matrix_out.cc precondition_block.cc precondition_block_ez.cc diff --git a/source/lac/block_matrix_array.cc b/source/lac/block_matrix_array.cc deleted file mode 100644 index a9f3f300b1..0000000000 --- a/source/lac/block_matrix_array.cc +++ /dev/null @@ -1,463 +0,0 @@ -// --------------------------------------------------------------------- -// -// Copyright (C) 2005 - 2018 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.md at -// the top level directory of deal.II. -// -// --------------------------------------------------------------------- - - -#include -#include -#include -#include -#include - -DEAL_II_NAMESPACE_OPEN - - -template -BlockMatrixArray::Entry::Entry(const Entry &e) - : row(e.row) - , col(e.col) - , prefix(e.prefix) - , transpose(e.transpose) - , matrix(e.matrix) -{ - Entry &e2 = const_cast(e); - e2.matrix = nullptr; -} - - - -template -BlockMatrixArray::Entry::~Entry() -{ - if (matrix) - delete matrix; -} - - -template -BlockMatrixArray::BlockMatrixArray() - : block_rows(0) - , block_cols(0) -{} - - - -template -BlockMatrixArray::BlockMatrixArray( - const unsigned int n_block_rows, - const unsigned int n_block_cols) - : block_rows(n_block_rows) - , block_cols(n_block_cols) -{} - - -template -void -BlockMatrixArray::initialize( - const unsigned int n_block_rows, - const unsigned int n_block_cols) -{ - block_rows = n_block_rows; - block_cols = n_block_cols; -} - - - -template -void -BlockMatrixArray::reinit( - const unsigned int n_block_rows, - const unsigned int n_block_cols) -{ - clear(); - block_rows = n_block_rows; - block_cols = n_block_cols; -} - - - -template -void -BlockMatrixArray::clear() -{ - entries.clear(); -} - - -template -void -BlockMatrixArray::vmult_add( - BlockVectorType & dst, - const BlockVectorType &src) const -{ - GrowingVectorMemory mem; - Assert(dst.n_blocks() == block_rows, - ExcDimensionMismatch(dst.n_blocks(), block_rows)); - Assert(src.n_blocks() == block_cols, - ExcDimensionMismatch(src.n_blocks(), block_cols)); - - typename VectorMemory::Pointer p_aux( - mem); - typename BlockVectorType::BlockType &aux = *p_aux; - - typename std::vector::const_iterator m = entries.begin(); - typename std::vector::const_iterator end = entries.end(); - - for (; m != end; ++m) - { - aux.reinit(dst.block(m->row)); - if (m->transpose) - m->matrix->Tvmult(aux, src.block(m->col)); - else - m->matrix->vmult(aux, src.block(m->col)); - dst.block(m->row).add(m->prefix, aux); - } -} - - - -template -void -BlockMatrixArray::vmult( - BlockVectorType & dst, - const BlockVectorType &src) const -{ - dst = 0.; - vmult_add(dst, src); -} - - - -template -void -BlockMatrixArray::Tvmult_add( - BlockVectorType & dst, - const BlockVectorType &src) const -{ - GrowingVectorMemory mem; - Assert(dst.n_blocks() == block_cols, - ExcDimensionMismatch(dst.n_blocks(), block_cols)); - Assert(src.n_blocks() == block_rows, - ExcDimensionMismatch(src.n_blocks(), block_rows)); - - typename std::vector::const_iterator m = entries.begin(); - typename std::vector::const_iterator end = entries.end(); - - typename VectorMemory::Pointer p_aux( - mem); - typename BlockVectorType::BlockType &aux = *p_aux; - - for (; m != end; ++m) - { - aux.reinit(dst.block(m->col)); - if (m->transpose) - m->matrix->vmult(aux, src.block(m->row)); - else - m->matrix->Tvmult(aux, src.block(m->row)); - dst.block(m->col).add(m->prefix, aux); - } -} - - - -template -void -BlockMatrixArray::Tvmult( - BlockVectorType & dst, - const BlockVectorType &src) const -{ - dst = 0.; - Tvmult_add(dst, src); -} - - - -template -number -BlockMatrixArray::matrix_scalar_product( - const BlockVectorType &u, - const BlockVectorType &v) const -{ - GrowingVectorMemory mem; - Assert(u.n_blocks() == block_rows, - ExcDimensionMismatch(u.n_blocks(), block_rows)); - Assert(v.n_blocks() == block_cols, - ExcDimensionMismatch(v.n_blocks(), block_cols)); - - typename VectorMemory::Pointer p_aux( - mem); - typename BlockVectorType::BlockType &aux = *p_aux; - - typename std::vector::const_iterator m; - typename std::vector::const_iterator end = entries.end(); - - number result = 0.; - - for (unsigned int i = 0; i < block_rows; ++i) - { - aux.reinit(u.block(i)); - for (m = entries.begin(); m != end; ++m) - { - if (m->row != i) - continue; - if (m->transpose) - m->matrix->Tvmult_add(aux, v.block(m->col)); - else - m->matrix->vmult(aux, v.block(m->col)); - } - result += u.block(i) * aux; - } - - return result; -} - - - -template -number -BlockMatrixArray::matrix_norm_square( - const BlockVectorType &u) const -{ - return matrix_scalar_product(u, u); -} - - - -template -unsigned int -BlockMatrixArray::n_block_rows() const -{ - return block_rows; -} - - - -template -unsigned int -BlockMatrixArray::n_block_cols() const -{ - return block_cols; -} - - - -//--------------------------------------------------------------------------- - -template -BlockTrianglePrecondition::BlockTrianglePrecondition() - : BlockMatrixArray() - , backward(false) -{} - - -template -BlockTrianglePrecondition::BlockTrianglePrecondition( - const unsigned int block_rows) - : BlockMatrixArray(block_rows, block_rows) - , backward(false) -{} - - -template -void -BlockTrianglePrecondition::reinit(const unsigned int n) -{ - BlockMatrixArray::reinit(n, n); -} - - -template -void -BlockTrianglePrecondition::do_row( - BlockVectorType &dst, - size_type row_num) const -{ - GrowingVectorMemory mem; - typename std::vector< - typename BlockMatrixArray::Entry>::const_iterator - m = this->entries.begin(); - typename std::vector< - typename BlockMatrixArray::Entry>::const_iterator - end = this->entries.end(); - std::vector::Entry>::const_iterator> - diagonals; - - typename VectorMemory::Pointer p_aux( - mem); - typename BlockVectorType::BlockType &aux = *p_aux; - - aux.reinit(dst.block(row_num), true); - - // Loop over all entries, since - // they are not ordered by rows. - for (; m != end; ++m) - { - const size_type i = m->row; - // Ignore everything not in - // this row - if (i != row_num) - continue; - const size_type j = m->col; - // Only use the lower (upper) - // triangle for forward - // (backward) substitution - if (((j > i) && !backward) || ((j < i) && backward)) - continue; - if (j == i) - { - diagonals.push_back(m); - } - else - { - if (m->transpose) - m->matrix->Tvmult(aux, dst.block(j)); - else - m->matrix->vmult(aux, dst.block(j)); - dst.block(i).add(-1.0 * m->prefix, aux); - } - } - Assert(diagonals.size() != 0, ExcNoDiagonal(row_num)); - - // Inverting the diagonal block is - // simple, if there is only one - // matrix - if (diagonals.size() == 1) - { - if (diagonals[0]->transpose) - diagonals[0]->matrix->Tvmult(aux, dst.block(row_num)); - else - diagonals[0]->matrix->vmult(aux, dst.block(row_num)); - dst.block(row_num).equ(diagonals[0]->prefix, aux); - } - else - { - aux = 0.; - for (size_type i = 0; i < diagonals.size(); ++i) - { - m = diagonals[i]; - // First, divide by the current - // factor, such that we can - // multiply by it later. - aux /= m->prefix; - if (m->transpose) - m->matrix->Tvmult_add(aux, dst.block(row_num)); - else - m->matrix->vmult_add(aux, dst.block(row_num)); - aux *= m->prefix; - } - dst.block(row_num) = aux; - } -} - - - -template -void -BlockTrianglePrecondition::vmult_add( - BlockVectorType & dst, - const BlockVectorType &src) const -{ - Assert(dst.n_blocks() == n_block_rows(), - ExcDimensionMismatch(dst.n_blocks(), n_block_rows())); - Assert(src.n_blocks() == n_block_cols(), - ExcDimensionMismatch(src.n_blocks(), n_block_cols())); - - BlockVectorType aux; - aux.reinit(dst); - vmult(aux, src); - dst += aux; -} - - - -template -void -BlockTrianglePrecondition::vmult( - BlockVectorType & dst, - const BlockVectorType &src) const -{ - Assert(dst.n_blocks() == n_block_rows(), - ExcDimensionMismatch(dst.n_blocks(), n_block_rows())); - Assert(src.n_blocks() == n_block_cols(), - ExcDimensionMismatch(src.n_blocks(), n_block_cols())); - - dst.equ(1., src); - - if (backward) - { - for (unsigned int i = n_block_rows(); i > 0;) - do_row(dst, --i); - } - else - { - for (unsigned int i = 0; i < n_block_rows(); ++i) - do_row(dst, i); - } -} - -template -void -BlockTrianglePrecondition::Tvmult( - BlockVectorType &, - const BlockVectorType &) const -{ - Assert(false, ExcNotImplemented()); -} - - -template -void -BlockTrianglePrecondition::Tvmult_add( - BlockVectorType &, - const BlockVectorType &) const -{ - Assert(false, ExcNotImplemented()); -} - -template class BlockMatrixArray; -template class BlockMatrixArray; -template class BlockTrianglePrecondition; -template class BlockTrianglePrecondition; - -#ifdef DEAL_II_WITH_TRILINOS -template class BlockMatrixArray; -template class BlockMatrixArray; -template class BlockTrianglePrecondition; -template class BlockTrianglePrecondition; -#endif - -#ifdef DEAL_II_WITH_PETSC -# ifdef PETSC_USE_COMPLEX -template class BlockMatrixArray, - PETScWrappers::MPI::BlockVector>; -template class BlockMatrixArray, - PETScWrappers::MPI::BlockVector>; -template class BlockTrianglePrecondition, - PETScWrappers::MPI::BlockVector>; -template class BlockTrianglePrecondition, - PETScWrappers::MPI::BlockVector>; -# else -template class BlockMatrixArray; -template class BlockMatrixArray; -template class BlockTrianglePrecondition; -template class BlockTrianglePrecondition; -# endif -#endif - -DEAL_II_NAMESPACE_CLOSE diff --git a/source/lac/matrix_lib.cc b/source/lac/matrix_lib.cc deleted file mode 100644 index 45cde9e336..0000000000 --- a/source/lac/matrix_lib.cc +++ /dev/null @@ -1,54 +0,0 @@ -// --------------------------------------------------------------------- -// -// Copyright (C) 2002 - 2018 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.md at -// the top level directory of deal.II. -// -// --------------------------------------------------------------------- - -#include -#include - -DEAL_II_NAMESPACE_OPEN - -MeanValueFilter::MeanValueFilter(size_type component) - : component(component) -{} - -template void -MeanValueFilter::filter(Vector &) const; -template void -MeanValueFilter::filter(Vector &) const; -template void -MeanValueFilter::filter(BlockVector &) const; -template void -MeanValueFilter::filter(BlockVector &) const; -template void -MeanValueFilter::vmult(Vector &, const Vector &) const; -template void -MeanValueFilter::vmult(Vector &, const Vector &) const; -template void -MeanValueFilter::vmult(BlockVector &, const BlockVector &) const; -template void -MeanValueFilter::vmult(BlockVector &, - const BlockVector &) const; - -template void -MeanValueFilter::vmult_add(Vector &, const Vector &) const; -template void -MeanValueFilter::vmult_add(Vector &, const Vector &) const; -template void -MeanValueFilter::vmult_add(BlockVector &, - const BlockVector &) const; -template void -MeanValueFilter::vmult_add(BlockVector &, - const BlockVector &) const; - -DEAL_II_NAMESPACE_CLOSE diff --git a/tests/bits/block_matrix_array_01.cc b/tests/bits/block_matrix_array_01.cc deleted file mode 100644 index 6c7ac546aa..0000000000 --- a/tests/bits/block_matrix_array_01.cc +++ /dev/null @@ -1,48 +0,0 @@ -// --------------------------------------------------------------------- -// -// Copyright (C) 2005 - 2018 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.md at -// the top level directory of deal.II. -// -// --------------------------------------------------------------------- - - - -// This tests the construction of a BlockMatrixArray and outputs the -// entered blocks using print_latex. - -#include -#include - -#include "../tests.h" - - -int -main() -{ - initlog(); - - FullMatrix A1(4, 4); - FullMatrix A2(4, 4); - FullMatrix B(4, 3); - FullMatrix C(3, 3); - - BlockMatrixArray block(2, 2); - - block.enter(A1, 0, 0); - block.enter(A2, 0, 0, 2, true); - block.enter(B, 0, 1, -3.); - block.enter(B, 0, 1, -3., true); - block.enter(C, 1, 1, 1., true); - - block.print_latex(deallog); - - return 0; -} diff --git a/tests/bits/block_matrix_array_01.output b/tests/bits/block_matrix_array_01.output deleted file mode 100644 index 084edb26a4..0000000000 --- a/tests/bits/block_matrix_array_01.output +++ /dev/null @@ -1,5 +0,0 @@ - -DEAL::\begin{array}{cc} -DEAL:: M0+2xM1^T & -3xM2-3xM3^T\\ -DEAL:: & M4^T -DEAL::\end{array} diff --git a/tests/bits/block_sparse_matrix_1.cc b/tests/bits/block_sparse_matrix_1.cc index b7782368f3..dca7043918 100644 --- a/tests/bits/block_sparse_matrix_1.cc +++ b/tests/bits/block_sparse_matrix_1.cc @@ -27,7 +27,6 @@ #include #include -#include #include #include #include diff --git a/tests/lac/filtered_matrix.cc b/tests/lac/filtered_matrix.cc index 3cac4a8632..1b91f02736 100644 --- a/tests/lac/filtered_matrix.cc +++ b/tests/lac/filtered_matrix.cc @@ -18,7 +18,6 @@ #include #include -#include #include #include diff --git a/tests/multigrid/smoother_block.cc b/tests/multigrid/smoother_block.cc deleted file mode 100644 index 01939cd915..0000000000 --- a/tests/multigrid/smoother_block.cc +++ /dev/null @@ -1,222 +0,0 @@ -// --------------------------------------------------------------------- -// -// Copyright (C) 2000 - 2018 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.md at -// the top level directory of deal.II. -// -// --------------------------------------------------------------------- - - - -#include - -#include -#include -#include - -#include - -#include - -#include "../tests.h" - -using namespace std; - -template -class ScalingMatrix : public Subscriptor -{ -public: - /** - * Constructor setting the - * scaling factor. Default is - * constructing the identity - * matrix. - */ - ScalingMatrix(number scaling_factor = 1.); - /** - * Apply preconditioner. - */ - template - void - vmult(VectorType &, const VectorType &) const; - - /** - * Apply transpose - * preconditioner. Since this is - * the identity, this function is - * the same as - * vmult(). - */ - template - void - Tvmult(VectorType &, const VectorType &) const; - /** - * Apply preconditioner, adding to the previous value. - */ - template - void - vmult_add(VectorType &, const VectorType &) const; - - /** - * Apply transpose - * preconditioner, adding. Since this is - * the identity, this function is - * the same as - * vmult_add(). - */ - template - void - Tvmult_add(VectorType &, const VectorType &) const; - -private: - number factor; -}; - - -//----------------------------------------------------------------------// - -template -ScalingMatrix::ScalingMatrix(number factor) - : factor(factor) -{} - - -template -template -inline void -ScalingMatrix::vmult(VectorType &dst, const VectorType &src) const -{ - dst.equ(factor, src); -} - -template -template -inline void -ScalingMatrix::Tvmult(VectorType &dst, const VectorType &src) const -{ - dst.equ(factor, src); -} - -template -template -inline void -ScalingMatrix::vmult_add(VectorType &dst, const VectorType &src) const -{ - dst.add(factor, src); -} - - - -template -template -inline void -ScalingMatrix::Tvmult_add(VectorType &dst, const VectorType &src) const -{ - dst.add(factor, src); -} - -//----------------------------------------------------------------------// - -template -void -check_smoother(const MGLevelObject &m, - const MGLevelObject & r) -{ - GrowingVectorMemory> mem; - MGSmootherBlock smoother; - - smoother.initialize(m, r); - - for (unsigned int l = m.min_level(); l <= m.max_level(); ++l) - { - deallog << "Level " << l << std::endl; - - BlockVector &u = *mem.alloc(); - BlockVector &f = *mem.alloc(); - u.reinit(m[l].n_block_rows(), 3); - f.reinit(u); - for (unsigned int b = 0; b < f.n_blocks(); ++b) - for (unsigned int i = 0; i < f.block(b).size(); ++i) - f.block(b)(i) = (b + 1) * (i + l); - - deallog << "First step" << std::endl; - smoother.set_steps(1); - smoother.smooth(l, u, f); - - for (unsigned int b = 0; b < u.n_blocks(); ++b) - { - for (unsigned int i = 0; i < u.block(b).size(); ++i) - deallog << '\t' << (int)(u.block(b)(i) + .5); - deallog << std::endl; - } - - deallog << "Second step" << std::endl; - smoother.smooth(l, u, f); - - for (unsigned int b = 0; b < u.n_blocks(); ++b) - { - for (unsigned int i = 0; i < u.block(b).size(); ++i) - deallog << '\t' << (int)(u.block(b)(i) + .5); - deallog << std::endl; - } - - deallog << "Two steps" << std::endl; - u = 0.; - smoother.set_steps(2); - smoother.smooth(l, u, f); - - for (unsigned int b = 0; b < u.n_blocks(); ++b) - { - for (unsigned int i = 0; i < u.block(b).size(); ++i) - deallog << '\t' << (int)(u.block(b)(i) + .5); - deallog << std::endl; - } - - mem.free(&u); - mem.free(&f); - } -} - -void -check() -{ - ScalingMatrix s1(-1.); - ScalingMatrix s2(2.); - ScalingMatrix s8(8.); - - GrowingVectorMemory> mem; - MGLevelObject> A(2, 4); - MGLevelObject> P(2, 4); - - for (unsigned int l = A.min_level(); l <= A.max_level(); ++l) - { - A[l].initialize(3, 3); - P[l].reinit(3); - for (unsigned int b = 0; b < A[l].n_block_rows(); ++b) - { - P[l].enter(s2, b, b, A[l].n_block_rows() - b); - A[l].enter(s8, b, b, 1); - for (unsigned int b2 = 0; b2 < A[l].n_block_rows(); ++b2) - A[l].enter(s1, b, b2, 1.); - } - } - - check_smoother(A, P); -} - - -int -main() -{ - initlog(); - deallog << std::setprecision(3); - - check(); -} diff --git a/tests/multigrid/smoother_block.output b/tests/multigrid/smoother_block.output deleted file mode 100644 index 2e2858b820..0000000000 --- a/tests/multigrid/smoother_block.output +++ /dev/null @@ -1,40 +0,0 @@ - -DEAL::Level 2 -DEAL::First step -DEAL:: 12 18 24 -DEAL:: 16 24 32 -DEAL:: 12 18 24 -DEAL::Second step -DEAL:: -311 -467 -623 -DEAL:: -319 -479 -639 -DEAL:: -87 -131 -175 -DEAL::Two steps -DEAL:: -311 -467 -623 -DEAL:: -319 -479 -639 -DEAL:: -87 -131 -175 -DEAL::Level 3 -DEAL::First step -DEAL:: 18 24 30 -DEAL:: 24 32 40 -DEAL:: 18 24 30 -DEAL::Second step -DEAL:: -467 -623 -779 -DEAL:: -479 -639 -799 -DEAL:: -131 -175 -219 -DEAL::Two steps -DEAL:: -467 -623 -779 -DEAL:: -479 -639 -799 -DEAL:: -131 -175 -219 -DEAL::Level 4 -DEAL::First step -DEAL:: 24 30 36 -DEAL:: 32 40 48 -DEAL:: 24 30 36 -DEAL::Second step -DEAL:: -623 -779 -935 -DEAL:: -639 -799 -959 -DEAL:: -175 -219 -263 -DEAL::Two steps -DEAL:: -623 -779 -935 -DEAL:: -639 -799 -959 -DEAL:: -175 -219 -263 -- 2.39.5