libMesh
gnuplot_io.C
Go to the documentation of this file.
1 // The libMesh Finite Element Library.
2 // Copyright (C) 2002-2017 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 <fstream>
20 #include <sstream>
21 #include <map>
22 
23 // Local includes
24 #include "libmesh/elem.h"
25 #include "libmesh/libmesh_logging.h"
26 #include "libmesh/mesh_base.h"
27 #include "libmesh/gnuplot_io.h"
28 
29 namespace libMesh
30 {
31 
33  const std::string & title,
34  int mesh_properties)
35  :
36  MeshOutput<MeshBase> (mesh_in),
37  _title(title)
38 {
39  _grid = (mesh_properties & GRID_ON);
40  _png_output = (mesh_properties & PNG_OUTPUT);
41 }
42 
43 void GnuPlotIO::write(const std::string & fname)
44 {
45  this->write_solution(fname);
46 }
47 
48 void GnuPlotIO::write_nodal_data (const std::string & fname,
49  const std::vector<Number> & soln,
50  const std::vector<std::string> & names)
51 {
52  LOG_SCOPE("write_nodal_data()", "GnuPlotIO");
53  this->write_solution(fname, &soln, &names);
54 }
55 
56 
57 
58 
59 void GnuPlotIO::write_solution(const std::string & fname,
60  const std::vector<Number> * soln,
61  const std::vector<std::string> * names)
62 {
63  // Even when writing on a serialized DistributedMesh, we expect
64  // non-proc-0 help with calls like n_active_elem
65  // libmesh_assert_equal_to (this->mesh().processor_id(), 0);
66 
67  const MeshBase & the_mesh = MeshOutput<MeshBase>::mesh();
68 
69  dof_id_type n_active_elem = the_mesh.n_active_elem();
70 
71  if (this->mesh().processor_id() == 0)
72  {
73  std::stringstream data_stream_name;
74  data_stream_name << fname << "_data";
75  const std::string data_file_name = data_stream_name.str();
76 
77  // This class is designed only for use with 1D meshes
78  libmesh_assert_equal_to (the_mesh.mesh_dimension(), 1);
79 
80  // Make sure we have a solution to plot
81  libmesh_assert ((names != libmesh_nullptr) && (soln != libmesh_nullptr));
82 
83  // Create an output stream for script file
84  std::ofstream out_stream(fname.c_str());
85 
86  // Make sure it opened correctly
87  if (!out_stream.good())
88  libmesh_file_error(fname.c_str());
89 
90  // The number of variables in the equation system
91  const unsigned int n_vars =
92  cast_int<unsigned int>(names->size());
93 
94  // Write header to stream
95  out_stream << "# This file was generated by gnuplot_io.C\n"
96  << "# Stores 1D solution data in GNUplot format\n"
97  << "# Execute this by loading gnuplot and typing "
98  << "\"call '" << fname << "'\"\n"
99  << "reset\n"
100  << "set title \"" << _title << "\"\n"
101  << "set xlabel \"x\"\n"
102  << "set xtics nomirror\n";
103 
104  // Loop over the elements to find the minimum and maximum x values,
105  // and also to find the element boundaries to write out as xtics
106  // if requested.
107  Real x_min=0., x_max=0.;
108 
109  // construct string for xtic positions at element edges
110  std::stringstream xtics_stream;
111 
112  unsigned int count = 0;
113 
114  for (const auto & el : the_mesh.active_element_ptr_range())
115  {
116  // if el is the left edge of the mesh, print its left node position
117  if (el->neighbor_ptr(0) == libmesh_nullptr)
118  {
119  x_min = (el->point(0))(0);
120  xtics_stream << "\"\" " << x_min << ", \\\n";
121  }
122  if (el->neighbor_ptr(1) == libmesh_nullptr)
123  {
124  x_max = (el->point(1))(0);
125  }
126  xtics_stream << "\"\" " << (el->point(1))(0);
127 
128  if (count+1 != n_active_elem)
129  {
130  xtics_stream << ", \\\n";
131  }
132  count++;
133  }
134 
135  out_stream << "set xrange [" << x_min << ":" << x_max << "]\n";
136 
137  if (_grid)
138  out_stream << "set x2tics (" << xtics_stream.str() << ")\nset grid noxtics noytics x2tics\n";
139 
140  if (_png_output)
141  {
142  out_stream << "set terminal png\n";
143  out_stream << "set output \"" << fname << ".png\"\n";
144  }
145 
146  out_stream << "plot "
147  << axes_limits
148  << " \"" << data_file_name << "\" using 1:2 title \"" << (*names)[0]
149  << "\" with lines";
150  if (n_vars > 1)
151  {
152  for (unsigned int i=1; i<n_vars; i++)
153  {
154  out_stream << ", \\\n\"" << data_file_name << "\" using 1:" << i+2
155  << " title \"" << (*names)[i] << "\" with lines";
156  }
157  }
158 
159  out_stream.close();
160 
161 
162  // Create an output stream for data file
163  std::ofstream data(data_file_name.c_str());
164 
165  if (!data.good())
166  libmesh_error_msg("ERROR: opening output data file " << data_file_name);
167 
168  // get ordered nodal data using a map
169  typedef std::pair<Real, std::vector<Number>> key_value_pair;
170  typedef std::map<Real, std::vector<Number>> map_type;
171  typedef map_type::iterator map_iterator;
172 
173  map_type node_map;
174 
175  for (const auto & elem : the_mesh.active_element_ptr_range())
176  for (unsigned int i=0; i<elem->n_nodes(); i++)
177  {
178  std::vector<Number> values;
179 
180  // Get the global id of the node
181  dof_id_type global_id = elem->node_id(i);
182 
183  for (unsigned int c=0; c<n_vars; c++)
184  values.push_back( (*soln)[global_id*n_vars + c] );
185 
186  node_map[ the_mesh.point(global_id)(0) ] = values;
187  }
188 
189  map_iterator map_it = node_map.begin();
190  const map_iterator end_map_it = node_map.end();
191 
192  for ( ; map_it != end_map_it; ++map_it)
193  {
194  key_value_pair kvp = *map_it;
195  std::vector<Number> values = kvp.second;
196 
197  data << kvp.first << "\t";
198 
199  for (std::size_t i=0; i<values.size(); i++)
200  data << values[i] << "\t";
201 
202  data << "\n";
203  }
204 
205  data.close();
206  }
207 }
208 
209 } // namespace libMesh
std::string _title
Definition: gnuplot_io.h:119
virtual dof_id_type n_active_elem() const =0
virtual const Point & point(const dof_id_type i) const =0
const class libmesh_nullptr_t libmesh_nullptr
virtual SimpleRange< element_iterator > active_element_ptr_range()=0
const MT & mesh() const
Definition: mesh_output.h:216
This class defines an abstract interface for Mesh output.
Definition: mesh_output.h:53
const unsigned int n_vars
Definition: tecplot_io.C:68
The libMesh namespace provides an interface to certain functionality in the library.
virtual void write(const std::string &) libmesh_override
Write the mesh to the specified file.
Definition: gnuplot_io.C:43
This is the MeshBase class.
Definition: mesh_base.h:68
libmesh_assert(j)
std::string axes_limits
GNUplot automatically adjusts the x and y-axes of 2D plots to "zoom in" on the data.
Definition: gnuplot_io.h:107
GnuPlotIO(const MeshBase &, const std::string &=std::string("FE 1D Solution"), int properties=0)
Constructor.
Definition: gnuplot_io.C:32
virtual void write_nodal_data(const std::string &, const std::vector< Number > &, const std::vector< std::string > &) libmesh_override
This method implements writing a mesh with nodal data to a specified file where the nodal data and va...
Definition: gnuplot_io.C:48
void write_solution(const std::string &, const std::vector< Number > *=libmesh_nullptr, const std::vector< std::string > *=libmesh_nullptr)
This method implements writing a mesh with nodal data to a specified file where the nodal data and va...
Definition: gnuplot_io.C:59
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real
unsigned int mesh_dimension() const
Definition: mesh_base.C:148
IterBase * data
Ideally this private member data should have protected access.
processor_id_type processor_id()
Definition: libmesh_base.h:96
uint8_t dof_id_type
Definition: id_types.h:64