]> https://gitweb.dealii.org/ - dealii-svn.git/commitdiff
Add new test.
authorwolf <wolf@0785d39b-7218-0410-832d-ea1e28bc413d>
Wed, 22 Oct 2003 22:32:46 +0000 (22:32 +0000)
committerwolf <wolf@0785d39b-7218-0410-832d-ea1e28bc413d>
Wed, 22 Oct 2003 22:32:46 +0000 (22:32 +0000)
git-svn-id: https://svn.dealii.org/trunk@8130 0785d39b-7218-0410-832d-ea1e28bc413d

tests/bits/Makefile
tests/bits/solution_transfer_1.cc [new file with mode: 0644]

index 312ecc61c710dea913b1478ff7626bc632798b8a..74135457dc6cdae005a0ca282e2f0847affbffd4 100644 (file)
@@ -134,7 +134,7 @@ volume_1.exe            : volume_1.g.$(OBJEXT)             $(libraries)
 volume_2.exe            : volume_2.g.$(OBJEXT)             $(libraries)
 volume_3.exe            : volume_3.g.$(OBJEXT)             $(libraries)
 volume_4.exe            : volume_4.g.$(OBJEXT)             $(libraries)
-
+solution_transfer_1.exe : solution_transfer_1.g.$(OBJEXT)  $(libraries)
 
 tests = anna_1 anna_2 anna_3 anna_4 anna_5 anna_6 \
        geometry_info_1 point_inside_1 point_inside_2 \
@@ -172,7 +172,7 @@ tests = anna_1 anna_2 anna_3 anna_4 anna_5 anna_6 \
        volume_1 volume_2 volume_3 volume_4 \
        mapping_cartesian_1 \
        mapping_q4_3d \
-       q_points
+       q_points solution_transfer
 
 ############################################################
 
diff --git a/tests/bits/solution_transfer_1.cc b/tests/bits/solution_transfer_1.cc
new file mode 100644 (file)
index 0000000..170ed2a
--- /dev/null
@@ -0,0 +1,74 @@
+//----------------------------  solution_transfer_1.cc  ---------------------------
+//    $Id$
+//    Version: $Name$
+//
+//    Copyright (C) 2002, 2003 by the deal.II authors and Anna Schneebeli
+//
+//    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.
+//
+//----------------------------  solution_transfer_1.cc  ---------------------------
+
+
+// something went wrong with the SolutionTransfer class: when we had a
+// vector or vectors as input and an empty vector of vectors as
+// output, then the output was a zero vector even if the input vector
+// was nonzero
+
+#include <base/logstream.h>
+#include <grid/tria.h>
+#include <dofs/dof_handler.h>
+#include <grid/grid_generator.h>
+#include <grid/tria_iterator.h>
+#include <dofs/dof_accessor.h>
+#include <dofs/dof_tools.h>
+#include <fe/fe_q.h>
+#include <numerics/solution_transfer.h>
+
+#include <fstream>
+    
+
+int main () 
+{
+  std::ofstream logfile("solution_transfer_1.output");
+  deallog.attach(logfile);
+  deallog.depth_console(0);
+
+  Triangulation<2> tria;
+  GridGenerator::hyper_cube (tria);
+  tria.refine_global (2);
+
+  FE_Q<2> fe(1);
+  
+  DoFHandler<2> dof_handler (tria);
+  dof_handler.distribute_dofs (fe);
+  
+  SolutionTransfer<2, double> soltrans (dof_handler);
+  std::vector<Vector<double> > solution_in, solution_out;
+  Vector<double> solution (dof_handler.n_dofs());
+  solution(0) = 1;
+  solution_in.push_back (solution);
+
+  tria.begin_active()->set_refine_flag();
+  tria.prepare_coarsening_and_refinement ();
+  soltrans.prepare_for_coarsening_and_refinement(solution_in);
+
+  tria.execute_coarsening_and_refinement ();
+  dof_handler.distribute_dofs (fe);
+
+  soltrans.interpolate(solution_in, solution_out);
+
+  Assert (solution_out.size() == solution_in.size(),
+         ExcInternalError());
+  Assert (solution_out[0].size() == dof_handler.n_dofs(),
+         ExcInternalError());
+
+                                  // this is the assertion that used
+                                  // to fail: the input vector is
+                                  // non-zero, so the output vector
+                                  // should not be either!
+  Assert (solution_out[0].l2_norm() > 0,
+         ExcInternalError());
+}

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.