tecplot_io.C
Go to the documentation of this file.
1 // The libMesh Finite Element Library.
2 // Copyright (C) 2002-2014 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 
19 
20 // C++ includes
21 #include <fstream>
22 #include <iomanip>
23 #include <sstream>
24 
25 // Local includes
26 #include "libmesh/libmesh_config.h"
28 #include "libmesh/tecplot_io.h"
29 #include "libmesh/mesh_base.h"
30 #include "libmesh/elem.h"
31 #include "libmesh/parallel.h"
32 
33 #ifdef LIBMESH_HAVE_TECPLOT_API
34 extern "C" {
35 # include <TECIO.h>
36 }
37 #endif
38 
39 
40 namespace libMesh
41 {
42 
43 
44 //--------------------------------------------------------
45 // Macros for handling Tecplot API data
46 
47 #ifdef LIBMESH_HAVE_TECPLOT_API
48 
49 namespace
50 {
51  class TecplotMacros
52  {
53  public:
54  TecplotMacros(const unsigned int n_nodes,
55  const unsigned int n_vars,
56  const unsigned int n_cells,
57  const unsigned int n_vert);
58  float & nd(const unsigned int i, const unsigned int j);
59  int & cd(const unsigned int i, const unsigned int j);
60  std::vector<float> nodalData;
61  std::vector<int> connData;
62  //float* nodalData;
63  //int* connData;
64 
65  void set_n_cells (const unsigned int nc);
66 
67  const unsigned int n_nodes;
68  const unsigned int n_vars;
69  unsigned int n_cells;
70  const unsigned int n_vert;
71  };
72 }
73 
74 
75 
76 inline
77 TecplotMacros::TecplotMacros(const unsigned int nn,
78  const unsigned int nvar,
79  const unsigned int nc,
80  const unsigned int nvrt) :
81  n_nodes(nn),
82  n_vars(nvar),
83  n_cells(nc),
84  n_vert(nvrt)
85 {
86  nodalData.resize(n_nodes*n_vars);
87  connData.resize(n_cells*n_vert);
88 }
89 
90 
91 
92 inline
93 float & TecplotMacros::nd(const unsigned int i, const unsigned int j)
94 {
95  return nodalData[(i)*(n_nodes) + (j)];
96 }
97 
98 
99 
100 inline
101 int & TecplotMacros::cd(const unsigned int i, const unsigned int j)
102 {
103  return connData[(i) + (j)*(n_vert)];
104 }
105 
106 
107 inline
108 void TecplotMacros::set_n_cells (const unsigned int nc)
109 {
110  n_cells = nc;
111  connData.resize(n_cells*n_vert);
112 }
113 #endif
114 //--------------------------------------------------------
115 
116 
117 
118 // ------------------------------------------------------------
119 // TecplotIO members
120 TecplotIO::TecplotIO (const MeshBase& mesh_in,
121  const bool binary_in,
122  const double time_in,
123  const int strand_offset_in) :
124  MeshOutput<MeshBase> (mesh_in),
125  _binary (binary_in),
126  _time (time_in),
127  _strand_offset (strand_offset_in),
128  _zone_title ("zone")
129 {
130  // Gather a list of subdomain ids in the mesh.
131  // We must do this now, while we have every
132  // processor's attention
133  // (some of the write methods only execute on processor 0).
134  mesh_in.subdomain_ids (_subdomain_ids);
135 }
136 
137 
138 
139 void TecplotIO::write (const std::string& fname)
140 {
141  if (this->mesh().processor_id() == 0)
142  {
143  if (this->binary())
144  this->write_binary (fname);
145  else
146  this->write_ascii (fname);
147  }
148 }
149 
150 
151 
152 void TecplotIO::write_nodal_data (const std::string& fname,
153  const std::vector<Number>& soln,
154  const std::vector<std::string>& names)
155 {
156  START_LOG("write_nodal_data()", "TecplotIO");
157 
158  if (this->mesh().processor_id() == 0)
159  {
160  if (this->binary())
161  this->write_binary (fname, &soln, &names);
162  else
163  this->write_ascii (fname, &soln, &names);
164  }
165 
166  STOP_LOG("write_nodal_data()", "TecplotIO");
167 }
168 
169 
170 
171 void TecplotIO::write_ascii (const std::string& fname,
172  const std::vector<Number>* v,
173  const std::vector<std::string>* solution_names)
174 {
175  // Should only do this on processor 0!
176  libmesh_assert_equal_to (this->mesh().processor_id(), 0);
177 
178  // Create an output stream
179  std::ofstream out_stream(fname.c_str());
180 
181  // Make sure it opened correctly
182  if (!out_stream.good())
183  libmesh_file_error(fname.c_str());
184 
185  // Get a constant reference to the mesh.
186  const MeshBase& the_mesh = MeshOutput<MeshBase>::mesh();
187 
188  // Write header to stream
189  {
190  {
191  // TODO: We used to print out the SVN revision here when we did keyword expansions...
192  out_stream << "# For a description of the Tecplot format see the Tecplot User's guide.\n"
193  << "#\n";
194  }
195 
196  out_stream << "Variables=x,y,z";
197 
198  if (solution_names != NULL)
199  for (unsigned int n=0; n<solution_names->size(); n++)
200  {
201 #ifdef LIBMESH_USE_REAL_NUMBERS
202 
203  // Write variable names for real variables
204  out_stream << "," << (*solution_names)[n];
205 
206 #else
207 
208  // Write variable names for complex variables
209  out_stream << "," << "r_" << (*solution_names)[n]
210  << "," << "i_" << (*solution_names)[n]
211  << "," << "a_" << (*solution_names)[n];
212 
213 #endif
214  }
215 
216  out_stream << '\n';
217 
218  out_stream << "Zone f=fepoint, n=" << the_mesh.n_nodes() << ", e=" << the_mesh.n_active_sub_elem();
219 
220  if (the_mesh.mesh_dimension() == 1)
221  out_stream << ", et=lineseg";
222  else if (the_mesh.mesh_dimension() == 2)
223  out_stream << ", et=quadrilateral";
224  else if (the_mesh.mesh_dimension() == 3)
225  out_stream << ", et=brick";
226  else
227  {
228  // Dimension other than 1, 2, or 3?
229  libmesh_error();
230  }
231 
232  // Use default mesh color = black
233  out_stream << ", c=black\n";
234 
235  } // finished writing header
236 
237  for (unsigned int i=0; i<the_mesh.n_nodes(); i++)
238  {
239  // Print the point without a newline
240  the_mesh.point(i).write_unformatted(out_stream, false);
241 
242  if ((v != NULL) && (solution_names != NULL))
243  {
244  const std::size_t n_vars = solution_names->size();
245 
246 
247  for (std::size_t c=0; c<n_vars; c++)
248  {
249 #ifdef LIBMESH_USE_REAL_NUMBERS
250  // Write real data
251  out_stream << std::setprecision(this->ascii_precision())
252  << (*v)[i*n_vars + c] << " ";
253 
254 #else
255  // Write complex data
256  out_stream << std::setprecision(this->ascii_precision())
257  << (*v)[i*n_vars + c].real() << " "
258  << (*v)[i*n_vars + c].imag() << " "
259  << std::abs((*v)[i*n_vars + c]) << " ";
260 
261 #endif
262  }
263  }
264 
265  // Write a new line after the data for this node
266  out_stream << '\n';
267  }
268 
271 
272  for ( ; it != end; ++it)
273  (*it)->write_connectivity(out_stream, TECPLOT);
274 }
275 
276 
277 
278 void TecplotIO::write_binary (const std::string& fname,
279  const std::vector<Number>* vec,
280  const std::vector<std::string>* solution_names)
281 {
282  //-----------------------------------------------------------
283  // Call the ASCII output function if configure did not detect
284  // the Tecplot binary API
285 #ifndef LIBMESH_HAVE_TECPLOT_API
286 
287  libMesh::err << "WARNING: Tecplot Binary files require the Tecplot API." << std::endl
288  << "Continuing with ASCII output."
289  << std::endl;
290 
291  if (this->mesh().processor_id() == 0)
292  this->write_ascii (fname, vec, solution_names);
293  return;
294 
295 
296 
297  //------------------------------------------------------------
298  // New binary formats, time aware and whatnot
299 #elif defined(LIBMESH_HAVE_TECPLOT_API_112)
300 
301  // Get a constant reference to the mesh.
302  const MeshBase& the_mesh = MeshOutput<MeshBase>::mesh();
303 
304  // Required variables
305  std::string tecplot_variable_names;
306  int
307  ierr = 0,
308  file_type = 0, // full
309  is_double = 0,
310 #ifdef DEBUG
311  tec_debug = 1,
312 #else
313  tec_debug = 0,
314 #endif
315  cell_type = -1,
316  nn_per_elem = -1;
317 
318  switch (the_mesh.mesh_dimension())
319  {
320  case 1:
321  cell_type = 1; // FELINESEG
322  nn_per_elem = 2;
323  break;
324 
325  case 2:
326  cell_type = 3; // FEQUADRILATERAL
327  nn_per_elem = 4;
328  break;
329 
330  case 3:
331  cell_type = 5; // FEBRICK
332  nn_per_elem = 8;
333  break;
334 
335  default:
336  libmesh_error();
337  }
338 
339  // Build a string containing all the variable names to pass to Tecplot
340  {
341  tecplot_variable_names += "x, y, z";
342 
343  if (solution_names != NULL)
344  {
345  for (unsigned int name=0; name<solution_names->size(); name++)
346  {
347 #ifdef LIBMESH_USE_REAL_NUMBERS
348 
349  tecplot_variable_names += ", ";
350  tecplot_variable_names += (*solution_names)[name];
351 
352 #else
353 
354  tecplot_variable_names += ", ";
355  tecplot_variable_names += "r_";
356  tecplot_variable_names += (*solution_names)[name];
357  tecplot_variable_names += ", ";
358  tecplot_variable_names += "i_";
359  tecplot_variable_names += (*solution_names)[name];
360  tecplot_variable_names += ", ";
361  tecplot_variable_names += "a_";
362  tecplot_variable_names += (*solution_names)[name];
363 
364 #endif
365  }
366  }
367  }
368 
369  // Instantiate a TecplotMacros interface. In 2D the most nodes per
370  // face should be 4, in 3D it's 8.
371 
372 
373  TecplotMacros tm(the_mesh.n_nodes(),
374 #ifdef LIBMESH_USE_REAL_NUMBERS
375  (3 + ((solution_names == NULL) ? 0 : solution_names->size())),
376 #else
377  (3 + 3*((solution_names == NULL) ? 0 : solution_names->size())),
378 #endif
379  the_mesh.n_active_sub_elem(),
380  nn_per_elem
381  );
382 
383 
384  // Copy the nodes and data to the TecplotMacros class. Note that we store
385  // everything as a float here since the eye doesn't require a double to
386  // understand what is going on
387  for (unsigned int v=0; v<the_mesh.n_nodes(); v++)
388  {
389  tm.nd(0,v) = static_cast<float>(the_mesh.point(v)(0));
390  tm.nd(1,v) = static_cast<float>(the_mesh.point(v)(1));
391  tm.nd(2,v) = static_cast<float>(the_mesh.point(v)(2));
392 
393  if ((vec != NULL) &&
394  (solution_names != NULL))
395  {
396  const unsigned int n_vars = solution_names->size();
397 
398  for (unsigned int c=0; c<n_vars; c++)
399  {
400 #ifdef LIBMESH_USE_REAL_NUMBERS
401 
402  tm.nd((3+c),v) = static_cast<float>((*vec)[v*n_vars + c]);
403 #else
404  tm.nd((3+3*c),v) = static_cast<float>((*vec)[v*n_vars + c].real());
405  tm.nd((3+3*c+1),v) = static_cast<float>((*vec)[v*n_vars + c].imag());
406  tm.nd((3+3*c+2),v) = static_cast<float>(std::abs((*vec)[v*n_vars + c]));
407 #endif
408  }
409  }
410  }
411 
412 
413  // Initialize the file
414  ierr = TECINI112 (NULL,
415  (char*) tecplot_variable_names.c_str(),
416  (char*) fname.c_str(),
417  (char*) ".",
418  &file_type,
419  &tec_debug,
420  &is_double);
421 
422  libmesh_assert_equal_to (ierr, 0);
423 
424  // A zone for each subdomain
425  bool firstzone=true;
426  for (std::set<subdomain_id_type>::const_iterator sbd_it=_subdomain_ids.begin();
427  sbd_it!=_subdomain_ids.end(); ++sbd_it)
428  {
429  // Copy the connectivity for this subdomain
430  {
433 
434  unsigned int n_subcells_in_subdomain=0;
435 
436  for (; it != end; ++it)
437  n_subcells_in_subdomain += (*it)->n_sub_elem();
438 
439  // update the connectivty array to include only the elements in this subdomain
440  tm.set_n_cells (n_subcells_in_subdomain);
441 
442  unsigned int te = 0;
443 
444  for (it = the_mesh.active_subdomain_elements_begin (*sbd_it);
445  it != end; ++it)
446  {
447  std::vector<dof_id_type> conn;
448  for (unsigned int se=0; se<(*it)->n_sub_elem(); se++)
449  {
450  (*it)->connectivity(se, TECPLOT, conn);
451 
452  for (unsigned int node=0; node<conn.size(); node++)
453  tm.cd(node,te) = conn[node];
454 
455  te++;
456  }
457  }
458  }
459 
460 
461  // Ready to call the Tecplot API for this subdomain
462  {
463  int
464  num_nodes = static_cast<int>(the_mesh.n_nodes()),
465  num_cells = static_cast<int>(tm.n_cells),
466  num_faces = 0,
467  i_cell_max = 0,
468  j_cell_max = 0,
469  k_cell_max = 0,
470  strand_id = std::max(*sbd_it,static_cast<subdomain_id_type>(1)) + this->strand_offset(),
471  parent_zone = 0,
472  is_block = 1,
473  num_face_connect = 0,
474  face_neighbor_mode = 0,
475  tot_num_face_nodes = 0,
476  num_connect_boundary_faces = 0,
477  tot_num_boundary_connect = 0,
478  share_connect_from_zone=0;
479 
480  std::vector<int>
481  passive_var_list (tm.n_vars, 0),
482  share_var_from_zone (tm.n_vars, 1); // We only write data for the first zone, all other
483  // zones will share from this one.
484 
485  // get the subdomain name from libMesh, if there is one.
486  std::string subdomain_name = the_mesh.subdomain_name(*sbd_it);
487  std::ostringstream zone_name;
488  zone_name << this->zone_title();
489 
490  // We will title this
491  // "{zone_title()}_{subdomain_name}", or
492  // "{zone_title()}_{subdomain_id}", or
493  // "{zone_title()}"
494  if (subdomain_name.size())
495  {
496  zone_name << "_";
497  zone_name << subdomain_name;
498  }
499  else if (_subdomain_ids.size() > 1)
500  {
501  zone_name << "_";
502  zone_name << *sbd_it;
503  }
504 
505  ierr = TECZNE112 ((char*) zone_name.str().c_str(),
506  &cell_type,
507  &num_nodes,
508  &num_cells,
509  &num_faces,
510  &i_cell_max,
511  &j_cell_max,
512  &k_cell_max,
513  &_time,
514  &strand_id,
515  &parent_zone,
516  &is_block,
517  &num_face_connect,
518  &face_neighbor_mode,
519  &tot_num_face_nodes,
520  &num_connect_boundary_faces,
521  &tot_num_boundary_connect,
522  &passive_var_list[0],
523  NULL, // = all are node centered
524  (firstzone) ? NULL : &share_var_from_zone[0],
525  &share_connect_from_zone);
526 
527  libmesh_assert_equal_to (ierr, 0);
528 
529  // Write *all* the data for the first zone, then share it with the others
530  if (firstzone)
531  {
532  int total =
533 #ifdef LIBMESH_USE_REAL_NUMBERS
534  ((3 + ((solution_names == NULL) ? 0 : solution_names->size()))*num_nodes);
535 #else
536  ((3 + 3*((solution_names == NULL) ? 0 : solution_names->size()))*num_nodes);
537 #endif
538 
539 
540  ierr = TECDAT112 (&total,
541  &tm.nodalData[0],
542  &is_double);
543 
544  libmesh_assert_equal_to (ierr, 0);
545  }
546 
547  // Write the connectivity
548  ierr = TECNOD112 (&tm.connData[0]);
549 
550  libmesh_assert_equal_to (ierr, 0);
551  }
552 
553  firstzone = false;
554  }
555 
556  // Done, close the file.
557  ierr = TECEND112 ();
558 
559  libmesh_assert_equal_to (ierr, 0);
560 
561 
562 
563 
564  //------------------------------------------------------------
565  // Legacy binary format
566 #else
567 
568  // Get a constant reference to the mesh.
569  const MeshBase& the_mesh = MeshOutput<MeshBase>::mesh();
570 
571  // Tecplot binary output only good for dim=2,3
572  if (the_mesh.mesh_dimension() == 1)
573  {
574  this->write_ascii (fname, vec, solution_names);
575 
576  return;
577  }
578 
579  // Required variables
580  std::string tecplot_variable_names;
581  int is_double = 0,
582  tec_debug = 0,
583  cell_type = ((the_mesh.mesh_dimension()==2) ? (1) : (3));
584 
585  // Build a string containing all the variable names to pass to Tecplot
586  {
587  tecplot_variable_names += "x, y, z";
588 
589  if (solution_names != NULL)
590  {
591  for (unsigned int name=0; name<solution_names->size(); name++)
592  {
593 #ifdef LIBMESH_USE_REAL_NUMBERS
594 
595  tecplot_variable_names += ", ";
596  tecplot_variable_names += (*solution_names)[name];
597 
598 #else
599 
600  tecplot_variable_names += ", ";
601  tecplot_variable_names += "r_";
602  tecplot_variable_names += (*solution_names)[name];
603  tecplot_variable_names += ", ";
604  tecplot_variable_names += "i_";
605  tecplot_variable_names += (*solution_names)[name];
606  tecplot_variable_names += ", ";
607  tecplot_variable_names += "a_";
608  tecplot_variable_names += (*solution_names)[name];
609 
610 #endif
611  }
612  }
613  }
614 
615  // Instantiate a TecplotMacros interface. In 2D the most nodes per
616  // face should be 4, in 3D it's 8.
617 
618 
619  TecplotMacros tm(the_mesh.n_nodes(),
620 #ifdef LIBMESH_USE_REAL_NUMBERS
621  (3 + ((solution_names == NULL) ? 0 : solution_names->size())),
622 #else
623  (3 + 3*((solution_names == NULL) ? 0 : solution_names->size())),
624 #endif
625  the_mesh.n_active_sub_elem(),
626  ((the_mesh.mesh_dimension() == 2) ? 4 : 8)
627  );
628 
629 
630  // Copy the nodes and data to the TecplotMacros class. Note that we store
631  // everything as a float here since the eye doesn't require a double to
632  // understand what is going on
633  for (unsigned int v=0; v<the_mesh.n_nodes(); v++)
634  {
635  tm.nd(0,v) = static_cast<float>(the_mesh.point(v)(0));
636  tm.nd(1,v) = static_cast<float>(the_mesh.point(v)(1));
637  tm.nd(2,v) = static_cast<float>(the_mesh.point(v)(2));
638 
639  if ((vec != NULL) &&
640  (solution_names != NULL))
641  {
642  const unsigned int n_vars = solution_names->size();
643 
644  for (unsigned int c=0; c<n_vars; c++)
645  {
646 #ifdef LIBMESH_USE_REAL_NUMBERS
647 
648  tm.nd((3+c),v) = static_cast<float>((*vec)[v*n_vars + c]);
649 #else
650  tm.nd((3+3*c),v) = static_cast<float>((*vec)[v*n_vars + c].real());
651  tm.nd((3+3*c+1),v) = static_cast<float>((*vec)[v*n_vars + c].imag());
652  tm.nd((3+3*c+2),v) = static_cast<float>(std::abs((*vec)[v*n_vars + c]));
653 #endif
654  }
655  }
656  }
657 
658 
659  // Copy the connectivity
660  {
661  unsigned int te = 0;
662 
665 
666  for ( ; it != end; ++it)
667  {
668  std::vector<dof_id_type> conn;
669  for (unsigned int se=0; se<(*it)->n_sub_elem(); se++)
670  {
671  (*it)->connectivity(se, TECPLOT, conn);
672 
673  for (unsigned int node=0; node<conn.size(); node++)
674  tm.cd(node,te) = conn[node];
675 
676  te++;
677  }
678  }
679  }
680 
681 
682  // Ready to call the Tecplot API
683  {
684  int ierr = 0,
685  num_nodes = static_cast<int>(the_mesh.n_nodes()),
686  num_cells = static_cast<int>(the_mesh.n_active_sub_elem());
687 
688 
689  ierr = TECINI (NULL,
690  (char*) tecplot_variable_names.c_str(),
691  (char*) fname.c_str(),
692  (char*) ".",
693  &tec_debug,
694  &is_double);
695 
696  libmesh_assert_equal_to (ierr, 0);
697 
698  ierr = TECZNE (NULL,
699  &num_nodes,
700  &num_cells,
701  &cell_type,
702  (char*) "FEBLOCK",
703  NULL);
704 
705  libmesh_assert_equal_to (ierr, 0);
706 
707 
708  int total =
709 #ifdef LIBMESH_USE_REAL_NUMBERS
710  ((3 + ((solution_names == NULL) ? 0 : solution_names->size()))*num_nodes);
711 #else
712  ((3 + 3*((solution_names == NULL) ? 0 : solution_names->size()))*num_nodes);
713 #endif
714 
715 
716  ierr = TECDAT (&total,
717  &tm.nodalData[0],
718  &is_double);
719 
720  libmesh_assert_equal_to (ierr, 0);
721 
722  ierr = TECNOD (&tm.connData[0]);
723 
724  libmesh_assert_equal_to (ierr, 0);
725 
726  ierr = TECEND ();
727 
728  libmesh_assert_equal_to (ierr, 0);
729  }
730 
731 #endif
732 }
733 
734 } // namespace libMesh

Site Created By: libMesh Developers
Last modified: February 07 2014 16:57:07 UTC

Hosted By:
SourceForge.net Logo