From: Wolfgang Bangerth Date: Mon, 25 Oct 2021 03:49:10 +0000 (-0600) Subject: Do not preset Patch::reference_cell. X-Git-Tag: v9.4.0-rc1~885^2 X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=3b51eae73cb130cb832e37b041320693b2305fc7;p=dealii.git Do not preset Patch::reference_cell. --- diff --git a/include/deal.II/base/bounding_box_data_out.h b/include/deal.II/base/bounding_box_data_out.h index 2dd9738de4..7469889710 100644 --- a/include/deal.II/base/bounding_box_data_out.h +++ b/include/deal.II/base/bounding_box_data_out.h @@ -141,9 +141,10 @@ BoundingBoxDataOut::build_patches( boost::geometry::convert(getter(*value), box); for (unsigned int v = 0; v < GeometryInfo::vertices_per_cell; ++v) { - patches[i].vertices[v] = box.vertex(v); - patches[i].patch_index = i; - patches[i].n_subdivisions = 1; + patches[i].vertices[v] = box.vertex(v); + patches[i].patch_index = i; + patches[i].n_subdivisions = 1; + patches[i].reference_cell = ReferenceCells::get_hypercube(); patches[i].points_are_available = false; } ++i; diff --git a/include/deal.II/lac/matrix_out.h b/include/deal.II/lac/matrix_out.h index 21319b292e..96866fb679 100644 --- a/include/deal.II/lac/matrix_out.h +++ b/include/deal.II/lac/matrix_out.h @@ -349,12 +349,15 @@ MatrixOut::build_patches(const Matrix & matrix, for (size_type i = 0; i < gridpoints_y; ++i) for (size_type j = 0; j < gridpoints_x; ++j, ++index) { + patches[index].n_subdivisions = 1; + patches[index].reference_cell = ReferenceCells::Quadrilateral; + // within each patch, order the points in such a way that if some // graphical output program (such as gnuplot) plots the quadrilaterals // as two triangles, then the diagonal of the quadrilateral which cuts // it into the two printed triangles is parallel to the diagonal of the // matrix, rather than perpendicular to it. this has the advantage that, - // for example, the unit matrix is plotted as a straight rim, rather + // for example, the unit matrix is plotted as a straight ridge, rather // than as a series of bumps and valleys along the diagonal patches[index].vertices[0](0) = j; patches[index].vertices[0](1) = -static_cast(i); diff --git a/source/base/data_out_base.cc b/source/base/data_out_base.cc index b0c68aa430..cfdf5227d8 100644 --- a/source/base/data_out_base.cc +++ b/source/base/data_out_base.cc @@ -959,6 +959,12 @@ namespace n_cells = 0; for (const auto &patch : patches) { + Assert(patch.reference_cell != ReferenceCells::Invalid, + ExcMessage( + "The reference cell for this patch is set to 'Invalid', " + "but that is clearly not a valid choice. Did you forget " + "to set the reference cell for the patch?")); + // The following formula doesn't hold for non-tensor products. if (patch.reference_cell == ReferenceCells::get_hypercube()) { @@ -1964,7 +1970,7 @@ namespace DataOutBase : patch_index(no_neighbor) , n_subdivisions(1) , points_are_available(false) - , reference_cell(ReferenceCells::get_hypercube()) + , reference_cell(ReferenceCells::Invalid) // all the other data has a constructor of its own, except for the "neighbors" // field, which we set to invalid values. { diff --git a/source/numerics/data_out_rotation.cc b/source/numerics/data_out_rotation.cc index 87fe4d7269..4b12a279ca 100644 --- a/source/numerics/data_out_rotation.cc +++ b/source/numerics/data_out_rotation.cc @@ -524,6 +524,8 @@ DataOutRotation::build_patches( for (unsigned int i = 0; i < new_patches.size(); ++i) { new_patches[i].n_subdivisions = n_subdivisions; + new_patches[i].reference_cell = ReferenceCells::get_hypercube(); + new_patches[i].data.reinit( n_datasets, Utilities::fixed_power(n_subdivisions + 1)); } diff --git a/source/numerics/data_out_stack.cc b/source/numerics/data_out_stack.cc index f2011bb21f..68a871c286 100644 --- a/source/numerics/data_out_stack.cc +++ b/source/numerics/data_out_stack.cc @@ -310,6 +310,7 @@ DataOutStack::build_patches( // the time direction) points dealii::DataOutBase::Patch default_patch; default_patch.n_subdivisions = n_subdivisions; + default_patch.reference_cell = ReferenceCells::get_hypercube(); default_patch.data.reinit(n_datasets, n_q_points * (n_subdivisions + 1)); patches.insert(patches.end(), n_patches, default_patch); diff --git a/tests/base/functions_04.cc b/tests/base/functions_04.cc index 2b2b1edd46..87ee67901d 100644 --- a/tests/base/functions_04.cc +++ b/tests/base/functions_04.cc @@ -61,6 +61,7 @@ check_function(const Functions::FlowFunction &f, patches[0].neighbors[i] = numbers::invalid_unsigned_int; patches[0].patch_index = 0; patches[0].n_subdivisions = sub; + patches[0].reference_cell = ReferenceCells::get_hypercube(); patches[0].points_are_available = false; vertex_number = 1; diff --git a/tests/data_out/data_out_base.cc b/tests/data_out/data_out_base.cc index ff7e7fa1f5..3a6dfa1274 100644 --- a/tests/data_out/data_out_base.cc +++ b/tests/data_out/data_out_base.cc @@ -161,6 +161,7 @@ create_patches(std::vector> &patches) DataOutBase::Patch &p = patches[c]; p.patch_index = c; p.n_subdivisions = nsub; + p.reference_cell = ReferenceCells::get_hypercube(); for (unsigned int i = 0; i < ncells; ++i) for (unsigned int j = 0; j < spacedim; ++j) diff --git a/tests/data_out/data_out_hdf5_02.cc b/tests/data_out/data_out_hdf5_02.cc index 2fea33737f..4573f9099e 100644 --- a/tests/data_out/data_out_hdf5_02.cc +++ b/tests/data_out/data_out_hdf5_02.cc @@ -41,6 +41,7 @@ create_patches(std::vector> &patches) const unsigned int nsubp = nsub + 1; patch.n_subdivisions = nsub; + patch.reference_cell = ReferenceCells::get_hypercube(); for (const unsigned int v : GeometryInfo::vertex_indices()) for (unsigned int d = 0; d < spacedim; ++d) patch.vertices[v](d) = diff --git a/tests/data_out/patches.h b/tests/data_out/patches.h index d3ec419cd6..0f6269e12f 100644 --- a/tests/data_out/patches.h +++ b/tests/data_out/patches.h @@ -35,6 +35,14 @@ create_patches(std::vector> &patches) const unsigned int nsubp = nsub + 1; patch.n_subdivisions = nsub; +#if DEAL_II_HAVE_CXX17 + if constexpr (dim > 0) + patch.reference_cell = ReferenceCells::get_hypercube(); +#else + if (dim > 0) + const_cast(patch.reference_cell) = + ReferenceCells::get_hypercube(); +#endif for (const unsigned int v : GeometryInfo::vertex_indices()) for (unsigned int d = 0; d < spacedim; ++d) patch.vertices[v](d) = @@ -98,6 +106,7 @@ create_continuous_patches(std::vector> &patches, { DataOutBase::Patch patch; patch.n_subdivisions = n_sub; + patch.reference_cell = ReferenceCells::get_hypercube(); for (unsigned int k = 0; k < trapez.size(); ++k) { Point p = trapez.point(k);