MeshKit
1.0
|
00001 00007 // ia_main.cpp 00008 // test IntervalMatching interface IAInterface for MeshKit 00009 00010 #include "TestUtil.hpp" 00011 #include "meshkit/MKCore.hpp" 00012 00013 #include "meshkit/IAInterface.hpp" 00014 #include "meshkit/IAVariable.hpp" 00015 #include "meshkit/TFIMapping.hpp" 00016 00017 #include <stdio.h> 00018 #include <iostream> 00019 #ifdef HAVE_ACIS 00020 #define INTASSIGN_TEST_FILE "quadface.sat" 00021 #elif defined(HAVE_OCC) 00022 #define INTASSIGN_TEST_FILE "quadface.stp" 00023 #endif 00024 00025 MeshKit::MKCore *mk; 00026 00027 MeshKit::IAInterface *new_ia_interface() 00028 { 00029 return 00030 (MeshKit::IAInterface*) mk->construct_meshop("IntervalAssignment"); 00031 } 00032 00033 void delete_ia_interface(MeshKit::IAInterface *) 00034 { 00035 //nada, mk takes care of it 00036 ; 00037 } 00038 00039 bool check_solution_correctness( MeshKit::IAInterface *ia_interface, 00040 std::vector< std::pair<int,int> > &correct_solution) 00041 { 00042 const bool verbose_output = true; 00043 const bool debug = false; 00044 bool all_correct = true; 00045 MeshKit::IAInterface::VariableVec::const_iterator b = ia_interface->variables_begin(); 00046 MeshKit::IAInterface::VariableVec::const_iterator e = ia_interface->variables_end(); 00047 MeshKit::IAInterface::VariableVec::const_iterator i = b; 00048 unsigned int c = 0; 00049 if (debug) 00050 std::cout << "Checking Solution Correctness" << std::endl; 00051 for ( ; i != e; ++i, ++c ) 00052 { 00053 const MeshKit::IAVariable *v = *i; 00054 assert(v); 00055 const int x = v->get_solution(); 00056 assert(c < correct_solution.size() ); 00057 const int lo = correct_solution[c].first; 00058 const int hi = correct_solution[c].second; 00059 if (debug) 00060 std::cout << "Checking variable " << c << " solution " << x << " in " 00061 << "[" << lo << "," << hi << "]?" << std::endl; 00062 if (x < lo) 00063 { 00064 if (verbose_output) 00065 std::cout << "ERROR: Variable " << c << " solution " << x << " BELOW " 00066 << "[" << lo << "," << hi << "]" << std::endl; 00067 all_correct = false; 00068 } 00069 if (x > hi) 00070 { 00071 if (verbose_output) 00072 std::cout << "ERROR: Variable " << c << " solution " << x << " ABOVE " 00073 << "[" << lo << "," << hi << "]" << std::endl; 00074 all_correct = false; 00075 } 00076 } 00077 if (debug) 00078 std::cout << "done checking solution correctness." << std::endl; 00079 return all_correct; 00080 } 00081 00082 void set_decoupled_pairs(MeshKit::IAInterface *ia_interface, 00083 int num_pairs, double goal1, double goal2, 00084 std::vector< std::pair<int,int> > &correct_solution) 00085 { 00086 // trivial 2-sided mapping problem 00087 // we can make multiple pairs, each pair is independent, 00088 // and pair i (in 0..num_pairs-1) has sides with one curve each with goals 00089 // i+goal1 and i+goal2, 00090 // 00091 // test scalability, relaxed nlp: 100,000 constraints in 1 second. milp: 40 variables in 1 second, grows exponentially! 00092 for (int i = 0; i<num_pairs; ++i) 00093 { 00094 // goals x_{2i} = 2, x_{2i+1} = 2 00095 // x_{2i}, goal: i + goal1 00096 const double g1 = i + goal1; 00097 const double g2 = i + goal2; 00098 MeshKit::IAVariable *v1 = ia_interface->create_variable( NULL, MeshKit::SOFT, g1); 00099 MeshKit::IAVariable *v2 = ia_interface->create_variable( NULL, MeshKit::SOFT, g2); 00100 const double compromise = sqrt( g1 * g2 ); 00101 double lo = floor(compromise); 00102 if ( ( compromise - lo ) < 0.1 ) 00103 --lo; 00104 if ( lo < 1. ) 00105 lo = 1.; 00106 double hi = ceil(compromise); 00107 if ( (hi - compromise) < 0.1 ) 00108 ++hi; 00109 correct_solution.push_back( std::make_pair( lo, hi ) ); 00110 correct_solution.push_back( std::make_pair( lo, hi ) ); 00111 00112 // constrain x_{2i} - x_{2i+1} = 0 00113 MeshKit::IAInterface::IAVariableVec side1, side2; 00114 side1.push_back(v1); 00115 side2.push_back(v2); 00116 ia_interface->constrain_sum_equal(side1, side2); 00117 } 00118 } 00119 00120 00121 void set_mapping_chain( MeshKit::IAInterface *ia_interface, const int num_sides, 00122 const bool grow_goal_by_i, 00123 const int goal_m1, const int goal_m2, 00124 const int num_curve_min, const int num_curve_max ) 00125 { 00126 // test problem 3, sides with more than one variable, with random goals 00127 printf("constructing coupled test problem - mapping chain\n"); 00128 srand(10234); 00129 MeshKit::IAInterface::IAVariableVec side1, side2; 00130 int num_vars = 0; 00131 for (int i = 0; i<num_sides; ++i) 00132 { 00133 // move side2 to side1 00134 side2.swap( side1 ); 00135 00136 // create new side2 00137 side2.clear(); 00138 assert( num_curve_min > 0 ); 00139 int num_curves = num_curve_min; 00140 if ( num_curve_max > num_curve_min ) 00141 num_curves += (rand() % (1 + num_curve_max - num_curve_min) ); 00142 for (int j = 0; j < num_curves; j++) 00143 { 00144 int goal_intervals = (1 + (rand() % goal_m1)) * (1 + (rand() % goal_m2)); 00145 if (grow_goal_by_i) 00146 goal_intervals += num_vars; 00147 MeshKit::IAVariable *v = ia_interface->create_variable( NULL, MeshKit::SOFT, goal_intervals); 00148 side2.push_back(v); 00149 } 00150 00151 // if we have two non-trivial opposite sides, then constrain them to be equal 00152 if (side1.size() && side2.size()) 00153 { 00154 ia_interface->constrain_sum_equal(side1, side2); 00155 } 00156 00157 // add a sum-even constraint 00158 if (i==0) 00159 { 00160 // printf("sum-even side: %d", i); 00161 assert( side2.size() ); 00162 ia_interface->constrain_sum_even(side2); 00163 } 00164 00165 // todo: try some hard-sets and non-trivial rhs 00166 } 00167 } 00168 00169 // sum-even constraints test problems 00170 /* 00171 // test problem 4, a simple sum-even constraint 00172 int num_surfaces = 12; // 12 00173 int num_curves_per_surface = 4; // 4 00174 int num_shared_curves = 1; // 2 00175 00176 int num_curves = 0; 00177 for (int i = 0; i < num_surfaces; ++i) 00178 { 00179 // gather the indices for the sum-even constraint 00180 int start_curve = num_curves - num_shared_curves; 00181 if (start_curve < 0) 00182 start_curve = 0; 00183 std::vector<int>curve_indices; 00184 if (debugging) 00185 printf("%d sum-even:",i); 00186 for (int j = 0; j < num_curves_per_surface; ++j) 00187 { 00188 curve_indices.push_back(start_curve+j); 00189 if (debugging) 00190 printf(" %d",start_curve+j); 00191 } 00192 num_curves = start_curve + num_curves_per_surface; 00193 const int rhs = 0; // test 0, -1 00194 constrain_sum_even(curve_indices,rhs); 00195 if (debugging) 00196 printf(" =%d\n",rhs); 00197 } 00198 // assign random goals to the curves 00199 for (int i = (int) I.size(); i < num_curves; ++i ) 00200 { 00201 double goal = 1 + ((double) (rand() % 59)) / 10.; // 1 to 6.9 00202 // force an odd sum for testing purposes 00203 //if (i==0) 00204 // goal += 1.; 00205 I.push_back(goal); 00206 } 00207 */ 00208 00209 void test_one_pair() 00210 { 00211 MeshKit::IAInterface *ia_interface = new_ia_interface(); 00212 ia_interface->destroy_data(); 00213 00214 std::vector< std::pair<int,int> > correct_solution; 00215 set_decoupled_pairs(ia_interface, 1, 1, 3, correct_solution); 00216 // set_decoupled_pairs(ia_interface, 1, 3.2, 12.1, correct_solution); 00217 ia_interface->execute_this(); 00218 bool solution_correct = check_solution_correctness( ia_interface, correct_solution ); 00219 CHECK( solution_correct ); 00220 } 00221 00222 void test_many_pairs() 00223 { 00224 MeshKit::IAInterface *ia_interface = new_ia_interface(); 00225 ia_interface->destroy_data(); 00226 00227 std::vector< std::pair<int,int> > correct_solution; 00228 set_decoupled_pairs(ia_interface, 8, 3.2, 12.1, correct_solution); 00229 set_decoupled_pairs(ia_interface, 1, 3.2, 12.1, correct_solution); 00230 set_decoupled_pairs(ia_interface, 8, 7.7, 4.2, correct_solution); 00231 set_decoupled_pairs(ia_interface, 40, 1.1, 5.2, correct_solution); 00232 set_decoupled_pairs(ia_interface, 40, 1.6, 4.5, correct_solution); 00233 set_decoupled_pairs(ia_interface, 4, 1.5, 1.5, correct_solution); 00234 set_decoupled_pairs(ia_interface, 4, 1, 1, correct_solution); 00235 00236 ia_interface->execute_this(); 00237 00238 bool solution_correct = check_solution_correctness( ia_interface, correct_solution ); 00239 CHECK( solution_correct ); 00240 00241 delete_ia_interface( ia_interface ); 00242 } 00243 00244 void test_long_chain() 00245 { 00246 MeshKit::IAInterface *ia_interface = new_ia_interface(); 00247 ia_interface->destroy_data(); 00248 00249 // test scalability: 20000 gives 20,000 constraints, 100,000 variables in 1 second relaxed solution 00250 set_mapping_chain(ia_interface, 16000, false, 3, 15, 2, 11); 00251 // goal distribution is gaussian in [1, 32] 00252 00253 ia_interface->execute_this(); 00254 00255 // bool solution_defined = check_solution( ia_interface ); 00256 00257 delete_ia_interface( ia_interface ); 00258 } 00259 00260 00261 void test_growing_chain() 00262 { 00263 // test problem 2 00264 // printf("constructing growing chain, coupled test problem\n"); 00265 MeshKit::IAInterface *ia_interface = new_ia_interface(); 00266 ia_interface->destroy_data(); 00267 00268 // goals are 1, 2, 3, 4, ... 16 00269 // one curve per side 00270 set_mapping_chain(ia_interface, 16, true, 1, 1, 1, 1); 00271 00272 ia_interface->execute_this(); 00273 00274 // bool solution_defined = check_solution( ia_interface ); 00275 00276 delete_ia_interface( ia_interface ); 00277 } 00278 00279 void mapping_test() 00280 { 00281 MeshKit::IAInterface *ia_interface = new_ia_interface(); 00282 ia_interface->destroy_data(); 00283 00284 std::string file_name = TestDir + "/" + INTASSIGN_TEST_FILE; 00285 printf("opening %s\n", file_name.c_str()); 00286 mk->load_geometry_mesh(file_name.c_str(), NULL); 00287 00288 //check the number of geometrical edges 00289 MeshKit::MEntVector surfs, loops; 00290 mk->get_entities_by_dimension(2, surfs); 00291 MeshKit::ModelEnt *this_surf = (*surfs.rbegin()); 00292 00293 // request a specific size 00294 //mk->sizing_function(0.1, true); 00295 MeshKit::MEntVector curves; 00296 std::vector<int> senses, loop_sizes; 00297 this_surf->boundary(1, curves, &senses, &loop_sizes); 00298 CHECK_EQUAL(4, (int)curves.size()); 00299 00300 // MeshKit::SizingFunction esize(mk, 3, 0.01); 00301 MeshKit::SizingFunction esize(mk, -1, -1); 00302 surfs[0]->sizing_function_index(esize.core_index()); 00303 00304 MeshKit::MEntVector side1, side2; 00305 side1.push_back(curves[0]); side2.push_back(curves[2]); 00306 ia_interface->constrain_sum_equal(ia_interface->make_constraint_group(side1), 00307 ia_interface->make_constraint_group(side2)); 00308 side1.clear(); side2.clear(); 00309 side1.push_back(curves[1]); side2.push_back(curves[3]); 00310 ia_interface->constrain_sum_equal(ia_interface->make_constraint_group(side1), 00311 ia_interface->make_constraint_group(side2)); 00312 00313 // if there are loops, and the loops have strictly less than 4 curves, then 00314 // ia_interface->constrain_sum_even( ia_interface->make_constraint_group(curves in loop) ); 00315 00316 //now, do the TFIMapping 00317 (MeshKit::TFIMapping*) mk->construct_meshop("TFIMapping", surfs); 00318 mk->setup_and_execute(); 00319 mk->save_mesh("intassign.exo"); 00320 // 00321 delete_ia_interface( ia_interface ); 00322 } 00323 00324 int main(int argv, char* argc[]) 00325 { 00326 // currently unable to create more than one mk called IntervalAssignment 00327 mk = new MeshKit::MKCore(); 00328 00329 int one_pair = RUN_TEST(test_one_pair); 00330 // run same test twice to check data clearing integrity 00331 int one_pair2 = RUN_TEST(test_one_pair); 00332 int many_pairs = RUN_TEST(test_many_pairs); 00333 int map_test = RUN_TEST(mapping_test); 00334 int growing_chain = RUN_TEST(test_growing_chain); 00335 // int long_chain = RUN_TEST(test_long_chain); 00336 00337 // int expt = RUN_TEST(test_exception); 00338 // int succ = RUN_TEST(test_success); 00339 00340 delete mk; 00341 int success = one_pair + one_pair2 + many_pairs + map_test + growing_chain; // + !abrt + !expt + succ; 00342 return success; 00343 }