--- /dev/null
+############################################################
+# 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
--- /dev/null
+// $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();
+ }
+}
--- /dev/null
+// $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;
--- /dev/null
+// $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;
+};