]> https://gitweb.dealii.org/ - dealii-svn.git/commitdiff
Add another test.
authorbangerth <bangerth@0785d39b-7218-0410-832d-ea1e28bc413d>
Thu, 28 Oct 2010 23:44:37 +0000 (23:44 +0000)
committerbangerth <bangerth@0785d39b-7218-0410-832d-ea1e28bc413d>
Thu, 28 Oct 2010 23:44:37 +0000 (23:44 +0000)
git-svn-id: https://svn.dealii.org/trunk@22543 0785d39b-7218-0410-832d-ea1e28bc413d

tests/codim_one/data_out_02.cc [new file with mode: 0644]
tests/codim_one/data_out_02/cmp/generic [new file with mode: 0644]
tests/codim_one/data_out_02/square.msh [new file with mode: 0644]

diff --git a/tests/codim_one/data_out_02.cc b/tests/codim_one/data_out_02.cc
new file mode 100644 (file)
index 0000000..002d3ef
--- /dev/null
@@ -0,0 +1,70 @@
+
+//----------------------------  template.cc  ---------------------------
+//    $Id$
+//    Version: $Name$
+//
+//    Copyright (C) 2010 by the deal.II authors and Sebastian Pauletti
+//
+//    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.
+//
+//----------------------------  template.cc  ---------------------------
+
+
+// DataOut::build_patches appeared to have a problem with outputting
+// lines in 2d where nodes were numbered differently when writing data
+// vectors as opposed to writing node locations. in the end this
+// turned out to be a feature: the mesh was a circle of lines, so
+// there are equally many cells as their were nodes, and consequently
+// DataOut assumed that it had cell_data, rather than
+// dof_data. passing the correct argument fixed the problem, but it
+// won't hurt to have this test anyway.
+
+#include "../tests.h"
+#include <fstream>
+#include <base/logstream.h>
+
+#include <grid/tria.h>
+#include <dofs/dof_handler.h>
+#include <fe/fe_q.h>
+#include <grid/grid_in.h>
+#include <lac/vector.h>
+#include <numerics/data_out.h>
+
+std::ofstream logfile("data_out_02/output");
+
+
+int main ()
+{
+  deallog.attach(logfile);
+  deallog.depth_console(0);
+
+  const unsigned int dim = 1;
+
+  Triangulation<dim,dim+1>    triangulation;
+  FE_Q<dim,dim+1>             fe(1);
+  DoFHandler<dim,dim+1>       dof_handler(triangulation);
+  Vector<double> soln;
+
+  GridIn<dim,dim+1> grid_in;
+  grid_in.attach_triangulation (triangulation);
+  std::ifstream fname("data_out_02/square.msh");
+  grid_in.read_msh (fname);
+
+  dof_handler.distribute_dofs (fe);
+  soln.reinit (dof_handler.n_dofs());
+  soln = 0;
+  for (unsigned int i=0; i<soln.size(); ++i)
+    soln(i) = i;
+  DataOut<dim, DoFHandler<dim, dim+1> > data_out;
+  data_out.attach_dof_handler (dof_handler);
+
+  data_out.add_data_vector (soln, "scalar_data",
+                           DataOut<dim,DoFHandler<dim, dim+1> >::type_dof_data);
+  data_out.build_patches ();
+  data_out.write_vtk (deallog.get_file_stream());
+
+  return 0;
+}
diff --git a/tests/codim_one/data_out_02/cmp/generic b/tests/codim_one/data_out_02/cmp/generic
new file mode 100644 (file)
index 0000000..06306fd
--- /dev/null
@@ -0,0 +1,28 @@
+
+# vtk DataFile Version 3.0
+#This file was generated 
+ASCII
+DATASET UNSTRUCTURED_GRID
+
+POINTS 8 double
+0.00000 0.00000 0
+1.00000 0.00000 0
+1.00000 0.00000 0
+1.00000 1.00000 0
+1.00000 1.00000 0
+0.00000 1.00000 0
+0.00000 1.00000 0
+0.00000 0.00000 0
+
+CELLS 4 12
+2      0       1
+2      2       3
+2      4       5
+2      6       7
+
+CELL_TYPES 4
+ 3 3 3 3
+POINT_DATA 8
+SCALARS scalar_data double 1
+LOOKUP_TABLE default
+0.00000 1.00000 1.00000 2.00000 2.00000 3.00000 3.00000 0.00000 
diff --git a/tests/codim_one/data_out_02/square.msh b/tests/codim_one/data_out_02/square.msh
new file mode 100644 (file)
index 0000000..997b625
--- /dev/null
@@ -0,0 +1,17 @@
+$MeshFormat
+2.1 0 8
+$EndMeshFormat
+$Nodes
+4
+1 0 0 0
+2 1 0 0
+3 1 1 0
+4 0 1 0
+$EndNodes
+$Elements
+4
+1 1 3 0 1 0 1 2
+2 1 3 0 1 0 2 3
+3 1 3 0 1 0 3 4
+4 1 3 0 1 0 4 1
+$EndElements

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.