]> https://gitweb.dealii.org/ - dealii.git/commitdiff
example for use of BlockMatrixArray and BlockTrianglePrecondition
authorGuido Kanschat <dr.guido.kanschat@gmail.com>
Sun, 6 Mar 2005 02:39:23 +0000 (02:39 +0000)
committerGuido Kanschat <dr.guido.kanschat@gmail.com>
Sun, 6 Mar 2005 02:39:23 +0000 (02:39 +0000)
git-svn-id: https://svn.dealii.org/trunk@10008 0785d39b-7218-0410-832d-ea1e28bc413d

deal.II/examples/doxygen/Makefile [new file with mode: 0644]
deal.II/examples/doxygen/block_matrix_array.cc [new file with mode: 0644]

diff --git a/deal.II/examples/doxygen/Makefile b/deal.II/examples/doxygen/Makefile
new file mode 100644 (file)
index 0000000..0a8ac81
--- /dev/null
@@ -0,0 +1,40 @@
+###########################################################################
+#    $Id$
+#    Version: $Name$
+#
+#    Copyright (C) 1998 - 2005 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.
+#
+###########################################################################
+include ../../common/Make.global_options
+
+all: block_matrix_array.exe
+
+######################################################################
+# Compilation of source code
+######################################################################
+
+%.o : %.cc
+       @echo =====debug======$(MT)== $<
+       @$(CXX) $(CXXFLAGS.g) $(CXXFLAGS) -c $< -o $@
+
+######################################################################
+# Specific targets
+######################################################################
+
+block_matrix_array.exe: block_matrix_array.o $(lib-base.g) $(lib-lac.g)
+       @echo =====linking====$(MT)== $<
+       @$(CXX) $(CXXFLAGS.g) $(LDFLAGS) -o $@ $^
+
+
+
+######################################################################
+# Pseudo target for cleaning directory
+######################################################################
+
+clean:
+       rm *~ *.o *.exe
diff --git a/deal.II/examples/doxygen/block_matrix_array.cc b/deal.II/examples/doxygen/block_matrix_array.cc
new file mode 100644 (file)
index 0000000..a69b774
--- /dev/null
@@ -0,0 +1,122 @@
+//---------------------------------------------------------------------------
+//    $Id$
+//    Version: $Name$ 
+//
+//    Copyright (C) 2005 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.
+//
+//---------------------------------------------------------------------------
+
+// See documentation of BlockMatrixArray for documentation of this example
+
+#include <base/logstream.h>
+#include <lac/block_matrix_array.h>
+#include <lac/full_matrix.h>
+#include <lac/vector.h>
+#include <lac/block_vector.h>
+#include <lac/vector_memory.h>
+#include <lac/precondition.h>
+#include <lac/solver_cg.h>
+#include <lac/solver_gmres.h>
+
+#include <iostream>
+#include <fstream>
+
+double Adata[] =
+{
+      4., .5, .1, 0.,
+      .5, 4., .5, .1,
+      .1, .5, 4., .5,
+      0., .1, .5, 4.
+};
+
+double B1data[] =
+{
+      .5, .1,
+      .4, .2,
+      .3, .3,
+      .2, .4
+};
+
+double B2data[] =
+{
+      .3, 0., -.3, 0.,
+      -.3, 0., .3, 0.
+};
+
+double Cdata[] =
+{
+      8., 1.,
+      1., 8.
+};
+
+int main () 
+{
+  FullMatrix<float> A(4,4);
+  FullMatrix<float> B1(4,2);
+  FullMatrix<float> B2(2,4);
+  FullMatrix<float> C(2,2);
+
+  A.fill(Adata);
+  B1.fill(B1data);
+  B2.fill(B2data);
+  C.fill(Cdata);
+  
+  BlockMatrixArray<FullMatrix<float> > matrix(2,2);
+  
+  matrix.enter(A,0,0,2.);
+  matrix.enter(B1,0,1,-1.);
+  matrix.enter(B2,0,1,1., true);
+  matrix.enter(B2,1,0,1.);
+  matrix.enter(B1,1,0,-1., true);
+  matrix.enter(C,1,1);
+  matrix.print_latex(deallog);
+  
+  std::vector<unsigned int> block_sizes(2);
+  block_sizes[0] = 4;
+  block_sizes[1] = 2;
+  
+  BlockVector<double> result(block_sizes);
+  BlockVector<double> x(block_sizes);
+  BlockVector<double> y(block_sizes);
+  for (unsigned int i=0;i<result.size();++i)
+    result(i) = i;
+
+  matrix.vmult(y, result);
+
+  SolverControl control(100,1.e-10);
+  GrowingVectorMemory<BlockVector<double> > mem;
+  PreconditionIdentity id;
+
+  SolverCG<BlockVector<double> > cg(control, mem);
+  cg.solve(matrix, x, y, id);
+  x.add(-1., result);
+  deallog << "Error " << x.l2_norm() << std::endl;
+  
+  FullMatrix<float> Ainv(4,4);
+  Ainv.invert(A);
+  FullMatrix<float> Cinv(2,2);
+  Cinv.invert(C);
+  
+  BlockTrianglePrecondition<FullMatrix<float> > precondition(2,2);
+  precondition.enter(Ainv,0,0,.5);
+  precondition.enter(Cinv,1,1);
+
+  cg.solve(matrix, x, y, precondition);
+  x.add(-1., result);
+  deallog << "Error " << x.l2_norm() << std::endl;
+  
+  precondition.enter(B1,1,0,-1., true);
+  precondition.enter(B2,1,0,1.);
+  
+  SolverGMRES<BlockVector<double> > gmres(control, mem);
+  gmres.solve(matrix, x, y, precondition);
+  x.add(-1., result);
+  deallog << "Error " << x.l2_norm() << std::endl;
+  
+  return 0;
+}

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.