MeshKit
1.0
|
00001 00007 #include "meshkit/MKCore.hpp" 00008 #include "meshkit/MeshOp.hpp" 00009 #include "meshkit/EBMesher.hpp" 00010 #include "meshkit/SCDMesh.hpp" 00011 #include "meshkit/ModelEnt.hpp" 00012 00013 //#include "CGMApp.hpp" 00014 using namespace MeshKit; 00015 00016 #include "TestUtil.hpp" 00017 00018 #ifdef HAVE_ACIS 00019 #define DEFAULT_TEST_FILE "sphere.sat" 00020 #elif defined(HAVE_OCC) 00021 #define DEFAULT_TEST_FILE "sphere.stp" 00022 #endif 00023 00024 const bool debug_ebmesher = false; 00025 00026 int load_and_mesh(const char *input_filename, 00027 const char *output_filename, 00028 int whole_geom, int* n_intervals, int mesh_based_geom, 00029 double box_increase, int vol_frac_res); 00030 00031 int main(int argc, char **argv) 00032 { 00033 // check command line arg 00034 std::string input_filename; 00035 const char *output_filename = NULL; 00036 int whole_geom = 1; 00037 int n_interval[3] = {10, 10, 10}; 00038 int mesh_based_geom = 0; 00039 double box_increase = .03; 00040 int vol_frac_res = 0; 00041 00042 if (argc > 2 && argc < 11) { 00043 input_filename += argv[1]; 00044 if (argc > 2) whole_geom = atoi(argv[2]); 00045 if (argc > 3) n_interval[0] = atoi(argv[3]); 00046 if (argc > 4) n_interval[1] = atoi(argv[4]); 00047 if (argc > 5) n_interval[2] = atoi(argv[5]); 00048 if (argc > 6) mesh_based_geom = atoi(argv[6]); 00049 if (argc > 7) output_filename = argv[7]; 00050 if (argc > 8) box_increase = atof(argv[8]); 00051 if (argc > 9) vol_frac_res = atoi(argv[9]); 00052 } 00053 else { 00054 std::cout << "Usage: " << argv[0] << "<input_geom_filename> <whole_geom> {x: # of intervals} {y: # of intervals} {z: # of intervals} {mesh_based_geom} {output_mesh_filename} {box_size_increase} {vol_frac_res}" << std::endl; 00055 std::cout << "<input_geom_filename> : input geometry file name" << std::endl; 00056 std::cout << "<whole_geom> : make mesh for whole geom or individually(1/0), default whole geom(1)" << std::endl; 00057 std::cout << "{x/y/z: # ofintervals} : optional argument. # of intervals. if it is not set, set to 10." << std::endl; 00058 std::cout << "<mesh_based_geom> : use mesh based geometry(1/0), default not-use(0)" << std::endl; 00059 std::cout << "{output_mesh_filename} : optional argument. if it is not set, dosn't export. can output mesh file (e.g. output.vtk.)" << std::endl; 00060 std::cout << "{box size increase} : optional argument. Cartesian mesh box increase form geometry. default 0.03" << std::endl; 00061 std::cout << "{vol_frac_res} : optional argument, volume fraction resolution of boundary cells for each material, you can specify it as # of divisions (e.g. 4)." << std::endl; 00062 std::cout << std::endl; 00063 if (argc != 1) return 1; 00064 std::cout << "No file specified. Defaulting to: " << DEFAULT_TEST_FILE << std::endl; 00065 std::string file_name = TestDir + "/" + DEFAULT_TEST_FILE; 00066 input_filename += TestDir; 00067 input_filename += "/"; 00068 input_filename += DEFAULT_TEST_FILE; 00069 } 00070 00071 if (load_and_mesh(input_filename.c_str(), output_filename, 00072 whole_geom, n_interval, mesh_based_geom, box_increase, vol_frac_res)) return 1; 00073 00074 return 0; 00075 } 00076 00077 int load_and_mesh(const char *input_filename, 00078 const char *output_filename, 00079 int whole_geom, int* n_interval, int mesh_based_geom, 00080 double box_increase, int vol_frac_res) 00081 { 00082 bool result; 00083 time_t start_time, load_time, mesh_time, vol_frac_time, 00084 export_time, query_time_techX, query_time; 00085 00086 // start up MK and load the geometry 00087 MKCore mk; 00088 time(&start_time); 00089 // CGMApp::instance()->attrib_manager()->auto_flag( CUBIT_TRUE ); 00090 mk.load_mesh(input_filename, NULL, 0, 0, 0, true); 00091 time(&load_time); 00092 00093 if (debug_ebmesher) { 00094 mk.save_mesh("input.vtk"); 00095 } 00096 00097 // get the volumes 00098 MEntVector vols; 00099 mk.get_entities_by_dimension(3, vols); 00100 00101 // make EBMesher 00102 EBMesher *ebm = (EBMesher*) mk.construct_meshop("EBMesher", vols); 00103 ebm->use_whole_geom(whole_geom); 00104 ebm->use_mesh_geometry(mesh_based_geom); 00105 ebm->set_num_interval(n_interval); 00106 ebm->increase_box(box_increase); 00107 if (mesh_based_geom) ebm->set_obb_tree_box_dimension(); 00108 00109 // mesh embedded boundary mesh, by calling execute 00110 mk.setup_and_execute(); 00111 time(&mesh_time); 00112 00113 // caculate volume fraction, only for geometry input 00114 if (vol_frac_res > 0) { 00115 result = ebm->get_volume_fraction(vol_frac_res); 00116 if (!result) { 00117 std::cerr << "Couldn't get volume fraction." << std::endl; 00118 return 1; 00119 } 00120 } 00121 time(&vol_frac_time); 00122 00123 // export mesh 00124 if (output_filename != NULL) { 00125 ebm->export_mesh(output_filename); 00126 } 00127 time(&export_time); 00128 00129 if (whole_geom && debug_ebmesher) { 00130 // techX query function test 00131 double boxMin[3], boxMax[3]; 00132 int nDiv[3]; 00133 std::map< CutCellSurfEdgeKey, std::vector<double>, LessThan > mdCutCellSurfEdge; 00134 std::vector<int> vnInsideCellTechX; 00135 00136 ebm->get_grid_and_edges_techX(boxMin, boxMax, nDiv, 00137 mdCutCellSurfEdge, vnInsideCellTechX); 00138 time(&query_time_techX); 00139 00140 // multiple intersection fraction query test 00141 std::map< CutCellSurfEdgeKey, std::vector<double>, LessThan > mdCutCellEdge; 00142 std::vector<int> vnInsideCell; 00143 result = ebm->get_grid_and_edges(boxMin, boxMax, nDiv, 00144 mdCutCellEdge, vnInsideCell); 00145 if (!result) { 00146 std::cerr << "Couldn't get mesh information." << std::endl; 00147 return 1; 00148 } 00149 time(&query_time); 00150 std::cout << "# of TechX cut-cell surfaces: " << mdCutCellSurfEdge.size() 00151 << ", # of nInsideCell: " << vnInsideCell.size()/3 << std::endl; 00152 } 00153 00154 std::cout << "EBMesh is succesfully finished." << std::endl; 00155 std::cout << "Time including loading: " 00156 << difftime(mesh_time, start_time) 00157 << " secs, Time excluding loading: " 00158 << difftime(mesh_time, load_time) 00159 << " secs, Time volume fraction: " 00160 << difftime(vol_frac_time, mesh_time) << " secs"; 00161 00162 if (output_filename != NULL) { 00163 std::cout << ", Time export mesh: " 00164 << difftime(export_time, vol_frac_time) << " secs"; 00165 } 00166 00167 if (whole_geom && debug_ebmesher) { 00168 std::cout << ", TechX query time: " 00169 << difftime(query_time_techX, export_time) 00170 << " secs, multiple intersection fraction query (elems, edge-cut fractions): " 00171 << difftime(query_time, query_time_techX) << " secs."; 00172 } 00173 00174 std::cout << std::endl; 00175 mk.clear_graph(); 00176 00177 return 0; 00178 }