From: Ralf Hartmann Date: Fri, 5 Jun 2009 04:18:11 +0000 (+0000) Subject: Output of curved patches for DataOutFaces if a Mapping is provided to the build_patch... X-Git-Tag: v8.0.0~7638 X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=7e8109994682eb5919e8a5d1c4b58af25246e24b;p=dealii.git Output of curved patches for DataOutFaces if a Mapping is provided to the build_patches() call (5560). git-svn-id: https://svn.dealii.org/trunk@18909 0785d39b-7218-0410-832d-ea1e28bc413d --- diff --git a/deal.II/deal.II/include/numerics/data_out_faces.h b/deal.II/deal.II/include/numerics/data_out_faces.h index c889a67ed7..a4cdc00d91 100644 --- a/deal.II/deal.II/include/numerics/data_out_faces.h +++ b/deal.II/deal.II/include/numerics/data_out_faces.h @@ -45,10 +45,12 @@ namespace internal const unsigned int n_datasets, const unsigned int n_subdivisions, const std::vector &n_postprocessor_outputs, + const Mapping &mapping, const FE &finite_elements, const UpdateFlags update_flags); const dealii::hp::QCollection q_collection; + const dealii::hp::MappingCollection mapping_collection; dealii::hp::FEFaceValues x_fe_values; std::vector > patch_normals; @@ -146,6 +148,50 @@ class DataOutFaces : public DataOut_DoFDatan_subdivisions>1, the + * points interior of subdivided + * patches which originate from + * cells at the boundary of the + * domain can be computed using the + * mapping, i.e. a higher order + * mapping leads to a + * representation of a curved + * boundary by using more + * subdivisions. + * + * Even for non-curved cells the + * mapping argument can be used + * for the Eulerian mappings (see + * class MappingQ1Eulerian) where + * a mapping is used not only to + * determine the position of + * points interior to a cell, but + * also of the vertices. It + * offers an opportunity to watch + * the solution on a deformed + * triangulation on which the + * computation was actually + * carried out, even if the mesh + * is internally stored in its + * undeformed configuration and + * the deformation is only + * tracked by an additional + * vector that holds the + * deformation of each vertex. + * + * @todo The @p mapping argument should be + * replaced by a hp::MappingCollection in + * case of a hp::DoFHandler. + */ + virtual void build_patches (const Mapping &mapping, + const unsigned int n_subdivisions = 0); + /** * Declare a way to describe a * face which we would like to diff --git a/deal.II/deal.II/source/numerics/data_out_faces.cc b/deal.II/deal.II/source/numerics/data_out_faces.cc index 399a7af9b8..9ae1ef49f5 100644 --- a/deal.II/deal.II/source/numerics/data_out_faces.cc +++ b/deal.II/deal.II/source/numerics/data_out_faces.cc @@ -41,6 +41,7 @@ namespace internal const unsigned int n_datasets, const unsigned int n_subdivisions, const std::vector &n_postprocessor_outputs, + const Mapping &mapping, const FE &finite_elements, const UpdateFlags update_flags) : @@ -52,7 +53,9 @@ namespace internal n_postprocessor_outputs, finite_elements), q_collection (quadrature), - x_fe_values (this->fe_collection, + mapping_collection (mapping), + x_fe_values (this->mapping_collection, + this->fe_collection, q_collection, update_flags) {} @@ -84,8 +87,24 @@ build_one_patch (const FaceDescriptor *cell_and_face, internal::DataOutFaces::ParallelData &data, DataOutBase::Patch &patch) { + // we use the mapping to transform the + // vertices. However, the mapping works on + // cells, not faces, so transform the face + // vertex to a cell vertex, that to a unit + // cell vertex and then, finally, that to + // the mapped vertex. In most cases this + // complicated procedure will be the + // identity. for (unsigned int vertex=0; vertex::vertices_per_cell; ++vertex) - patch.vertices[vertex] = cell_and_face->first->face(cell_and_face->second)->vertex(vertex); + patch.vertices[vertex] = data.mapping_collection[0].transform_unit_to_real_cell + (cell_and_face->first, + GeometryInfo::unit_cell_vertex + (GeometryInfo::face_to_cell_vertices + (cell_and_face->second, + vertex, + cell_and_face->first->face_orientation(cell_and_face->second), + cell_and_face->first->face_flip(cell_and_face->second), + cell_and_face->first->face_rotation(cell_and_face->second)))); if (data.n_datasets > 0) { @@ -94,8 +113,25 @@ build_one_patch (const FaceDescriptor *cell_and_face, = data.x_fe_values.get_present_fe_values (); const unsigned int n_q_points = fe_patch_values.n_quadrature_points; - - + + // store the intermediate points + Assert(patch.space_dim==DH::dimension, ExcInternalError()); + const std::vector > & q_points=fe_patch_values.get_quadrature_points(); + // resize the patch.data member + // in order to have enough memory + // for the quadrature points as + // well + patch.data.reinit(patch.data.size(0)+DH::dimension, + patch.data.size(1)); + // set the flag indicating that + // for this cell the points are + // explicitly given + patch.points_are_available=true; + // copy points to patch.data + for (unsigned int i=0; i void DataOutFaces::build_patches (const unsigned int n_subdivisions_) +{ + build_patches (StaticMappingQ1::mapping, n_subdivisions_); +} + + + +template +void DataOutFaces::build_patches (const Mapping &mapping, + const unsigned int n_subdivisions_) { // Check consistency of redundant // template parameter @@ -286,11 +331,14 @@ void DataOutFaces::build_patches (const unsigned int n_subdivisions_) for (unsigned int i=0; idof_data.size(); ++i) if (this->dof_data[i]->postprocessor) update_flags |= this->dof_data[i]->postprocessor->get_needed_update_flags(); + update_flags |= update_quadrature_points; internal::DataOutFaces::ParallelData thread_data (patch_points, n_components, n_datasets, n_subdivisions, - n_postprocessor_outputs, this->dofs->get_fe(), + n_postprocessor_outputs, + mapping, + this->dofs->get_fe(), update_flags); DataOutBase::Patch sample_patch; sample_patch.n_subdivisions = n_subdivisions; @@ -388,7 +436,7 @@ DataOutFaces::next_face (const FaceDescriptor &old_face) // explicit instantiations -// don't instantiate anything for the 1d and 2d cases +// don't instantiate anything for the 1d #if deal_II_dimension >=2 template class DataOutFaces >; template class DataOutFaces >;