]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Fix a bug in DataOut::build_patches.
authorWolfgang Bangerth <bangerth@math.tamu.edu>
Wed, 19 Jan 2011 01:09:35 +0000 (01:09 +0000)
committerWolfgang Bangerth <bangerth@math.tamu.edu>
Wed, 19 Jan 2011 01:09:35 +0000 (01:09 +0000)
git-svn-id: https://svn.dealii.org/trunk@23215 0785d39b-7218-0410-832d-ea1e28bc413d

deal.II/doc/news/changes.h
deal.II/source/numerics/data_out.cc
tests/codim_one/data_out_03.cc [new file with mode: 0644]
tests/codim_one/data_out_03/cmp/generic [new file with mode: 0644]
tests/codim_one/data_out_03/square.msh [new file with mode: 0644]

index 895d9de98634668a5614c53efacc9fcbbb710b56..35efa57dd17a8a1ce8015f0fdb1eecaab18c5ed5 100644 (file)
@@ -47,7 +47,14 @@ should be fixed now.
 <h3>Specific improvements</h3>
 
 <ol>
-<li> Restructured the internals of PETScWrappers::Precondition* to allow a
+<li> Fixed: When calling DataOut::build_patches with a mapping, requesting more
+than one subdivision, and when <code>dim@<spacedim</code>, then some cells
+were not properly mapped. This is now fixed.
+<br>
+(Wolfgang Bangerth, 2011/01/18)
+</li>
+
+<li> New: Restructured the internals of PETScWrappers::Precondition* to allow a
 PETSc PC object to exist without a solver. New: use Precondition*::vmult() to
 apply the preconditioner once. Preconditioners now have a default constructor
 and an initialize() function and are no longer initialized in the solver call,
index b1467ebd1c219bb66bd1cd86405e34dafe9dd8c6..90d85c1af2d7d83e3e7deea602d27e9eeaba5e72 100644 (file)
@@ -2,7 +2,7 @@
 //    $Id$
 //    Version: $Name$
 //
-//    Copyright (C) 1999, 2000, 2001, 2002, 2003, 2004, 2005, 2006, 2007, 2008, 2009, 2010 by the deal.II authors
+//    Copyright (C) 1999, 2000, 2001, 2002, 2003, 2004, 2005, 2006, 2007, 2008, 2009, 2010, 2011 by the deal.II authors
 //
 //    This file is subject to QPL and may not be  distributed
 //    without copyright and license information. Please refer
@@ -708,14 +708,22 @@ build_one_patch (const std::pair<cell_iterator, unsigned int> *cell_and_index,
                                       // of curved cells, if necessary
                                       // append the quadrature points to
                                       // the last rows of the patch.data
-                                      // member. THis is the case if we
+                                      // member. This is the case if we
                                       // want to produce curved cells at
                                       // the boundary and this cell
                                       // actually is at the boundary, or
                                       // else if we want to produce curved
                                       // cells everywhere
-      if (curved_cell_region==curved_inner_cells ||
-         (curved_cell_region==curved_boundary && cell_and_index->first->at_boundary()))
+                                      //
+                                      // note: a cell is *always* at
+                                      // the boundary if dim<spacedim
+      if (curved_cell_region==curved_inner_cells
+         ||
+         (curved_cell_region==curved_boundary
+          &&
+          (cell_and_index->first->at_boundary()
+           ||
+           (DH::dimension != DH::space_dimension))))
        {
          Assert(patch.space_dim==DH::space_dimension, ExcInternalError());
          const std::vector<Point<DH::space_dimension> > & q_points=fe_patch_values.get_quadrature_points();
@@ -768,11 +776,11 @@ build_one_patch (const std::pair<cell_iterator, unsigned int> *cell_and_index,
                  if (update_flags & update_hessians)
                    this->dof_data[dataset]->get_function_hessians (fe_patch_values,
                                                                    data.patch_hessians);
-                                                                   
+
                  if (update_flags & update_quadrature_points)
                    data.patch_evaluation_points = fe_patch_values.get_quadrature_points();
-                   
-                   
+
+
                  std::vector<Point<DH::space_dimension> > dummy_normals;
                  postprocessor->
                    compute_derived_quantities_scalar(data.patch_values,
@@ -797,12 +805,12 @@ build_one_patch (const std::pair<cell_iterator, unsigned int> *cell_and_index,
                  if (update_flags & update_hessians)
                    this->dof_data[dataset]->get_function_hessians (fe_patch_values,
                                                                    data.patch_hessians_system);
-                 
+
                  if (update_flags & update_quadrature_points)
                    data.patch_evaluation_points = fe_patch_values.get_quadrature_points();
-                 
+
                  std::vector<Point<DH::space_dimension> > dummy_normals;
-                 
+
                  postprocessor->
                    compute_derived_quantities_vector(data.patch_values_system,
                                                      data.patch_gradients_system,
diff --git a/tests/codim_one/data_out_03.cc b/tests/codim_one/data_out_03.cc
new file mode 100644 (file)
index 0000000..5040bce
--- /dev/null
@@ -0,0 +1,129 @@
+
+//----------------------------  template.cc  ---------------------------
+//    $Id$
+//    Version: $Name$
+//
+//    Copyright (C) 2010, 2011 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 did not compute mapped coordinates for all cells
+// if dim<spacedim even if a mapping was explicitly requested.
+
+#include "../tests.h"
+#include <fstream>
+#include <base/logstream.h>
+#include <base/function.h>
+
+#include <grid/tria.h>
+#include <grid/tria_boundary_lib.h>
+#include <grid/grid_generator.h>
+#include <grid/grid_tools.h>
+#include <dofs/dof_handler.h>
+#include <fe/fe_q.h>
+#include <fe/fe_system.h>
+#include <grid/grid_in.h>
+#include <lac/vector.h>
+#include <numerics/data_out.h>
+#include <numerics/vectors.h>
+#include <fe/mapping_q.h>
+
+std::ofstream logfile("data_out_03/output");
+
+template <int dim>
+class Identity : public Function<dim>
+{
+  public:
+    Identity ()
+                   :
+                   Function<dim>(dim)
+      {
+      }
+
+
+    virtual double value (const Point<dim> &p,
+                         const unsigned int component) const
+      {
+       return p(component);
+      }
+
+    virtual void vector_value(const Point<dim> &p,
+                             Vector< double > & values) const
+      {
+       for (unsigned int i=0;i<dim;i++)
+         values(i) = p(i);
+      }
+};
+
+
+int main ()
+{
+  deallog.attach(logfile);
+  deallog.depth_console(0);
+
+  int fe_degree =2;
+  int mapping_degree = 2;
+
+
+  Triangulation<2,3> tria;
+
+  std::map< Triangulation<2,3>::cell_iterator,
+    Triangulation<3,3>::face_iterator> surface_to_volume_mapping;
+
+  HyperBallBoundary<3> boundary_description;
+  Triangulation<3> volume_mesh;
+  GridGenerator::half_hyper_ball(volume_mesh);
+
+
+  volume_mesh.set_boundary (1, boundary_description);
+  volume_mesh.set_boundary (0, boundary_description);
+
+  static HyperBallBoundary<2,3> surface_description;
+  tria.set_boundary (1, surface_description);
+  tria.set_boundary (0, surface_description);
+
+  std::set<unsigned char> boundary_ids;
+  boundary_ids.insert(0);
+
+  GridTools::extract_boundary_mesh (volume_mesh, tria,
+                                   boundary_ids);
+
+  // test for the position
+  MappingQ<2,3>   mapping(mapping_degree,true);
+
+  FESystem<2,3> fe_test (FE_Q<2,3>
+                                                           (fe_degree),3);
+  DoFHandler<2,3> dh_test(tria);
+  dh_test.distribute_dofs(fe_test);
+
+  Vector<double> position(dh_test.n_dofs());
+  VectorTools::interpolate(mapping,dh_test,Identity<3>(),position);
+
+  std::vector<DataComponentInterpretation::DataComponentInterpretation>
+    data_component_interpretation
+    (3, DataComponentInterpretation::component_is_part_of_vector);
+
+  std::vector<std::string> solution_names (3, "position");
+
+  DataOut<2, DoFHandler<2,3> > data_out;
+  data_out.attach_dof_handler (dh_test);
+  data_out.add_data_vector (position, solution_names,
+                           DataOut<2,DoFHandler<2,
+                           3> >::type_dof_data,
+                           data_component_interpretation);
+  data_out.build_patches (mapping, 2);
+
+  data_out.write_gnuplot (deallog.get_file_stream());
+
+
+  dh_test.clear();
+  tria.set_boundary(0);
+
+  return 0;
+}
diff --git a/tests/codim_one/data_out_03/cmp/generic b/tests/codim_one/data_out_03/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_03/square.msh b/tests/codim_one/data_out_03/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.