LehrFEM++ 1.0.0
A simple Finite Element Library for teaching
Loading...
Searching...
No Matches
Public Types | Public Member Functions | Private Member Functions | Static Private Member Functions | Private Attributes | List of all members
lf::io::VtkWriter Class Reference

Write a mesh along with mesh data into a vtk file. More...

#include <lf/io/vtk_writer.h>

Public Types

using dim_t = base::dim_t
 

Public Member Functions

 VtkWriter (const VtkWriter &)=delete
 
 VtkWriter (VtkWriter &&)=delete
 
VtkWriteroperator= (const VtkWriter &)=delete
 
VtkWriteroperator= (VtkWriter &&)=delete
 
 VtkWriter (std::shared_ptr< const mesh::Mesh > mesh, std::string filename, dim_t codim=0, unsigned char order=1)
 Construct a new VtkWriter.
 
void setBinary (bool binary)
 Determines whether the Vtk file is written in binary or ASCII mode (default).
 
void WritePointData (const std::string &name, const mesh::utils::MeshDataSet< unsigned char > &mds, unsigned char undefined_value=0)
 Add a new unsigned char attribute dataset that attaches data to the points/nodes of the mesh.
 
void WritePointData (const std::string &name, const mesh::utils::MeshDataSet< char > &mds, char undefined_value=0)
 Add a new char attribute dataset that attaches data to the points/nodes of the mesh.
 
void WritePointData (const std::string &name, const mesh::utils::MeshDataSet< unsigned int > &mds, unsigned int undefined_value=0)
 Add a new unsigned int attribute dataset that attaches data to the points/nodes of the mesh.
 
void WritePointData (const std::string &name, const mesh::utils::MeshDataSet< int > &mds, int undefined_value=0)
 Add a new unsigned int attribute dataset that attaches data to the points/nodes of the mesh.
 
void WritePointData (const std::string &name, const mesh::utils::MeshDataSet< float > &mds, float undefined_value=0.F)
 Add a new float attribute dataset that attaches data to the points/nodes of the mesh.
 
void WritePointData (const std::string &name, const mesh::utils::MeshDataSet< double > &mds, double undefined_value=0.)
 Add a new double attribute dataset that attaches data to the points/nodes of the mesh.
 
void WritePointData (const std::string &name, const mesh::utils::MeshDataSet< Eigen::Vector2d > &mds, const Eigen::Vector2d &undefined_value={0, 0})
 Add a new vector attribute dataset that attaches vectors to the points/nodes of the mesh.
 
void WritePointData (const std::string &name, const mesh::utils::MeshDataSet< Eigen::Vector2f > &mds, const Eigen::Vector2f &undefined_value={0, 0})
 Add a new vector attribute dataset that attaches vectors to the points/nodes of the mesh.
 
void WritePointData (const std::string &name, const mesh::utils::MeshDataSet< Eigen::Vector3d > &mds, const Eigen::Vector3d &undefined_value={0, 0, 0})
 Add a new vector attribute dataset that attaches vectors to the points/nodes of the mesh.
 
void WritePointData (const std::string &name, const mesh::utils::MeshDataSet< Eigen::Vector3f > &mds, const Eigen::Vector3f &undefined_value={0, 0, 0})
 Add a new vector attribute dataset that attaches vectors to the points/nodes of the mesh.
 
void WritePointData (const std::string &name, const mesh::utils::MeshDataSet< Eigen::VectorXd > &mds, const Eigen::VectorXd &undefined_value=Eigen::Vector3d(0, 0, 0))
 Add a new vector attribute dataset that attaches vectors to the points/nodes of the mesh.
 
void WritePointData (const std::string &name, const mesh::utils::MeshDataSet< Eigen::VectorXf > &mds, const Eigen::VectorXf &undefined_value=Eigen::Vector3f(0, 0, 0))
 Add a new vector attribute dataset that attaches vectors to the points/nodes of the mesh.
 
template<mesh::utils::MeshFunction MESH_FUNCTION>
void WritePointData (const std::string &name, const MESH_FUNCTION &mesh_function)
 Sample a MeshFunction at points of the mesh and possible edges/interior (for order>1).
 
void WriteCellData (const std::string &name, const mesh::utils::MeshDataSet< unsigned char > &mds, unsigned char undefined_value=0)
 Add a new unsigned char attribute dataset that attaches data to the cells of the mesh (i.e. to entities with codim = "codim that was specified in constructor of VtkWriter")
 
void WriteCellData (const std::string &name, const mesh::utils::MeshDataSet< char > &mds, char undefined_value=0)
 Add a new char attribute dataset that attaches data to the cells of the mesh (i.e. to entities with codim = "codim that was specified in constructor of VtkWriter")
 
void WriteCellData (const std::string &name, const mesh::utils::MeshDataSet< unsigned int > &mds, unsigned int undefined_value=0)
 Add a new unsigned int attribute dataset that attaches data to the cells of the mesh (i.e. to entities with codim = "codim that was specified in constructor of VtkWriter")
 
void WriteCellData (const std::string &name, const mesh::utils::MeshDataSet< int > &mds, int undefined_value=0)
 Add a new int attribute dataset that attaches data to the cells of the mesh (i.e. to entities with codim = "codim that was specified in constructor of VtkWriter")
 
void WriteCellData (const std::string &name, const mesh::utils::MeshDataSet< float > &mds, float undefined_value=0)
 Add a new float attribute dataset that attaches data to the cells of the mesh (i.e. to entities with codim = "codim that was specified in constructor of VtkWriter")
 
void WriteCellData (const std::string &name, const mesh::utils::MeshDataSet< double > &mds, double undefined_value=0)
 Add a new double attribute dataset that attaches data to the cells of the mesh (i.e. to entities with codim = "codim that was specified in constructor of VtkWriter")
 
void WriteCellData (const std::string &name, const mesh::utils::MeshDataSet< Eigen::Vector2d > &mds, const Eigen::Vector2d &undefined_value={0, 0})
 Add a new vector attribute dataset that attaches vectors to the cells of the mesh.
 
void WriteCellData (const std::string &name, const mesh::utils::MeshDataSet< Eigen::Vector2f > &mds, const Eigen::Vector2f &undefined_value={0, 0})
 Add a new vector attribute dataset that attaches vectors to the cells of the mesh.
 
void WriteCellData (const std::string &name, const mesh::utils::MeshDataSet< Eigen::Vector3d > &mds, const Eigen::Vector3d &undefined_value={0, 0, 0})
 Add a new vector attribute dataset that attaches vectors to the cells of the mesh.
 
void WriteCellData (const std::string &name, const mesh::utils::MeshDataSet< Eigen::Vector3f > &mds, const Eigen::Vector3f &undefined_value={0, 0, 0})
 Add a new vector attribute dataset that attaches vectors to the cells of the mesh.
 
void WriteCellData (const std::string &name, const mesh::utils::MeshDataSet< Eigen::VectorXd > &mds, const Eigen::VectorXd &undefined_value=Eigen::Vector3d(0, 0, 0))
 Add a new vector attribute dataset that attaches vectors to the cells of the mesh.
 
void WriteCellData (const std::string &name, const mesh::utils::MeshDataSet< Eigen::VectorXf > &mds, const Eigen::VectorXf &undefined_value=Eigen::Vector3f(0, 0, 0))
 Add a new vector attribute dataset that attaches vectors to the cells of the mesh.
 
template<mesh::utils::MeshFunction MESH_FUNCTION>
void WriteCellData (const std::string &name, const MESH_FUNCTION &mesh_function)
 Sample a MeshFunction at the barycenter of the cell and visualize it as cell data in the Vtk File.
 
void WriteGlobalData (const std::string &name, std::vector< int > data)
 Write global data into the vtk file that is not related to the mesh at all.
 
void WriteGlobalData (const std::string &name, std::vector< float > data)
 Write global data into the vtk file that is not related to the mesh at all.
 
void WriteGlobalData (const std::string &name, std::vector< double > data)
 Write global data into the vtk file that is not related to the mesh at all.
 
 ~VtkWriter ()
 Destructor, writes everything into the file and closes it.
 

Private Member Functions

template<class T >
void WriteScalarPointData (const std::string &name, const mesh::utils::MeshDataSet< T > &mds, T undefined_value)
 
template<int ROWS, class T >
void WriteVectorPointData (const std::string &name, const mesh::utils::MeshDataSet< Eigen::Matrix< T, ROWS, 1 > > &mds, const Eigen::Matrix< T, ROWS, 1 > &undefined_value)
 
template<class T >
void WriteScalarCellData (const std::string &name, const mesh::utils::MeshDataSet< T > &mds, T undefined_value)
 
template<int ROWS, class T >
void WriteVectorCellData (const std::string &name, const mesh::utils::MeshDataSet< Eigen::Matrix< T, ROWS, 1 > > &mds, const Eigen::Matrix< T, ROWS, 1 > &undefined_value)
 
template<class T >
void WriteFieldData (const std::string &name, std::vector< T > data)
 
template<class DATA >
void CheckAttributeSetName (const DATA &data, const std::string &name)
 

Static Private Member Functions

template<int ROWS, class T >
static void PadWithZeros (Eigen::Matrix< T, 3, 1 > &out, const Eigen::Matrix< T, ROWS, 1 > &in)
 

Private Attributes

std::shared_ptr< const mesh::Meshmesh_
 
VtkFile vtk_file_
 
std::string filename_
 
dim_t codim_
 
unsigned char order_
 
std::array< Eigen::MatrixXd, 5 > aux_nodes_
 
std::array< mesh::utils::CodimMeshDataSet< unsigned int >, 2 > aux_node_offset_
 

Detailed Description

Write a mesh along with mesh data into a vtk file.

Sample usage:

std::shared_ptr<mesh::Mesh> mesh; // initialize mesh somehow
// data stored with codim=0 entities
std::shared_ptr<mesh::utils::MeshDataSet<double>> cell_data;
// vector data stored with points of the mesh
std::shared_ptr<mesh::utils::MeshDataSet<Eigen::Vector2d>> point_data;
// A Mesh function representing the function sin(x)*cos(y)
auto mf = mesh::utils::MeshFunctionGlobal(
[](const Eigen::Vector2d& x) { return std::sin(x[0]) * std::cos(x[1]); });
io::VtkWriter vtk_writer(mesh, "filename.vtk");
vtk_writer.WriteCellData("cellData", *cell_data);
vtk_writer.WritePointData("pointData", *point_data);
vtk_writer.WritePointData("mfPoint", mf);

Mesh

VtkWriter can work with any mesh implementing the mesh::Mesh interface. By default, it will translate codim=0 mesh entities into VTK cells. If necessary, VtkWriter can also write out all codim=1 mesh entities (skeleton), see documentation of the Constructor.

Data

VtkWriter can attach different types of data to either cells (codim=0) or points (codim=dimMesh) of the mesh (see the different overloads of WritePointData() and WriteCellData()). Each dataset is translated into a VTK AttributeSet which can then be visualized in ParaView.

In addition to point- and cell-data VtkWriter can also write out global data (in VTK terms this is called "field attribute sets"). This is useful e.g. to write along the current time of the simulation and visualize it with a label in ParaView.

Note
The actual Vtk file is written when the VtkWriter is destructed! It is very important that this destructor is called.
By default, VtkWriter outputs files in text-format. However, for larger meshes you can switch to the more efficient binary mode (setBinary()).

How to view the output in Paraview

Once a *.vtk file has been written with the VtkWriter, this file can be opened in Paraview as follows:

  1. Open the file via File->Open (here we open the file smiley.vtk which is produced by the example code in examples/io/vtk_gmsh_demo)
  2. Note that Paraview added the Vtk-File to the "Pipeline Browser". You can activate the file and show it in the view by selecting the Vtk file in the Pipeline Browser and clicking on "Apply" below:
  1. Paraview will now automatically select a set of data that was written with WritePointData() or WriteCellData() and visualize it on the mesh. You can select a different dataset or change the representation mode (e.g. to show the underlying mesh):
  1. There are many ways to visualize vector datasets (i.e. if the MeshDataSet is Eigen::Vector2d or Eigen::Vector3d valued) in Paraview. If you select a vector dataset, you can additionally select how the mesh should be colored: Magnitude of the vector, or the x-, y- or z-component of the vector.
  1. If you want to visualize the vectors using arrows, you can use a "Glyph" filter:
    1. Click on "Filters->Common->Glyph"
    2. Note that Paraview added the filter to the "Pipeline Browser"
    3. Select the new filter and adjust the properties in the "Properties Window" below:
      Subsection Property Name Description
      Glyph Source Glyph Type Arrow is the most common choice, 2D Glyph is also interesting
      Active Attributes Scalars The (scalar) dataset that is used to color the arrows.
      Active Attributes Vector The vector dataset that defines the vector direction. If Scalars is a cell-based dataset, this must also be cell-based. If it Scalars is a point dataset, this must be a vector point dataset
      Scaling Scale Mode Set this to vector if the length of the arrows should be proportional to the magnitude of the vectors
      Scaling Scale Factor If vectors are too small/too big adjust this factor
    4. Click on "Apply"
  2. Finally, if you want to take a look at the raw data you can e.g. click on "Split Vertical" and then select "Spreadsheet View". You can then select rows and see to what cell/point they correspond in the mesh:

Higher order cells

Recent versions of Paraview support higher order Lagrange cells of order up to 10. VtkWriter can write such higher order cells so that e.g. higher-order geometry mappings or higher order BVP solutions can be visualized appropriately.

Example usage:

The following code will output a mesh of the unit circle ( \( \{\vec{x} \in \mathbb{R}^2 | \norm{x}<1 \}\)) consisting of second order quadrilaterals and the function \( \sin(\pi x) \cos(\pi x)\) to vtk:

// read a 2nd order gmsh mesh:
GmshReader reader(std::make_unique<mesh::hybrid2d::MeshFactory>(2),
"unit_circle.msh");
std::shared_ptr<mesh::Mesh> mesh = reader.mesh();
// define mesh function of the form sin(\pi*x)*cos(\pi*y) (in global
// coordinates)
auto mfTrig = mesh::utils::MeshFunctionGlobal([](const Eigen::Vector2d& x) {
return std::sin(base::kPi * x[0]) * std::cos(base::kPi * x[1]);
});
// output the mesh + mesh function using 1st order cells:
VtkWriter vtk1(mesh, "1storder.vtk");
vtk1.WritePointData("trig", mfTrig);
// output the mesh + mesh function using 2nd order cells:
VtkWriter vtk2(mesh, "2ndorder.vtk", 0, 2);
vtk2.WritePointData("trig", mfTrig);

The two vtk files generated by this code can be visualized in Paraview:

We can see that the first order mesh approximates the boundary of the circle with line segments and that the mesh function is approximated by bilinear basis functions. In the present example this is not very satisfactory, but luckily second order cells remedy this situation.

Because it is rather hard to visualize second or higher order basis functions directly, Paraview interpolates such basis functions linearly. By default, Paraview approximates second order quadrilaterals by dividing them into four first order quadrilaterals. This is what is shown in the middle column above where we can see that the boundary of the circle is still not very smooth. In order to get a better representation of the second order basis functions, we can instruct Paraview to subdivide the cells even more with the tesselate filter (see below for instructions). This gives the desired result (on the right).

Note
you need a recent version of Paraview to visualize higher-order cells. Version 5.6.0 is known to work.
By using higher order cells, not only the geometry is better approximated, the Point datasets are also much better approximated. This is unfortunately not true for cell datasets.
If you construct a VtkWriter with order greater than one, you cannot call WritePointData() methods anymore that accept a mesh::utils::MeshDataSet. If you want to write PointData, you need to call WritePointData(const std::string &name, const MESH_FUNCTION &mesh_function) which accepts a MeshFunction MeshFunction. This is because a mesh::utils::MeshDataSet can only provide values at nodes of the mesh, but not at other (Lagrange) points that lie e.g. on the middle of an edge or inside of a cell.
The current version of Paraview (5.6.0) has some issues with cells of order 7 or higher and it just crashes. Also there are some problems with the visualization of higher order Lagrange curves. This is only important when VtkWriter is constructed with codim>0.

How to apply the tesselate filter

As pointed out above, Paraview doesn't approximate higher-order Lagrange cells nicely by default. This can be resolved by applying the tesselate filter in Paraview:

1) After you have loaded a VtkFile with order>1 cells, click on "Filters->Search" in the menu bar. Then search for the "tesselate" filter. 2) Make the following settings for the tesselate filter in the "properties" window:

3) click on "Apply"

Definition at line 272 of file vtk_writer.h.

Member Typedef Documentation

◆ dim_t

Definition at line 274 of file vtk_writer.h.

Constructor & Destructor Documentation

◆ VtkWriter() [1/3]

lf::io::VtkWriter::VtkWriter ( const VtkWriter & )
delete

◆ VtkWriter() [2/3]

lf::io::VtkWriter::VtkWriter ( VtkWriter && )
delete

◆ VtkWriter() [3/3]

lf::io::VtkWriter::VtkWriter ( std::shared_ptr< const mesh::Mesh > mesh,
std::string filename,
dim_t codim = 0,
unsigned char order = 1 )

Construct a new VtkWriter.

Parameters
meshThe underlying mesh that should be written into the VtkFile.
filenameThe filename of the Vtk File
codim(Optional) the codimension of the cells, default is 0, i.e. the codim=0 entities of the mesh are written into the vtk files and data can be attached to these entities. If you set codim=1, the sekelton (only codim=1 entities) will be written into the vtk file and you can visualize data on the skeleton
order(Optional) the order of the cells written to Paraview. By default first order cells are used. Should be in [0,10], see also the documentation of the VtkWriter class for more info.

Definition at line 578 of file vtk_writer.cc.

References mesh_.

◆ ~VtkWriter()

lf::io::VtkWriter::~VtkWriter ( )
inline

Destructor, writes everything into the file and closes it.

Definition at line 784 of file vtk_writer.h.

References filename_, vtk_file_, and lf::io::WriteToFile().

Member Function Documentation

◆ CheckAttributeSetName()

template<class DATA >
void lf::io::VtkWriter::CheckAttributeSetName ( const DATA & data,
const std::string & name )
private

◆ operator=() [1/2]

VtkWriter & lf::io::VtkWriter::operator= ( const VtkWriter & )
delete

◆ operator=() [2/2]

VtkWriter & lf::io::VtkWriter::operator= ( VtkWriter && )
delete

◆ PadWithZeros()

template<int ROWS, class T >
void lf::io::VtkWriter::PadWithZeros ( Eigen::Matrix< T, 3, 1 > & out,
const Eigen::Matrix< T, ROWS, 1 > & in )
staticprivate

Definition at line 992 of file vtk_writer.h.

◆ setBinary()

void lf::io::VtkWriter::setBinary ( bool binary)
inline

Determines whether the Vtk file is written in binary or ASCII mode (default).

Parameters
binarytrue if binary mode should be used, otherwise false

Definition at line 303 of file vtk_writer.h.

References lf::io::VtkFile::ASCII, lf::io::VtkFile::BINARY, lf::io::VtkFile::format, and vtk_file_.

◆ WriteCellData() [1/13]

void lf::io::VtkWriter::WriteCellData ( const std::string & name,
const mesh::utils::MeshDataSet< char > & mds,
char undefined_value = 0 )

Add a new char attribute dataset that attaches data to the cells of the mesh (i.e. to entities with codim = "codim that was specified in constructor of VtkWriter")

Parameters
nameThe name of the attribute set.
mdsThe mesh dataset that that attaches the data to the cells the mesh.
undefined_valueThe value that should be written for a cell to which mds does not attach data (i.e. if mds.DefinedOn() == false)

Definition at line 868 of file vtk_writer.cc.

References WriteScalarCellData().

◆ WriteCellData() [2/13]

void lf::io::VtkWriter::WriteCellData ( const std::string & name,
const mesh::utils::MeshDataSet< double > & mds,
double undefined_value = 0 )

Add a new double attribute dataset that attaches data to the cells of the mesh (i.e. to entities with codim = "codim that was specified in constructor of VtkWriter")

Parameters
nameThe name of the attribute set.
mdsThe mesh dataset that that attaches the data to the cells the mesh.
undefined_valueThe value that should be written for a cell to which mds does not attach data (i.e. if mds.DefinedOn() == false)

Definition at line 892 of file vtk_writer.cc.

References WriteScalarCellData().

◆ WriteCellData() [3/13]

void lf::io::VtkWriter::WriteCellData ( const std::string & name,
const mesh::utils::MeshDataSet< Eigen::Vector2d > & mds,
const Eigen::Vector2d & undefined_value = {0, 0} )

Add a new vector attribute dataset that attaches vectors to the cells of the mesh.

Parameters
nameThe name of the attribute set.
mdsThe mesh dataset that that attaches the data to the cells of the mesh.
undefined_valueThe value that should be written for a cell to which mds does not attach data (i.e. if mds.DefinedOn() == false)

Definition at line 898 of file vtk_writer.cc.

◆ WriteCellData() [4/13]

void lf::io::VtkWriter::WriteCellData ( const std::string & name,
const mesh::utils::MeshDataSet< Eigen::Vector2f > & mds,
const Eigen::Vector2f & undefined_value = {0, 0} )

Add a new vector attribute dataset that attaches vectors to the cells of the mesh.

Parameters
nameThe name of the attribute set.
mdsThe mesh dataset that that attaches the data to the cells of the mesh.
undefined_valueThe value that should be written for a cell to which mds does not attach data (i.e. if mds.DefinedOn() == false)

Definition at line 905 of file vtk_writer.cc.

◆ WriteCellData() [5/13]

void lf::io::VtkWriter::WriteCellData ( const std::string & name,
const mesh::utils::MeshDataSet< Eigen::Vector3d > & mds,
const Eigen::Vector3d & undefined_value = {0, 0, 0} )

Add a new vector attribute dataset that attaches vectors to the cells of the mesh.

Parameters
nameThe name of the attribute set.
mdsThe mesh dataset that that attaches the data to the cells of the mesh.
undefined_valueThe value that should be written for a cell to which mds does not attach data (i.e. if mds.DefinedOn() == false)

Definition at line 912 of file vtk_writer.cc.

◆ WriteCellData() [6/13]

void lf::io::VtkWriter::WriteCellData ( const std::string & name,
const mesh::utils::MeshDataSet< Eigen::Vector3f > & mds,
const Eigen::Vector3f & undefined_value = {0, 0, 0} )

Add a new vector attribute dataset that attaches vectors to the cells of the mesh.

Parameters
nameThe name of the attribute set.
mdsThe mesh dataset that that attaches the data to the cells of the mesh.
undefined_valueThe value that should be written for a cell to which mds does not attach data (i.e. if mds.DefinedOn() == false)

Definition at line 919 of file vtk_writer.cc.

◆ WriteCellData() [7/13]

void lf::io::VtkWriter::WriteCellData ( const std::string & name,
const mesh::utils::MeshDataSet< Eigen::VectorXd > & mds,
const Eigen::VectorXd & undefined_value = Eigen::Vector3d(0, 0, 0) )

Add a new vector attribute dataset that attaches vectors to the cells of the mesh.

Parameters
nameThe name of the attribute set.
mdsThe mesh dataset that that attaches the data to the cells of the mesh.
undefined_valueThe value that should be written for a cell to which mds does not attach data (i.e. if mds.DefinedOn() == false)
Note
This version accepts in principle arbitrary size double Vectors, however only the first three components are visualized!

Definition at line 926 of file vtk_writer.cc.

◆ WriteCellData() [8/13]

void lf::io::VtkWriter::WriteCellData ( const std::string & name,
const mesh::utils::MeshDataSet< Eigen::VectorXf > & mds,
const Eigen::VectorXf & undefined_value = Eigen::Vector3f(0, 0, 0) )

Add a new vector attribute dataset that attaches vectors to the cells of the mesh.

Parameters
nameThe name of the attribute set.
mdsThe mesh dataset that that attaches the data to the cells of the mesh.
undefined_valueThe value that should be written for a cell to which mds does not attach data (i.e. if mds.DefinedOn() == false)
Note
This version accepts in principle arbitrary size float Vectors, however only the first three components are visualized!

Definition at line 933 of file vtk_writer.cc.

◆ WriteCellData() [9/13]

void lf::io::VtkWriter::WriteCellData ( const std::string & name,
const mesh::utils::MeshDataSet< float > & mds,
float undefined_value = 0 )

Add a new float attribute dataset that attaches data to the cells of the mesh (i.e. to entities with codim = "codim that was specified in constructor of VtkWriter")

Parameters
nameThe name of the attribute set.
mdsThe mesh dataset that that attaches the data to the cells the mesh.
undefined_valueThe value that should be written for a cell to which mds does not attach data (i.e. if mds.DefinedOn() == false)

Definition at line 886 of file vtk_writer.cc.

References WriteScalarCellData().

◆ WriteCellData() [10/13]

void lf::io::VtkWriter::WriteCellData ( const std::string & name,
const mesh::utils::MeshDataSet< int > & mds,
int undefined_value = 0 )

Add a new int attribute dataset that attaches data to the cells of the mesh (i.e. to entities with codim = "codim that was specified in constructor of VtkWriter")

Parameters
nameThe name of the attribute set.
mdsThe mesh dataset that that attaches the data to the cells the mesh.
undefined_valueThe value that should be written for a cell to which mds does not attach data (i.e. if mds.DefinedOn() == false)

Definition at line 880 of file vtk_writer.cc.

References WriteScalarCellData().

◆ WriteCellData() [11/13]

void lf::io::VtkWriter::WriteCellData ( const std::string & name,
const mesh::utils::MeshDataSet< unsigned char > & mds,
unsigned char undefined_value = 0 )

Add a new unsigned char attribute dataset that attaches data to the cells of the mesh (i.e. to entities with codim = "codim that was specified in constructor of VtkWriter")

Parameters
nameThe name of the attribute set.
mdsThe mesh dataset that that attaches the data to the cells the mesh.
undefined_valueThe value that should be written for a cell to which mds does not attach data (i.e. if mds.DefinedOn() == false)

Definition at line 862 of file vtk_writer.cc.

References WriteScalarCellData().

◆ WriteCellData() [12/13]

void lf::io::VtkWriter::WriteCellData ( const std::string & name,
const mesh::utils::MeshDataSet< unsigned int > & mds,
unsigned int undefined_value = 0 )

Add a new unsigned int attribute dataset that attaches data to the cells of the mesh (i.e. to entities with codim = "codim that was specified in constructor of VtkWriter")

Parameters
nameThe name of the attribute set.
mdsThe mesh dataset that that attaches the data to the cells the mesh.
undefined_valueThe value that should be written for a cell to which mds does not attach data (i.e. if mds.DefinedOn() == false)

Definition at line 874 of file vtk_writer.cc.

References WriteScalarCellData().

◆ WriteCellData() [13/13]

template<mesh::utils::MeshFunction MESH_FUNCTION>
void lf::io::VtkWriter::WriteCellData ( const std::string & name,
const MESH_FUNCTION & mesh_function )

Sample a MeshFunction at the barycenter of the cell and visualize it as cell data in the Vtk File.

Template Parameters
MESH_FUNCTIONAn object fulfilling the MeshFunction concept. The mesh::utils::MeshFunctionReturnType should be one of
  • unsigned char
  • char
  • unsigned int
  • int
  • float
  • double
  • Eigen::Vector2d
  • Eigen::Vector2f
  • Eigen::Vector3d
  • Eigen::Vector3f
  • Eigen::VectorXd
  • Eigen::VectorXf
Parameters
nameThe name of the dataset, shouldn't contain any spaces!
mesh_functionThe Mesh Function to be sampled.
Note
this function will evaluate the MeshFunction on barycenters of the cells of the mesh, i.e. at codim=0 entities. Some MeshFunction are not well defined on cells (e.g. boundary conditions) and cannot be visualized with this function. (An assert will fail)

Example usage

std::shared_ptr<mesh::Mesh> mesh; // initialize mesh somehow
// construct a first order lagrange fe space
auto fes = std::make_shared<uscalfe::FeSpaceLagrangeO1<double>>(mesh);
// A Mesh function representing the gradient of the first shape function
Eigen::VectorXd x(fes->LocGlobMap().NumDofs());
x[0] = 1.0;
auto mfShapeFun = fe::MeshFunctionGradFE(fes, x);
// write the gradient to vtk (cell based)
io::VtkWriter vtk_writer(mesh, "filename.vtk");
vtk_writer.WriteCellData("shapeFun", mfShapeFun);

Definition at line 923 of file vtk_writer.h.

References lf::io::VtkFile::cell_data, CheckAttributeSetName(), codim_, lf::io::VtkFile::ScalarData< T >::data, lf::io::VtkFile::VectorData< T >::data, lf::base::RefEl::Id(), lf::base::RefEl::kPoint(), lf::base::RefEl::kQuad(), lf::base::RefEl::kSegment(), lf::base::RefEl::kTria(), mesh_, lf::base::RefEl::NodeCoords(), and vtk_file_.

◆ WriteFieldData()

template<class T >
void lf::io::VtkWriter::WriteFieldData ( const std::string & name,
std::vector< T > data )
private

◆ WriteGlobalData() [1/3]

void lf::io::VtkWriter::WriteGlobalData ( const std::string & name,
std::vector< double > data )

Write global data into the vtk file that is not related to the mesh at all.

Parameters
nameName of the global dataset.
dataThe data series that should be written.

This function writes so-called "Field Data" into the vtk file which can be read and visualized by paraview. It is e.g. a convenient way to export the simulation time or the global mesh size.

Definition at line 950 of file vtk_writer.cc.

References WriteFieldData().

◆ WriteGlobalData() [2/3]

void lf::io::VtkWriter::WriteGlobalData ( const std::string & name,
std::vector< float > data )

Write global data into the vtk file that is not related to the mesh at all.

Parameters
nameName of the global dataset.
dataThe data series that should be written.

This function writes so-called "Field Data" into the vtk file which can be read and visualized by paraview. It is e.g. a convenient way to export the simulation time or the global mesh size.

Definition at line 945 of file vtk_writer.cc.

References WriteFieldData().

◆ WriteGlobalData() [3/3]

void lf::io::VtkWriter::WriteGlobalData ( const std::string & name,
std::vector< int > data )

Write global data into the vtk file that is not related to the mesh at all.

Parameters
nameName of the global dataset.
dataThe data series that should be written.

This function writes so-called "Field Data" into the vtk file which can be read and visualized by paraview. It is e.g. a convenient way to export the simulation time or the global mesh size.

Definition at line 940 of file vtk_writer.cc.

References WriteFieldData().

◆ WritePointData() [1/13]

void lf::io::VtkWriter::WritePointData ( const std::string & name,
const mesh::utils::MeshDataSet< char > & mds,
char undefined_value = 0 )

Add a new char attribute dataset that attaches data to the points/nodes of the mesh.

Parameters
nameThe name of the attribute set.
mdsThe mesh dataset that that attaches the data to the points of the mesh.
undefined_valueThe value that should be written for a point to which mds does not attach data (i.e. if mds.DefinedOn() == false)
Note
This overload is only supported if VtkWriter has been instantiated with order=1. For order>1 use WritePointData(const std::string& name, const MESH_FUNCTION& mesh_function).

Definition at line 790 of file vtk_writer.cc.

References WriteScalarPointData().

◆ WritePointData() [2/13]

void lf::io::VtkWriter::WritePointData ( const std::string & name,
const mesh::utils::MeshDataSet< double > & mds,
double undefined_value = 0. )

Add a new double attribute dataset that attaches data to the points/nodes of the mesh.

Parameters
nameThe name of the attribute set.
mdsThe mesh dataset that that attaches the data to the points of the mesh.
undefined_valueThe value that should be written for a point to which mds does not attach data (i.e. if mds.DefinedOn() == false)
Note
This overload is only supported if VtkWriter has been instantiated with order=1. For order>1 use WritePointData(const std::string& name, const MESH_FUNCTION& mesh_function).

Definition at line 814 of file vtk_writer.cc.

References WriteScalarPointData().

◆ WritePointData() [3/13]

void lf::io::VtkWriter::WritePointData ( const std::string & name,
const mesh::utils::MeshDataSet< Eigen::Vector2d > & mds,
const Eigen::Vector2d & undefined_value = {0, 0} )

Add a new vector attribute dataset that attaches vectors to the points/nodes of the mesh.

Parameters
nameThe name of the attribute set.
mdsThe mesh dataset that that attaches the data to the points of the mesh.
undefined_valueThe value that should be written for a point to which mds does not attach data (i.e. if mds.DefinedOn() == false)
Note
This overload is only supported if VtkWriter has been instantiated with order=1. For order>1 use WritePointData(const std::string& name, const MESH_FUNCTION& mesh_function).

Definition at line 820 of file vtk_writer.cc.

◆ WritePointData() [4/13]

void lf::io::VtkWriter::WritePointData ( const std::string & name,
const mesh::utils::MeshDataSet< Eigen::Vector2f > & mds,
const Eigen::Vector2f & undefined_value = {0, 0} )

Add a new vector attribute dataset that attaches vectors to the points/nodes of the mesh.

Parameters
nameThe name of the attribute set.
mdsThe mesh dataset that that attaches the data to the points of the mesh.
undefined_valueThe value that should be written for a point to which mds does not attach data (i.e. if mds.DefinedOn() == false)
Note
This overload is only supported if VtkWriter has been instantiated with order=1. For order>1 use WritePointData(const std::string& name, const MESH_FUNCTION& mesh_function).

Definition at line 827 of file vtk_writer.cc.

◆ WritePointData() [5/13]

void lf::io::VtkWriter::WritePointData ( const std::string & name,
const mesh::utils::MeshDataSet< Eigen::Vector3d > & mds,
const Eigen::Vector3d & undefined_value = {0, 0, 0} )

Add a new vector attribute dataset that attaches vectors to the points/nodes of the mesh.

Parameters
nameThe name of the attribute set.
mdsThe mesh dataset that that attaches the data to the points of the mesh.
undefined_valueThe value that should be written for a point to which mds does not attach data (i.e. if mds.DefinedOn() == false)
Note
This overload is only supported if VtkWriter has been instantiated with order=1. For order>1 use WritePointData(const std::string& name, const MESH_FUNCTION& mesh_function).

Definition at line 834 of file vtk_writer.cc.

◆ WritePointData() [6/13]

void lf::io::VtkWriter::WritePointData ( const std::string & name,
const mesh::utils::MeshDataSet< Eigen::Vector3f > & mds,
const Eigen::Vector3f & undefined_value = {0, 0, 0} )

Add a new vector attribute dataset that attaches vectors to the points/nodes of the mesh.

Parameters
nameThe name of the attribute set.
mdsThe mesh dataset that that attaches the data to the points of the mesh.
undefined_valueThe value that should be written for a point to which mds does not attach data (i.e. if mds.DefinedOn() == false)
Note
This overload is only supported if VtkWriter has been instantiated with order=1. For order>1 use WritePointData(const std::string& name, const MESH_FUNCTION& mesh_function).

Definition at line 841 of file vtk_writer.cc.

◆ WritePointData() [7/13]

void lf::io::VtkWriter::WritePointData ( const std::string & name,
const mesh::utils::MeshDataSet< Eigen::VectorXd > & mds,
const Eigen::VectorXd & undefined_value = Eigen::Vector3d(0, 0, 0) )

Add a new vector attribute dataset that attaches vectors to the points/nodes of the mesh.

Parameters
nameThe name of the attribute set.
mdsThe mesh dataset that that attaches the data to the points of the mesh.
undefined_valueThe value that should be written for a point to which mds does not attach data (i.e. if mds.DefinedOn() == false)
Note
This version accepts in principle arbitrary size double Vectors, however only the first three components are visualized!
This overload is only supported if VtkWriter has been instantiated with order=1. For order>1 use WritePointData(const std::string& name, const MESH_FUNCTION& mesh_function).

Definition at line 848 of file vtk_writer.cc.

◆ WritePointData() [8/13]

void lf::io::VtkWriter::WritePointData ( const std::string & name,
const mesh::utils::MeshDataSet< Eigen::VectorXf > & mds,
const Eigen::VectorXf & undefined_value = Eigen::Vector3f(0, 0, 0) )

Add a new vector attribute dataset that attaches vectors to the points/nodes of the mesh.

Parameters
nameThe name of the attribute set.
mdsThe mesh dataset that that attaches the data to the points of the mesh.
undefined_valueThe value that should be written for a point to which mds does not attach data (i.e. if mds.DefinedOn() == false)
Note
This version accepts in principle arbitrary size float Vectors, however only the first three components are visualized!
This overload is only supported if VtkWriter has been instantiated with order=1. For order>1 use WritePointData(const std::string& name, const MESH_FUNCTION& mesh_function).

Definition at line 855 of file vtk_writer.cc.

◆ WritePointData() [9/13]

void lf::io::VtkWriter::WritePointData ( const std::string & name,
const mesh::utils::MeshDataSet< float > & mds,
float undefined_value = 0.F )

Add a new float attribute dataset that attaches data to the points/nodes of the mesh.

Parameters
nameThe name of the attribute set.
mdsThe mesh dataset that that attaches the data to the points of the mesh.
undefined_valueThe value that should be written for a point to which mds does not attach data (i.e. if mds.DefinedOn() == false)
Note
This overload is only supported if VtkWriter has been instantiated with order=1. For order>1 use WritePointData(const std::string& name, const MESH_FUNCTION& mesh_function).

Definition at line 808 of file vtk_writer.cc.

References WriteScalarPointData().

◆ WritePointData() [10/13]

void lf::io::VtkWriter::WritePointData ( const std::string & name,
const mesh::utils::MeshDataSet< int > & mds,
int undefined_value = 0 )

Add a new unsigned int attribute dataset that attaches data to the points/nodes of the mesh.

Parameters
nameThe name of the attribute set.
mdsThe mesh dataset that that attaches the data to the points of the mesh.
undefined_valueThe value that should be written for a point to which mds does not attach data (i.e. if mds.DefinedOn() == false)
Note
This overload is only supported if VtkWriter has been instantiated with order=1. For order>1 use WritePointData(const std::string& name, const MESH_FUNCTION& mesh_function).

Definition at line 802 of file vtk_writer.cc.

References WriteScalarPointData().

◆ WritePointData() [11/13]

void lf::io::VtkWriter::WritePointData ( const std::string & name,
const mesh::utils::MeshDataSet< unsigned char > & mds,
unsigned char undefined_value = 0 )

Add a new unsigned char attribute dataset that attaches data to the points/nodes of the mesh.

Parameters
nameThe name of the attribute set.
mdsThe mesh dataset that that attaches the data to the points of the mesh.
undefined_valueThe value that should be written for a point to which mds does not attach data (i.e. if mds.DefinedOn() == false)
Note
This overload is only supported if VtkWriter has been instantiated with order=1. For order>1 use WritePointData(const std::string& name, const MESH_FUNCTION& mesh_function).

Definition at line 784 of file vtk_writer.cc.

References WriteScalarPointData().

◆ WritePointData() [12/13]

void lf::io::VtkWriter::WritePointData ( const std::string & name,
const mesh::utils::MeshDataSet< unsigned int > & mds,
unsigned int undefined_value = 0 )

Add a new unsigned int attribute dataset that attaches data to the points/nodes of the mesh.

Parameters
nameThe name of the attribute set.
mdsThe mesh dataset that that attaches the data to the points of the mesh.
undefined_valueThe value that should be written for a point to which mds does not attach data (i.e. if mds.DefinedOn() == false)
Note
This overload is only supported if VtkWriter has been instantiated with order=1. For order>1 use WritePointData(const std::string& name, const MESH_FUNCTION& mesh_function).

Definition at line 796 of file vtk_writer.cc.

References WriteScalarPointData().

◆ WritePointData() [13/13]

template<mesh::utils::MeshFunction MESH_FUNCTION>
void lf::io::VtkWriter::WritePointData ( const std::string & name,
const MESH_FUNCTION & mesh_function )

Sample a MeshFunction at points of the mesh and possible edges/interior (for order>1).

Template Parameters
MESH_FUNCTIONAn object fulfilling the MeshFunction concept. The mesh::utils::MeshFunctionReturnType should be one of
  • unsigned char
  • char
  • unsigned int
  • int
  • float
  • double
  • Eigen::Vector2d
  • Eigen::Vector2f
  • Eigen::Vector3d
  • Eigen::Vector3f
Parameters
nameThe name of the dataset, cannot contain spaces!
mesh_functionThe Mesh Function to be sampled.
Note
this function will evaluate the MeshFunction on the points of the mesh, i.e. at entities with codim=dimMesh. If the VtkWriter has been instantiated with order>1, the mesh function will also be evaluated on edges and possibly in the interior of the elements. Some MeshFunction are not well defined on points, e.g. the gradient of the solution of a BVP is not well defined on the points of the mesh. In this case WriteCellData(const std::string &name, const MESH_FUNCTION &mesh_function) may be more appropriate.

Example usage

std::shared_ptr<mesh::Mesh> mesh; // initialize mesh somehow
// construct a first order lagrange fe space
auto fes = std::make_shared<uscalfe::FeSpaceLagrangeO1<double>>(mesh);
// A Mesh function representing the first shape function of fes
Eigen::VectorXd x(fes->LocGlobMap().NumDofs());
x[0] = 1.0;
auto mfShapeFun = fe::MeshFunctionFE(fes, x);
// A Meshfunction representing the gradient of the first shape function:
auto mfGradShapeFun = fe::MeshFunctionGradFE(fes, x);
// write to VTK
io::VtkWriter vtk_writer(mesh, "filename.vtk");
vtk_writer.WritePointData("shapeFun", mfShapeFun);
// Trying to sample the gradient directly on the points of the mesh
// will fail an assert, see note above:
// vtk_writer.WritePointData("GradShapeFun", mfGradShapeFun);
// Sample the geometry and the mesh function with 3rd order lagrange basis
// functions (see documentation of VtkWriter):
io::VtkWriter vtk_writerOrder3(mesh, "order3.vtk", 0, 3);
auto mfTrig = mesh::utils::MeshFunctionGlobal(
[](const auto& x) { return std::sin(x[0]) * std::cos(x[1]); });
vtk_writer.WritePointData("mfTrig", mfTrig);

Definition at line 836 of file vtk_writer.h.

References aux_node_offset_, aux_nodes_, CheckAttributeSetName(), codim_, lf::io::VtkFile::ScalarData< T >::data, lf::io::VtkFile::VectorData< T >::data, lf::base::RefEl::kTria(), mesh_, order_, lf::io::VtkFile::point_data, lf::io::VtkFile::UnstructuredGrid::points, lf::io::VtkFile::unstructured_grid, and vtk_file_.

◆ WriteScalarCellData()

template<class T >
void lf::io::VtkWriter::WriteScalarCellData ( const std::string & name,
const mesh::utils::MeshDataSet< T > & mds,
T undefined_value )
private

◆ WriteScalarPointData()

template<class T >
void lf::io::VtkWriter::WriteScalarPointData ( const std::string & name,
const mesh::utils::MeshDataSet< T > & mds,
T undefined_value )
private

◆ WriteVectorCellData()

template<int ROWS, class T >
void lf::io::VtkWriter::WriteVectorCellData ( const std::string & name,
const mesh::utils::MeshDataSet< Eigen::Matrix< T, ROWS, 1 > > & mds,
const Eigen::Matrix< T, ROWS, 1 > & undefined_value )
private

◆ WriteVectorPointData()

template<int ROWS, class T >
void lf::io::VtkWriter::WriteVectorPointData ( const std::string & name,
const mesh::utils::MeshDataSet< Eigen::Matrix< T, ROWS, 1 > > & mds,
const Eigen::Matrix< T, ROWS, 1 > & undefined_value )
private

Member Data Documentation

◆ aux_node_offset_

std::array<mesh::utils::CodimMeshDataSet<unsigned int>, 2> lf::io::VtkWriter::aux_node_offset_
private

Definition at line 799 of file vtk_writer.h.

Referenced by WritePointData().

◆ aux_nodes_

std::array<Eigen::MatrixXd, 5> lf::io::VtkWriter::aux_nodes_
private

Definition at line 795 of file vtk_writer.h.

Referenced by WritePointData().

◆ codim_

dim_t lf::io::VtkWriter::codim_
private

◆ filename_

std::string lf::io::VtkWriter::filename_
private

Definition at line 789 of file vtk_writer.h.

Referenced by ~VtkWriter().

◆ mesh_

std::shared_ptr<const mesh::Mesh> lf::io::VtkWriter::mesh_
private

◆ order_

unsigned char lf::io::VtkWriter::order_
private

Definition at line 791 of file vtk_writer.h.

Referenced by WritePointData(), WriteScalarPointData(), and WriteVectorPointData().

◆ vtk_file_

VtkFile lf::io::VtkWriter::vtk_file_
private

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