From c958471a7776d988ac17f4b8113dfca737709313 Mon Sep 17 00:00:00 2001 From: wolf Date: Wed, 25 Feb 2004 19:46:34 +0000 Subject: [PATCH] Fold the implementations of testmatrix.cc into the header file. This way we can also use this file from the tests/bits directory. git-svn-id: https://svn.dealii.org/trunk@8545 0785d39b-7218-0410-832d-ea1e28bc413d --- tests/lac/Makefile | 16 +-- tests/lac/testmatrix.cc | 269 ---------------------------------------- tests/lac/testmatrix.h | 251 ++++++++++++++++++++++++++++++++++++- 3 files changed, 257 insertions(+), 279 deletions(-) delete mode 100644 tests/lac/testmatrix.cc diff --git a/tests/lac/Makefile b/tests/lac/Makefile index 75fb41e693..96235c7e50 100644 --- a/tests/lac/Makefile +++ b/tests/lac/Makefile @@ -1,6 +1,6 @@ ############################################################ # $Id$ -# Copyright (C) 2000, 2001, 2002, 2003 by the deal.II authors +# Copyright (C) 2000, 2001, 2002, 2003, 2004 by the deal.II authors ############################################################ ############################################################ @@ -25,15 +25,15 @@ block_vector.exe : block_vector.g.$(OBJEXT) $(librar block_vector_iterator.exe : block_vector_iterator.g.$(OBJEXT) $(libraries) full_matrix.exe : full_matrix.g.$(OBJEXT) $(libraries) matrix_out.exe : matrix_out.g.$(OBJEXT) $(libraries) -solver.exe : solver.g.$(OBJEXT) testmatrix.g.$(OBJEXT) $(libraries) -eigen.exe : eigen.g.$(OBJEXT) testmatrix.g.$(OBJEXT) $(libraries) -sparse_matrices.exe : sparse_matrices.g.$(OBJEXT) testmatrix.g.$(OBJEXT) $(libraries) -sparse_matrices.opt.exe : sparse_matrices.$(OBJEXT) testmatrix.$(OBJEXT) $(libraries_o) -sparsity_pattern.exe : sparsity_pattern.g.$(OBJEXT) testmatrix.g.$(OBJEXT) $(libraries) -sparse_ilu.exe : sparse_ilu.g.$(OBJEXT) testmatrix.g.$(OBJEXT) $(libraries) +solver.exe : solver.g.$(OBJEXT) $(libraries) +eigen.exe : eigen.g.$(OBJEXT) $(libraries) +sparse_matrices.exe : sparse_matrices.g.$(OBJEXT) $(libraries) +sparse_matrices.opt.exe : sparse_matrices.$(OBJEXT) $(libraries_o) +sparsity_pattern.exe : sparsity_pattern.g.$(OBJEXT) $(libraries) +sparse_ilu.exe : sparse_ilu.g.$(OBJEXT) $(libraries) vector-vector.exe : vector-vector.g.$(OBJEXT) $(libraries) #mgbase.exe : mgbase.g.$(OBJEXT) $(libraries) -#mg.exe : mg.g.$(OBJEXT) testmatrix.g.$(OBJEXT) $(libraries) +#mg.exe : mg.g.$(OBJEXT) $(libraries) tests = $(sort \ diff --git a/tests/lac/testmatrix.cc b/tests/lac/testmatrix.cc deleted file mode 100644 index a4ee0ceceb..0000000000 --- a/tests/lac/testmatrix.cc +++ /dev/null @@ -1,269 +0,0 @@ -//---------------------------- testmatrix.cc --------------------------- -// testmatrix.cc,v 1.22 2003/09/09 16:09:15 wolf Exp -// Version: -// -// Copyright (C) 1998, 1999, 2000, 2002, 2003 by the deal.II authors -// -// This file is subject to QPL and may not be distributed -// without copyright and license information. Please refer -// to the file deal.II/doc/license.html for the text and -// further information on this license. -// -//---------------------------- testmatrix.cc --------------------------- - - -#include "../tests.h" -#include "testmatrix.h" -#include -#include -#include -#include - -FDMatrix::FDMatrix(unsigned int nx, unsigned int ny) - : - nx(nx), ny(ny) -{} - -void -FDMatrix::five_point_structure(SparsityPattern& structure) const -{ - for(unsigned int i=0;i<=ny-2;i++) - { - for(unsigned int j=0;j<=nx-2; j++) - { - // Number of the row to be entered - unsigned int row = j+(nx-1)*i; - structure.add(row, row); - if (j>0) - { - structure.add(row-1, row); - structure.add(row, row-1); - } - if (j0) - { - structure.add(row-(nx-1), row); - structure.add(row, row-(nx-1)); - } - if (i0) - { - structure.add(row-1, row); - structure.add(row, row-1); - if (i>0) - { - structure.add(row-1, row-(nx-1)); - structure.add(row-(nx-1), row-1); - } - if (i0) - { - structure.add(row+1, row-(nx-1)); - structure.add(row-(nx-1), row+1); - } - if (i0) - { - structure.add(row-(nx-1), row); - structure.add(row, row-(nx-1)); - } - if (i -void -FDMatrix::nine_point(MATRIX& A, bool) const -{ - for(unsigned int i=0;i<=ny-2;i++) - { - for(unsigned int j=0;j<=nx-2; j++) - { - // Number of the row to be entered - unsigned int row = j+(nx-1)*i; - - A.set(row, row, 20.); - if (j>0) - { - A.set(row-1, row, -4.); - A.set(row, row-1, -4.); - if (i>0) - { - A.set(row-1, row-(nx-1), -1.); - A.set(row-(nx-1), row-1, -1.); - } - if (i0) - { - A.set(row+1, row-(nx-1), -1.); - A.set(row-(nx-1), row+1, -1.); - } - if (i0) - { - A.set(row-(nx-1), row, -4.); - A.set(row, row-(nx-1), -4.); - } - if (i -void -FDMatrix::five_point(MATRIX& A, bool nonsymmetric) const -{ - for(unsigned int i=0;i<=ny-2;i++) - { - for(unsigned int j=0;j<=nx-2; j++) - { - // Number of the row to be entered - unsigned int row = j+(nx-1)*i; - if (nonsymmetric) - A.set(row, row, 5.); - else - A.set(row, row, 4.); - if (j>0) - { - if (nonsymmetric) - A.set(row-1, row, -2.); - else - A.set(row-1, row, -1.); - A.set(row, row-1, -1.); - } - if (j0) - { - A.set(row-(nx-1), row, -1.); - A.set(row, row-(nx-1), -1.); - } - if (i -void -FDMatrix::upwind(MATRIX& A, bool back) const -{ - for(unsigned int i=0;i<=ny-2;i++) - { - for(unsigned int j=0;j<=nx-2; j++) - { - // Number of the row to be entered - unsigned int row = j+(nx-1)*i; - A.set(row, row, 3.); - - if (j>0 && !back) - A.set(row, row-1, -1.); - - if (j -void -FDMatrix::gnuplot_print(std::ostream& s, const Vector& V) const -{ - for(unsigned int i=0;i<=ny-2;i++) - { - for(unsigned int j=0;j<=nx-2; j++) - { - // Number of the row to be entered - unsigned int row = j+(nx-1)*i; - s << (j+1)/(float)nx << '\t' << (i+1)/(float)ny << '\t' << V(row) << std::endl; - } - s << std::endl; - } - s << std::endl; -} - - - -template void FDMatrix::five_point(SparseMatrix&, bool) const; -template void FDMatrix::five_point(SparseMatrix&, bool) const; -template void FDMatrix::five_point(BlockSparseMatrix&, bool) const; -template void FDMatrix::five_point(BlockSparseMatrix&, bool) const; -template void FDMatrix::five_point(SparseMatrixEZ&, bool) const; -template void FDMatrix::five_point(SparseMatrix&, bool) const; -template void FDMatrix::nine_point(SparseMatrix&, bool) const; -template void FDMatrix::nine_point(SparseMatrix&, bool) const; -template void FDMatrix::nine_point(SparseMatrixEZ&, bool) const; -template void FDMatrix::nine_point(SparseMatrix&, bool) const; -template void FDMatrix::upwind(SparseMatrix&, bool) const; -template void FDMatrix::upwind(SparseMatrix&, bool) const; -template void FDMatrix::upwind(SparseMatrixEZ&, bool) const; -template void FDMatrix::upwind(SparseMatrix&, bool) const; -template void FDMatrix::gnuplot_print(std::ostream& s, const Vector& V) const; -template void FDMatrix::gnuplot_print(std::ostream& s, const Vector& V) const; -template void FDMatrix::gnuplot_print(std::ostream& s, const Vector& V) const; - diff --git a/tests/lac/testmatrix.h b/tests/lac/testmatrix.h index 04ab4a5fd8..f7579d2969 100644 --- a/tests/lac/testmatrix.h +++ b/tests/lac/testmatrix.h @@ -3,9 +3,11 @@ #include "../tests.h" #include #include +#include +#include +#include +#include -template class Vector; -class SparsityPattern; /** * Finite difference matrix on uniform grid. @@ -62,3 +64,248 @@ class FDMatrix */ unsigned int ny; }; + + +// --------------- inline and template functions ----------------- + +inline +FDMatrix::FDMatrix(unsigned int nx, unsigned int ny) + : + nx(nx), ny(ny) +{} + + + +inline +void +FDMatrix::five_point_structure(SparsityPattern& structure) const +{ + for(unsigned int i=0;i<=ny-2;i++) + { + for(unsigned int j=0;j<=nx-2; j++) + { + // Number of the row to be entered + unsigned int row = j+(nx-1)*i; + structure.add(row, row); + if (j>0) + { + structure.add(row-1, row); + structure.add(row, row-1); + } + if (j0) + { + structure.add(row-(nx-1), row); + structure.add(row, row-(nx-1)); + } + if (i0) + { + structure.add(row-1, row); + structure.add(row, row-1); + if (i>0) + { + structure.add(row-1, row-(nx-1)); + structure.add(row-(nx-1), row-1); + } + if (i0) + { + structure.add(row+1, row-(nx-1)); + structure.add(row-(nx-1), row+1); + } + if (i0) + { + structure.add(row-(nx-1), row); + structure.add(row, row-(nx-1)); + } + if (i +void +FDMatrix::nine_point(MATRIX& A, bool) const +{ + for(unsigned int i=0;i<=ny-2;i++) + { + for(unsigned int j=0;j<=nx-2; j++) + { + // Number of the row to be entered + unsigned int row = j+(nx-1)*i; + + A.set(row, row, 20.); + if (j>0) + { + A.set(row-1, row, -4.); + A.set(row, row-1, -4.); + if (i>0) + { + A.set(row-1, row-(nx-1), -1.); + A.set(row-(nx-1), row-1, -1.); + } + if (i0) + { + A.set(row+1, row-(nx-1), -1.); + A.set(row-(nx-1), row+1, -1.); + } + if (i0) + { + A.set(row-(nx-1), row, -4.); + A.set(row, row-(nx-1), -4.); + } + if (i +inline +void +FDMatrix::five_point(MATRIX& A, bool nonsymmetric) const +{ + for(unsigned int i=0;i<=ny-2;i++) + { + for(unsigned int j=0;j<=nx-2; j++) + { + // Number of the row to be entered + unsigned int row = j+(nx-1)*i; + if (nonsymmetric) + A.set(row, row, 5.); + else + A.set(row, row, 4.); + if (j>0) + { + if (nonsymmetric) + A.set(row-1, row, -2.); + else + A.set(row-1, row, -1.); + A.set(row, row-1, -1.); + } + if (j0) + { + A.set(row-(nx-1), row, -1.); + A.set(row, row-(nx-1), -1.); + } + if (i +inline +void +FDMatrix::upwind(MATRIX& A, bool back) const +{ + for(unsigned int i=0;i<=ny-2;i++) + { + for(unsigned int j=0;j<=nx-2; j++) + { + // Number of the row to be entered + unsigned int row = j+(nx-1)*i; + A.set(row, row, 3.); + + if (j>0 && !back) + A.set(row, row-1, -1.); + + if (j +inline +void +FDMatrix::gnuplot_print(std::ostream& s, const Vector& V) const +{ + for(unsigned int i=0;i<=ny-2;i++) + { + for(unsigned int j=0;j<=nx-2; j++) + { + // Number of the row to be entered + unsigned int row = j+(nx-1)*i; + s << (j+1)/(float)nx << '\t' << (i+1)/(float)ny << '\t' << V(row) << std::endl; + } + s << std::endl; + } + s << std::endl; +} -- 2.39.5