Dear all,
I am trying to read a grid from a vtk file created by Paraview.
Unfortunately, this does not work.
Moreover, reading a vtk file that has been created using deal.ii also does
not work for me. It seems to me that deal.ii reads the file but does not
connect the faces to form a mesh. This is underlined by the attached
minimal toy example which creates a two dimensional cross, prints the
cross' interior and boundary faces, writes it to a vtk file, reads the mesh
from the file, and again prints the interior and boundary faces. Execution
of the program gives the following output:
Self Interface 1 0 -0.5
Self Interface 2 -0.5 0
Self Interface 3 0.5 0
Self Boundaryface 1 1 -0.5
Self Interface 4 0 0.5
Self Boundaryface 2 -0.5 1
Self Boundaryface 3 1 0.5
Self Boundaryface 4 0.5 1
Self Boundaryface 5 -1 -0.5
Self Boundaryface 6 -1.5 0
Self Boundaryface 7 -1 0.5
Self Boundaryface 8 1.5 0
Self Boundaryface 9 -0.5 -1
Self Boundaryface 10 0 -1.5
Self Boundaryface 11 0.5 -1
Self Boundaryface 12 0 1.5
Read Boundaryface 1 0 -0.5
Read Boundaryface 2 -0.5 0
Read Boundaryface 3 0.5 0
Read Boundaryface 4 0 0.5
Read Boundaryface 5 -1 -0.5
Read Boundaryface 6 -1.5 0
Read Boundaryface 7 -0.5 0
Read Boundaryface 8 -1 0.5
Read Boundaryface 9 1 -0.5
Read Boundaryface 10 0.5 0
Read Boundaryface 11 1.5 0
Read Boundaryface 12 1 0.5
Read Boundaryface 13 0 -1.5
Read Boundaryface 14 -0.5 -1
Read Boundaryface 15 0.5 -1
Read Boundaryface 16 0 -0.5
Read Boundaryface 17 0 0.5
Read Boundaryface 18 -0.5 1
Read Boundaryface 19 0.5 1
Read Boundaryface 20 0 1.5
Is there any quick fix to this issue (except for not using vtk)? I failed
with merging the mesh with itself to eliminate the double boundary faces
(e.g. 1 and 16).
Thank you very much!
Best regards,
Andreas
--
The deal.II project is located at http://www.dealii.org/
For mailing list/forum options, see
https://groups.google.com/d/forum/dealii?hl=en
---
You received this message because you are subscribed to the Google Groups
"deal.II User Group" group.
To unsubscribe from this group and stop receiving emails from it, send an email
to [email protected].
To view this discussion on the web visit
https://groups.google.com/d/msgid/dealii/bfb4e975-63ab-4be0-9e3e-f9bf0a131ce6%40googlegroups.com.
#include <iostream>
#include <fstream>
#include <deal.II/grid/grid_generator.h>
#include <deal.II/base/conditional_ostream.h>
#include <deal.II/distributed/tria.h>
#include <deal.II/grid/grid_in.h>
#include <deal.II/grid/grid_out.h>
namespace ElasticityDG
{
using namespace dealii;
template <int dim>
class InteriorPenaltyProblem
{
public:
InteriorPenaltyProblem();
private:
MPI_Comm mpi_communicator;
ConditionalOStream pcout;
parallel::distributed::Triangulation<dim> triangulationSelf, triangulationRead;
};
template <int dim>
InteriorPenaltyProblem<dim>::InteriorPenaltyProblem()
:
mpi_communicator(MPI_COMM_WORLD),
pcout(std::cout, (Utilities::MPI::this_mpi_process(mpi_communicator) == 0)),
triangulationSelf (mpi_communicator), triangulationRead(mpi_communicator)
{
GridGenerator::hyper_cross(triangulationSelf, std::vector<unsigned int>(2*dim, 1));
GridOut grid_writer;
std::ofstream output_file("./DealOutput.vtk");
grid_writer.write_vtk(triangulationSelf, output_file);
GridIn<dim> grid_reader;
grid_reader.attach_triangulation(triangulationRead);
std::ifstream input_file("./DealOutput.vtk");
grid_reader.read_vtk(input_file);
int interfacesSelf = 0, boundaryfacesSelf = 0;
int interfacesRead = 0, boundaryfacesRead = 0;
Point<dim> barycenter;
for(auto face = triangulationSelf.begin_active_face(); face != triangulationSelf.end_face(); ++face)
{
barycenter = face -> center();
if (face -> at_boundary()) pcout << "Self Boundaryface " << ++boundaryfacesSelf << " " << barycenter << std::endl;
if (!(face -> at_boundary())) pcout << "Self Interface " << ++interfacesSelf << " " << barycenter << std::endl;
}
pcout << std::endl;
for(auto face = triangulationRead.begin_active_face(); face != triangulationRead.end_face(); ++face)
{
barycenter = face -> center();
if (face -> at_boundary()) pcout << "Read Boundaryface " << ++boundaryfacesRead << " " << barycenter << std::endl;
if (!(face -> at_boundary())) pcout << "Read Interface " << ++interfacesRead << " " << barycenter << std::endl;
}
}
}
int main(int argc, char *argv[])
{
try
{
using namespace dealii;
using namespace ElasticityDG;
const int spatial_dimension = 2;
Utilities::MPI::MPI_InitFinalize mpi_initialization(argc, argv, numbers::invalid_unsigned_int);
InteriorPenaltyProblem<spatial_dimension> test1;
}
catch (std::exception &exc)
{
std::cerr << std::endl << std::endl
<< "----------------------------------------------------"
<< std::endl;
std::cerr << "Exception on processing: " << std::endl
<< exc.what() << std::endl
<< "Aborting!" << std::endl
<< "----------------------------------------------------"
<< std::endl;
return 1;
}
catch (...)
{
std::cerr << std::endl << std::endl
<< "----------------------------------------------------"
<< std::endl;
std::cerr << "Unknown exception!" << std::endl
<< "Aborting!" << std::endl
<< "----------------------------------------------------"
<< std::endl;
return 1;
}
return 0;
}
