From 126b25e6902b64042d81922a40977e8bb0f7ec6b Mon Sep 17 00:00:00 2001 From: logari81 Date: Sat, 14 Apr 2012 12:10:32 +0200 Subject: [PATCH] Sketcher: basic graph based system partitioning --- src/Mod/Sketcher/App/Sketch.cpp | 15 +-- src/Mod/Sketcher/App/Sketch.h | 3 +- src/Mod/Sketcher/App/freegcs/GCS.cpp | 142 ++++++++++++++------------- src/Mod/Sketcher/App/freegcs/GCS.h | 10 +- 4 files changed, 93 insertions(+), 77 deletions(-) diff --git a/src/Mod/Sketcher/App/Sketch.cpp b/src/Mod/Sketcher/App/Sketch.cpp index 97e26f1d8..608d9da0e 100644 --- a/src/Mod/Sketcher/App/Sketch.cpp +++ b/src/Mod/Sketcher/App/Sketch.cpp @@ -1563,13 +1563,14 @@ bool Sketch::updateGeometry() // solving ========================================================== -int Sketch::solve() +int Sketch::solve(void) { Base::TimeInfo start_time; if (!isInitMove) { // make sure we are in single subsystem mode GCSsys.clearByTag(-1); GCSsys.clearByTag(-2); + isFine = true; } int ret; @@ -1580,15 +1581,15 @@ int Sketch::solve() case 0: // solving with the default DogLeg solver // (or with SQP if we are in moving mode) solvername = isInitMove ? "SQP" : "DogLeg"; - ret = GCSsys.solve(true, GCS::DogLeg); + ret = GCSsys.solve(isFine, GCS::DogLeg); break; case 1: // solving with the LevenbergMarquardt solver solvername = "LevenbergMarquardt"; - ret = GCSsys.solve(true, GCS::LevenbergMarquardt); + ret = GCSsys.solve(isFine, GCS::LevenbergMarquardt); break; case 2: // solving with the BFGS solver solvername = "BFGS"; - ret = GCSsys.solve(true, GCS::BFGS); + ret = GCSsys.solve(isFine, GCS::BFGS); break; case 3: // last resort: augment the system with a second subsystem and use the SQP solver solvername = "SQP(augmented system)"; @@ -1601,7 +1602,7 @@ int Sketch::solve() GCSsys.addConstraintEqual(*it, &InitParameters[i], -2); } GCSsys.initSolution(Parameters); - ret = GCSsys.solve(true); + ret = GCSsys.solve(isFine); break; } @@ -1647,8 +1648,10 @@ int Sketch::solve() return ret; } -int Sketch::initMove(int geoId, PointPos pos) +int Sketch::initMove(int geoId, PointPos pos, bool fine) { + isFine = fine; + geoId = checkGeoId(geoId); GCSsys.clearByTag(-1); diff --git a/src/Mod/Sketcher/App/Sketch.h b/src/Mod/Sketcher/App/Sketch.h index 0d1543cb3..ec108421b 100644 --- a/src/Mod/Sketcher/App/Sketch.h +++ b/src/Mod/Sketcher/App/Sketch.h @@ -94,7 +94,7 @@ public: /** initializes a point (or curve) drag by setting the current * sketch status as a reference */ - int initMove(int geoId, PointPos pos); + int initMove(int geoId, PointPos pos, bool fine=true); /** move this point (or curve) to a new location and solve. * This will introduce some additional weak constraints expressing @@ -212,6 +212,7 @@ protected: std::vector Circles; bool isInitMove; + bool isFine; private: /// retrieves the index of a point diff --git a/src/Mod/Sketcher/App/freegcs/GCS.cpp b/src/Mod/Sketcher/App/freegcs/GCS.cpp index b0ee985fb..10b2492cd 100644 --- a/src/Mod/Sketcher/App/freegcs/GCS.cpp +++ b/src/Mod/Sketcher/App/freegcs/GCS.cpp @@ -27,9 +27,14 @@ #include "qp_eq.h" #include +#include +#include + namespace GCS { +typedef boost::adjacency_list Graph; + /////////////////////////////////////// // Solver /////////////////////////////////////// @@ -38,9 +43,7 @@ namespace GCS System::System() : clist(0), c2p(), p2c(), - subsys0(0), - subsys1(0), - subsys2(0), + subsyslist(0), reference(), init(false) { @@ -48,9 +51,7 @@ System::System() System::System(std::vector clist_) : c2p(), p2c(), - subsys0(0), - subsys1(0), - subsys2(0), + subsyslist(0), reference(), init(false) { @@ -553,6 +554,28 @@ void System::initSolution(VEC_pD ¶ms) // - Organizes the rest of constraints into two subsystems for // tag ids >=0 and < 0 respectively and applies the // system reduction specified in the previous step + MAP_pD_I params_index; + for (int i=0; i < int(params.size()); ++i) + params_index[params[i]] = i; + + Graph g; + for (int i=0; i < int(params.size() + clist.size()); i++) + boost::add_vertex(g); + + int cvtid = int(params.size()); + for (std::vector::const_iterator constr=clist.begin(); + constr != clist.end(); ++constr, cvtid++) { + VEC_pD &cparams = c2p[*constr]; + for (VEC_pD::const_iterator param=cparams.begin(); + param != cparams.end(); ++param) { + MAP_pD_I::const_iterator it = params_index.find(*param); + if (it != params_index.end()) + boost::add_edge(cvtid, it->second, g); + } + } + + VEC_I components(boost::num_vertices(g)); + int components_size = boost::connected_components(g, &components[0]); clearReference(); for (VEC_pD::const_iterator param=params.begin(); @@ -564,9 +587,6 @@ void System::initSolution(VEC_pD ¶ms) reductionmap.clear(); { VEC_pD reduced_params=params; - MAP_pD_I params_index; - for (int i=0; i < int(params.size()); ++i) - params_index[params[i]] = i; for (std::vector::const_iterator constr=clist.begin(); constr != clist.end(); ++constr) { @@ -589,27 +609,39 @@ void System::initSolution(VEC_pD ¶ms) reductionmap[params[i]] = reduced_params[i]; } - int i=0; - std::vector clist0, clist1, clist2; + std::vector< std::vector > clists0(components_size), + clists1(components_size), + clists2(components_size); + int i = int(params.size()); for (std::vector::const_iterator constr=clist.begin(); constr != clist.end(); ++constr, i++) { if (eliminated.count(*constr) == 0) { + int id = components[i]; if ((*constr)->getTag() >= 0) - clist0.push_back(*constr); + clists0[id].push_back(*constr); else if ((*constr)->getTag() == -1) // move constraints - clist1.push_back(*constr); + clists1[id].push_back(*constr); else // distance from reference constraints - clist2.push_back(*constr); + clists2[id].push_back(*constr); } } + std::vector< std::vector > plists(components_size); + for (int i=0; i < int(params.size()); ++i) { + int id = components[i]; + plists[id].push_back(params[i]); + } + clearSubSystems(); - if (clist0.size() > 0) - subsys0 = new SubSystem(clist0, params, reductionmap); - if (clist1.size() > 0) - subsys1 = new SubSystem(clist1, params, reductionmap); - if (clist2.size() > 0) - subsys2 = new SubSystem(clist2, params, reductionmap); + for (int cid=0; cid < components_size; cid++) { + subsyslist.push_back(std::vector(0)); + if (clists0[cid].size() > 0) + subsyslist[cid].push_back(new SubSystem(clists0[cid], plists[cid], reductionmap)); + if (clists1[cid].size() > 0) + subsyslist[cid].push_back(new SubSystem(clists1[cid], plists[cid], reductionmap)); + if (clists2[cid].size() > 0) + subsyslist[cid].push_back(new SubSystem(clists2[cid], plists[cid], reductionmap)); + } init = true; } @@ -634,31 +666,26 @@ int System::solve(VEC_pD ¶ms, bool isFine, Algorithm alg) int System::solve(bool isFine, Algorithm alg) { - if (subsys0) { - resetToReference(); - if (subsys2) { - int ret = solve(subsys0, subsys2, isFine); - if (subsys1) // give subsys1 higher priority than subsys2 - // in this case subsys2 acts like a preconditioner - return solve(subsys0, subsys1, isFine); - else - return ret; + bool isReset = false; + // return success by default in order to permit coincidence constraints to be applied + // even if no other system has to be solved + int res = Success; + for (int cid=0; cid < int(subsyslist.size()); cid++) { + if (subsyslist[cid].size() > 0 && !isReset) { + resetToReference(); + isReset = true; } - else if (subsys1) - return solve(subsys0, subsys1, isFine); - else - return solve(subsys0, isFine, alg); + if (subsyslist[cid].size() == 1) + res = std::max(res, solve(subsyslist[cid][0], isFine, alg)); + else if (subsyslist[cid].size() == 2) + res = std::max(res, solve(subsyslist[cid][0], subsyslist[cid][1], isFine)); + else if (subsyslist[cid].size() > 2) + // subsystem 1 has higher priority than subsystems 2,3,... + // these subsystems act like a preconditioner + for (int i=subsyslist[cid].size()-1; i > 0; i--) + res = std::max(res, solve(subsyslist[cid][0], subsyslist[cid][i], isFine)); } - else if (subsys1) { - resetToReference(); - if (subsys2) - return solve(subsys1, subsys2, isFine); - else - return solve(subsys1, isFine, alg); - } - else - // return success in order to permit coincidence constraints to be applied - return Success; + return res; } int System::solve(SubSystem *subsys, bool isFine, Algorithm alg) @@ -1193,25 +1220,11 @@ int System::solve(SubSystem *subsysA, SubSystem *subsysB, bool isFine) } -void System::getSubSystems(std::vector &subsysvec) -{ - subsysvec.clear(); - if (subsys0) - subsysvec.push_back(subsys0); - if (subsys1) - subsysvec.push_back(subsys1); - if (subsys2) - subsysvec.push_back(subsys2); -} - void System::applySolution() { - if (subsys2) - subsys2->applySolution(); - if (subsys1) - subsys1->applySolution(); - if (subsys0) - subsys0->applySolution(); + for (int cid=0; cid < int(subsyslist.size()); cid++) + for (int i=subsyslist[cid].size()-1; i >= 0; i--) + subsyslist[cid][i]->applySolution(); for (MAP_pD_pD::const_iterator it=reductionmap.begin(); it != reductionmap.end(); ++it) @@ -1303,12 +1316,9 @@ int System::diagnose(VEC_pD ¶ms, VEC_I &conflicting) void System::clearSubSystems() { init = false; - std::vector subsystems; - getSubSystems(subsystems); - free(subsystems); - subsys0 = NULL; - subsys1 = NULL; - subsys2 = NULL; + for (int i=0; i < int(subsyslist.size()); i++) + free(subsyslist[i]); + subsyslist.clear(); } double lineSearch(SubSystem *subsys, Eigen::VectorXd &xdir) diff --git a/src/Mod/Sketcher/App/freegcs/GCS.h b/src/Mod/Sketcher/App/freegcs/GCS.h index 844e9365e..f17f585c9 100644 --- a/src/Mod/Sketcher/App/freegcs/GCS.h +++ b/src/Mod/Sketcher/App/freegcs/GCS.h @@ -54,9 +54,12 @@ namespace GCS std::map c2p; // constraint to parameter adjacency list std::map > p2c; // parameter to constraint adjacency list - SubSystem *subsys0; // has the highest priority, always used as the primary subsystem - SubSystem *subsys1; // normally used as secondary subsystem, it is considered primary only if subsys0 is missing - SubSystem *subsys2; // has the lowest priority, always used as secondary system + // each row of subsyslist contains 3 subsystems. + // the first one has the highest priority, always used as the primary subsystem + // the second one is normally used as secondary subsystem, it is considered primary + // only if the first one is missing + // the rhird one has the lowest priority, always used as secondary system + std::vector< std::vector > subsyslist; void clearSubSystems(); MAP_pD_D reference; @@ -152,7 +155,6 @@ namespace GCS int solve(SubSystem *subsys, bool isFine=true, Algorithm alg=DogLeg); int solve(SubSystem *subsysA, SubSystem *subsysB, bool isFine=true); - void getSubSystems(std::vector &subsysvec); void applySolution(); void undoSolution();