Poster of Linux kernelThe best gift for a Linux geek
VTKIO

VTKIO

Section: C Library Functions (3) Updated: Thu Apr 7 2011
Local index Up
 

NAME

VTKIO -  

SYNOPSIS


#include <vtk_io.h>

Inherits MeshInput< MeshBase >, and MeshOutput< MeshBase >.  

Public Member Functions


VTKIO (MeshBase &mesh, MeshData *mesh_data=NULL)

VTKIO (const MeshBase &mesh, MeshData *mesh_data=NULL)

virtual void write_equation_systems (const std::string &fname, const EquationSystems &es)

virtual void read (const std::string &)

virtual void write (const std::string &)

vtkUnstructuredGrid * get_vtk_grid ()

virtual void write_nodal_data (const std::string &, const std::vector< Number > &, const std::vector< std::string > &)
 

Protected Member Functions


MeshBase & mesh ()

void skip_comment_lines (std::istream &in, const char comment_start)

const MeshBase & mesh () const
 

Private Member Functions


vtkPoints * nodes_to_vtk (const MeshBase &mesh)

vtkCellArray * cells_to_vtk (const MeshBase &mesh, int *&types)

void solution_to_vtk (const EquationSystems &es, vtkUnstructuredGrid *&grid)

void system_vectors_to_vtk (const EquationSystems &es, vtkUnstructuredGrid *&grid)
 

Private Attributes


vtkUnstructuredGrid * _vtk_grid

MeshData * _mesh_data
 

Detailed Description

This class implements reading and writing meshes in the VTK format. Format description: cf. VTK home page.

This class will not have any functionality unless VTK is detected during configure and hence LIBMESH_HAVE_VTK is defined.

Author:

Wout Ruijter, 2007 (Checked in to LibMesh by J.W. Peterson)

Definition at line 55 of file vtk_io.h.  

Constructor & Destructor Documentation

 

VTKIO::VTKIO (MeshBase &mesh, MeshData *mesh_data = NULL) [inline]Constructor. Takes a writeable reference to a mesh object. This is the constructor required to read a mesh.

Definition at line 144 of file vtk_io.h.

References _vtk_grid.

                                                 :
        MeshInput<MeshBase> (mesh),
        MeshOutput<MeshBase>(mesh),
        _mesh_data(mesh_data)
{
  _vtk_grid = NULL;     
  libmesh_experimental();
}
 

VTKIO::VTKIO (const MeshBase &mesh, MeshData *mesh_data = NULL) [inline]Constructor. Takes a read-only reference to a mesh object. This is the constructor required to write a mesh.

Definition at line 156 of file vtk_io.h.

References _vtk_grid.

                                                       :
        MeshOutput<MeshBase>(mesh),
        _mesh_data(mesh_data)
{
  _vtk_grid = NULL;     
  libmesh_experimental();
}
 

Member Function Documentation

 

vtkCellArray * VTKIO::cells_to_vtk (const MeshBase &mesh, int *&types) [private]write the cells from the mesh into a vtkUnstructuredGrid

Definition at line 138 of file vtk_io.C.

References Elem::connectivity(), libMeshEnums::EDGE2, libMeshEnums::EDGE3, libMeshEnums::EDGE4, MeshBase::elem(), libMeshEnums::HEX20, libMeshEnums::HEX27, libMeshEnums::HEX8, DofObject::id(), libMeshEnums::INFEDGE2, libMeshEnums::INFHEX16, libMeshEnums::INFHEX18, libMeshEnums::INFHEX8, libMeshEnums::INFPRISM12, libMeshEnums::INFPRISM6, libMeshEnums::INFQUAD4, libMeshEnums::INFQUAD6, libMeshEnums::INVALID_ELEM, MeshBase::n_active_elem(), Elem::n_nodes(), libMeshEnums::NODEELEM, libMeshEnums::PRISM15, libMeshEnums::PRISM18, libMeshEnums::PRISM6, libMesh::processor_id(), libMeshEnums::PYRAMID5, libMeshEnums::QUAD4, libMeshEnums::QUAD8, libMeshEnums::QUAD9, libMeshEnums::TET10, libMeshEnums::TET4, libMeshEnums::TRI3, libMeshEnums::TRI6, Elem::type(), and libMeshEnums::VTK.

Referenced by write(), and write_equation_systems().

                                                                  {
        //, vtkUnstructuredGrid*& grid){
        //  MeshBase::const_element_iterator       it  = mesh.active_local_elements_begin();
        //  const MeshBase::const_element_iterator end = mesh.active_local_elements_end(); 
        if (libMesh::processor_id() == 0)
        {

                vtkCellArray* cells = vtkCellArray::New();
                //     cells->Allocate(100000);
                vtkIdList *pts = vtkIdList::New();

                // for some reason the iterator variant of this code throws an error
                //     MeshBase::const_element_iterator       it  = mesh.active_elements_begin();
                //     MeshBase::const_element_iterator end = mesh.active_elements_end(); 

                //     for ( ; it != end; ++it)
                //     {
                for(unsigned int el_nr =0; el_nr< mesh.n_active_elem(); ++el_nr){
                        Elem *elem  = mesh.elem(el_nr);
                        //         if(elem->active()){ 
                        //       std::cout<<'sub elem '<<elem->has_children()<<std::endl;  
                        //       (*it);
                        vtkIdType celltype = VTK_EMPTY_CELL; // initialize to something to avoid compiler warning

                        switch(elem->type())
                        {
                                case EDGE2:                             
                                        celltype = VTK_LINE;
                                        break;
                                case EDGE3:      
                                        celltype = VTK_QUADRATIC_EDGE;
                                        break;// 1
                                case TRI3:       
                                        celltype = VTK_TRIANGLE;
                                        break;// 3
                                case TRI6:       
                                        celltype = VTK_QUADRATIC_TRIANGLE;
                                        break;// 4
                                case QUAD4:      
                                        celltype = VTK_QUAD;
                                        break;// 5
                                case QUAD8:      
                                        celltype = VTK_QUADRATIC_QUAD;
                                        break;// 6
                                case TET4:      
                                        celltype = VTK_TETRA;
                                        break;// 8
                                case TET10:      
                                        celltype = VTK_QUADRATIC_TETRA;
                                        break;// 9
                                case HEX8:    
                                        celltype = VTK_HEXAHEDRON;
                                        break;// 10
                                case HEX20:      
                                        celltype = VTK_QUADRATIC_HEXAHEDRON;
                                        break;// 12
                                case PRISM6:     
                                        celltype = VTK_WEDGE;
                                        break;// 13
                                case PRISM15:   
                                        celltype = VTK_HIGHER_ORDER_WEDGE;
                                        break;// 14
                                case PRISM18:    
                                        break;// 15
                                case PYRAMID5:
                                        celltype = VTK_PYRAMID;
                                        break;// 16
#if VTK_MAJOR_VERSION > 5 || (VTK_MAJOR_VERSION == 5 && VTK_MINOR_VERSION > 0)
                                case QUAD9:      
                                        celltype = VTK_BIQUADRATIC_QUAD;
                                        break;
#else
                                case QUAD9:
#endif 
                                case EDGE4:      
                                case HEX27:      
                                case INFEDGE2:   
                                case INFQUAD4:   
                                case INFQUAD6:   
                                case INFHEX8:    
                                case INFHEX16:   
                                case INFHEX18:   
                                case INFPRISM6:  
                                case INFPRISM12: 
                                case NODEELEM:   
                                case INVALID_ELEM:
                                default:
                                        {
                                                std::cerr<<'element type '<<elem->type()<<' not implemented'<<std::endl;
                                                libmesh_error();
                                        }
                        }
                        //        std::cout<<'elem '<<elem->n_nodes()<<' '<<elem->n_sub_elem()<<std::endl;  
                        //        std::cout<<'type '<<celltype<<' '<<VTK_BIQUADRATIC_QUAD<<std::endl;  
                        types[elem->id()] = celltype;
                        pts->SetNumberOfIds(elem->n_nodes());
                        // get the connectivity for this element
                        std::vector<unsigned int> conn;
                        elem->connectivity(0,VTK,conn);
                        //        std::cout<<'conn size '<<conn.size()<<std::endl;  
                        for(unsigned int i=0;i<conn.size();++i)
                        {
                                //           std::cout<<' sz '<<pts->GetNumberOfIds()<<std::endl;  
                                pts->InsertId(i,conn[i]);
                        } 
                        //        std::cout<<'set cell'<<std::endl;  
                        cells->InsertNextCell(pts);
                        //         }
                        //         std::cout<<'finish set cell'<<std::endl;  
                } // end loop over active elements
                pts->Delete();

                return cells;
        }
        return NULL;
}
 

vtkUnstructuredGrid* VTKIO::get_vtk_grid () [inline]Get a pointer to the VTK datastructure

Definition at line 101 of file vtk_io.h.

References _vtk_grid.

{return _vtk_grid;}
 

MeshBase & MeshInput< MeshBase >::mesh () [protected, inherited]Returns the object as a writeable reference.

Referenced by GMVIO::_read_materials(), GMVIO::_read_nodes(), GMVIO::_read_one_cell(), GMVIO::add_cell_centered_data(), GMVIO::copy_nodal_solution(), UNVIO::element_in(), TetGenIO::element_in(), UNVIO::element_out(), UNVIO::node_in(), TetGenIO::node_in(), UNVIO::node_out(), XdrIO::read(), read(), GMVIO::read(), ExodusII_IO::read(), LegacyXdrIO::read_ascii(), LegacyXdrIO::read_binary(), UCDIO::read_implementation(), LegacyXdrIO::read_mesh(), GmshIO::read_mesh(), XdrIO::read_serialized_bcs(), XdrIO::read_serialized_connectivity(), XdrIO::read_serialized_nodes(), OFFIO::read_stream(), MatlabIO::read_stream(), solution_to_vtk(), XdrIO::write(), write(), TetGenIO::write(), GMVIO::write_ascii_new_impl(), GMVIO::write_ascii_old_impl(), GMVIO::write_binary(), GMVIO::write_discontinuous_gmv(), UNVIO::write_implementation(), UCDIO::write_implementation(), LegacyXdrIO::write_mesh(), GmshIO::write_mesh(), GmshIO::write_post(), XdrIO::write_serialized_bcs(), XdrIO::write_serialized_connectivity(), XdrIO::write_serialized_nodes(), and LegacyXdrIO::write_soln().  

const MeshBase & MeshOutput< MeshBase >::mesh () const [protected, inherited]Returns the object as a read-only reference.

Referenced by PostscriptIO::write(), FroIO::write(), MEDITIO::write_ascii(), EnsightIO::write_geometry_ascii(), EnsightIO::write_scalar_ascii(), GnuPlotIO::write_solution(), DivaIO::write_stream(), and EnsightIO::write_vector_ascii().  

vtkPoints * VTKIO::nodes_to_vtk (const MeshBase &mesh) [private]write the nodes from the mesh into a vtkUnstructuredGrid

Definition at line 101 of file vtk_io.C.

References DofObject::id(), MeshBase::n_nodes(), MeshBase::nodes_begin(), MeshBase::nodes_end(), and libMesh::processor_id().

Referenced by write(), and write_equation_systems().

                                                  {
        if (libMesh::processor_id() == 0)
                {

                        vtkPoints* points = vtkPoints::New();
                        vtkDoubleArray* pcoords = vtkDoubleArray::New();
                        pcoords->SetNumberOfComponents(3);
                        pcoords->SetNumberOfTuples(mesh.n_nodes());

                        // FIXME change to local iterators when I've figured out how to write in parallel
                        //  MeshBase::const_node_iterator nd = mesh.local_nodes_begin();
                        //  MeshBase::const_node_iterator nd_end = mesh.local_nodes_end();
                        //  points->SetNumberOfPoints(mesh.n_local_nodes()); // it seems that it needs this to prevent a segfault
                        MeshBase::const_node_iterator nd = mesh.nodes_begin();
                        MeshBase::const_node_iterator nd_end = mesh.nodes_end();
                        // write out nodal data
                        for (;nd!=nd_end;++nd)
                        {
                                Node* node = (*nd);
                                float* pnt = new float[3];
                                for(unsigned int i=0;i<3;++i){
                                        pnt[i] = (*node)(i);
                                } 
                                //     pcoords->InsertNextTuple(pnt);
                                pcoords->SetTuple(node->id(),pnt);
                                //     std::cout<<'pnt '<<pnt[0]<<' '<<pnt[1]<<' '<<pcoords-><<std::endl;  
                                //     points->InsertPoint(node->id(),pnt);
                                //     std::cout<<'point '<<node->id()<<' '<<(*node)(0)<<' '<<(*node)(1)<<std::endl;  
                                delete [] pnt;
                        }
                        points->SetData(pcoords);
                        pcoords->Delete();
                        return points;
        }
        return NULL;
}
 

void VTKIO::read (const std::string &name) [virtual]This method implements reading a mesh from a specified file in VTK format.

Implements MeshInput< MeshBase >.

Definition at line 449 of file vtk_io.C.

References _mesh_data, _vtk_grid, MeshBase::add_elem(), MeshData::add_foreign_node_id(), MeshBase::add_point(), MeshBase::clear(), Elem::connectivity(), MeshInput< MeshBase >::mesh(), Elem::n_nodes(), MeshBase::node_ptr(), libMesh::processor_id(), DofObject::set_id(), Elem::set_node(), and libMeshEnums::VTK.

Referenced by UnstructuredMesh::read().

{
  // This is a serial-only process for now;
  // the Mesh should be read on processor 0 and
  // broadcast later
  libmesh_assert(libMesh::processor_id() == 0);

#ifndef LIBMESH_HAVE_VTK
  std::cerr << 'Cannot read VTK file: ' << name
            << 'ou must have VTK installed and correctly configured to read VTK meshes.'
            << std::endl;
  libmesh_error();

#else
//  std::cout<<'read '<<name <<std::endl;  
  vtkXMLUnstructuredGridReader *reader = vtkXMLUnstructuredGridReader::New();
  reader->SetFileName( name.c_str() );
  //std::cout<<'force read'<<std::endl;  
  // Force reading
  reader->Update();

  // read in the grid   
//  vtkUnstructuredGrid *grid = reader->GetOutput();
  _vtk_grid = reader->GetOutput();
  _vtk_grid->Update();
  reader->Delete();

  // Get a reference to the mesh
  MeshBase& mesh = MeshInput<MeshBase>::mesh();

  // Clear out any pre-existing data from the Mesh
  mesh.clear();
        
  // always numbered nicely??, so we can loop like this
  // I'm pretty sure it is numbered nicely
  for (unsigned int i=0; i < static_cast<unsigned int>(_vtk_grid->GetNumberOfPoints()); ++i)
    {
      // add to the id map
      // and add the actual point
      double * pnt = _vtk_grid->GetPoint(static_cast<vtkIdType>(i));
      Point xyz(pnt[0],pnt[1],pnt[2]);
      Node* newnode = mesh.add_point(xyz,i);

      // Add node to the nodes vector &
      // tell the MeshData object the foreign node id.
      if (this->_mesh_data != NULL)
        this->_mesh_data->add_foreign_node_id (newnode, i);
    }
        
  for (unsigned int i=0; i < static_cast<unsigned int>(_vtk_grid->GetNumberOfCells()); ++i)
    {
      vtkCell* cell = _vtk_grid->GetCell(i);
      Elem* elem = NULL;  // Initialize to avoid compiler warning
      switch(cell->GetCellType())
        {
        case VTK_TETRA:
          elem = new Tet4();
          break;
        case VTK_WEDGE:
          elem = new Prism6();
          break;
        case VTK_HEXAHEDRON:
          elem = new Hex8();
          break;
        case VTK_PYRAMID:       
          elem = new Pyramid5();                                        
          break;
        case VTK_QUADRATIC_HEXAHEDRON:
          elem = new Hex20();
          break;
        case VTK_QUADRATIC_TETRA:
          elem = new Tet10();
          break;
        default:
          std::cerr << 'element type not implemented in vtkinterface ' << cell->GetCellType() << std::endl;
          libmesh_error();
        }

  // get the straightforward numbering from the VTK cells
  for(unsigned int j=0;j<elem->n_nodes();++j){
          elem->set_node(j) = mesh.node_ptr(cell->GetPointId(j));
  } 
  // then get the connectivity 
  std::vector<unsigned int> conn;
  elem->connectivity(0,VTK,conn);
  // then reshuffle the nodes according to the connectivity, this
  // two-time-assign would evade the definition of the vtk_mapping
  for(unsigned int j=0;j<conn.size();++j){
          elem->set_node(j) = mesh.node_ptr(conn[j]);
  } 
  elem->set_id(i);
  mesh.add_elem(elem);
  } // end loop over VTK cells
  _vtk_grid->Delete();
#endif // LIBMESH_HAVE_VTK
}
 

void MeshInput< MeshBase >::skip_comment_lines (std::istream &in, const charcomment_start) [protected, inherited]Reads input from in, skipping all the lines that start with the character comment_start.

Referenced by TetGenIO::read(), and UCDIO::read_implementation().  

void VTKIO::solution_to_vtk (const EquationSystems &es, vtkUnstructuredGrid *&grid) [private]write the solution to a vtkUnstructuredGrid

Definition at line 300 of file vtk_io.C.

References System::current_solution(), DofObject::dof_number(), EquationSystems::get_mesh(), EquationSystems::get_system(), MeshInput< MeshBase >::mesh(), MeshBase::n_nodes(), MeshTools::n_nodes(), EquationSystems::n_systems(), System::n_vars(), Quality::name(), MeshBase::node(), libMesh::processor_id(), System::solution, and System::variable_name().

Referenced by write_equation_systems().

                                                                                {
        // write only single processor
        if(libMesh::processor_id()==0){

                const MeshBase& mesh = (MeshBase&)es.get_mesh();

                const unsigned int n_nodes = es.get_mesh().n_nodes();
                // loop over the systems
                for(unsigned int i=0;i<es.n_systems();++i){ // for all systems, regardless of whether they are active or not
                        
                        const System& sys = es.get_system(i);

                        const unsigned int n_vars = sys.n_vars();
                        // loop over variables
                        for(unsigned int j=0;j<n_vars;++j){

                                std::string name = sys.variable_name(j);

                                vtkFloatArray *data = vtkFloatArray::New(); 

                                data->SetName(name.c_str());

                                data->SetNumberOfValues(sys.solution->size());

                                for(unsigned int k=0;k<n_nodes;++k){

                                        const unsigned int dof_nr = mesh.node(k).dof_number(i,j,0);

#ifdef LIBMESH_USE_COMPLEX_NUMBERS
                                        data->SetValue(k,sys.current_solution(dof_nr).real());
#else
                                        data->SetValue(k,sys.current_solution(dof_nr));
#endif

                                } 
                                grid->GetPointData()->AddArray(data);
                                data->Delete();
                        } 
                        
                }
        }
}
 

void VTKIO::system_vectors_to_vtk (const EquationSystems &es, vtkUnstructuredGrid *&grid) [private]write the system vectors to vtk

Definition at line 347 of file vtk_io.C.

References EquationSystems::get_mesh(), EquationSystems::get_system(), MeshBase::n_nodes(), EquationSystems::n_systems(), libMesh::processor_id(), System::vectors_begin(), and System::vectors_end().

                                                                                     {
        if (libMesh::processor_id() == 0){

                std::map<std::string,std::vector<Number> > vecs; 
                for(unsigned int i=0;i<es.n_systems();++i){
                        const System& sys = es.get_system(i);
                        System::const_vectors_iterator v_end = sys.vectors_end();
                        System::const_vectors_iterator it = sys.vectors_begin();
                        for(;it!= v_end;++it){ // for all vectors on this system
                                std::vector<Number> values;     
                                std::cout<<'it '<<it->first<<std::endl;  

                                it->second->localize_to_one(values,0);
                                std::cout<<'finish localize'<<std::endl;  
                                vecs[it->first] = values;
                        }
                }


                std::map<std::string,std::vector<Number> >::iterator it = vecs.begin(); 

                for(;it!=vecs.end();++it){

                        vtkFloatArray *data = vtkFloatArray::New(); 

                        data->SetName(it->first.c_str());

                        libmesh_assert(it->second.size()==es.get_mesh().n_nodes());

                        data->SetNumberOfValues(it->second.size());

                        for(unsigned int i=0;i<it->second.size();++i){

#ifdef LIBMESH_USE_COMPLEX_NUMBERS
                                data->SetValue(i,it->second[i].real());
#else
                                data->SetValue(i,it->second[i]);
#endif

                        }

                        grid->GetPointData()->AddArray(data);
                        data->Delete();
                } 

        }

/*
        // write out the vectors added to the systems
        if (libMesh::processor_id() == 0)
                {
                        std::cout<<'write system vectors'<<std::endl;  
                const MeshBase& mesh = (MeshBase&)es.get_mesh();
                const unsigned int n_nodes = mesh.n_nodes();
                for(unsigned int i=0;i<es.n_systems();++i){ // for all systems, regardless of whether they are active or not
                        const System& sys = es.get_system(i);
                        std::cout<<'i '<<i<<std::endl;  
                        System::const_vectors_iterator v_end = sys.vectors_end();
                        System::const_vectors_iterator it = sys.vectors_begin();
                        for(;it!= v_end;++it){ // for all vectors on this system
                                std::cout<<'it- '<<it->second->size()<<' '<<n_nodes<<' '<<it->first<<std::endl;  
                                vtkFloatArray *data = vtkFloatArray::New(); 
                                data->SetName(it->first.c_str());
                                std::vector<Number> values;     
                                it->second->localize(values);
                                data->SetNumberOfValues(n_nodes);


        //         MeshBase::const_node_iterator it = mesh.active_nodes_begin();
        //         const MeshBase::const_node_iterator n_end = mesh.nodes_end();
        //         for(unsigned int count=0;it!=n_end;++it,++count){                    
                                for(unsigned int j=0;j<n_nodes;++j){
                                        std::cout<<'j '<<j<<' '<<(*it->second).size()<<std::endl;  
//               const unsigned int dof_nr = mesh.node(j).dof_number(i,0,0);
                                        data->SetValue(j,values[j]);
        //            it++;
                                } 
                                grid->GetPointData()->AddArray(data);           
                        } 
                } 
        }*/
}
 

void VTKIO::write (const std::string &name) [virtual]This method implements writing a mesh to a specified '.poly' file. '.poly' files defines so called Piecewise Linear Complex (PLC).

This method implements writing to a .vtu (VTK Unstructured Grid) file. This is one of the new style XML dataformats.

Implements MeshOutput< MeshBase >.

Definition at line 619 of file vtk_io.C.

References _vtk_grid, cells_to_vtk(), MeshInput< MeshBase >::mesh(), MeshBase::n_active_elem(), nodes_to_vtk(), and libMesh::processor_id().

Referenced by UnstructuredMesh::write().

{  
#ifndef LIBMESH_HAVE_VTK
  std::cerr << 'Cannot write VTK file: ' << name
            << 'ou must have VTK installed and correctly configured to write VTK meshes.'
            << std::endl;
  libmesh_error();

#else
  if (libMesh::processor_id() == 0)
  {

          MeshBase& mesh = MeshInput<MeshBase>::mesh();
          _vtk_grid = vtkUnstructuredGrid::New();
          vtkXMLUnstructuredGridWriter* writer = vtkXMLUnstructuredGridWriter::New();
          std::cout<<'write nodes '<<std::endl;  
          vtkPoints* pnts = nodes_to_vtk(mesh);
          _vtk_grid->SetPoints(pnts);

          std::cout<<'write elements '<<std::endl;  
          int * types = new int[mesh.n_active_elem()];
          vtkCellArray* cells = cells_to_vtk(mesh,types);
          _vtk_grid->SetCells(types,cells);
          //  , _vtk_grid);
          writer->SetInput(_vtk_grid);
          writer->SetDataModeToAscii();
          writer->SetFileName(name.c_str());
          writer->Write();
          writer->Delete();
          _vtk_grid->Delete();
  }
#endif // LIBMESH_HAVE_VTK
}
 

void VTKIO::write_equation_systems (const std::string &fname, const EquationSystems &es) [virtual]This method implements writing a mesh with nodal data to a specified file where the nodal data and variable names are provided. Overloads writing equation systems, this is done because when overloading write_nodal_data there would be no way to export cell centered data

This method writes out the equationsystems to a .pvtu file (VTK parallel unstructured grid).

Reimplemented from MeshOutput< MeshBase >.

Definition at line 554 of file vtk_io.C.

References _vtk_grid, cells_to_vtk(), EquationSystems::get_mesh(), MeshBase::n_active_elem(), EquationSystems::n_systems(), nodes_to_vtk(), libMesh::processor_id(), and solution_to_vtk().

{
#ifndef LIBMESH_HAVE_VTK

  // Do something with the es to avoid a compiler warning.
  es.n_systems();
  
  std::cerr << 'Cannot write VTK file: ' << fname
            << 'ou must have VTK installed and correctly configured to read VTK meshes.'
            << std::endl;
  libmesh_error();
  
#else
  if (libMesh::processor_id() == 0)
                {

          // check if the filename extension is pvtu
          libmesh_assert(fname.substr(fname.rfind('.'),fname.size())=='.pvtu');
          /*
                * we only use Unstructured grids
                */
          _vtk_grid = vtkUnstructuredGrid::New();
          vtkXMLPUnstructuredGridWriter* writer= vtkXMLPUnstructuredGridWriter::New();
          std::cout<<'get points '<<std::endl;  
          vtkPoints* pnts = nodes_to_vtk((const MeshBase&)es.get_mesh());
          _vtk_grid->SetPoints(pnts);
         
          int * types = new int[es.get_mesh().n_active_elem()];
          std::cout<<'get cells'<<std::endl;  
          vtkCellArray* cells = cells_to_vtk((const MeshBase&)es.get_mesh(), types);

          std::cout<<'set cells'<<std::endl;  
          _vtk_grid->SetCells(types,cells);
          
          // I'd like to write out meshdata, but this requires some coding, in
          // particular, non_const meshdata iterators
          //   const MeshData& md = es.get_mesh_data();
          //   if(es.has_mesh_data())
          //      meshdata_to_vtk(md,_vtk_grid);
          //   libmesh_assert (soln.size() ==mesh.n_nodes()*names.size());
          std::cout<<'write solution'<<std::endl;  
          solution_to_vtk(es,_vtk_grid);

//#ifdef DEBUG
//     if(true) // add some condition here, although maybe it is more sensible to give each vector a flag on whether it is to be written out or not
//       system_vectors_to_vtk(es,_vtk_grid);
//#endif
          writer->SetInput(_vtk_grid);
          writer->SetFileName(fname.c_str());
          writer->SetDataModeToAscii();
          writer->Write();

          delete [] types;
          pnts->Delete();
          cells->Delete();
          _vtk_grid->Delete();
          writer->Delete();
        }
#endif
}
 

virtual void MeshOutput< MeshBase >::write_nodal_data (const std::string &, const std::vector< Number > &, const std::vector< std::string > &) [inline, virtual, inherited]This method implements writing a mesh with nodal data to a specified file where the nodal data and variable names are provided.

Reimplemented in ExodusII_IO, GmshIO, GMVIO, GnuPlotIO, MEDITIO, and TecplotIO.

Definition at line 91 of file mesh_output.h.

  { libmesh_error(); }
 

Member Data Documentation

 

MeshData* VTKIO::_mesh_data [private]A pointer to the MeshData object you would like to use. with this VTKIO object. Can be NULL.

Definition at line 136 of file vtk_io.h.

Referenced by read().  

vtkUnstructuredGrid* VTKIO::_vtk_grid [private]pointer to the VTK grid

Definition at line 130 of file vtk_io.h.

Referenced by get_vtk_grid(), read(), VTKIO(), write(), and write_equation_systems().

 

Author

Generated automatically by Doxygen for libMesh from the source code.


 

Index

NAME
SYNOPSIS
Public Member Functions
Protected Member Functions
Private Member Functions
Private Attributes
Detailed Description
Constructor & Destructor Documentation
VTKIO::VTKIO (MeshBase &mesh, MeshData *mesh_data = NULL) [inline]Constructor. Takes a writeable reference to a mesh object. This is the constructor required to read a mesh.
VTKIO::VTKIO (const MeshBase &mesh, MeshData *mesh_data = NULL) [inline]Constructor. Takes a read-only reference to a mesh object. This is the constructor required to write a mesh.
Member Function Documentation
vtkCellArray * VTKIO::cells_to_vtk (const MeshBase &mesh, int *&types) [private]write the cells from the mesh into a vtkUnstructuredGrid
vtkUnstructuredGrid* VTKIO::get_vtk_grid () [inline]Get a pointer to the VTK datastructure
MeshBase & MeshInput< MeshBase >::mesh () [protected, inherited]Returns the object as a writeable reference.
const MeshBase & MeshOutput< MeshBase >::mesh () const [protected, inherited]Returns the object as a read-only reference.
vtkPoints * VTKIO::nodes_to_vtk (const MeshBase &mesh) [private]write the nodes from the mesh into a vtkUnstructuredGrid
void VTKIO::read (const std::string &name) [virtual]This method implements reading a mesh from a specified file in VTK format.
void MeshInput< MeshBase >::skip_comment_lines (std::istream &in, const charcomment_start) [protected, inherited]Reads input from in, skipping all the lines that start with the character comment_start.
void VTKIO::solution_to_vtk (const EquationSystems &es, vtkUnstructuredGrid *&grid) [private]write the solution to a vtkUnstructuredGrid
void VTKIO::system_vectors_to_vtk (const EquationSystems &es, vtkUnstructuredGrid *&grid) [private]write the system vectors to vtk
void VTKIO::write (const std::string &name) [virtual]This method implements writing a mesh to a specified '.poly' file. '.poly' files defines so called Piecewise Linear Complex (PLC).
void VTKIO::write_equation_systems (const std::string &fname, const EquationSystems &es) [virtual]This method implements writing a mesh with nodal data to a specified file where the nodal data and variable names are provided. Overloads writing equation systems, this is done because when overloading write_nodal_data there would be no way to export cell centered data
virtual void MeshOutput< MeshBase >::write_nodal_data (const std::string &, const std::vector< Number > &, const std::vector< std::string > &) [inline, virtual, inherited]This method implements writing a mesh with nodal data to a specified file where the nodal data and variable names are provided.
Member Data Documentation
MeshData* VTKIO::_mesh_data [private]A pointer to the MeshData object you would like to use. with this VTKIO object. Can be NULL.
vtkUnstructuredGrid* VTKIO::_vtk_grid [private]pointer to the VTK grid
Author

This document was created by man2html, using the manual pages.
Time: 21:58:09 GMT, April 16, 2011