From e8d1950b073d02795b89b2cb51c0f2d9ec2a8d7a Mon Sep 17 00:00:00 2001 From: Daniel Arndt Date: Mon, 30 Mar 2020 20:20:57 -0400 Subject: [PATCH] Remove deprecated class FilteredMatrix --- .../incompatibilities/20200401DanielArndt | 4 + include/deal.II/lac/filtered_matrix.h | 1008 ----------------- include/deal.II/numerics/matrix_tools.h | 4 +- .../deal.II/numerics/vector_tools.templates.h | 117 +- source/grid/grid_tools.cc | 42 +- tests/data_out/data_out_curved_cells.cc | 53 +- tests/grid/grid_transform_02.output | 52 +- tests/lac/filtered_matrix.cc | 99 -- tests/lac/filtered_matrix.output | 19 - tests/lac/filtered_matrix_01.cc | 66 -- tests/lac/filtered_matrix_01.output | 5 - tests/lac/filtered_matrix_02.cc | 66 -- tests/lac/filtered_matrix_02.output | 5 - tests/lac/filtered_matrix_03.cc | 70 -- tests/lac/filtered_matrix_03.output | 5 - tests/lac/filtered_matrix_04.cc | 70 -- tests/lac/filtered_matrix_04.output | 5 - tests/lac/filtered_matrix_05.cc | 64 -- tests/lac/filtered_matrix_05.output | 5 - tests/lac/filtered_matrix_06.cc | 224 ---- tests/lac/filtered_matrix_06.output | 480 -------- tests/mappings/mapping_q1_eulerian_01.cc | 1 - tests/mappings/mapping_q_eulerian_01.cc | 1 - tests/mappings/mapping_q_eulerian_02.cc | 1 - tests/mappings/mapping_q_eulerian_03.cc | 1 - tests/mappings/mapping_q_eulerian_07.cc | 1 - tests/mappings/mapping_q_eulerian_08.cc | 1 - tests/mappings/mapping_q_eulerian_11.cc | 1 - tests/matrix_free/interpolate_to_mg.cc | 1 - 29 files changed, 88 insertions(+), 2383 deletions(-) create mode 100644 doc/news/changes/incompatibilities/20200401DanielArndt delete mode 100644 include/deal.II/lac/filtered_matrix.h delete mode 100644 tests/lac/filtered_matrix.cc delete mode 100644 tests/lac/filtered_matrix.output delete mode 100644 tests/lac/filtered_matrix_01.cc delete mode 100644 tests/lac/filtered_matrix_01.output delete mode 100644 tests/lac/filtered_matrix_02.cc delete mode 100644 tests/lac/filtered_matrix_02.output delete mode 100644 tests/lac/filtered_matrix_03.cc delete mode 100644 tests/lac/filtered_matrix_03.output delete mode 100644 tests/lac/filtered_matrix_04.cc delete mode 100644 tests/lac/filtered_matrix_04.output delete mode 100644 tests/lac/filtered_matrix_05.cc delete mode 100644 tests/lac/filtered_matrix_05.output delete mode 100644 tests/lac/filtered_matrix_06.cc delete mode 100644 tests/lac/filtered_matrix_06.output diff --git a/doc/news/changes/incompatibilities/20200401DanielArndt b/doc/news/changes/incompatibilities/20200401DanielArndt new file mode 100644 index 0000000000..4ddea854ab --- /dev/null +++ b/doc/news/changes/incompatibilities/20200401DanielArndt @@ -0,0 +1,4 @@ +Removed: The deprecated FilteredMatrix class has been removed. +Use the LinearOperator class instead. +
+(Daniel Arndt, 2020/04/01) diff --git a/include/deal.II/lac/filtered_matrix.h b/include/deal.II/lac/filtered_matrix.h deleted file mode 100644 index 0fe20a7d5f..0000000000 --- a/include/deal.II/lac/filtered_matrix.h +++ /dev/null @@ -1,1008 +0,0 @@ -// --------------------------------------------------------------------- -// -// Copyright (C) 2001 - 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_filtered_matrix_h -# define dealii_filtered_matrix_h - - - -# include - -# include -# include -# include - -# include -# include - -# include -# include - -DEAL_II_NAMESPACE_OPEN - -// Forward declaration -# ifndef DOXYGEN -template -class FilteredMatrixBlock; -# endif - -/*! @addtogroup Matrix2 - *@{ - */ - - -/** - * This class is a wrapper for linear systems of equations with simple - * equality constraints fixing individual degrees of freedom to a certain - * value such as when using Dirichlet boundary values. - * - * In order to accomplish this, the vmult(), Tvmult(), vmult_add() and - * Tvmult_add functions modify the same function of the original matrix such - * as if all constrained entries of the source vector were zero. Additionally, - * all constrained entries of the destination vector are set to zero. - * - *

Usage

- * - * Usage is simple: create an object of this type, point it to a matrix that - * shall be used for $A$ above (either through the constructor, the copy - * constructor, or the set_referenced_matrix() function), specify the list of - * boundary values or other constraints (through the add_constraints() - * function), and then for each required solution modify the right hand side - * vector (through apply_constraints()) and use this object as matrix object - * in a linear solver. As linear solvers should only use vmult() and - * residual() functions of a matrix class, this class should be as good a - * matrix as any other for that purpose. - * - * Furthermore, also the precondition_Jacobi() function is provided (since the - * computation of diagonal elements of the filtered matrix $A_X$ is simple), - * so you can use this as a preconditioner. Some other functions useful for - * matrices are also available. - * - * A typical code snippet showing the above steps is as follows: - * @code - * // set up sparse matrix A and right hand side b somehow - * ... - * - * // initialize filtered matrix with matrix and boundary value constraints - * FilteredMatrix > filtered_A (A); - * filtered_A.add_constraints (boundary_values); - * - * // set up a linear solver - * SolverControl control (1000, 1.e-10, false, false); - * GrowingVectorMemory > mem; - * SolverCG > solver (control, mem); - * - * // set up a preconditioner object - * PreconditionJacobi > prec; - * prec.initialize (A, 1.2); - * FilteredMatrix > filtered_prec (prec); - * filtered_prec.add_constraints (boundary_values); - * - * // compute modification of right hand side - * filtered_A.apply_constraints (b, true); - * - * // solve for solution vector x - * solver.solve (filtered_A, x, b, filtered_prec); - * @endcode - * - * - *

Connection to other classes

- * - * The function MatrixTools::apply_boundary_values() does exactly the same - * that this class does, except for the fact that that function actually - * modifies the matrix. Consequently, it is only possible to solve with a - * matrix to which MatrixTools::apply_boundary_values() was applied for one - * right hand side and one set of boundary values since the modification of - * the right hand side depends on the original matrix. - * - * While this is a feasible method in cases where only one solution of the - * linear system is required, for example in solving linear stationary - * systems, one would often like to have the ability to solve multiple times - * with the same matrix in nonlinear problems (where one often does not want - * to update the Hessian between Newton steps, despite having different right - * hand sides in subsequent steps) or time dependent problems, without having - * to re-assemble the matrix or copy it to temporary matrices with which one - * then can work. For these cases, this class is meant. - * - * - *

Some background

Mathematically speaking, it is used to represent a - * system of linear equations $Ax=b$ with the constraint that $B_D x = g_D$, - * where $B_D$ is a rectangular matrix with exactly one $1$ in each row, and - * these $1$s in those columns representing constrained degrees of freedom - * (e.g. for Dirichlet boundary nodes, thus the index $D$) and zeroes for all - * other diagonal entries, and $g_D$ having the requested nodal values for - * these constrained nodes. Thus, the underdetermined equation $B_D x = g_D$ - * fixes only the constrained nodes and does not impose any condition on the - * others. We note that $B_D B_D^T = 1_D$, where $1_D$ is the identity matrix - * with dimension as large as the number of constrained degrees of freedom. - * Likewise, $B_D^T B_D$ is the diagonal matrix with diagonal entries $0$ or - * $1$ that, when applied to a vector, leaves all constrained nodes untouched - * and deletes all unconstrained ones. - * - * For solving such a system of equations, we first write down the Lagrangian - * $L=1/2 x^T A x - x^T b + l^T B_D x$, where $l$ is a Lagrange multiplier for - * the constraints. The stationarity condition then reads - * @code - * [ A B_D^T ] [x] = [b ] - * [ B_D 0 ] [l] = [g_D] - * @endcode - * - * The first equation then reads $B_D^T l = b-Ax$. On the other hand, if we - * left-multiply the first equation by $B_D^T B_D$, we obtain $B_D^T B_D A x + - * B_D^T l = B_D^T B_D b$ after equating $B_D B_D^T$ to the identity matrix. - * Inserting the previous equality, this yields $(A - B_D^T B_D A) x = (1 - - * B_D^T B_D)b$. Since $x=(1 - B_D^T B_D) x + B_D^T B_D x = (1 - B_D^T B_D) x - * + B_D^T g_D$, we can restate the linear system: $A_D x = (1 - B_D^T B_D)b - - * (1 - B_D^T B_D) A B^T g_D$, where $A_D = (1 - B_D^T B_D) A (1 - B_D^T B_D)$ - * is the matrix where all rows and columns corresponding to constrained nodes - * have been deleted. - * - * The last system of equation only defines the value of the unconstrained - * nodes, while the constrained ones are determined by the equation $B_D x = - * g_D$. We can combine these two linear systems by using the zeroed out rows - * of $A_D$: if we set the diagonal to $1$ and the corresponding zeroed out - * element of the right hand side to that of $g_D$, then this fixes the - * constrained elements as well. We can write this as follows: $A_X x = (1 - - * B_D^T B_D)b - (1 - B_D^T B_D) A B^T g_D + B_D^T g_D$, where $A_X = A_D + - * B_D^T B_D$. Note that the two parts of the latter matrix operate on - * disjoint subspaces (the first on the unconstrained nodes, the latter on the - * constrained ones). - * - * In iterative solvers, it is not actually necessary to compute $A_X$ - * explicitly, since only matrix-vector operations need to be performed. This - * can be done in a three-step procedure that first clears all elements in the - * incoming vector that belong to constrained nodes, then performs the product - * with the matrix $A$, then clears again. This class is a wrapper to this - * procedure, it takes a pointer to a matrix with which to perform matrix- - * vector products, and does the cleaning of constrained elements itself. This - * class therefore implements an overloaded @p vmult function that does the - * matrix-vector product, as well as @p Tvmult for transpose matrix-vector - * multiplication and @p residual for residual computation, and can thus be - * used as a matrix replacement in linear solvers. - * - * It also has the ability to generate the modification of the right hand - * side, through the apply_constraints() function. - * - * - * - *

Template arguments

- * - * This class takes as template arguments a matrix and a vector class. The - * former must provide @p vmult, @p vmult_add, @p Tvmult, and @p residual - * member function that operate on the vector type (the second template - * argument). The latter template parameter must provide access to individual - * elements through operator(), assignment through - * operator=. - * - * - *

Thread-safety

- * - * The functions that operate as a matrix and do not change the internal state - * of this object are synchronized and thus threadsafe. Consequently, you do - * not need to serialize calls to @p vmult or @p residual . - * - * @author Wolfgang Bangerth 2001, Luca Heltai 2006, Guido Kanschat 2007, 2008 - * - * @deprecated Use LinearOperator instead. See the documentation of - * constrained_linear_operator(). - */ -template -class DEAL_II_DEPRECATED FilteredMatrix : public Subscriptor -{ -public: - class const_iterator; - - /** - * Declare the type of container size. - */ - using size_type = types::global_dof_index; - - /** - * Accessor class for iterators - */ - class Accessor - { - /** - * Constructor. Since we use accessors only for read access, a const - * matrix pointer is sufficient. - */ - Accessor(const FilteredMatrix *matrix, const size_type index); - - public: - /** - * Row number of the element represented by this object. - */ - size_type - row() const; - - /** - * Column number of the element represented by this object. - */ - size_type - column() const; - - /** - * Value of the right hand side for this row. - */ - double - value() const; - - private: - /** - * Advance to next entry - */ - void - advance(); - - /** - * The matrix accessed. - */ - const FilteredMatrix *matrix; - - /** - * Current row number. - */ - size_type index; - - // Make enclosing class a friend. - friend class const_iterator; - }; - - /** - * Standard-conforming iterator. - */ - class const_iterator - { - public: - /** - * Constructor. - */ - const_iterator(const FilteredMatrix *matrix, - const size_type index); - - /** - * Prefix increment. - */ - const_iterator & - operator++(); - - /** - * Postfix increment. - */ - const_iterator & - operator++(int); - - /** - * Dereferencing operator. - */ - const Accessor &operator*() const; - - /** - * Dereferencing operator. - */ - const Accessor *operator->() const; - - /** - * Comparison. True, if both iterators point to the same matrix position. - */ - bool - operator==(const const_iterator &) const; - /** - * Inverse of ==. - */ - bool - operator!=(const const_iterator &) const; - - /** - * Comparison operator. Result is true if either the first row number is - * smaller or if the row numbers are equal and the first index is smaller. - */ - bool - operator<(const const_iterator &) const; - - /** - * Comparison operator. Compares just the other way around than the - * operator above. - */ - bool - operator>(const const_iterator &) const; - - private: - /** - * Store an object of the accessor class. - */ - Accessor accessor; - }; - - /** - * Typedef defining a type that represents a pair of degree of freedom index - * and the value it shall have. - */ - using IndexValuePair = std::pair; - - /** - * @name Constructors and initialization - */ - //@{ - /** - * Default constructor. You will have to set the matrix to be used later - * using initialize(). - */ - FilteredMatrix(); - - /** - * Copy constructor. Use the matrix and the constraints set in the given - * object for the present one as well. - */ - FilteredMatrix(const FilteredMatrix &fm); - - /** - * Constructor. Use the given matrix for future operations. - * - * @arg @p m: The matrix being used in multiplications. - * - * @arg @p expect_constrained_source: See documentation of - * #expect_constrained_source. - */ - template - FilteredMatrix(const MatrixType &matrix, - const bool expect_constrained_source = false); - - /** - * Copy operator. Take over matrix and constraints from the other object. - */ - FilteredMatrix & - operator=(const FilteredMatrix &fm); - - /** - * Set the matrix to be used further on. You will probably also want to call - * the clear_constraints() function if constraints were previously added. - * - * @arg @p m: The matrix being used in multiplications. - * - * @arg @p expect_constrained_source: See documentation of - * #expect_constrained_source. - */ - template - void - initialize(const MatrixType &m, const bool expect_constrained_source = false); - - /** - * Delete all constraints and the matrix pointer. - */ - void - clear(); - //@} - /** - * @name Managing constraints - */ - //@{ - /** - * Add the constraint that the value with index i should have the - * value v. - */ - void - add_constraint(const size_type i, const double v); - - /** - * Add a list of constraints to the ones already managed by this object. The - * actual data type of this list must be so that dereferenced iterators are - * pairs of indices and the corresponding values to be enforced on the - * respective solution vector's entry. Thus, the data type might be, for - * example, a @p std::list or @p std::vector of IndexValuePair objects, but - * also a std::map. - * - * The second component of these pairs will only be used in - * apply_constraints(). The first is used to set values to zero in matrix - * vector multiplications. - * - * It is an error if the argument contains an entry for a degree of freedom - * that has already been constrained previously. - */ - template - void - add_constraints(const ConstraintList &new_constraints); - - /** - * Delete the list of constraints presently in use. - */ - void - clear_constraints(); - //@} - /** - * Vector operations - */ - //@{ - /** - * Apply the constraints to a right hand side vector. This needs to be done - * before starting to solve with the filtered matrix. If the matrix is - * symmetric (i.e. the matrix itself, not only its sparsity pattern), set - * the second parameter to @p true to use a faster algorithm. Note: This - * method is deprecated as matrix_is_symmetric parameter is no longer used. - */ - DEAL_II_DEPRECATED - void - apply_constraints(VectorType &v, const bool matrix_is_symmetric) const; - /** - * Apply the constraints to a right hand side vector. This needs to be done - * before starting to solve with the filtered matrix. - */ - void - apply_constraints(VectorType &v) const; - - /** - * Matrix-vector multiplication: this operation performs pre_filter(), - * multiplication with the stored matrix and post_filter() in that order. - */ - void - vmult(VectorType &dst, const VectorType &src) const; - - /** - * Matrix-vector multiplication: this operation performs pre_filter(), - * transposed multiplication with the stored matrix and post_filter() in - * that order. - */ - void - Tvmult(VectorType &dst, const VectorType &src) const; - - /** - * Adding matrix-vector multiplication. - * - * @note The result vector of this multiplication will have the constraint - * entries set to zero, independent of the previous value of dst. - * We expect that in most cases this is the required behavior. - */ - void - vmult_add(VectorType &dst, const VectorType &src) const; - - /** - * Adding transpose matrix-vector multiplication: - * - * @note The result vector of this multiplication will have the constraint - * entries set to zero, independent of the previous value of dst. - * We expect that in most cases this is the required behavior. - */ - void - Tvmult_add(VectorType &dst, const VectorType &src) const; - //@} - - /** - * @name Iterators - */ - //@{ - /** - * Iterator to the first constraint. - */ - const_iterator - begin() const; - /** - * Final iterator. - */ - const_iterator - end() const; - //@} - - /** - * Determine an estimate for the memory consumption (in bytes) of this - * object. Since we are not the owner of the matrix referenced, its memory - * consumption is not included. - */ - std::size_t - memory_consumption() const; - -private: - /** - * Determine, whether multiplications can expect that the source vector has - * all constrained entries set to zero. - * - * If so, the auxiliary vector can be avoided and memory as well as time can - * be saved. - * - * We expect this for instance in Newton's method, where the residual - * already should be zero on constrained nodes. This is, because there is no - * test function in these nodes. - */ - bool expect_constrained_source; - - /** - * Declare an abbreviation for an iterator into the array constraint pairs, - * since that data type is so often used and is rather awkward to write out - * each time. - */ - using const_index_value_iterator = - typename std::vector::const_iterator; - - /** - * Helper class used to sort pairs of indices and values. Only the index is - * considered as sort key. - */ - struct PairComparison - { - /** - * Function comparing the pairs @p i1 and @p i2 for their keys. - */ - bool - operator()(const IndexValuePair &i1, const IndexValuePair &i2) const; - }; - - /** - * Pointer to the sparsity pattern used for this matrix. - */ - std::shared_ptr> matrix; - - /** - * Sorted list of pairs denoting the index of the variable and the value to - * which it shall be fixed. - */ - std::vector constraints; - - /** - * Do the pre-filtering step, i.e. zero out those components that belong to - * constrained degrees of freedom. - */ - void - pre_filter(VectorType &v) const; - - /** - * Do the postfiltering step, i.e. set constrained degrees of freedom to the - * value of the input vector, as the matrix contains only ones on the - * diagonal for these degrees of freedom. - */ - void - post_filter(const VectorType &in, VectorType &out) const; - - friend class Accessor; - // FilteredMatrixBlock accesses pre_filter() and post_filter(). - friend class FilteredMatrixBlock; -}; - -/*@}*/ -/*---------------------- Inline functions -----------------------------------*/ - - -//--------------------------------Iterators--------------------------------------// - -template -inline FilteredMatrix::Accessor::Accessor( - const FilteredMatrix *matrix, - const size_type index) - : matrix(matrix) - , index(index) -{ - AssertIndexRange(index, matrix->constraints.size() + 1); -} - - - -template -inline types::global_dof_index -FilteredMatrix::Accessor::row() const -{ - return matrix->constraints[index].first; -} - - - -template -inline types::global_dof_index -FilteredMatrix::Accessor::column() const -{ - return matrix->constraints[index].first; -} - - - -template -inline double -FilteredMatrix::Accessor::value() const -{ - return matrix->constraints[index].second; -} - - - -template -inline void -FilteredMatrix::Accessor::advance() -{ - Assert(index < matrix->constraints.size(), ExcIteratorPastEnd()); - ++index; -} - - - -template -inline FilteredMatrix::const_iterator::const_iterator( - const FilteredMatrix *matrix, - const size_type index) - : accessor(matrix, index) -{} - - - -template -inline typename FilteredMatrix::const_iterator & -FilteredMatrix::const_iterator::operator++() -{ - accessor.advance(); - return *this; -} - - -template -inline const typename FilteredMatrix::Accessor & - FilteredMatrix::const_iterator::operator*() const -{ - return accessor; -} - - -template -inline const typename FilteredMatrix::Accessor * - FilteredMatrix::const_iterator::operator->() const -{ - return &accessor; -} - - -template -inline bool -FilteredMatrix::const_iterator:: -operator==(const const_iterator &other) const -{ - return (accessor.index == other.accessor.index && - accessor.matrix == other.accessor.matrix); -} - - -template -inline bool -FilteredMatrix::const_iterator:: -operator!=(const const_iterator &other) const -{ - return !(*this == other); -} - - - -//------------------------- FilteredMatrix ----------------------------------// - -template -inline typename FilteredMatrix::const_iterator -FilteredMatrix::begin() const -{ - return const_iterator(this, 0); -} - - -template -inline typename FilteredMatrix::const_iterator -FilteredMatrix::end() const -{ - return const_iterator(this, constraints.size()); -} - - -template -inline bool -FilteredMatrix::PairComparison:: -operator()(const IndexValuePair &i1, const IndexValuePair &i2) const -{ - return (i1.first < i2.first); -} - - - -template -template -inline void -FilteredMatrix::initialize(const MatrixType &m, bool ecs) -{ - matrix.reset(new_pointer_matrix_base(m, VectorType())); - - expect_constrained_source = ecs; -} - - - -template -inline FilteredMatrix::FilteredMatrix() - : expect_constrained_source(false) -{} - - - -template -inline FilteredMatrix::FilteredMatrix(const FilteredMatrix &fm) - : Subscriptor() - , expect_constrained_source(fm.expect_constrained_source) - , matrix(fm.matrix) - , constraints(fm.constraints) -{} - - - -template -template -inline FilteredMatrix::FilteredMatrix(const MatrixType &m, - const bool ecs) - : expect_constrained_source(false) -{ - initialize(m, ecs); -} - - - -template -inline FilteredMatrix & -FilteredMatrix::operator=(const FilteredMatrix &fm) -{ - matrix = fm.matrix; - expect_constrained_source = fm.expect_constrained_source; - constraints = fm.constraints; - return *this; -} - - - -template -inline void -FilteredMatrix::add_constraint(const size_type index, - const double value) -{ - // add new constraint to end - constraints.push_back(IndexValuePair(index, value)); -} - - - -template -template -inline void -FilteredMatrix::add_constraints( - const ConstraintList &new_constraints) -{ - // add new constraints to end - const size_type old_size = constraints.size(); - constraints.reserve(old_size + new_constraints.size()); - constraints.insert(constraints.end(), - new_constraints.begin(), - new_constraints.end()); - // then merge the two arrays to - // form one sorted one - std::inplace_merge(constraints.begin(), - constraints.begin() + old_size, - constraints.end(), - PairComparison()); -} - - - -template -inline void -FilteredMatrix::clear_constraints() -{ - // swap vectors to release memory - std::vector empty; - constraints.swap(empty); -} - - - -template -inline void -FilteredMatrix::clear() -{ - clear_constraints(); - matrix.reset(); -} - - - -template -inline void -FilteredMatrix::apply_constraints( - VectorType &v, - const bool /* matrix_is_symmetric */) const -{ - apply_constraints(v); -} - - -template -inline void -FilteredMatrix::apply_constraints(VectorType &v) const -{ - GrowingVectorMemory mem; - typename VectorMemory::Pointer tmp_vector(mem); - tmp_vector->reinit(v); - const_index_value_iterator i = constraints.begin(); - const const_index_value_iterator e = constraints.end(); - for (; i != e; ++i) - { - AssertIsFinite(i->second); - (*tmp_vector)(i->first) = -i->second; - } - - // This vmult is without bc, to get - // the rhs correction in a correct - // way. - matrix->vmult_add(v, *tmp_vector); - // finally set constrained - // entries themselves - for (i = constraints.begin(); i != e; ++i) - { - AssertIsFinite(i->second); - v(i->first) = i->second; - } -} - - -template -inline void -FilteredMatrix::pre_filter(VectorType &v) const -{ - // iterate over all constraints and - // zero out value - const_index_value_iterator i = constraints.begin(); - const const_index_value_iterator e = constraints.end(); - for (; i != e; ++i) - v(i->first) = 0; -} - - - -template -inline void -FilteredMatrix::post_filter(const VectorType &in, - VectorType & out) const -{ - // iterate over all constraints and - // set value correctly - const_index_value_iterator i = constraints.begin(); - const const_index_value_iterator e = constraints.end(); - for (; i != e; ++i) - { - AssertIsFinite(in(i->first)); - out(i->first) = in(i->first); - } -} - - - -template -inline void -FilteredMatrix::vmult(VectorType &dst, const VectorType &src) const -{ - if (!expect_constrained_source) - { - GrowingVectorMemory mem; - typename VectorMemory::Pointer tmp_vector(mem); - // first copy over src vector and - // pre-filter - tmp_vector->reinit(src, true); - *tmp_vector = src; - pre_filter(*tmp_vector); - // then let matrix do its work - matrix->vmult(dst, *tmp_vector); - } - else - { - matrix->vmult(dst, src); - } - - // finally do post-filtering - post_filter(src, dst); -} - - - -template -inline void -FilteredMatrix::Tvmult(VectorType &dst, const VectorType &src) const -{ - if (!expect_constrained_source) - { - GrowingVectorMemory mem; - typename VectorMemory::Pointer tmp_vector(mem); - // first copy over src vector and - // pre-filter - tmp_vector->reinit(src, true); - *tmp_vector = src; - pre_filter(*tmp_vector); - // then let matrix do its work - matrix->Tvmult(dst, *tmp_vector); - } - else - { - matrix->Tvmult(dst, src); - } - - // finally do post-filtering - post_filter(src, dst); -} - - - -template -inline void -FilteredMatrix::vmult_add(VectorType & dst, - const VectorType &src) const -{ - if (!expect_constrained_source) - { - GrowingVectorMemory mem; - typename VectorMemory::Pointer tmp_vector(mem); - // first copy over src vector and - // pre-filter - tmp_vector->reinit(src, true); - *tmp_vector = src; - pre_filter(*tmp_vector); - // then let matrix do its work - matrix->vmult_add(dst, *tmp_vector); - } - else - { - matrix->vmult_add(dst, src); - } - - // finally do post-filtering - post_filter(src, dst); -} - - - -template -inline void -FilteredMatrix::Tvmult_add(VectorType & dst, - const VectorType &src) const -{ - if (!expect_constrained_source) - { - GrowingVectorMemory mem; - typename VectorMemory::Pointer tmp_vector(mem); - // first copy over src vector and - // pre-filter - tmp_vector->reinit(src, true); - *tmp_vector = src; - pre_filter(*tmp_vector); - // then let matrix do its work - matrix->Tvmult_add(dst, *tmp_vector); - } - else - { - matrix->Tvmult_add(dst, src); - } - - // finally do post-filtering - post_filter(src, dst); -} - - - -template -inline std::size_t -FilteredMatrix::memory_consumption() const -{ - return (MemoryConsumption::memory_consumption(matrix) + - MemoryConsumption::memory_consumption(constraints)); -} - - - -DEAL_II_NAMESPACE_CLOSE - -#endif -/*---------------------------- filtered_matrix.h ---------------------------*/ diff --git a/include/deal.II/numerics/matrix_tools.h b/include/deal.II/numerics/matrix_tools.h index 588d5962da..14a8bc5599 100644 --- a/include/deal.II/numerics/matrix_tools.h +++ b/include/deal.II/numerics/matrix_tools.h @@ -735,8 +735,8 @@ namespace MatrixCreator * depends on the original matrix, this is not possible without storing the * original matrix somewhere and applying the @p apply_boundary_conditions * function to a copy of it each time we want to solve. In that case, you can - * use the FilteredMatrix class in the @p LAC sublibrary. There you can also - * find a formal (mathematical) description of the process of modifying the + * use the constrained_linear_operator() function. In its documentation you can + * also find a formal (mathematical) description of the process of modifying the * matrix and right hand side vectors for boundary values. * * diff --git a/include/deal.II/numerics/vector_tools.templates.h b/include/deal.II/numerics/vector_tools.templates.h index 366263d6a0..8bce666efa 100644 --- a/include/deal.II/numerics/vector_tools.templates.h +++ b/include/deal.II/numerics/vector_tools.templates.h @@ -60,7 +60,6 @@ #include #include #include -#include #include #include #include @@ -3367,77 +3366,6 @@ namespace VectorTools return std::numeric_limits::min(); } - - - template - void - invert_mass_matrix( - const SparseMatrix & mass_matrix, - const FilteredMatrix> &filtered_mass_matrix, - FilteredMatrix> & filtered_preconditioner, - const Vector & rhs, - Vector & boundary_projection) - { - // Allow for a maximum of 5*n steps to reduce the residual by 10^-12. n - // steps may not be sufficient, since roundoff errors may accumulate for - // badly conditioned matrices. This behavior can be observed, e.g. for - // FE_Q_Hierarchical for degree higher than three. - ReductionControl control(5 * rhs.size(), 0., 1.e-12, false, false); - GrowingVectorMemory> memory; - SolverCG> cg(control, memory); - - PreconditionSSOR> prec; - prec.initialize(mass_matrix, 1.2); - filtered_preconditioner.initialize(prec, true); - // solve - cg.solve(filtered_mass_matrix, - boundary_projection, - rhs, - filtered_preconditioner); - filtered_preconditioner.apply_constraints(boundary_projection, true); - filtered_preconditioner.clear(); - } - - - - template - void - invert_mass_matrix( - const SparseMatrix & mass_matrix, - const FilteredMatrix> &filtered_mass_matrix, - FilteredMatrix> & filtered_preconditioner, - const Vector> & rhs, - Vector> & boundary_projection) - { - auto solve_for_one_component = [&](const bool real_part) { - // copy the real or imaginary part out of the rhs vector - Vector rhs_part(rhs.size()); - for (unsigned int i = 0; i < rhs.size(); ++i) - rhs_part(i) = (real_part ? rhs(i).real() : rhs(i).imag()); - - // then solve the linear system for this part - Vector boundary_projection_part(boundary_projection.size()); - invert_mass_matrix(mass_matrix, - filtered_mass_matrix, - filtered_preconditioner, - rhs_part, - boundary_projection_part); - - // finally copy the real or imaginary part of the - // solution back into the global solution vector - for (unsigned int i = 0; i < boundary_projection.size(); ++i) - if (real_part == true) - boundary_projection(i).real(boundary_projection_part(i)); - else - boundary_projection(i).imag(boundary_projection_part(i)); - }; - - // solve for real and imaginary parts of the solution separately - solve_for_one_component(true); - solve_for_one_component(false); - } - - template class DoFHandlerType, @@ -3572,42 +3500,6 @@ namespace VectorTools static_cast *>(nullptr), component_mapping); - // For certain weird elements, - // there might be degrees of - // freedom on the boundary, but - // their shape functions do not - // have support there. Let's - // eliminate them here. - - // The Bogner-Fox-Schmidt element - // is an example for those. - - // TODO: Maybe we should figure out if the element really needs this - - FilteredMatrix> filtered_mass_matrix(mass_matrix, true); - FilteredMatrix> filtered_precondition; - std::vector excluded_dofs(mass_matrix.m(), false); - - // we assemble mass matrix with unit weight, - // thus it will be real-valued irrespectively of the underlying algebra - // with positive elements on diagonal. - // Thus in order to extend this filtering to complex-algebra simply take - // the real-part of element. - number max_element = 0.; - for (unsigned int i = 0; i < mass_matrix.m(); ++i) - if (real_part_bigger_than(mass_matrix.diag_element(i), max_element)) - max_element = mass_matrix.diag_element(i); - - for (unsigned int i = 0; i < mass_matrix.m(); ++i) - if (real_part_bigger_than(1.e-8 * max_element, - mass_matrix.diag_element(i))) - { - filtered_mass_matrix.add_constraint(i, 0.); - filtered_precondition.add_constraint(i, 0.); - mass_matrix.diag_element(i) = 1.; - excluded_dofs[i] = true; - } - Vector boundary_projection(rhs.size()); // cannot reduce residual in a useful way if we are close to the square @@ -3616,16 +3508,11 @@ namespace VectorTools boundary_projection = 0; else { - invert_mass_matrix(mass_matrix, - filtered_mass_matrix, - filtered_precondition, - rhs, - boundary_projection); + invert_mass_matrix(mass_matrix, rhs, boundary_projection); } // fill in boundary values for (unsigned int i = 0; i < dof_to_boundary_mapping.size(); ++i) - if (dof_to_boundary_mapping[i] != numbers::invalid_dof_index && - !excluded_dofs[dof_to_boundary_mapping[i]]) + if (dof_to_boundary_mapping[i] != numbers::invalid_dof_index) { AssertIsFinite(boundary_projection(dof_to_boundary_mapping[i])); diff --git a/source/grid/grid_tools.cc b/source/grid/grid_tools.cc index 4a4cee06be..2fd1c5b98c 100644 --- a/source/grid/grid_tools.cc +++ b/source/grid/grid_tools.cc @@ -41,8 +41,8 @@ #include #include +#include #include -#include #include #include #include @@ -1104,15 +1104,15 @@ namespace GridTools * in order to allow parallel execution. */ void - laplace_solve(const SparseMatrix & S, - const std::map &fixed_dofs, - Vector & u) + laplace_solve(const SparseMatrix & S, + const AffineConstraints &constraints, + Vector & u) { - const unsigned int n_dofs = S.n(); - FilteredMatrix> SF(S); + const unsigned int n_dofs = S.n(); + const auto op = linear_operator(S); + const auto SF = constrained_linear_operator(constraints, op); PreconditionJacobi> prec; prec.initialize(S, 1.2); - FilteredMatrix> PF(prec); SolverControl control(n_dofs, 1.e-10, false, false); GrowingVectorMemory> mem; @@ -1120,9 +1120,11 @@ namespace GridTools Vector f(n_dofs); - SF.add_constraints(fixed_dofs); - SF.apply_constraints(f, true); - solver.solve(SF, u, f, PF); + const auto constrained_rhs = + constrained_right_hand_side(constraints, op, f); + solver.solve(SF, u, constrained_rhs, prec); + + constraints.distribute(u); } } // namespace @@ -1171,7 +1173,7 @@ namespace GridTools StaticMappingQ1::mapping, dof_handler, quadrature, S, coefficient); // set up the boundary values for the laplace problem - std::map fixed_dofs[dim]; + std::array, dim> constraints; typename std::map>::const_iterator map_end = new_points.end(); @@ -1190,14 +1192,20 @@ namespace GridTools if (map_iter != map_end) for (unsigned int i = 0; i < dim; ++i) - fixed_dofs[i].insert(std::pair( - cell->vertex_dof_index(vertex_no, 0), - (solve_for_absolute_positions ? - map_iter->second(i) : - map_iter->second(i) - vertex_point[i]))); + { + constraints[i].add_line(cell->vertex_dof_index(vertex_no, 0)); + constraints[i].set_inhomogeneity( + cell->vertex_dof_index(vertex_no, 0), + (solve_for_absolute_positions ? + map_iter->second(i) : + map_iter->second(i) - vertex_point[i])); + } } } + for (unsigned int i = 0; i < dim; ++i) + constraints[i].close(); + // solve the dim problems with different right hand sides. Vector us[dim]; for (unsigned int i = 0; i < dim; ++i) @@ -1206,7 +1214,7 @@ namespace GridTools // solve linear systems in parallel Threads::TaskGroup<> tasks; for (unsigned int i = 0; i < dim; ++i) - tasks += Threads::new_task(&laplace_solve, S, fixed_dofs[i], us[i]); + tasks += Threads::new_task(&laplace_solve, S, constraints[i], us[i]); tasks.join_all(); // change the coordinates of the points of the triangulation diff --git a/tests/data_out/data_out_curved_cells.cc b/tests/data_out/data_out_curved_cells.cc index 1e590aad09..f40bf99020 100644 --- a/tests/data_out/data_out_curved_cells.cc +++ b/tests/data_out/data_out_curved_cells.cc @@ -38,9 +38,11 @@ #include #include -#include +#include #include #include +#include +#include #include #include #include @@ -56,17 +58,14 @@ #include "../tests.h" -// this is copied from GridGenerator void -laplace_solve(const SparseMatrix & S, - const std::map &m, - Vector & u) +laplace_solve(const SparseMatrix & S, + const AffineConstraints &constraints, + Vector & u) { - const unsigned int n_dofs = S.n(); - FilteredMatrix> SF(S); - PreconditionJacobi> prec; - prec.initialize(S, 1.2); - FilteredMatrix> PF(prec); + const unsigned int n_dofs = S.n(); + const auto op = linear_operator(S); + const auto constrained_op = constrained_linear_operator(constraints, op); SolverControl control(10000, 1.e-10, false, false); GrowingVectorMemory> mem; @@ -74,13 +73,13 @@ laplace_solve(const SparseMatrix & S, Vector f(n_dofs); - SF.add_constraints(m); - SF.apply_constraints(f, true); - solver.solve(SF, u, f, PF); + const auto constrained_rhs = constrained_right_hand_side(constraints, op, f); + solver.solve(constrained_op, u, constrained_rhs, PreconditionIdentity()); + constraints.distribute(u); } -// create a rinf grid and compute a MappingQEuler to represent the inner +// create a ring grid and compute a MappingQEuler to represent the inner // boundary void curved_grid(std::ostream &out) @@ -151,27 +150,25 @@ curved_grid(std::ostream &out) MatrixCreator::create_laplace_matrix(mapping_q1, dof_handler, quadrature, S); // set up the boundary values for // the laplace problem - std::vector> m(2); - // fill these maps: on the inner boundary try + std::vector> m(2); + // fill these constraints: on the inner boundary try // to approximate a circle, on the outer // boundary use straight lines - DoFHandler<2>::cell_iterator cell = dof_handler.begin_active(), - endc = dof_handler.end(); DoFHandler<2>::face_iterator face; - for (; cell != endc; ++cell) + for (const auto &cell : dof_handler.active_cell_iterators()) { // fix all vertices for (const unsigned int vertex_no : GeometryInfo<2>::vertex_indices()) { for (unsigned int i = 0; i < 2; ++i) - m[i][cell->vertex_dof_index(vertex_no, 0)] = 0.; + m[i].add_line(cell->vertex_dof_index(vertex_no, 0)); } if (cell->at_boundary()) for (const unsigned int face_no : GeometryInfo<2>::face_indices()) { face = cell->face(face_no); - // insert a modifiued value for + // insert a modified value for // the middle point of boundary // faces if (face->at_boundary()) @@ -179,16 +176,24 @@ curved_grid(std::ostream &out) const double eps = 1e-4; if (std::fabs(face->vertex(1).norm() - r_i) < eps) for (unsigned int i = 0; i < 2; ++i) - m[i][face->dof_index(0)] = - (face->center() * (r_i / face->center().norm() - 1))(i); + { + m[i].add_line(face->dof_index(0)); + m[i].set_inhomogeneity(face->dof_index(0), + (face->center() * + (r_i / face->center().norm() - + 1))(i)); + } else if (std::fabs(face->vertex(1).norm() - r_a) < eps) for (unsigned int i = 0; i < 2; ++i) - m[i][face->dof_index(0)] = 0.; + m[i].add_line(face->dof_index(0)); else Assert(false, ExcInternalError()); } } } + m[0].close(); + m[1].close(); + // solve the 2 problems with // different right hand sides. Vector us[2]; diff --git a/tests/grid/grid_transform_02.output b/tests/grid/grid_transform_02.output index b35bbeabbd..88d3c9a87e 100644 --- a/tests/grid/grid_transform_02.output +++ b/tests/grid/grid_transform_02.output @@ -2,7 +2,7 @@ %%Title: deal.II Output %%Creator: the deal.II library -%%BoundingBox: 0 0 301 300 +%%BoundingBox: 0 0 301 301 /m {moveto} bind def /x {lineto stroke} bind def /b {0 0 0 setrgbcolor} def @@ -10,37 +10,37 @@ %%EndProlog 0.5 setlinewidth -b 1.86517e-13 0 m 0 37.5 x -b 37.5 2.66454e-14 m 38.4292 37.5 x -b 1.86517e-13 0 m 37.5 2.66454e-14 x +b 0 0 m 0 37.5 x +b 37.5 0 m 38.4292 37.5 x +b 0 0 m 37.5 0 x b 0 37.5 m 38.4292 37.5 x -b 37.5 2.66454e-14 m 38.4292 37.5 x -b 75 2.66454e-14 m 76.8146 37.5 x -b 37.5 2.66454e-14 m 75 2.66454e-14 x +b 37.5 0 m 38.4292 37.5 x +b 75 0 m 76.8146 37.5 x +b 37.5 0 m 75 0 x b 38.4292 37.5 m 76.8146 37.5 x -b 75 2.66454e-14 m 76.8146 37.5 x -b 112.5 2.66454e-14 m 115.008 37.5 x -b 75 2.66454e-14 m 112.5 2.66454e-14 x +b 75 0 m 76.8146 37.5 x +b 112.5 0 m 115.008 37.5 x +b 75 0 m 112.5 0 x b 76.8146 37.5 m 115.008 37.5 x -b 112.5 2.66454e-14 m 115.008 37.5 x -b 150 2.66454e-14 m 152.77 37.5 x -b 112.5 2.66454e-14 m 150 2.66454e-14 x +b 112.5 0 m 115.008 37.5 x +b 150 0 m 152.77 37.5 x +b 112.5 0 m 150 0 x b 115.008 37.5 m 152.77 37.5 x -b 150 2.66454e-14 m 152.77 37.5 x -b 187.5 2.66454e-14 m 190.008 37.5 x -b 150 2.66454e-14 m 187.5 2.66454e-14 x +b 150 0 m 152.77 37.5 x +b 187.5 0 m 190.008 37.5 x +b 150 0 m 187.5 0 x b 152.77 37.5 m 190.008 37.5 x -b 187.5 2.66454e-14 m 190.008 37.5 x -b 225 2.66454e-14 m 226.815 37.5 x -b 187.5 2.66454e-14 m 225 2.66454e-14 x +b 187.5 0 m 190.008 37.5 x +b 225 0 m 226.815 37.5 x +b 187.5 0 m 225 0 x b 190.008 37.5 m 226.815 37.5 x -b 225 2.66454e-14 m 226.815 37.5 x -b 262.5 2.66454e-14 m 263.429 37.5 x -b 225 2.66454e-14 m 262.5 2.66454e-14 x +b 225 0 m 226.815 37.5 x +b 262.5 0 m 263.429 37.5 x +b 225 0 m 262.5 0 x b 226.815 37.5 m 263.429 37.5 x -b 262.5 2.66454e-14 m 263.429 37.5 x +b 262.5 0 m 263.429 37.5 x b 300 0 m 300 37.5 x -b 262.5 2.66454e-14 m 300 0 x +b 262.5 0 m 300 0 x b 263.429 37.5 m 300 37.5 x b 0 37.5 m 0 75 x b 38.4292 37.5 m 39.3146 75 x @@ -234,10 +234,10 @@ b 264.315 225 m 263.429 262.5 x b 300 225 m 300 262.5 x b 264.315 225 m 300 225 x b 263.429 262.5 m 300 262.5 x -b 0 262.5 m 1.59872e-13 300 x +b 0 262.5 m 0 300 x b 38.4292 262.5 m 37.5 300 x b 0 262.5 m 38.4292 262.5 x -b 1.59872e-13 300 m 37.5 300 x +b 0 300 m 37.5 300 x b 38.4292 262.5 m 37.5 300 x b 76.8146 262.5 m 75 300 x b 38.4292 262.5 m 76.8146 262.5 x diff --git a/tests/lac/filtered_matrix.cc b/tests/lac/filtered_matrix.cc deleted file mode 100644 index 1b91f02736..0000000000 --- a/tests/lac/filtered_matrix.cc +++ /dev/null @@ -1,99 +0,0 @@ -// --------------------------------------------------------------------- -// -// Copyright (C) 2007 - 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. -// -// --------------------------------------------------------------------- - - -// Test properties of FilteredMatrix and iterators - -#include -#include -#include -#include - -#include "../tests.h" - - -template -void -test(const FilteredMatrix &M) -{ - deallog << "Iterator"; - - unsigned int max = 0; - for (typename FilteredMatrix::const_iterator i = M.begin(); - i != M.end(); - ++i) - { - Assert(i->row() == i->column(), ExcInternalError()); - deallog << ' ' << i->row() << ':' << i->value(); - max = i->row(); - } - VectorType v(max + 1); - VectorType w(max + 1); - - for (unsigned int i = 0; i < v.size(); ++i) - v(i) = 31 + i; - - deallog << std::endl << "vmult "; - - w = 0.; - M.vmult(w, v); - for (unsigned int i = 0; i < v.size(); ++i) - deallog << ' ' << w(i); - - deallog << std::endl << "Tvmult"; - - w = 0.; - M.Tvmult(w, v); - for (unsigned int i = 0; i < v.size(); ++i) - deallog << ' ' << w(i); - - deallog << std::endl << "vmult_add"; - - M.vmult_add(w, v); - for (unsigned int i = 0; i < v.size(); ++i) - deallog << ' ' << w(i); - - deallog << std::endl << "boundary"; - - M.apply_constraints(w, true); - for (unsigned int i = 0; i < v.size(); ++i) - deallog << ' ' << w(i); - - deallog << std::endl << "boundary"; - - M.apply_constraints(w, false); - for (unsigned int i = 0; i < v.size(); ++i) - deallog << ' ' << w(i); - - deallog << std::endl; -} - -int -main() -{ - initlog(); - - PreconditionIdentity identity; - FilteredMatrix> f; - f.initialize(identity, false); - test(f); - - for (unsigned int i = 0; i < 5; ++i) - f.add_constraint(i * i, i / 2.); - test(f); - - f.initialize(identity, true); - test(f); -} diff --git a/tests/lac/filtered_matrix.output b/tests/lac/filtered_matrix.output deleted file mode 100644 index c1fc45f14f..0000000000 --- a/tests/lac/filtered_matrix.output +++ /dev/null @@ -1,19 +0,0 @@ - -DEAL::Iterator -DEAL::vmult 31.0000 -DEAL::Tvmult 31.0000 -DEAL::vmult_add 62.0000 -DEAL::boundary 62.0000 -DEAL::boundary 62.0000 -DEAL::Iterator 0:0.00000 1:0.500000 4:1.00000 9:1.50000 16:2.00000 -DEAL::vmult 31.0000 32.0000 33.0000 34.0000 35.0000 36.0000 37.0000 38.0000 39.0000 40.0000 41.0000 42.0000 43.0000 44.0000 45.0000 46.0000 47.0000 -DEAL::Tvmult 31.0000 32.0000 33.0000 34.0000 35.0000 36.0000 37.0000 38.0000 39.0000 40.0000 41.0000 42.0000 43.0000 44.0000 45.0000 46.0000 47.0000 -DEAL::vmult_add 31.0000 32.0000 66.0000 68.0000 35.0000 72.0000 74.0000 76.0000 78.0000 40.0000 82.0000 84.0000 86.0000 88.0000 90.0000 92.0000 47.0000 -DEAL::boundary 0.00000 0.500000 66.0000 68.0000 1.00000 72.0000 74.0000 76.0000 78.0000 1.50000 82.0000 84.0000 86.0000 88.0000 90.0000 92.0000 2.00000 -DEAL::boundary 0.00000 0.500000 66.0000 68.0000 1.00000 72.0000 74.0000 76.0000 78.0000 1.50000 82.0000 84.0000 86.0000 88.0000 90.0000 92.0000 2.00000 -DEAL::Iterator 0:0.00000 1:0.500000 4:1.00000 9:1.50000 16:2.00000 -DEAL::vmult 31.0000 32.0000 33.0000 34.0000 35.0000 36.0000 37.0000 38.0000 39.0000 40.0000 41.0000 42.0000 43.0000 44.0000 45.0000 46.0000 47.0000 -DEAL::Tvmult 31.0000 32.0000 33.0000 34.0000 35.0000 36.0000 37.0000 38.0000 39.0000 40.0000 41.0000 42.0000 43.0000 44.0000 45.0000 46.0000 47.0000 -DEAL::vmult_add 31.0000 32.0000 66.0000 68.0000 35.0000 72.0000 74.0000 76.0000 78.0000 40.0000 82.0000 84.0000 86.0000 88.0000 90.0000 92.0000 47.0000 -DEAL::boundary 0.00000 0.500000 66.0000 68.0000 1.00000 72.0000 74.0000 76.0000 78.0000 1.50000 82.0000 84.0000 86.0000 88.0000 90.0000 92.0000 2.00000 -DEAL::boundary 0.00000 0.500000 66.0000 68.0000 1.00000 72.0000 74.0000 76.0000 78.0000 1.50000 82.0000 84.0000 86.0000 88.0000 90.0000 92.0000 2.00000 diff --git a/tests/lac/filtered_matrix_01.cc b/tests/lac/filtered_matrix_01.cc deleted file mode 100644 index b3958eedcb..0000000000 --- a/tests/lac/filtered_matrix_01.cc +++ /dev/null @@ -1,66 +0,0 @@ -// --------------------------------------------------------------------- -// -// Copyright (C) 2007 - 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 "../tests.h" - -template -void -checkVmult(FullMatrix &A, - Vector & V, - bool expect_constrained_source = false) -{ - deallog << "vmult" << std::endl; - - FilteredMatrix> F; - F.initialize(A, expect_constrained_source); - F.add_constraint(0, 1); - - Vector O(A.m()); - - F.vmult(O, V); - - for (unsigned int i = 0; i < O.size(); ++i) - deallog << O(i) << '\t'; - deallog << std::endl; -} - -int -main() -{ - std::ofstream logfile("output"); - deallog << std::fixed; - deallog << std::setprecision(4); - deallog.attach(logfile); - - const double Adata[] = {1, 2, 3, 4, 5, 6, 7, 8, 9}; - - FullMatrix A(3, 3); - - A.fill(Adata); - - Vector V1(3); - Vector V2(3); - - V1(0) = V2(0) = 1; - V1(1) = V2(1) = 2; - V1(2) = V2(2) = 3; - - checkVmult(A, V1, false); - checkVmult(A, V2, true); -} diff --git a/tests/lac/filtered_matrix_01.output b/tests/lac/filtered_matrix_01.output deleted file mode 100644 index e893ffef71..0000000000 --- a/tests/lac/filtered_matrix_01.output +++ /dev/null @@ -1,5 +0,0 @@ - -DEAL::vmult -DEAL::1.0000 28.0000 43.0000 -DEAL::vmult -DEAL::1.0000 32.0000 50.0000 diff --git a/tests/lac/filtered_matrix_02.cc b/tests/lac/filtered_matrix_02.cc deleted file mode 100644 index 86e6869ed4..0000000000 --- a/tests/lac/filtered_matrix_02.cc +++ /dev/null @@ -1,66 +0,0 @@ -// --------------------------------------------------------------------- -// -// Copyright (C) 2007 - 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 "../tests.h" - -template -void -checkTvmult(FullMatrix &A, - Vector & V, - bool expect_constrained_source = false) -{ - deallog << "Tvmult" << std::endl; - - FilteredMatrix> F; - F.initialize(A, expect_constrained_source); - F.add_constraint(0, 1); - - Vector O(A.n()); - - F.Tvmult(O, V); - - for (unsigned int i = 0; i < O.size(); ++i) - deallog << O(i) << '\t'; - deallog << std::endl; -} - -int -main() -{ - std::ofstream logfile("output"); - deallog << std::fixed; - deallog << std::setprecision(4); - deallog.attach(logfile); - - const double Adata[] = {1, 2, 3, 4, 5, 6, 7, 8, 9}; - - FullMatrix A(3, 3); - - A.fill(Adata); - - Vector V1(3); - Vector V2(3); - - V1(0) = V2(0) = 1; - V1(1) = V2(1) = 2; - V1(2) = V2(2) = 3; - - checkTvmult(A, V1, false); - checkTvmult(A, V2, true); -} diff --git a/tests/lac/filtered_matrix_02.output b/tests/lac/filtered_matrix_02.output deleted file mode 100644 index 5c3f2bee78..0000000000 --- a/tests/lac/filtered_matrix_02.output +++ /dev/null @@ -1,5 +0,0 @@ - -DEAL::Tvmult -DEAL::1.0000 34.0000 39.0000 -DEAL::Tvmult -DEAL::1.0000 36.0000 42.0000 diff --git a/tests/lac/filtered_matrix_03.cc b/tests/lac/filtered_matrix_03.cc deleted file mode 100644 index 1c334e1790..0000000000 --- a/tests/lac/filtered_matrix_03.cc +++ /dev/null @@ -1,70 +0,0 @@ -// --------------------------------------------------------------------- -// -// Copyright (C) 2007 - 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 "../tests.h" - -template -void -checkVmult_Add(FullMatrix &A, - Vector & V, - bool expect_constrained_source = false) -{ - deallog << "vmult_add" << std::endl; - - FilteredMatrix> F; - F.initialize(A, expect_constrained_source); - F.add_constraint(0, 1); - - Vector O(A.m()); - for (unsigned int i = 0; i < O.size(); ++i) - { - O(i) = 1; - } - - F.vmult_add(O, V); - - for (unsigned int i = 0; i < O.size(); ++i) - deallog << O(i) << '\t'; - deallog << std::endl; -} - -int -main() -{ - std::ofstream logfile("output"); - deallog << std::fixed; - deallog << std::setprecision(4); - deallog.attach(logfile); - - const double Adata[] = {1, 2, 3, 4, 5, 6, 7, 8, 9}; - - FullMatrix A(3, 3); - - A.fill(Adata); - - Vector V1(3); - Vector V2(3); - - V1(0) = V2(0) = 1; - V1(1) = V2(1) = 2; - V1(2) = V2(2) = 3; - - checkVmult_Add(A, V1, false); - checkVmult_Add(A, V2, true); -} diff --git a/tests/lac/filtered_matrix_03.output b/tests/lac/filtered_matrix_03.output deleted file mode 100644 index 9ea8103a31..0000000000 --- a/tests/lac/filtered_matrix_03.output +++ /dev/null @@ -1,5 +0,0 @@ - -DEAL::vmult_add -DEAL::1.0000 29.0000 44.0000 -DEAL::vmult_add -DEAL::1.0000 33.0000 51.0000 diff --git a/tests/lac/filtered_matrix_04.cc b/tests/lac/filtered_matrix_04.cc deleted file mode 100644 index 6ad73c901d..0000000000 --- a/tests/lac/filtered_matrix_04.cc +++ /dev/null @@ -1,70 +0,0 @@ -// --------------------------------------------------------------------- -// -// Copyright (C) 2007 - 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 "../tests.h" - -template -void -checkTvmult_Add(FullMatrix &A, - Vector & V, - bool expect_constrained_source = false) -{ - deallog << "Tvmult_add" << std::endl; - - FilteredMatrix> F; - F.initialize(A, expect_constrained_source); - F.add_constraint(0, 1); - - Vector O(A.n()); - for (unsigned int i = 0; i < O.size(); ++i) - { - O(i) = 1; - } - - F.Tvmult_add(O, V); - - for (unsigned int i = 0; i < O.size(); ++i) - deallog << O(i) << '\t'; - deallog << std::endl; -} - -int -main() -{ - std::ofstream logfile("output"); - deallog << std::fixed; - deallog << std::setprecision(4); - deallog.attach(logfile); - - const double Adata[] = {1, 2, 3, 4, 5, 6, 7, 8, 9}; - - FullMatrix A(3, 3); - - A.fill(Adata); - - Vector V1(3); - Vector V2(3); - - V1(0) = V2(0) = 1; - V1(1) = V2(1) = 2; - V1(2) = V2(2) = 3; - - checkTvmult_Add(A, V1, false); - checkTvmult_Add(A, V2, true); -} diff --git a/tests/lac/filtered_matrix_04.output b/tests/lac/filtered_matrix_04.output deleted file mode 100644 index b455ab3529..0000000000 --- a/tests/lac/filtered_matrix_04.output +++ /dev/null @@ -1,5 +0,0 @@ - -DEAL::Tvmult_add -DEAL::1.0000 35.0000 40.0000 -DEAL::Tvmult_add -DEAL::1.0000 37.0000 43.0000 diff --git a/tests/lac/filtered_matrix_05.cc b/tests/lac/filtered_matrix_05.cc deleted file mode 100644 index 4e99667888..0000000000 --- a/tests/lac/filtered_matrix_05.cc +++ /dev/null @@ -1,64 +0,0 @@ -// --------------------------------------------------------------------- -// -// Copyright (C) 2007 - 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 "../tests.h" - -template -void -checkApply_Constraints(FullMatrix &A, - Vector & V, - bool matrix_is_symmetric = false) -{ - deallog << "apply_constraints" << std::endl; - - FilteredMatrix> F; - F.initialize(A); - F.add_constraint(0, 1); - - F.apply_constraints(V, matrix_is_symmetric); - - for (unsigned int i = 0; i < V.size(); ++i) - deallog << V(i) << '\t'; - deallog << std::endl; -} - -int -main() -{ - std::ofstream logfile("output"); - deallog << std::fixed; - deallog << std::setprecision(4); - deallog.attach(logfile); - - const double Adata[] = {1, 2, 3, 4, 5, 6, 7, 8, 9}; - - FullMatrix A(3, 3); - - A.fill(Adata); - - Vector V1(3); - Vector V2(3); - - V1(0) = V2(0) = 1; - V1(1) = V2(1) = 2; - V1(2) = V2(2) = 3; - - checkApply_Constraints(A, V1, false); - checkApply_Constraints(A, V2, true); -} diff --git a/tests/lac/filtered_matrix_05.output b/tests/lac/filtered_matrix_05.output deleted file mode 100644 index 664c0a38fb..0000000000 --- a/tests/lac/filtered_matrix_05.output +++ /dev/null @@ -1,5 +0,0 @@ - -DEAL::apply_constraints -DEAL::1.0000 -2.0000 -4.0000 -DEAL::apply_constraints -DEAL::1.0000 -2.0000 -4.0000 diff --git a/tests/lac/filtered_matrix_06.cc b/tests/lac/filtered_matrix_06.cc deleted file mode 100644 index 4df8737b31..0000000000 --- a/tests/lac/filtered_matrix_06.cc +++ /dev/null @@ -1,224 +0,0 @@ -// --------------------------------------------------------------------- -// -// Copyright (C) 2013 - 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 - -#include -#include -#include - -#include -#include -#include -#include -#include -#include - -#include -#include - -#include -#include - -#include "../tests.h" - - - -void -solve_filtered(std::map &bv, - SparseMatrix & A, - Vector & u, - Vector & f) -{ - FilteredMatrix> A1(A); - A1.add_constraints(bv); - - SolverControl control(1000, 1.e-10, false, false); - PrimitiveVectorMemory> mem; - SolverCG> solver(control, mem); - PreconditionJacobi> prec; - FilteredMatrix> fprec; - prec.initialize(A, 1.2); - fprec.initialize(prec); - - Vector f1(f.size()); - f1 = f; - A1.apply_constraints(f1, true); - - solver.solve(A1, u, f1, fprec); - - for (std::map::const_iterator i = bv.begin(); - i != bv.end(); - ++i) - AssertThrow(std::fabs(u(i->first) - i->second) < 1e-8, ExcInternalError()); -} - - - -template -void -solve_eliminated(std::map &bv, - SparseMatrix & A, - Vector & u, - Vector & f) -{ - MatrixTools::apply_boundary_values(bv, A, u, f); - - SolverControl control(1000, 1.e-10, false, false); - PrimitiveVectorMemory> mem; - SolverCG> solver(control, mem); - PreconditionJacobi<> prec; - prec.initialize(A, 1.2); - - solver.solve(A, u, f, prec); -} - - - -template -void -check() -{ - Triangulation tr; - - Functions::CosineFunction cosine; - - if (dim == 2) - GridGenerator::hyper_ball(tr, Point(), 1); - else - GridGenerator::hyper_cube(tr, -1, 1); - tr.reset_manifold(0); - - tr.refine_global(5 - dim); - - MappingQ mapping(2); - FE_Q element(1); - QGauss quadrature(4); - - DoFHandler dof(tr); - dof.distribute_dofs(element); - - FEValues fe(mapping, - element, - quadrature, - update_values | update_gradients | update_quadrature_points | - update_JxW_values); - - std::vector global_dofs(element.dofs_per_cell); - std::vector function(quadrature.size()); - - Vector f(dof.n_dofs()); - - SparsityPattern A_pattern(dof.n_dofs(), - dof.n_dofs(), - dof.max_couplings_between_dofs()); - DoFTools::make_sparsity_pattern(dof, A_pattern); - A_pattern.compress(); - - SparseMatrix A(A_pattern); - - typename DoFHandler::active_cell_iterator cell = dof.begin_active(); - const typename DoFHandler::cell_iterator end = dof.end(); - - for (; cell != end; ++cell) - { - fe.reinit(cell); - cell->get_dof_indices(global_dofs); - cosine.value_list(fe.get_quadrature_points(), function); - - for (unsigned int k = 0; k < quadrature.size(); ++k) - { - double dx = fe.JxW(k); - - for (unsigned int i = 0; i < element.dofs_per_cell; ++i) - { - const double v = fe.shape_value(i, k); - const Tensor<1, dim> grad_v = fe.shape_grad(i, k); - - double rhs = dx * v * (function[k]); - - f(global_dofs[i]) += rhs; - for (unsigned int j = 0; j < element.dofs_per_cell; ++j) - { - const Tensor<1, dim> grad_u = fe.shape_grad(j, k); - double el = dx * (grad_u * grad_v); - A.add(global_dofs[i], global_dofs[j], el); - } - } - } - } - - // interpolate boundary values - std::map bv; - VectorTools::interpolate_boundary_values( - mapping, dof, 0, cosine, bv, std::vector()); - // the cosine has too many zero - // values on the boundary of the - // domain, so reset the elements to - // some other value - for (typename std::map::iterator i = - bv.begin(); - i != bv.end(); - ++i) - i->second = std::sin(i->second + 0.5) + 1.0; - - // first solve filtered. this does - // not change the matrix - Vector u_filtered(dof.n_dofs()); - solve_filtered(bv, A, u_filtered, f); - - // then solve by eliminating in the - // matrix. since this changes the - // matrix, this call must come - // second - Vector u_eliminated(dof.n_dofs()); - solve_eliminated(bv, A, u_eliminated, f); - - // output and check - for (unsigned int i = 0; i < dof.n_dofs(); ++i) - { - deallog << u_filtered(i) << std::endl; - Assert(std::fabs(u_filtered(i) - u_eliminated(i)) < 1e-8, - ExcInternalError()); - }; -} - - -int -main() -{ - initlog(); - deallog << std::setprecision(2) << std::fixed; - - deallog.push("1d"); - check<1>(); - deallog.pop(); - deallog.push("2d"); - check<2>(); - deallog.pop(); - deallog.push("3d"); - check<3>(); - deallog.pop(); -} diff --git a/tests/lac/filtered_matrix_06.output b/tests/lac/filtered_matrix_06.output deleted file mode 100644 index 20dfc60ff1..0000000000 --- a/tests/lac/filtered_matrix_06.output +++ /dev/null @@ -1,480 +0,0 @@ - -DEAL:1d::1.48 -DEAL:1d::1.64 -DEAL:1d::1.79 -DEAL:1d::1.94 -DEAL:1d::2.08 -DEAL:1d::2.21 -DEAL:1d::2.33 -DEAL:1d::2.43 -DEAL:1d::2.52 -DEAL:1d::2.59 -DEAL:1d::2.65 -DEAL:1d::2.69 -DEAL:1d::2.72 -DEAL:1d::2.74 -DEAL:1d::2.75 -DEAL:1d::2.75 -DEAL:1d::2.75 -DEAL:2d::1.64 -DEAL:2d::1.72 -DEAL:2d::1.68 -DEAL:2d::1.74 -DEAL:2d::1.77 -DEAL:2d::1.78 -DEAL:2d::1.72 -DEAL:2d::1.77 -DEAL:2d::1.80 -DEAL:2d::1.80 -DEAL:2d::1.81 -DEAL:2d::1.81 -DEAL:2d::1.82 -DEAL:2d::1.82 -DEAL:2d::1.83 -DEAL:2d::1.75 -DEAL:2d::1.79 -DEAL:2d::1.81 -DEAL:2d::1.78 -DEAL:2d::1.81 -DEAL:2d::1.83 -DEAL:2d::1.83 -DEAL:2d::1.84 -DEAL:2d::1.84 -DEAL:2d::1.85 -DEAL:2d::1.80 -DEAL:2d::1.81 -DEAL:2d::1.77 -DEAL:2d::1.78 -DEAL:2d::1.82 -DEAL:2d::1.80 -DEAL:2d::1.72 -DEAL:2d::1.74 -DEAL:2d::1.64 -DEAL:2d::1.68 -DEAL:2d::1.77 -DEAL:2d::1.72 -DEAL:2d::1.83 -DEAL:2d::1.81 -DEAL:2d::1.84 -DEAL:2d::1.83 -DEAL:2d::1.79 -DEAL:2d::1.75 -DEAL:2d::1.81 -DEAL:2d::1.78 -DEAL:2d::1.80 -DEAL:2d::1.83 -DEAL:2d::1.84 -DEAL:2d::1.82 -DEAL:2d::1.84 -DEAL:2d::1.85 -DEAL:2d::1.85 -DEAL:2d::1.86 -DEAL:2d::1.86 -DEAL:2d::1.87 -DEAL:2d::1.84 -DEAL:2d::1.86 -DEAL:2d::1.87 -DEAL:2d::1.86 -DEAL:2d::1.87 -DEAL:2d::1.87 -DEAL:2d::1.87 -DEAL:2d::1.87 -DEAL:2d::1.88 -DEAL:2d::1.88 -DEAL:2d::1.85 -DEAL:2d::1.84 -DEAL:2d::1.86 -DEAL:2d::1.85 -DEAL:2d::1.83 -DEAL:2d::1.80 -DEAL:2d::1.84 -DEAL:2d::1.82 -DEAL:2d::1.87 -DEAL:2d::1.87 -DEAL:2d::1.88 -DEAL:2d::1.87 -DEAL:2d::1.86 -DEAL:2d::1.84 -DEAL:2d::1.87 -DEAL:2d::1.86 -DEAL:2d::1.72 -DEAL:2d::1.74 -DEAL:2d::1.77 -DEAL:2d::1.77 -DEAL:2d::1.78 -DEAL:2d::1.80 -DEAL:2d::1.79 -DEAL:2d::1.81 -DEAL:2d::1.81 -DEAL:2d::1.83 -DEAL:2d::1.80 -DEAL:2d::1.81 -DEAL:2d::1.82 -DEAL:2d::1.81 -DEAL:2d::1.82 -DEAL:2d::1.83 -DEAL:2d::1.83 -DEAL:2d::1.84 -DEAL:2d::1.84 -DEAL:2d::1.85 -DEAL:2d::1.83 -DEAL:2d::1.84 -DEAL:2d::1.84 -DEAL:2d::1.85 -DEAL:2d::1.86 -DEAL:2d::1.87 -DEAL:2d::1.87 -DEAL:2d::1.87 -DEAL:2d::1.85 -DEAL:2d::1.86 -DEAL:2d::1.86 -DEAL:2d::1.87 -DEAL:2d::1.87 -DEAL:2d::1.88 -DEAL:2d::1.87 -DEAL:2d::1.88 -DEAL:2d::1.80 -DEAL:2d::1.81 -DEAL:2d::1.82 -DEAL:2d::1.77 -DEAL:2d::1.78 -DEAL:2d::1.80 -DEAL:2d::1.83 -DEAL:2d::1.84 -DEAL:2d::1.81 -DEAL:2d::1.83 -DEAL:2d::1.72 -DEAL:2d::1.74 -DEAL:2d::1.77 -DEAL:2d::1.64 -DEAL:2d::1.68 -DEAL:2d::1.72 -DEAL:2d::1.79 -DEAL:2d::1.81 -DEAL:2d::1.75 -DEAL:2d::1.78 -DEAL:2d::1.85 -DEAL:2d::1.86 -DEAL:2d::1.84 -DEAL:2d::1.85 -DEAL:2d::1.87 -DEAL:2d::1.88 -DEAL:2d::1.87 -DEAL:2d::1.87 -DEAL:2d::1.83 -DEAL:2d::1.84 -DEAL:2d::1.80 -DEAL:2d::1.82 -DEAL:2d::1.86 -DEAL:2d::1.87 -DEAL:2d::1.84 -DEAL:2d::1.86 -DEAL:2d::1.88 -DEAL:2d::1.88 -DEAL:2d::1.88 -DEAL:2d::1.89 -DEAL:2d::1.89 -DEAL:2d::1.89 -DEAL:2d::1.89 -DEAL:2d::1.90 -DEAL:2d::1.89 -DEAL:2d::1.89 -DEAL:2d::1.89 -DEAL:2d::1.90 -DEAL:2d::1.90 -DEAL:2d::1.90 -DEAL:2d::1.90 -DEAL:2d::1.90 -DEAL:2d::1.89 -DEAL:2d::1.88 -DEAL:2d::1.89 -DEAL:2d::1.89 -DEAL:2d::1.88 -DEAL:2d::1.87 -DEAL:2d::1.88 -DEAL:2d::1.87 -DEAL:2d::1.90 -DEAL:2d::1.89 -DEAL:2d::1.90 -DEAL:2d::1.90 -DEAL:2d::1.89 -DEAL:2d::1.88 -DEAL:2d::1.89 -DEAL:2d::1.88 -DEAL:2d::1.89 -DEAL:2d::1.89 -DEAL:2d::1.88 -DEAL:2d::1.89 -DEAL:2d::1.90 -DEAL:2d::1.90 -DEAL:2d::1.89 -DEAL:2d::1.90 -DEAL:2d::1.88 -DEAL:2d::1.88 -DEAL:2d::1.87 -DEAL:2d::1.87 -DEAL:2d::1.89 -DEAL:2d::1.89 -DEAL:2d::1.88 -DEAL:2d::1.88 -DEAL:2d::1.90 -DEAL:2d::1.89 -DEAL:2d::1.89 -DEAL:2d::1.89 -DEAL:2d::1.89 -DEAL:2d::1.88 -DEAL:2d::1.88 -DEAL:2d::1.87 -DEAL:2d::1.89 -DEAL:2d::1.88 -DEAL:2d::1.88 -DEAL:2d::1.87 -DEAL:2d::1.88 -DEAL:2d::1.87 -DEAL:2d::1.87 -DEAL:2d::1.86 -DEAL:2d::1.72 -DEAL:2d::1.74 -DEAL:2d::1.77 -DEAL:2d::1.78 -DEAL:2d::1.77 -DEAL:2d::1.80 -DEAL:2d::1.80 -DEAL:2d::1.81 -DEAL:2d::1.81 -DEAL:2d::1.82 -DEAL:2d::1.82 -DEAL:2d::1.83 -DEAL:2d::1.79 -DEAL:2d::1.81 -DEAL:2d::1.81 -DEAL:2d::1.83 -DEAL:2d::1.83 -DEAL:2d::1.84 -DEAL:2d::1.84 -DEAL:2d::1.85 -DEAL:2d::1.80 -DEAL:2d::1.81 -DEAL:2d::1.77 -DEAL:2d::1.78 -DEAL:2d::1.82 -DEAL:2d::1.80 -DEAL:2d::1.72 -DEAL:2d::1.74 -DEAL:2d::1.64 -DEAL:2d::1.68 -DEAL:2d::1.77 -DEAL:2d::1.72 -DEAL:2d::1.83 -DEAL:2d::1.81 -DEAL:2d::1.84 -DEAL:2d::1.83 -DEAL:2d::1.79 -DEAL:2d::1.75 -DEAL:2d::1.81 -DEAL:2d::1.78 -DEAL:2d::1.83 -DEAL:2d::1.84 -DEAL:2d::1.84 -DEAL:2d::1.85 -DEAL:2d::1.85 -DEAL:2d::1.86 -DEAL:2d::1.86 -DEAL:2d::1.87 -DEAL:2d::1.86 -DEAL:2d::1.87 -DEAL:2d::1.87 -DEAL:2d::1.87 -DEAL:2d::1.85 -DEAL:2d::1.84 -DEAL:2d::1.86 -DEAL:2d::1.85 -DEAL:2d::1.83 -DEAL:2d::1.80 -DEAL:2d::1.84 -DEAL:2d::1.82 -DEAL:2d::1.87 -DEAL:2d::1.87 -DEAL:2d::1.86 -DEAL:2d::1.84 -DEAL:2d::1.72 -DEAL:2d::1.74 -DEAL:2d::1.77 -DEAL:2d::1.77 -DEAL:2d::1.78 -DEAL:2d::1.80 -DEAL:2d::1.79 -DEAL:2d::1.81 -DEAL:2d::1.81 -DEAL:2d::1.83 -DEAL:2d::1.80 -DEAL:2d::1.81 -DEAL:2d::1.82 -DEAL:2d::1.81 -DEAL:2d::1.82 -DEAL:2d::1.83 -DEAL:2d::1.83 -DEAL:2d::1.84 -DEAL:2d::1.84 -DEAL:2d::1.85 -DEAL:2d::1.83 -DEAL:2d::1.84 -DEAL:2d::1.84 -DEAL:2d::1.85 -DEAL:2d::1.86 -DEAL:2d::1.87 -DEAL:2d::1.85 -DEAL:2d::1.86 -DEAL:2d::1.86 -DEAL:2d::1.87 -DEAL:2d::1.87 -DEAL:2d::1.87 -DEAL:2d::1.80 -DEAL:2d::1.81 -DEAL:2d::1.82 -DEAL:2d::1.77 -DEAL:2d::1.78 -DEAL:2d::1.80 -DEAL:2d::1.83 -DEAL:2d::1.84 -DEAL:2d::1.81 -DEAL:2d::1.83 -DEAL:2d::1.72 -DEAL:2d::1.74 -DEAL:2d::1.77 -DEAL:2d::1.79 -DEAL:2d::1.81 -DEAL:2d::1.85 -DEAL:2d::1.86 -DEAL:2d::1.84 -DEAL:2d::1.85 -DEAL:2d::1.87 -DEAL:2d::1.87 -DEAL:2d::1.83 -DEAL:2d::1.84 -DEAL:2d::1.86 -DEAL:3d::1.48 -DEAL:3d::1.48 -DEAL:3d::1.48 -DEAL:3d::1.48 -DEAL:3d::1.48 -DEAL:3d::1.48 -DEAL:3d::1.48 -DEAL:3d::1.53 -DEAL:3d::1.48 -DEAL:3d::1.48 -DEAL:3d::1.48 -DEAL:3d::1.55 -DEAL:3d::1.48 -DEAL:3d::1.48 -DEAL:3d::1.48 -DEAL:3d::1.55 -DEAL:3d::1.48 -DEAL:3d::1.59 -DEAL:3d::1.48 -DEAL:3d::1.48 -DEAL:3d::1.48 -DEAL:3d::1.55 -DEAL:3d::1.48 -DEAL:3d::1.59 -DEAL:3d::1.48 -DEAL:3d::1.59 -DEAL:3d::1.63 -DEAL:3d::1.48 -DEAL:3d::1.48 -DEAL:3d::1.48 -DEAL:3d::1.53 -DEAL:3d::1.48 -DEAL:3d::1.48 -DEAL:3d::1.48 -DEAL:3d::1.48 -DEAL:3d::1.48 -DEAL:3d::1.55 -DEAL:3d::1.48 -DEAL:3d::1.48 -DEAL:3d::1.48 -DEAL:3d::1.55 -DEAL:3d::1.48 -DEAL:3d::1.48 -DEAL:3d::1.59 -DEAL:3d::1.48 -DEAL:3d::1.48 -DEAL:3d::1.48 -DEAL:3d::1.48 -DEAL:3d::1.53 -DEAL:3d::1.48 -DEAL:3d::1.55 -DEAL:3d::1.48 -DEAL:3d::1.48 -DEAL:3d::1.48 -DEAL:3d::1.48 -DEAL:3d::1.48 -DEAL:3d::1.48 -DEAL:3d::1.48 -DEAL:3d::1.55 -DEAL:3d::1.59 -DEAL:3d::1.48 -DEAL:3d::1.48 -DEAL:3d::1.48 -DEAL:3d::1.48 -DEAL:3d::1.53 -DEAL:3d::1.48 -DEAL:3d::1.48 -DEAL:3d::1.48 -DEAL:3d::1.48 -DEAL:3d::1.48 -DEAL:3d::1.48 -DEAL:3d::1.55 -DEAL:3d::1.48 -DEAL:3d::1.48 -DEAL:3d::1.48 -DEAL:3d::1.48 -DEAL:3d::1.48 -DEAL:3d::1.48 -DEAL:3d::1.53 -DEAL:3d::1.48 -DEAL:3d::1.55 -DEAL:3d::1.48 -DEAL:3d::1.55 -DEAL:3d::1.59 -DEAL:3d::1.48 -DEAL:3d::1.48 -DEAL:3d::1.48 -DEAL:3d::1.48 -DEAL:3d::1.48 -DEAL:3d::1.48 -DEAL:3d::1.48 -DEAL:3d::1.48 -DEAL:3d::1.48 -DEAL:3d::1.48 -DEAL:3d::1.53 -DEAL:3d::1.48 -DEAL:3d::1.48 -DEAL:3d::1.55 -DEAL:3d::1.48 -DEAL:3d::1.48 -DEAL:3d::1.48 -DEAL:3d::1.48 -DEAL:3d::1.48 -DEAL:3d::1.48 -DEAL:3d::1.48 -DEAL:3d::1.48 -DEAL:3d::1.53 -DEAL:3d::1.55 -DEAL:3d::1.48 -DEAL:3d::1.48 -DEAL:3d::1.48 -DEAL:3d::1.48 -DEAL:3d::1.48 -DEAL:3d::1.48 -DEAL:3d::1.48 -DEAL:3d::1.48 -DEAL:3d::1.48 -DEAL:3d::1.53 -DEAL:3d::1.48 -DEAL:3d::1.48 -DEAL:3d::1.48 -DEAL:3d::1.48 -DEAL:3d::1.48 -DEAL:3d::1.48 -DEAL:3d::1.48 diff --git a/tests/mappings/mapping_q1_eulerian_01.cc b/tests/mappings/mapping_q1_eulerian_01.cc index 31571942ff..d06bb86dfc 100644 --- a/tests/mappings/mapping_q1_eulerian_01.cc +++ b/tests/mappings/mapping_q1_eulerian_01.cc @@ -41,7 +41,6 @@ #include #include -#include #include #include #include diff --git a/tests/mappings/mapping_q_eulerian_01.cc b/tests/mappings/mapping_q_eulerian_01.cc index a171851b92..4728db96af 100644 --- a/tests/mappings/mapping_q_eulerian_01.cc +++ b/tests/mappings/mapping_q_eulerian_01.cc @@ -39,7 +39,6 @@ #include #include -#include #include #include #include diff --git a/tests/mappings/mapping_q_eulerian_02.cc b/tests/mappings/mapping_q_eulerian_02.cc index 9002e1fa69..41362d5819 100644 --- a/tests/mappings/mapping_q_eulerian_02.cc +++ b/tests/mappings/mapping_q_eulerian_02.cc @@ -41,7 +41,6 @@ #include #include -#include #include #include #include diff --git a/tests/mappings/mapping_q_eulerian_03.cc b/tests/mappings/mapping_q_eulerian_03.cc index 953cb01b11..95dff81694 100644 --- a/tests/mappings/mapping_q_eulerian_03.cc +++ b/tests/mappings/mapping_q_eulerian_03.cc @@ -42,7 +42,6 @@ #include #include -#include #include #include #include diff --git a/tests/mappings/mapping_q_eulerian_07.cc b/tests/mappings/mapping_q_eulerian_07.cc index 598478ea1e..539b182591 100644 --- a/tests/mappings/mapping_q_eulerian_07.cc +++ b/tests/mappings/mapping_q_eulerian_07.cc @@ -42,7 +42,6 @@ #include #include -#include #include #include #include diff --git a/tests/mappings/mapping_q_eulerian_08.cc b/tests/mappings/mapping_q_eulerian_08.cc index c5a2cd74d2..870be2c220 100644 --- a/tests/mappings/mapping_q_eulerian_08.cc +++ b/tests/mappings/mapping_q_eulerian_08.cc @@ -44,7 +44,6 @@ #include #include -#include #include #include #include diff --git a/tests/mappings/mapping_q_eulerian_11.cc b/tests/mappings/mapping_q_eulerian_11.cc index 78dd5191f9..5c1a6e7497 100644 --- a/tests/mappings/mapping_q_eulerian_11.cc +++ b/tests/mappings/mapping_q_eulerian_11.cc @@ -41,7 +41,6 @@ #include #include -#include #include #include #include diff --git a/tests/matrix_free/interpolate_to_mg.cc b/tests/matrix_free/interpolate_to_mg.cc index a844db8c0a..6ddcccc3ee 100644 --- a/tests/matrix_free/interpolate_to_mg.cc +++ b/tests/matrix_free/interpolate_to_mg.cc @@ -42,7 +42,6 @@ #include #include -#include #include #include #include -- 2.39.5