Commit f30759ed authored by Jakub Klinkovský's avatar Jakub Klinkovský
Browse files

used the VTU writer from TNL for output in cgal-mesher

closes #2
parent 3131ceae
Loading
Loading
Loading
Loading
Loading
+8 −0
Original line number Diff line number Diff line
@@ -26,6 +26,7 @@ endif()

# bundled packages
include_directories(libs/lyra/include)
include_directories(libs/tnl/src libs/tnl/src/3rdparty)

cmake_policy(SET CMP0074 NEW)
find_package(CGAL REQUIRED OPTIONAL_COMPONENTS Qt5)
@@ -48,3 +49,10 @@ if( CGAL_ACTIVATE_CONCURRENT_MESH_3 AND TBB_FOUND )
    target_compile_definitions(cgal-mesher PRIVATE -DCGAL_CONCURRENT_MESH_3)
    target_link_libraries(cgal-mesher PUBLIC CGAL::TBB_support)
endif()

find_package(ZLIB)
if( ZLIB_FOUND )
    target_compile_definitions(cgal-mesher PRIVATE "-DHAVE_ZLIB")
    target_include_directories(cgal-mesher PUBLIC ${ZLIB_INCLUDE_DIRS})
    target_link_libraries(cgal-mesher PUBLIC ${ZLIB_LIBRARIES})
endif()

cgal-mesh-to-tnl.h

0 → 100644
+88 −0
Original line number Diff line number Diff line
#pragma once

#include <TNL/Meshes/MeshBuilder.h>
#include <TNL/Meshes/Mesh.h>
#include <TNL/Meshes/DefaultConfig.h>
#include <TNL/Meshes/Topologies/Tetrahedron.h>

template< typename Cell,
          int SpaceDimension = Cell::dimension,
          typename Real = double,
          typename GlobalIndex = int,
          typename LocalIndex = GlobalIndex >
struct MinimalConfig
: public TNL::Meshes::DefaultConfig< Cell, SpaceDimension, Real, GlobalIndex, LocalIndex >
{
	static constexpr bool subentityStorage( int entityDimension, int subentityDimension )
	{
	   return subentityDimension == 0 && entityDimension == Cell::dimension;
	}

	static constexpr bool superentityStorage( int entityDimension, int superentityDimension )
	{
	   return false;
	}

	static constexpr bool entityTagsStorage( int entityDimension )
	{
	   return false;
	}

	static constexpr bool dualGraphStorage()
	{
	   return false;
	}
};

template< typename Real >
using TetrahedralMesh = TNL::Meshes::Mesh< MinimalConfig<
                        TNL::Meshes::Topologies::Tetrahedron,
                        3,
                        Real,
                        std::size_t,
                        short int > >;

template< typename C3T3 >
TetrahedralMesh< typename C3T3::Triangulation::Geom_traits::FT >
cgal_mesh_to_tnl( const C3T3& c3t3 )
{
	using MeshType = TetrahedralMesh< typename C3T3::Triangulation::Geom_traits::FT >;
	using PointType = typename MeshType::PointType;

	// initialize TNL mesh builder
	TNL::Meshes::MeshBuilder< MeshType > builder;
	// NOTE: c3t3's triangulation contains even points outside the domain (on
	// the bounding ball), which is why this function must work with c3t3 even
	// though working with just the triangulation would be easier
	builder.setPointsCount(c3t3.triangulation().number_of_vertices() - c3t3.number_of_far_points());
	builder.setCellsCount(c3t3.number_of_cells());

	// index mapping for vertices in the complex
	std::unordered_map<typename C3T3::Triangulation::Vertex_handle, std::size_t> V;

	// set points
	std::size_t idx = 0;
	for( auto vit : c3t3.triangulation().finite_vertex_handles() )
	{
		if(vit->in_dimension() <= -1) continue;
		V[vit] = idx;
		builder.setPoint(idx++, {vit->point()[0], vit->point()[1], vit->point()[2]});
	}

	// set cell seeds
	idx = 0;
	for( auto cit = c3t3.cells_in_complex_begin();
		cit != c3t3.cells_in_complex_end();
		cit++ )
	{
		auto& seed = builder.getCellSeed(idx++);
		for (int i=0; i<4; i++)
			seed.setCornerId(i, V[cit->vertex(i)]);
	}

	// build the mesh
	MeshType mesh;
	if (!builder.build(mesh))
		throw std::runtime_error("TNL's MeshBuilder failed");
	return mesh;
}
+15 −24
Original line number Diff line number Diff line
@@ -9,12 +9,13 @@
#include <CGAL/Mesh_complex_3_in_triangulation_3.h>
#include <CGAL/Mesh_criteria_3.h>
#include <CGAL/Polyhedron_3.h>
//#include <CGAL/Polyhedral_mesh_domain_3.h>
#include <CGAL/Polyhedral_mesh_domain_with_features_3.h>
#include <CGAL/make_mesh_3.h>
#include <CGAL/refine_mesh_3.h>
#include <CGAL/tetrahedral_remeshing.h>
#include <CGAL/IO/output_to_vtu.h>

#include <TNL/Meshes/Writers/VTUWriter.h>
#include "cgal-mesh-to-tnl.h"

#include "lyra/lyra.hpp"

@@ -480,33 +481,23 @@ void polyDomain(const std::filesystem::path & input_file,
	// Remeshing
	if (args.remesh_iterations > 0) {
		std::cout << "Remeshing..." << std::endl;
		// Extract triangulation from c3t3 and discard the c3t3 object
		using Triangulation_3 = CGAL::Triangulation_3<typename Tr::Geom_traits, typename Tr::Triangulation_data_structure>;
		Triangulation_3 t3 = CGAL::convert_to_triangulation_3(std::move(c3t3));

		// Remesh
		CGAL::tetrahedral_isotropic_remeshing(t3,
		// NOTE: passing c3t3 directly to the remeshing function is an undocumented feature of CGAL
		CGAL::tetrahedral_isotropic_remeshing(c3t3,
			args.remesh_target_edge_length,
			CGAL::parameters::number_of_iterations(args.remesh_iterations)
			.remesh_boundaries(args.remesh_boundaries)
		);

		// Create another c3t3 for output (but of different type)
		CGAL::Mesh_complex_3_in_triangulation_3<Triangulation_3, int, int> c3t3;
		c3t3.triangulation() = std::move(t3);
		c3t3.rescan_after_load_of_triangulation();

		// Output
		std::ofstream out(output_file);
		output_to_vtu(out, c3t3);
		out.close();
	}
	else {
		// Output

	// Output to VTU
	// (we use TNL instead of CGAL because the latter uses the "appended" data
	// format which is difficult to parse for libraries except VTK itself)
	using MeshType = TetrahedralMesh<typename C3T3::Triangulation::Geom_traits::FT>;
	const MeshType tnlMesh = cgal_mesh_to_tnl(c3t3);
	std::ofstream out(output_file);
		output_to_vtu(out, c3t3);
		out.close();
	}
	using MeshWriter = TNL::Meshes::Writers::VTUWriter<MeshType>;
	MeshWriter writer(out);
	writer.template writeEntities<3>(tnlMesh);
}

int main(int argc, char * argv[])