From: David Wells Date: Sun, 29 Apr 2018 18:08:50 +0000 (-0400) Subject: Avoid writing colinear points in GridOut::write_gnuplot. X-Git-Tag: v9.1.0-rc1~1105^2~3 X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=a2d7de113d0f57398e167878892090d00db884fe;p=dealii.git Avoid writing colinear points in GridOut::write_gnuplot. GridOut::write_gnuplot has an option where, for plotting curved grids, it can plot multiple line segments per face to give the illusion of a cell with curved faces. This procedure creates massive output files since every face of every corresponds to many individual line segments. This commit implements a compression scheme where we do not output intermediate points that are colinear. This lowers the number of points in a polar grid output file by about 50% (i.e., the radial cell faces are not curved and can be plotted with two points instead of perhaps five or six; same for the interior cell faces). --- diff --git a/source/grid/grid_out.cc b/source/grid/grid_out.cc index 42b8ad118a..3f09caab8d 100644 --- a/source/grid/grid_out.cc +++ b/source/grid/grid_out.cc @@ -3157,6 +3157,34 @@ namespace internal { namespace { + /** + * GNUPlot output can, optionally, output multiple line segments on each + * grid line to make curved lines look curved. However, this is very + * wasteful when the line itself is straight since we do not need multiple + * line segments to draw a straight line. This function tries to identify + * whether or not the collection of points corresponds to a straight line: + * if it does then all but the first and last points can be removed. + */ + template + void + remove_colinear_points(std::vector> &points) + { + while (points.size() > 2) + { + Tensor<1, spacedim> first_difference = points[1] - points[0]; + first_difference /= first_difference.norm(); + Tensor<1, spacedim> second_difference = points[2] - points[1]; + second_difference /= second_difference.norm(); + // If the three points are colinear then remove the middle one. + if ((first_difference - second_difference).norm() < 1e-10) + points.erase(points.begin() + 1); + else + break; + } + } + + + template void write_gnuplot (const dealii::Triangulation<1,spacedim> &tria, std::ostream &out, @@ -3263,12 +3291,24 @@ namespace internal face = cell->face(face_no); if (face->at_boundary() || gnuplot_flags.curved_inner_cells) { + // Save the points on each face to a vector and then try + // to remove colinear points that won't show up in the + // generated plot. + std::vector> line_points; // compute offset of quadrature points within set of // projected points const unsigned int offset=face_no*n_points; - for (unsigned int i=0; itransform_unit_to_real_cell - (cell, q_projector->point(offset+i))) + // we don't need to transform the vertices: they are + // already in physical space. + line_points.push_back(face->vertex(0)); + for (unsigned int i=1; itransform_unit_to_real_cell + (cell, q_projector->point(offset+i))); + line_points.push_back(face->vertex(1)); + internal::remove_colinear_points(line_points); + + for (const Point &point : line_points) + out << point << ' ' << cell->level() << ' ' << static_cast(cell->material_id()) << '\n'; @@ -3462,19 +3502,30 @@ namespace internal const typename dealii::Triangulation::line_iterator line=face->line(l); - const Point &v0=line->vertex(0), - &v1=line->vertex(1); + const Point &v0=line->vertex(0), &v1=line->vertex(1); if (line->at_boundary() || gnuplot_flags.curved_inner_cells) { + // Save the points on each face to a vector and + // then try to remove colinear points that won't + // show up in the generated plot. + std::vector> line_points; // transform_real_to_unit_cell could be replaced // by using QProjector::project_to_line // which is not yet implemented const Point u0=mapping->transform_real_to_unit_cell(cell, v0), u1=mapping->transform_real_to_unit_cell(cell, v1); - for (unsigned int i=0; itransform_unit_to_real_cell - (cell, (1-boundary_points[i][0])*u0+boundary_points[i][0]*u1)) + const Point center; + // we don't need to transform the vertices: they + // are already in physical space. + line_points.push_back(v0); + for (unsigned int i=1; itransform_unit_to_real_cell + (cell, (1-boundary_points[i][0])*u0 + boundary_points[i][0]*u1)); + line_points.push_back(v1); + internal::remove_colinear_points(line_points); + for (const Point &point : line_points) + out << point << ' ' << cell->level() << ' ' << static_cast(cell->material_id()) << '\n'; }