libMesh::TecplotIO Class Reference

#include <tecplot_io.h>

Inheritance diagram for libMesh::TecplotIO:

List of all members.

Public Member Functions

 TecplotIO (const MeshBase &, const bool binary=false, const double time=0., const int strand_offset=0)
virtual void write (const std::string &)
virtual void write_nodal_data (const std::string &, const std::vector< Number > &, const std::vector< std::string > &)
bool & binary ()
double & time ()
int & strand_offset ()
std::string & zone_title ()
virtual void write_equation_systems (const std::string &, const EquationSystems &, const std::set< std::string > *system_names=NULL)
unsigned int & ascii_precision ()

Protected Member Functions

const MeshBasemesh () const

Protected Attributes

const bool _is_parallel_format

Private Member Functions

void write_ascii (const std::string &, const std::vector< Number > *=NULL, const std::vector< std::string > *=NULL)
void write_binary (const std::string &, const std::vector< Number > *=NULL, const std::vector< std::string > *=NULL)

Private Attributes

bool _binary
double _time
int _strand_offset
std::string _zone_title
std::set< subdomain_id_type_subdomain_ids

Detailed Description

This class implements writing meshes in the Tecplot format.

Author:
Benjamin S. Kirk, 2004

Definition at line 47 of file tecplot_io.h.


Constructor & Destructor Documentation

libMesh::TecplotIO::TecplotIO ( const MeshBase mesh_in,
const bool  binary = false,
const double  time = 0.,
const int  strand_offset = 0 
) [explicit]

Constructor. Takes a reference to a constant mesh object. This constructor will only allow us to write the mesh. The optional parameter binary can be used to switch between ASCII (false, the default) or binary (true) output files.

Definition at line 120 of file tecplot_io.C.

References _subdomain_ids, and libMesh::MeshBase::subdomain_ids().

00123                                                   :
00124   MeshOutput<MeshBase> (mesh_in),
00125   _binary (binary_in),
00126   _time (time_in),
00127   _strand_offset (strand_offset_in),
00128   _zone_title ("zone")
00129 {
00130   // Gather a list of subdomain ids in the mesh.
00131   // We must do this now, while we have every
00132   // processor's attention
00133   // (some of the write methods only execute on processor 0).
00134   mesh_in.subdomain_ids (_subdomain_ids);
00135 }


Member Function Documentation

unsigned int& libMesh::MeshOutput< MeshBase >::ascii_precision (  )  [inherited]

Return/set the precision to use when writing ASCII files.

By default we use numeric_limits<Real>::digits10 + 2, which should be enough to write out to ASCII and get the exact same Real back when reading in.

Referenced by write_ascii(), libMesh::GMVIO::write_ascii_new_impl(), and libMesh::GMVIO::write_ascii_old_impl().

bool & libMesh::TecplotIO::binary (  )  [inline]

Flag indicating whether or not to write a binary file (if the tecio.a library was found by configure).

Definition at line 156 of file tecplot_io.h.

References _binary.

Referenced by write(), and write_nodal_data().

00157 {
00158   return _binary;
00159 }

int & libMesh::TecplotIO::strand_offset (  )  [inline]

Strand offset for this file. Each mesh block will be written to (strand_id=block_id+1+strand_offset). Written to newer binary formats that are time-aware, defaults to 0.

Definition at line 172 of file tecplot_io.h.

References _strand_offset.

Referenced by write_binary().

00173 {
00174   return _strand_offset;
00175 }

double & libMesh::TecplotIO::time (  )  [inline]

Solution time for transient data. Written to newer binary formats that are time-aware.

Definition at line 164 of file tecplot_io.h.

References _time.

00165 {
00166   return _time;
00167 }

void libMesh::TecplotIO::write ( const std::string &  fname  )  [virtual]

This method implements writing a mesh to a specified file.

Implements libMesh::MeshOutput< MeshBase >.

Definition at line 139 of file tecplot_io.C.

References binary(), libMesh::processor_id(), write_ascii(), and write_binary().

00140 {
00141   if (libMesh::processor_id() == 0)
00142     {
00143       if (this->binary())
00144         this->write_binary (fname);
00145       else
00146         this->write_ascii  (fname);
00147     }
00148 }

void libMesh::TecplotIO::write_ascii ( const std::string &  fname,
const std::vector< Number > *  v = NULL,
const std::vector< std::string > *  solution_names = NULL 
) [private]

This method implements writing a mesh with nodal data to a specified file where the nodal data and variable names are optionally provided. This will write an ASCII file.

Definition at line 171 of file tecplot_io.C.

References std::abs(), libMesh::MeshBase::active_elements_begin(), libMesh::MeshBase::active_elements_end(), libMesh::MeshOutput< MeshBase >::ascii_precision(), end, libMesh::MeshOutput< MeshBase >::mesh(), libMesh::MeshBase::mesh_dimension(), libMesh::MeshBase::n_active_sub_elem(), libMesh::MeshBase::n_nodes(), n_vars, libMesh::MeshBase::point(), libMesh::processor_id(), libMeshEnums::TECPLOT, and libMesh::TypeVector< T >::write_unformatted().

Referenced by write(), write_binary(), and write_nodal_data().

00174 {
00175   // Should only do this on processor 0!
00176   libmesh_assert_equal_to (libMesh::processor_id(), 0);
00177 
00178   // Create an output stream
00179   std::ofstream out_stream(fname.c_str());
00180 
00181   // Make sure it opened correctly
00182   if (!out_stream.good())
00183     libmesh_file_error(fname.c_str());
00184 
00185   // Get a constant reference to the mesh.
00186   const MeshBase& the_mesh = MeshOutput<MeshBase>::mesh();
00187 
00188   // Write header to stream
00189   {
00190     {
00191       // TODO: We used to print out the SVN revision here when we did keyword expansions...
00192       out_stream << "# For a description of the Tecplot format see the Tecplot User's guide.\n"
00193           << "#\n";
00194     }
00195 
00196     out_stream << "Variables=x,y,z";
00197 
00198     if (solution_names != NULL)
00199       for (unsigned int n=0; n<solution_names->size(); n++)
00200         {
00201 #ifdef LIBMESH_USE_REAL_NUMBERS
00202 
00203           // Write variable names for real variables
00204           out_stream << "," << (*solution_names)[n];
00205 
00206 #else
00207 
00208           // Write variable names for complex variables
00209           out_stream << "," << "r_"   << (*solution_names)[n]
00210               << "," << "i_"   << (*solution_names)[n]
00211               << "," << "a_"   << (*solution_names)[n];
00212 
00213 #endif
00214         }
00215 
00216     out_stream << '\n';
00217 
00218     out_stream << "Zone f=fepoint, n=" << the_mesh.n_nodes() << ", e=" << the_mesh.n_active_sub_elem();
00219 
00220     if (the_mesh.mesh_dimension() == 1)
00221       out_stream << ", et=lineseg";
00222     else if (the_mesh.mesh_dimension() == 2)
00223       out_stream << ", et=quadrilateral";
00224     else if (the_mesh.mesh_dimension() == 3)
00225       out_stream << ", et=brick";
00226     else
00227       {
00228         // Dimension other than 1, 2, or 3?
00229         libmesh_error();
00230       }
00231 
00232     // Use default mesh color = black
00233     out_stream << ", c=black\n";
00234 
00235   } // finished writing header
00236 
00237   for (unsigned int i=0; i<the_mesh.n_nodes(); i++)
00238     {
00239       // Print the point without a newline
00240       the_mesh.point(i).write_unformatted(out_stream, false);
00241 
00242       if ((v != NULL) && (solution_names != NULL))
00243         {
00244           const std::size_t n_vars = solution_names->size();
00245 
00246 
00247           for (std::size_t c=0; c<n_vars; c++)
00248             {
00249 #ifdef LIBMESH_USE_REAL_NUMBERS
00250               // Write real data
00251               out_stream << std::setprecision(this->ascii_precision())
00252                   << (*v)[i*n_vars + c] << " ";
00253 
00254 #else
00255               // Write complex data
00256               out_stream << std::setprecision(this->ascii_precision())
00257                   << (*v)[i*n_vars + c].real() << " "
00258                   << (*v)[i*n_vars + c].imag() << " "
00259                   << std::abs((*v)[i*n_vars + c]) << " ";
00260 
00261 #endif
00262             }
00263         }
00264 
00265       // Write a new line after the data for this node
00266       out_stream << '\n';
00267     }
00268 
00269   MeshBase::const_element_iterator       it  = the_mesh.active_elements_begin();
00270   const MeshBase::const_element_iterator end = the_mesh.active_elements_end();
00271 
00272   for ( ; it != end; ++it)
00273     (*it)->write_connectivity(out_stream, TECPLOT);
00274 }

void libMesh::TecplotIO::write_binary ( const std::string &  fname,
const std::vector< Number > *  vec = NULL,
const std::vector< std::string > *  solution_names = NULL 
) [private]

This method implements writing a mesh with nodal data to a specified file where the nodal data and variable names are optionally provided. This will write a binary file if the tecio.a library was found at compile time, otherwise a warning message will be printed and an ASCII file will be created.

Definition at line 278 of file tecplot_io.C.

References _subdomain_ids, _time, std::abs(), libMesh::MeshBase::active_elements_begin(), libMesh::MeshBase::active_elements_end(), libMesh::MeshBase::active_subdomain_elements_begin(), libMesh::MeshBase::active_subdomain_elements_end(), end, libMesh::err, std::max(), libMesh::MeshOutput< MeshBase >::mesh(), libMesh::MeshBase::mesh_dimension(), libMesh::MeshBase::n_active_sub_elem(), libMesh::MeshBase::n_nodes(), n_vars, libMesh::Quality::name(), libMesh::MeshBase::point(), libMesh::processor_id(), strand_offset(), libMesh::MeshBase::subdomain_name(), libMeshEnums::TECPLOT, write_ascii(), and zone_title().

Referenced by write(), and write_nodal_data().

00281 {
00282   //-----------------------------------------------------------
00283   // Call the ASCII output function if configure did not detect
00284   // the Tecplot binary API
00285 #ifndef LIBMESH_HAVE_TECPLOT_API
00286 
00287   libMesh::err << "WARNING: Tecplot Binary files require the Tecplot API." << std::endl
00288                << "Continuing with ASCII output."
00289                << std::endl;
00290   
00291   if (libMesh::processor_id() == 0)
00292     this->write_ascii (fname, vec, solution_names);
00293   return;
00294 
00295 
00296   
00297   //------------------------------------------------------------
00298   // New binary formats, time aware and whatnot
00299 #elif defined(LIBMESH_HAVE_TECPLOT_API_112)
00300   
00301   // Get a constant reference to the mesh.
00302   const MeshBase& the_mesh = MeshOutput<MeshBase>::mesh();
00303     
00304   // Required variables
00305   std::string tecplot_variable_names;
00306   int
00307     ierr      =  0,
00308     file_type =  0, // full
00309     is_double =  0,
00310 #ifdef DEBUG
00311     tec_debug =  1,
00312 #else
00313     tec_debug =  0,
00314 #endif
00315     cell_type   = -1,
00316     nn_per_elem = -1;
00317 
00318   switch (the_mesh.mesh_dimension())
00319     {
00320     case 1:
00321       cell_type   = 1;  // FELINESEG
00322       nn_per_elem = 2;
00323       break;
00324       
00325     case 2:
00326       cell_type   = 3; // FEQUADRILATERAL
00327       nn_per_elem = 4;
00328       break;
00329       
00330     case 3:
00331       cell_type   = 5; // FEBRICK
00332       nn_per_elem = 8;
00333       break;
00334 
00335     default:
00336       libmesh_error();
00337     }
00338 
00339   // Build a string containing all the variable names to pass to Tecplot
00340   {
00341     tecplot_variable_names += "x, y, z";
00342 
00343     if (solution_names != NULL)
00344       {
00345         for (unsigned int name=0; name<solution_names->size(); name++)
00346           {
00347 #ifdef LIBMESH_USE_REAL_NUMBERS
00348 
00349             tecplot_variable_names += ", ";
00350             tecplot_variable_names += (*solution_names)[name];
00351 
00352 #else
00353 
00354             tecplot_variable_names += ", ";
00355             tecplot_variable_names += "r_";
00356             tecplot_variable_names += (*solution_names)[name];
00357             tecplot_variable_names += ", ";
00358             tecplot_variable_names += "i_";
00359             tecplot_variable_names += (*solution_names)[name];
00360             tecplot_variable_names += ", ";
00361             tecplot_variable_names += "a_";
00362             tecplot_variable_names += (*solution_names)[name];
00363 
00364 #endif
00365           }
00366       }
00367   }
00368 
00369   // Instantiate a TecplotMacros interface.  In 2D the most nodes per
00370   // face should be 4, in 3D it's 8.
00371 
00372 
00373   TecplotMacros tm(the_mesh.n_nodes(),
00374 #ifdef LIBMESH_USE_REAL_NUMBERS
00375                    (3 + ((solution_names == NULL) ? 0 : solution_names->size())),
00376 #else
00377                    (3 + 3*((solution_names == NULL) ? 0 : solution_names->size())),
00378 #endif
00379                    the_mesh.n_active_sub_elem(),
00380                    nn_per_elem
00381                    );
00382 
00383 
00384   // Copy the nodes and data to the TecplotMacros class. Note that we store
00385   // everything as a float here since the eye doesn't require a double to
00386   // understand what is going on
00387   for (unsigned int v=0; v<the_mesh.n_nodes(); v++)
00388     {
00389       tm.nd(0,v) = static_cast<float>(the_mesh.point(v)(0));
00390       tm.nd(1,v) = static_cast<float>(the_mesh.point(v)(1));
00391       tm.nd(2,v) = static_cast<float>(the_mesh.point(v)(2));
00392 
00393       if ((vec != NULL) &&
00394           (solution_names != NULL))
00395         {
00396           const unsigned int n_vars = solution_names->size();
00397 
00398           for (unsigned int c=0; c<n_vars; c++)
00399             {
00400 #ifdef LIBMESH_USE_REAL_NUMBERS
00401 
00402               tm.nd((3+c),v)     = static_cast<float>((*vec)[v*n_vars + c]);
00403 #else
00404               tm.nd((3+3*c),v)   = static_cast<float>((*vec)[v*n_vars + c].real());
00405               tm.nd((3+3*c+1),v) = static_cast<float>((*vec)[v*n_vars + c].imag());
00406               tm.nd((3+3*c+2),v) = static_cast<float>(std::abs((*vec)[v*n_vars + c]));
00407 #endif
00408             }
00409         }
00410     }
00411 
00412 
00413   // Initialize the file
00414   ierr = TECINI112 (NULL,
00415                     (char*) tecplot_variable_names.c_str(),
00416                     (char*) fname.c_str(),
00417                     (char*) ".",
00418                     &file_type,
00419                     &tec_debug,
00420                     &is_double);
00421   
00422   libmesh_assert_equal_to (ierr, 0);
00423 
00424   // A zone for each subdomain
00425   bool firstzone=true;
00426   for (std::set<subdomain_id_type>::const_iterator sbd_it=_subdomain_ids.begin();
00427        sbd_it!=_subdomain_ids.end(); ++sbd_it)
00428     {      
00429       // Copy the connectivity for this subdomain
00430       {
00431         MeshBase::const_element_iterator       it  = the_mesh.active_subdomain_elements_begin (*sbd_it);
00432         const MeshBase::const_element_iterator end = the_mesh.active_subdomain_elements_end   (*sbd_it);
00433 
00434         unsigned int n_subcells_in_subdomain=0;
00435 
00436         for (; it != end; ++it)
00437           n_subcells_in_subdomain += (*it)->n_sub_elem();
00438         
00439         // update the connectivty array to include only the elements in this subdomain
00440         tm.set_n_cells (n_subcells_in_subdomain);
00441 
00442         unsigned int te = 0;
00443           
00444         for (it  = the_mesh.active_subdomain_elements_begin (*sbd_it);
00445              it != end; ++it)
00446           {
00447             std::vector<dof_id_type> conn;
00448             for (unsigned int se=0; se<(*it)->n_sub_elem(); se++)
00449               {
00450                 (*it)->connectivity(se, TECPLOT, conn);
00451                   
00452                 for (unsigned int node=0; node<conn.size(); node++)
00453                   tm.cd(node,te) = conn[node];
00454                   
00455                 te++;
00456               }
00457           }
00458       }
00459         
00460         
00461       // Ready to call the Tecplot API for this subdomain
00462       { 
00463         int
00464           num_nodes   = static_cast<int>(the_mesh.n_nodes()),
00465           num_cells   = static_cast<int>(tm.n_cells),
00466           num_faces   = 0,
00467           i_cell_max  = 0,
00468           j_cell_max  = 0,
00469           k_cell_max  = 0,
00470           strand_id   = std::max(*sbd_it,static_cast<subdomain_id_type>(1)) + this->strand_offset(),
00471           parent_zone = 0,
00472           is_block    = 1,
00473           num_face_connect   = 0,
00474           face_neighbor_mode = 0,
00475           tot_num_face_nodes = 0,
00476           num_connect_boundary_faces = 0,
00477           tot_num_boundary_connect   = 0,
00478           share_connect_from_zone=0;
00479           
00480         std::vector<int>
00481           passive_var_list    (tm.n_vars, 0),
00482           share_var_from_zone (tm.n_vars, 1); // We only write data for the first zone, all other
00483                                               // zones will share from this one.
00484 
00485         // get the subdomain name from libMesh, if there is one.
00486         std::string subdomain_name = the_mesh.subdomain_name(*sbd_it);
00487         std::ostringstream zone_name;
00488         zone_name << this->zone_title();
00489         
00490         // We will title this
00491         // "{zone_title()}_{subdomain_name}", or
00492         // "{zone_title()}_{subdomain_id}", or
00493         // "{zone_title()}"
00494         if (subdomain_name.size())      
00495           {
00496             zone_name << "_";
00497             zone_name << subdomain_name;
00498           }
00499         else if (_subdomain_ids.size() > 1)
00500           {
00501             zone_name << "_";
00502             zone_name << *sbd_it;
00503           }
00504             
00505         ierr = TECZNE112 ((char*) zone_name.str().c_str(),
00506                           &cell_type,
00507                           &num_nodes,
00508                           &num_cells,
00509                           &num_faces,
00510                           &i_cell_max,
00511                           &j_cell_max,
00512                           &k_cell_max,
00513                           &_time,
00514                           &strand_id,
00515                           &parent_zone,
00516                           &is_block,
00517                           &num_face_connect,
00518                           &face_neighbor_mode,
00519                           &tot_num_face_nodes,
00520                           &num_connect_boundary_faces,
00521                           &tot_num_boundary_connect,
00522                           &passive_var_list[0],
00523                           NULL, // = all are node centered
00524                           (firstzone) ? NULL : &share_var_from_zone[0],
00525                           &share_connect_from_zone);
00526           
00527         libmesh_assert_equal_to (ierr, 0);
00528                   
00529         // Write *all* the data for the first zone, then share it with the others
00530         if (firstzone)
00531           {
00532             int total =
00533 #ifdef LIBMESH_USE_REAL_NUMBERS
00534               ((3 + ((solution_names == NULL) ? 0 : solution_names->size()))*num_nodes);
00535 #else
00536               ((3 + 3*((solution_names == NULL) ? 0 : solution_names->size()))*num_nodes);
00537 #endif
00538             
00539           
00540             ierr = TECDAT112 (&total,
00541                               &tm.nodalData[0],
00542                               &is_double);
00543               
00544             libmesh_assert_equal_to (ierr, 0);
00545           }
00546 
00547         // Write the connectivity
00548         ierr = TECNOD112 (&tm.connData[0]);
00549           
00550         libmesh_assert_equal_to (ierr, 0);
00551       }
00552 
00553       firstzone = false;
00554     }
00555 
00556   // Done, close the file.
00557   ierr = TECEND112 ();
00558     
00559   libmesh_assert_equal_to (ierr, 0);
00560   
00561 
00562 
00563   
00564   //------------------------------------------------------------
00565   // Legacy binary format
00566 #else
00567 
00568   // Get a constant reference to the mesh.
00569   const MeshBase& the_mesh = MeshOutput<MeshBase>::mesh();
00570 
00571   // Tecplot binary output only good for dim=2,3
00572   if (the_mesh.mesh_dimension() == 1)
00573     {
00574       this->write_ascii (fname, vec, solution_names);
00575 
00576       return;
00577     }
00578 
00579   // Required variables
00580   std::string tecplot_variable_names;
00581   int is_double =  0,
00582     tec_debug =  0,
00583     cell_type = ((the_mesh.mesh_dimension()==2) ? (1) : (3));
00584 
00585   // Build a string containing all the variable names to pass to Tecplot
00586   {
00587     tecplot_variable_names += "x, y, z";
00588 
00589     if (solution_names != NULL)
00590       {
00591         for (unsigned int name=0; name<solution_names->size(); name++)
00592           {
00593 #ifdef LIBMESH_USE_REAL_NUMBERS
00594 
00595             tecplot_variable_names += ", ";
00596             tecplot_variable_names += (*solution_names)[name];
00597 
00598 #else
00599 
00600             tecplot_variable_names += ", ";
00601             tecplot_variable_names += "r_";
00602             tecplot_variable_names += (*solution_names)[name];
00603             tecplot_variable_names += ", ";
00604             tecplot_variable_names += "i_";
00605             tecplot_variable_names += (*solution_names)[name];
00606             tecplot_variable_names += ", ";
00607             tecplot_variable_names += "a_";
00608             tecplot_variable_names += (*solution_names)[name];
00609 
00610 #endif
00611           }
00612       }
00613   }
00614 
00615   // Instantiate a TecplotMacros interface.  In 2D the most nodes per
00616   // face should be 4, in 3D it's 8.
00617 
00618 
00619   TecplotMacros tm(the_mesh.n_nodes(),
00620 #ifdef LIBMESH_USE_REAL_NUMBERS
00621                    (3 + ((solution_names == NULL) ? 0 : solution_names->size())),
00622 #else
00623                    (3 + 3*((solution_names == NULL) ? 0 : solution_names->size())),
00624 #endif
00625                    the_mesh.n_active_sub_elem(),
00626                    ((the_mesh.mesh_dimension() == 2) ? 4 : 8)
00627                    );
00628 
00629 
00630   // Copy the nodes and data to the TecplotMacros class. Note that we store
00631   // everything as a float here since the eye doesn't require a double to
00632   // understand what is going on
00633   for (unsigned int v=0; v<the_mesh.n_nodes(); v++)
00634     {
00635       tm.nd(0,v) = static_cast<float>(the_mesh.point(v)(0));
00636       tm.nd(1,v) = static_cast<float>(the_mesh.point(v)(1));
00637       tm.nd(2,v) = static_cast<float>(the_mesh.point(v)(2));
00638 
00639       if ((vec != NULL) &&
00640           (solution_names != NULL))
00641         {
00642           const unsigned int n_vars = solution_names->size();
00643 
00644           for (unsigned int c=0; c<n_vars; c++)
00645             {
00646 #ifdef LIBMESH_USE_REAL_NUMBERS
00647 
00648               tm.nd((3+c),v)     = static_cast<float>((*vec)[v*n_vars + c]);
00649 #else
00650               tm.nd((3+3*c),v)   = static_cast<float>((*vec)[v*n_vars + c].real());
00651               tm.nd((3+3*c+1),v) = static_cast<float>((*vec)[v*n_vars + c].imag());
00652               tm.nd((3+3*c+2),v) = static_cast<float>(std::abs((*vec)[v*n_vars + c]));
00653 #endif
00654             }
00655         }
00656     }
00657 
00658 
00659   // Copy the connectivity
00660   {
00661     unsigned int te = 0;
00662 
00663     MeshBase::const_element_iterator       it  = the_mesh.active_elements_begin();
00664     const MeshBase::const_element_iterator end = the_mesh.active_elements_end();
00665 
00666     for ( ; it != end; ++it)
00667       {
00668         std::vector<dof_id_type> conn;
00669         for (unsigned int se=0; se<(*it)->n_sub_elem(); se++)
00670           {
00671             (*it)->connectivity(se, TECPLOT, conn);
00672 
00673             for (unsigned int node=0; node<conn.size(); node++)
00674               tm.cd(node,te) = conn[node];
00675 
00676             te++;
00677           }
00678       }
00679   }
00680 
00681 
00682   // Ready to call the Tecplot API
00683   {
00684     int ierr = 0,
00685       num_nodes = static_cast<int>(the_mesh.n_nodes()),
00686       num_cells = static_cast<int>(the_mesh.n_active_sub_elem());
00687 
00688 
00689     ierr = TECINI (NULL,
00690                    (char*) tecplot_variable_names.c_str(),
00691                    (char*) fname.c_str(),
00692                    (char*) ".",
00693                    &tec_debug,
00694                    &is_double);
00695 
00696     libmesh_assert_equal_to (ierr, 0);
00697 
00698     ierr = TECZNE (NULL,
00699                    &num_nodes,
00700                    &num_cells,
00701                    &cell_type,
00702                    (char*) "FEBLOCK",
00703                    NULL);
00704 
00705     libmesh_assert_equal_to (ierr, 0);
00706 
00707 
00708     int total =
00709 #ifdef LIBMESH_USE_REAL_NUMBERS
00710       ((3 + ((solution_names == NULL) ? 0 : solution_names->size()))*num_nodes);
00711 #else
00712       ((3 + 3*((solution_names == NULL) ? 0 : solution_names->size()))*num_nodes);
00713 #endif
00714 
00715 
00716     ierr = TECDAT (&total,
00717                    &tm.nodalData[0],
00718                    &is_double);
00719 
00720     libmesh_assert_equal_to (ierr, 0);
00721 
00722     ierr = TECNOD (&tm.connData[0]);
00723 
00724     libmesh_assert_equal_to (ierr, 0);
00725 
00726     ierr = TECEND ();
00727 
00728     libmesh_assert_equal_to (ierr, 0);
00729   }
00730 
00731 #endif
00732 }

virtual void libMesh::MeshOutput< MeshBase >::write_equation_systems ( const std::string &  ,
const EquationSystems ,
const std::set< std::string > *  system_names = NULL 
) [virtual, inherited]

This method implements writing a mesh with data to a specified file where the data is taken from the EquationSystems object.

Referenced by libMesh::Nemesis_IO::write_timestep().

void libMesh::TecplotIO::write_nodal_data ( const std::string &  fname,
const std::vector< Number > &  soln,
const std::vector< std::string > &  names 
) [virtual]

This method implements writing a mesh with nodal data to a specified file where the nodal data and variable names are provided.

Reimplemented from libMesh::MeshOutput< MeshBase >.

Definition at line 152 of file tecplot_io.C.

References binary(), libMesh::processor_id(), write_ascii(), and write_binary().

00155 {
00156   START_LOG("write_nodal_data()", "TecplotIO");
00157 
00158   if (libMesh::processor_id() == 0)
00159     {
00160       if (this->binary())
00161         this->write_binary (fname, &soln, &names);
00162       else
00163         this->write_ascii  (fname, &soln, &names);
00164     }
00165 
00166   STOP_LOG("write_nodal_data()", "TecplotIO");
00167 }

std::string & libMesh::TecplotIO::zone_title (  )  [inline]

The zone title to write.

Definition at line 180 of file tecplot_io.h.

References _zone_title.

Referenced by write_binary().

00181 {
00182   return _zone_title;
00183 }


Member Data Documentation

Flag to write binary data.

Definition at line 128 of file tecplot_io.h.

Referenced by binary().

const bool libMesh::MeshOutput< MeshBase >::_is_parallel_format [protected, inherited]

Flag specifying whether this format is parallel-capable. If this is false (default) I/O is only permitted when the mesh has been serialized.

Definition at line 126 of file mesh_output.h.

Referenced by libMesh::PostscriptIO::write(), libMesh::FroIO::write(), libMesh::EnsightIO::write(), and libMesh::DivaIO::write().

Offset for Tecplot's STRANDID.

Definition at line 138 of file tecplot_io.h.

Referenced by strand_offset().

The subdomains in the mesh.

Definition at line 148 of file tecplot_io.h.

Referenced by TecplotIO(), and write_binary().

double libMesh::TecplotIO::_time [private]

Solution time.

Definition at line 133 of file tecplot_io.h.

Referenced by time(), and write_binary().

std::string libMesh::TecplotIO::_zone_title [private]

The zone title to write.

Definition at line 143 of file tecplot_io.h.

Referenced by zone_title().


The documentation for this class was generated from the following files:

Site Created By: libMesh Developers
Last modified: February 05 2013 19:55:40 UTC

Hosted By:
SourceForge.net Logo