intx_in_plane_test.cpp

Go to the documentation of this file.

/*
 * intx_in_plane_test.cpp
 *
 * Created on: Oct 4, 2012
 */
#include <iostream>
#include <sstream>
#include <ctime>
#include <cstdlib>
#include <cstdio>
#include <cstring>
#include "moab/Core.hpp"
#include "moab/Interface.hpp"
#include "moab/IntxMesh/Intx2MeshInPlane.hpp"
#include "moab/IntxMesh/IntxUtils.hpp"
#include "TestUtil.hpp"
#include <cmath>
 
using namespace moab;
 
int main( int argc, char* argv[] )
{
 // check command line arg
 const char* filename_mesh1 = STRINGIFY( MESHDIR ) "/mbcslam/m1.vtk";
 const char* filename_mesh2 = STRINGIFY( MESHDIR ) "/mbcslam/m2.vtk";
 const char* newFile = "intx1.vtk";
 const char* edgesFile = "polyWithEdges.vtk";
 if( argc == 4 )
 {
 filename_mesh1 = argv[1];
 filename_mesh2 = argv[2];
 newFile = argv[3];
 }
 else
 {
 printf( "Usage: %s <mesh_filename1> <mesh_filename2> <newFile>\n", argv[0] );
 if( argc != 1 ) return 1;
 printf( "No files specified. Defaulting to: %s %s %s\n", filename_mesh1, filename_mesh2, newFile );
 }
 
 // read meshes in 2 file sets
 Core moab;
 Interface* mb = &moab; // global
 EntityHandle sf1, sf2;
 ErrorCode rval = mb->create_meshset( MESHSET_SET, sf1 );
 if( MB_SUCCESS != rval ) return 1;
 rval = mb->create_meshset( MESHSET_SET, sf2 );
 if( MB_SUCCESS != rval ) return 1;
 rval = mb->load_file( filename_mesh1, &sf1 );
 if( MB_SUCCESS != rval ) return 1;
 rval = mb->load_file( filename_mesh2, &sf2 );
 if( MB_SUCCESS != rval ) return 1;
 
 EntityHandle outputSet;
 rval = mb->create_meshset( MESHSET_SET, outputSet );
 if( MB_SUCCESS != rval ) return 1;
 
 Intx2MeshInPlane worker( mb );
 IntxAreaUtils areaAdaptor;
 rval = areaAdaptor.positive_orientation( mb, sf1, -1 );MB_CHK_ERR( rval );
 rval = areaAdaptor.positive_orientation( mb, sf2, -1 );MB_CHK_ERR( rval );
 
 worker.set_error_tolerance( 1.e-5 );
 rval = worker.FindMaxEdges( sf1, sf2 );MB_CHK_ERR( rval );
 // worker.enable_debug();
 rval = worker.intersect_meshes( sf1, sf2, outputSet );
 if( MB_SUCCESS != rval ) return 1;
 rval = mb->write_mesh( newFile, &outputSet, 1 );
 if( MB_SUCCESS != rval ) return 1;
 
 // retrieve polygons and get adjacent edges
 Range polygons;
 rval = mb->get_entities_by_type( outputSet, MBPOLYGON, polygons );
 if( MB_SUCCESS != rval ) return 1;
 
 Range edges;
 rval = mb->get_adjacencies( polygons, 1, true, edges, Interface::UNION );
 if( MB_SUCCESS != rval ) return 1;
 
 std::cout << "number of edges:" << edges.size() << "\n";
 // add edges to the output set
 rval = mb->add_entities( outputSet, edges );
 if( MB_SUCCESS != rval ) return 1;
 rval = mb->write_mesh( edgesFile, &outputSet, 1 );
 if( MB_SUCCESS != rval ) return 1;
 return 0;
}