From: Wolfgang Bangerth Date: Wed, 19 Jan 2011 01:09:35 +0000 (+0000) Subject: Fix a bug in DataOut::build_patches. X-Git-Tag: v8.0.0~4472 X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=628810ad713573b025a8a4ec19de547ff3dabff1;p=dealii.git Fix a bug in DataOut::build_patches. git-svn-id: https://svn.dealii.org/trunk@23215 0785d39b-7218-0410-832d-ea1e28bc413d --- diff --git a/deal.II/doc/news/changes.h b/deal.II/doc/news/changes.h index 895d9de986..35efa57dd1 100644 --- a/deal.II/doc/news/changes.h +++ b/deal.II/doc/news/changes.h @@ -47,7 +47,14 @@ should be fixed now.

Specific improvements

    -
  1. Restructured the internals of PETScWrappers::Precondition* to allow a +
  2. Fixed: When calling DataOut::build_patches with a mapping, requesting more +than one subdivision, and when dim@, then some cells +were not properly mapped. This is now fixed. +
    +(Wolfgang Bangerth, 2011/01/18) +
  3. + +
  4. 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, diff --git a/deal.II/source/numerics/data_out.cc b/deal.II/source/numerics/data_out.cc index b1467ebd1c..90d85c1af2 100644 --- a/deal.II/source/numerics/data_out.cc +++ b/deal.II/source/numerics/data_out.cc @@ -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_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 dimfirst->at_boundary() + || + (DH::dimension != DH::space_dimension)))) { Assert(patch.space_dim==DH::space_dimension, ExcInternalError()); const std::vector > & q_points=fe_patch_values.get_quadrature_points(); @@ -768,11 +776,11 @@ build_one_patch (const std::pair *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 > dummy_normals; postprocessor-> compute_derived_quantities_scalar(data.patch_values, @@ -797,12 +805,12 @@ build_one_patch (const std::pair *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 > 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 index 0000000000..5040bce743 --- /dev/null +++ b/tests/codim_one/data_out_03.cc @@ -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 +#include +#include + +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include + +std::ofstream logfile("data_out_03/output"); + +template +class Identity : public Function +{ + public: + Identity () + : + Function(dim) + { + } + + + virtual double value (const Point &p, + const unsigned int component) const + { + return p(component); + } + + virtual void vector_value(const Point &p, + Vector< double > & values) const + { + for (unsigned int i=0;i 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 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 position(dh_test.n_dofs()); + VectorTools::interpolate(mapping,dh_test,Identity<3>(),position); + + std::vector + data_component_interpretation + (3, DataComponentInterpretation::component_is_part_of_vector); + + std::vector 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 index 0000000000..06306fd562 --- /dev/null +++ b/tests/codim_one/data_out_03/cmp/generic @@ -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 index 0000000000..997b625492 --- /dev/null +++ b/tests/codim_one/data_out_03/square.msh @@ -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