]> https://gitweb.dealii.org/ - dealii-svn.git/commitdiff
Test fuer iterative Loeser
authorguido <guido@0785d39b-7218-0410-832d-ea1e28bc413d>
Mon, 26 Apr 1999 20:04:34 +0000 (20:04 +0000)
committerguido <guido@0785d39b-7218-0410-832d-ea1e28bc413d>
Mon, 26 Apr 1999 20:04:34 +0000 (20:04 +0000)
git-svn-id: https://svn.dealii.org/trunk@1210 0785d39b-7218-0410-832d-ea1e28bc413d

tests/lac/Makefile [new file with mode: 0644]
tests/lac/solver.cc [new file with mode: 0644]
tests/lac/testmatrix.cc [new file with mode: 0644]
tests/lac/testmatrix.h [new file with mode: 0644]

diff --git a/tests/lac/Makefile b/tests/lac/Makefile
new file mode 100644 (file)
index 0000000..3d1a1cf
--- /dev/null
@@ -0,0 +1,109 @@
+############################################################
+# Directory with the DEAL distribution
+############################################################
+
+D = ../..
+
+
+############################################################
+# Include general settings for including DEAL libraries
+############################################################
+
+include $D/common/Make.global_options
+
+
+############################################################
+# Set debug-mode as a default
+############################################################
+
+debug-mode = on
+
+
+############################################################
+# Define library names
+############################################################
+
+libs   = -llac -lbase
+libs.g = -llac.g -lbase.g
+
+
+############################################################
+# Select compiler flags according to debug-mode
+############################################################
+
+ifeq ($(debug-mode),on)
+libraries = $(libs.g)
+flags     = $(CXXFLAGS.g) $(CXXFLAGS)
+endif
+
+ifeq ($(debug-mode),off)
+libraries = $(libs)
+flags     = $(CXXFLAGS.o) $(CXXFLAGS)
+endif
+
+all: solver
+
+############################################################
+# Typical block for building a running program
+#
+# 1. provide a list of source files in ...-cc-files
+#
+# 2. generate the list of object files according to debug-mode
+#
+# 3. make executable
+#
+# 4. Explicit dependencies of object files (will be automatic soon)
+#
+############################################################
+
+solver-cc-files = solver.cc testmatrix.cc
+
+ifeq ($(debug-mode),on)
+solver-o-files = $(solver-cc-files:.cc=.go)
+else
+solver-o-files = $(solver-cc-files:.cc=.o)
+endif
+
+solver: $(solver-o-files) $(libraries)
+       $(CXX) $(flags) -o $@ $^
+
+############################################################
+# Continue with other targets if needed
+############################################################
+
+
+target1-cc-files = t1.cc t2.cc t3.cc
+
+ifeq ($(debug-mode),on)
+target1-o-files = $(target1-cc-files:.cc=.go)
+else
+target1-o-files = $(target1-cc-files:.cc=.o)
+endif
+
+target1: $(target1-o-files) $(libraries)
+       $(CXX) $(flags) -o $@ $^
+
+
+
+############################################################
+# Cleanup targets
+############################################################
+
+clean:
+       rm -f *.o *.go
+
+veryclean: clean
+       rm -f solver *.inp *.gpl *.eps *.gnuplot
+
+############################################################
+# Automatic generation of dependencies
+############################################################
+
+all-cc-files = $(solver-cc-files) # $(target2-cc-files) ...
+
+Make.depend: $(all-cc-files)
+       @echo =====Dependecies=== Make.depend
+       @$(CXX) $(flags) $^ -M > $@
+       @perl -pi~ -e 's/(^[^.]+)\.o:/\1.o \1.go:/;' $@
+
+include Make.depend
diff --git a/tests/lac/solver.cc b/tests/lac/solver.cc
new file mode 100644 (file)
index 0000000..0ca3ec9
--- /dev/null
@@ -0,0 +1,87 @@
+// $Id$
+//
+// Test program for linear solvers.
+
+#include <cmath>
+#include "testmatrix.h"
+#include <base/logstream.h>
+#include <lac/sparsematrix.h>
+#include <lac/vector.h>
+#include <lac/vector_memory.h>
+#include <lac/solver_control.h>
+#include <lac/solver_pcg.h>
+#include <lac/solver_pgmres.h>
+#include <lac/solver_bicgstab.h>
+#include <lac/solver_richardson.h>
+#include <lac/precondition.h>
+
+#define MATRIX SparseMatrix<float> 
+#define VECTOR Vector<double> 
+
+main()
+{
+  PrimitiveVectorMemory<VECTOR > mem;
+  SolverControl control(100, 1.e-5, true);
+  SolverPCG<MATRIX, VECTOR > cg(control, mem);
+  SolverGMRES<MATRIX, VECTOR > gmres(control, mem,20);
+  SolverBicgstab<MATRIX, VECTOR > bicgstab(control, mem);
+  SolverRichardson<MATRIX, VECTOR > rich(control, mem);
+
+  for (unsigned int size=10; size <= 10; size *= 10)
+    {
+      deallog << "Size " << size << endl;
+      
+      unsigned int dim = (size+1)*(size+1);
+
+                                      // Make matrix
+      FDMatrix testproblem(size, size);
+      SparseMatrixStruct structure(dim, dim, 5);
+      testproblem.build_structure(structure);
+      structure.compress();
+      MATRIX A(structure);
+      testproblem.laplacian(A);
+
+      PreconditionIdentity<VECTOR >
+       prec_no;
+      PreconditionRelaxation<MATRIX, VECTOR>
+       prec_ssor(A, &MATRIX::template precondition_SSOR<double>, 1.2);
+      
+      VECTOR f(dim);
+      VECTOR u(dim);
+      
+      deallog.push("no");
+
+      f = 1.;
+      u = 0.;
+      cg.solve(A,u,f,prec_no);
+
+      f = 1.;
+      u = 0.;
+      bicgstab.solve(A,u,f,prec_no);
+
+      f = 1.;
+      u = 0.;
+      gmres.solve(A,u,f,prec_no);
+
+      deallog.pop();
+      deallog.push("ssor");      
+
+      f = 1.;
+      u = 0.;
+      rich.solve(A,u,f,prec_ssor);
+
+      f = 1.;
+      u = 0.;
+      cg.solve(A,u,f,prec_ssor);
+
+      f = 1.;
+      u = 0.;
+      bicgstab.solve(A,u,f,prec_ssor);
+
+      f = 1.;
+      u = 0.;
+      gmres.solve(A,u,f,prec_ssor);
+
+      deallog.pop();
+    }
+}
diff --git a/tests/lac/testmatrix.cc b/tests/lac/testmatrix.cc
new file mode 100644 (file)
index 0000000..656cade
--- /dev/null
@@ -0,0 +1,82 @@
+// $Id$
+
+#include "testmatrix.h"
+#include <lac/sparsematrix.h>
+
+FDMatrix::FDMatrix(unsigned int nx, unsigned int ny)
+               :
+               nx(nx), ny(ny)
+{}
+
+void
+FDMatrix::build_structure(SparseMatrixStruct& structure) const
+{
+  for(unsigned int i=0;i<=ny;i++)
+    {
+      for(unsigned int j=0;j<=nx; 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 (j<nx)
+           {
+             structure.add(row+1, row);
+             structure.add(row, row+1);
+           }
+         if (i>0)
+           {
+             structure.add(row-nx, row);
+             structure.add(row, row-nx);
+           }
+         if (i<ny)
+           {
+             structure.add(row+nx, row);
+             structure.add(row, row+nx);
+           }
+       }
+    }
+}
+
+template<typename number>
+void
+FDMatrix::laplacian(SparseMatrix<number>& A) const
+{
+   for(unsigned int i=0;i<=ny;i++)
+    {
+      for(unsigned int j=0;j<=nx; j++)
+       {
+                                          // Number of the row to be entered
+         unsigned int row = j+(nx+1)*i;
+         
+         A.set(row, row, 4.);
+         if (j>0)
+           {
+             A.set(row-1, row, -1.);
+             A.set(row, row-1, -1.);
+           }
+         if (j<nx)
+           {
+             A.set(row+1, row, -1.);
+             A.set(row, row+1, -1.);
+           }
+         if (i>0)
+           {
+             A.set(row-nx, row, -1.);
+             A.set(row, row-nx, -1.);
+           }
+         if (i<ny)
+           {
+             A.set(row+nx, row, -1.);
+             A.set(row, row+nx, -1.);
+           }
+       }
+    } 
+}
+
+template void FDMatrix::laplacian(SparseMatrix<double>&) const;
+template void FDMatrix::laplacian(SparseMatrix<float>&) const;
diff --git a/tests/lac/testmatrix.h b/tests/lac/testmatrix.h
new file mode 100644 (file)
index 0000000..ffb32d8
--- /dev/null
@@ -0,0 +1,39 @@
+// $Id$
+
+#include <lac/forward-declarations.h>
+
+/**
+ * Finite difference matrix on uniform grid.
+ * Generator for simple 5-point discretization of Laplace problem.
+ */
+
+class FDMatrix
+{
+  public:
+                                    /**
+                                     * Constructor specifying grid resolution.
+                                     */
+    FDMatrix(unsigned int nx, unsigned int ny);
+    
+                                    /**
+                                     * Generate the matrix structure.
+                                     */
+    void build_structure(SparseMatrixStruct& structure) const;
+    
+                                    /**
+                                     * Fill the matrix with values.
+                                     */
+    template <typename number>
+    void laplacian(SparseMatrix<number>&) const;
+
+  private:
+                                    /**
+                                     * Number of gridpoints in x-direction.
+                                     */
+    unsigned int nx;
+
+                                    /**
+                                     * Number of gridpoints in y-direction.
+                                     */
+    unsigned int ny;
+};

In the beginning the Universe was created. This has made a lot of people very angry and has been widely regarded as a bad move.

Douglas Adams


Typeset in Trocchi and Trocchi Bold Sans Serif.