Line data Source code
1 : // The libMesh Finite Element Library.
2 : // Copyright (C) 2002-2025 Benjamin S. Kirk, John W. Peterson, Roy H. Stogner
3 :
4 : // This library is free software; you can redistribute it and/or
5 : // modify it under the terms of the GNU Lesser General Public
6 : // License as published by the Free Software Foundation; either
7 : // version 2.1 of the License, or (at your option) any later version.
8 :
9 : // This library is distributed in the hope that it will be useful,
10 : // but WITHOUT ANY WARRANTY; without even the implied warranty of
11 : // MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
12 : // Lesser General Public License for more details.
13 :
14 : // You should have received a copy of the GNU Lesser General Public
15 : // License along with this library; if not, write to the Free Software
16 : // Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA
17 :
18 : // C++ includes
19 : #include <ctime>
20 : #include <iomanip>
21 : #include <iostream>
22 : #include <sstream>
23 :
24 : // Local includes
25 : #include "libmesh/postscript_io.h"
26 : #include "libmesh/mesh_tools.h"
27 : #include "libmesh/elem.h"
28 : #include "libmesh/int_range.h"
29 :
30 : namespace libMesh
31 : {
32 :
33 :
34 : // Transformation map between monomial (physical space) and Bezier bases.
35 : const float PostscriptIO::_bezier_transform[3][3] =
36 : {
37 : {-1.f/6.f, 1.f/6.f, 1.},
38 : {-1.f/6.f, 0.5, 1.f/6.f},
39 : {0., 1., 0.}
40 : };
41 :
42 :
43 0 : PostscriptIO::PostscriptIO (const MeshBase & mesh_in) :
44 : MeshOutput<MeshBase> (mesh_in),
45 0 : shade_value(0.0),
46 0 : line_width(0.5),
47 : //_M(3,3),
48 : _offset(0., 0.),
49 0 : _scale(1.0),
50 0 : _current_point(0., 0.)
51 : {
52 : // This code is still undergoing some development.
53 : libmesh_experimental();
54 :
55 : // Entries of transformation matrix from physical to Bezier coords.
56 : // _M(0,0) = -1./6.; _M(0,1) = 1./6.; _M(0,2) = 1.;
57 : // _M(1,0) = -1./6.; _M(1,1) = 0.5 ; _M(1,2) = 1./6.;
58 : // _M(2,0) = 0. ; _M(2,1) = 1. ; _M(2,2) = 0.;
59 :
60 : // Make sure there is enough room to store Bezier coefficients.
61 0 : _bezier_coeffs.resize(3);
62 0 : }
63 :
64 :
65 :
66 0 : PostscriptIO::~PostscriptIO () = default;
67 :
68 :
69 :
70 0 : void PostscriptIO::write (const std::string & fname)
71 : {
72 : // We may need to gather a DistributedMesh to output it, making that
73 : // const qualifier in our constructor a dirty lie
74 0 : MeshSerializer serialize(const_cast<MeshBase &>(this->mesh()), !_is_parallel_format);
75 :
76 0 : if (this->mesh().processor_id() == 0)
77 : {
78 : // Get a constant reference to the mesh.
79 0 : const MeshBase & the_mesh = MeshOutput<MeshBase>::mesh();
80 :
81 : // Only works in 2D
82 0 : libmesh_assert_equal_to (the_mesh.mesh_dimension(), 2);
83 :
84 : // Create output file stream.
85 : // _out is now a private member of the class.
86 0 : _out.open(fname.c_str());
87 :
88 : // Make sure it opened correctly
89 0 : if (!_out.good())
90 0 : libmesh_file_error(fname.c_str());
91 :
92 : // The mesh bounding box gives us info about what the
93 : // Postscript bounding box should be.
94 0 : BoundingBox bbox = MeshTools::create_bounding_box(the_mesh);
95 :
96 : // Add a little extra padding to the "true" bounding box so
97 : // that we can still see the boundary
98 0 : const Real percent_padding = 0.01;
99 0 : const Real dx=bbox.second(0)-bbox.first(0); libmesh_assert_greater (dx, 0.0);
100 0 : const Real dy=bbox.second(1)-bbox.first(1); libmesh_assert_greater (dy, 0.0);
101 :
102 0 : const Real x_min = bbox.first(0) - percent_padding*dx;
103 0 : const Real y_min = bbox.first(1) - percent_padding*dy;
104 0 : const Real x_max = bbox.second(0) + percent_padding*dx;
105 0 : const Real y_max = bbox.second(1) + percent_padding*dy;
106 :
107 : // Width of the output as given in postscript units.
108 : // This usually is given by the strange unit 1/72 inch.
109 : // A width of 300 represents a size of roughly 10 cm.
110 0 : const Real width = 300;
111 0 : _scale = width / (x_max-x_min);
112 0 : _offset(0) = x_min;
113 0 : _offset(1) = y_min;
114 :
115 : // Header writing stuff stolen from Deal.II
116 0 : std::time_t time1= std::time (0);
117 0 : std::tm * time = std::localtime(&time1);
118 : _out << "%!PS-Adobe-2.0 EPSF-1.2" << '\n'
119 : //<< "%!PS-Adobe-1.0" << '\n' // Lars' PS version
120 : << "%%Filename: " << fname << '\n'
121 : << "%%Title: LibMesh Output" << '\n'
122 : << "%%Creator: LibMesh: A C++ finite element library" << '\n'
123 0 : << "%%Creation Date: "
124 0 : << time->tm_year+1900 << "/"
125 0 : << time->tm_mon+1 << "/"
126 0 : << time->tm_mday << " - "
127 : << time->tm_hour << ":"
128 0 : << std::setw(2) << time->tm_min << ":"
129 0 : << std::setw(2) << time->tm_sec << '\n'
130 : << "%%BoundingBox: "
131 : // lower left corner
132 0 : << "0 0 "
133 : // upper right corner
134 0 : << static_cast<unsigned int>( rint(double((x_max-x_min) * _scale )))
135 0 : << ' '
136 0 : << static_cast<unsigned int>( rint(double((y_max-y_min) * _scale )))
137 0 : << '\n';
138 :
139 : // define some abbreviations to keep
140 : // the output small:
141 : // m=move turtle to
142 : // l=define a line
143 : // s=set rgb color
144 : // sg=set gray value
145 : // lx=close the line and plot the line
146 : // lf=close the line and fill the interior
147 : _out << "/m {moveto} bind def" << '\n'
148 : << "/l {lineto} bind def" << '\n'
149 : << "/s {setrgbcolor} bind def" << '\n'
150 : << "/sg {setgray} bind def" << '\n'
151 : << "/cs {curveto stroke} bind def" << '\n'
152 : << "/lx {lineto closepath stroke} bind def" << '\n'
153 0 : << "/lf {lineto closepath fill} bind def" << '\n';
154 :
155 0 : _out << "%%EndProlog" << '\n';
156 : // << '\n';
157 :
158 : // Set line width in the postscript file.
159 0 : _out << line_width << " setlinewidth" << '\n';
160 :
161 : // Set line cap and join options
162 0 : _out << "1 setlinecap" << '\n';
163 0 : _out << "1 setlinejoin" << '\n';
164 :
165 : // allow only five digits for output (instead of the default
166 : // six); this should suffice even for fine grids, but reduces
167 : // the file size significantly
168 0 : _out << std::setprecision (5);
169 :
170 : // Loop over the active elements, draw lines for the edges. We
171 : // draw even quadratic elements with straight sides, i.e. a straight
172 : // line sits between each pair of vertices. Also we draw every edge
173 : // for an element regardless of the fact that it may overlap with
174 : // another. This would probably be a useful optimization...
175 0 : for (const auto & elem : the_mesh.active_element_ptr_range())
176 0 : this->plot_linear_elem(elem);
177 :
178 : // Issue the showpage command, and we're done.
179 0 : _out << "showpage" << std::endl;
180 :
181 : } // end if (this->mesh().processor_id() == 0)
182 0 : }
183 :
184 :
185 :
186 :
187 :
188 :
189 0 : void PostscriptIO::plot_linear_elem(const Elem * elem)
190 : {
191 : // Clear the string contents. Yes, this really is how you do that...
192 0 : _cell_string.str("");
193 :
194 : // The general strategy is:
195 : // 1.) Use m := {moveto} to go to vertex 0.
196 : // 2.) Use l := {lineto} commands to draw lines to vertex 1, 2, ... N-1.
197 : // 3a.) Use lx := {lineto closepath stroke} command at vertex N to draw the last line.
198 : // 3b.) lf := {lineto closepath fill} command to shade the cell just drawn
199 : // All of our 2D elements' vertices are numbered in counterclockwise order,
200 : // so we can just draw them in the same order.
201 :
202 : // 1.)
203 0 : _current_point = (elem->point(0) - _offset) * _scale;
204 0 : _cell_string << _current_point(0) << " " << _current_point(1) << " "; // write x y
205 0 : _cell_string << "m ";
206 :
207 : // 2.)
208 0 : const unsigned int nv=elem->n_vertices();
209 0 : for (unsigned int v=1; v<nv-1; ++v)
210 : {
211 0 : _current_point = (elem->point(v) - _offset) * _scale;
212 0 : _cell_string << _current_point(0) << " " << _current_point(1) << " "; // write x y
213 0 : _cell_string << "l ";
214 : }
215 :
216 : // 3.)
217 0 : _current_point = (elem->point(nv-1) - _offset) * _scale;
218 0 : _cell_string << _current_point(0) << " " << _current_point(1) << " "; // write x y
219 :
220 : // We draw the shaded (interior) parts first, if applicable.
221 0 : if (shade_value > 0.0)
222 0 : _out << shade_value << " sg " << _cell_string.str() << "lf\n";
223 :
224 : // Draw the black lines (I guess we will always do this)
225 0 : _out << "0 sg " << _cell_string.str() << "lx\n";
226 0 : }
227 :
228 :
229 :
230 :
231 0 : void PostscriptIO::plot_quadratic_elem(const Elem * elem)
232 : {
233 0 : for (auto ns : elem->side_index_range())
234 : {
235 : // Build the quadratic side
236 0 : std::unique_ptr<const Elem> side = elem->build_side_ptr(ns);
237 :
238 : // Be sure it's quadratic (Edge2). Eventually we could
239 : // handle cubic elements as well...
240 0 : libmesh_assert_equal_to ( side->type(), EDGE3 );
241 :
242 0 : _out << "0 sg ";
243 :
244 : // Move to the first point on this side.
245 0 : _current_point = (side->point(0) - _offset) * _scale;
246 0 : _out << _current_point(0) << " " << _current_point(1) << " "; // write x y
247 0 : _out << "m ";
248 :
249 : // Compute _bezier_coeffs for this edge. This fills up
250 : // the _bezier_coeffs vector.
251 0 : this->_compute_edge_bezier_coeffs(side.get());
252 :
253 : // Print curveto path to file
254 0 : for (auto i : index_range(_bezier_coeffs))
255 0 : _out << _bezier_coeffs[i](0) << " " << _bezier_coeffs[i](1) << " ";
256 0 : _out << " cs\n";
257 0 : }
258 0 : }
259 :
260 :
261 :
262 :
263 0 : void PostscriptIO::_compute_edge_bezier_coeffs(const Elem * elem)
264 : {
265 : // I only know how to do this for an Edge3!
266 0 : libmesh_assert_equal_to (elem->type(), EDGE3);
267 :
268 : // Get x-coordinates into an array, transform them,
269 : // and repeat for y.
270 0 : float phys_coords[3] = {0., 0., 0.};
271 0 : float bez_coords[3] = {0., 0., 0.};
272 :
273 0 : for (unsigned int i=0; i<2; ++i)
274 : {
275 : // Initialize vectors. Physical coordinates are initialized
276 : // by their postscript-scaled values.
277 0 : for (unsigned int j=0; j<3; ++j)
278 : {
279 0 : phys_coords[j] = static_cast<float>
280 0 : ((elem->point(j)(i) - _offset(i)) * _scale);
281 0 : bez_coords[j] = 0.; // zero out result vector
282 : }
283 :
284 : // Multiply matrix times vector
285 0 : for (unsigned int j=0; j<3; ++j)
286 0 : for (unsigned int k=0; k<3; ++k)
287 0 : bez_coords[j] += _bezier_transform[j][k]*phys_coords[k];
288 :
289 : // Store result in _bezier_coeffs
290 0 : for (unsigned int j=0; j<3; ++j)
291 0 : _bezier_coeffs[j](i) = phys_coords[j];
292 : }
293 0 : }
294 :
295 : } // namespace libMesh
|