#ifndef _GRADIENT_PROJECTION_H #define _GRADIENT_PROJECTION_H #include <libvpsc/solve_VPSC.h> #include <libvpsc/variable.h> #include <libvpsc/constraint.h> #include <libvpsc/generate-constraints.h> #include <vector> #include <iostream> #include <math.h> using namespace std; typedef vector<vpsc::Constraint*> Constraints; typedef vector<vpsc::Variable*> Variables; typedef vector<pair<unsigned,double> > OffsetList; 00018 class SimpleConstraint { public: SimpleConstraint(unsigned l, unsigned r, double g) : left(l), right(r), gap(g) {} unsigned left; unsigned right; double gap; }; typedef vector<SimpleConstraint*> SimpleConstraints; 00027 class AlignmentConstraint { friend class GradientProjection; public: AlignmentConstraint(double pos) : position(pos), variable(NULL) {} void updatePosition() { position = variable->position(); } OffsetList offsets; void* guide; double position; private: vpsc::Variable* variable; }; typedef vector<AlignmentConstraint*> AlignmentConstraints; 00042 class PageBoundaryConstraints { public: PageBoundaryConstraints(double lm, double rm, double w) : leftMargin(lm), rightMargin(rm), weight(w) { } void createVarsAndConstraints(Variables &vs, Constraints &cs) { vpsc::Variable* vl, * vr; // create 2 dummy vars, based on the dimension we are in vs.push_back(vl=new vpsc::Variable(vs.size(), leftMargin, weight)); vs.push_back(vr=new vpsc::Variable(vs.size(), rightMargin, weight)); // for each of the "real" variables, create a constraint that puts that var // between our two new dummy vars, depending on the dimension. for(OffsetList::iterator o=offsets.begin(); o!=offsets.end(); ++o) { cs.push_back(new vpsc::Constraint(vl, vs[o->first], o->second)); cs.push_back(new vpsc::Constraint(vs[o->first], vr, o->second)); } } OffsetList offsets; private: double leftMargin; double rightMargin; double weight; }; typedef vector<pair<unsigned,double> > CList; /** * A DummyVarPair is a pair of variables with an ideal distance between them and which have no * other interaction with other variables apart from through constraints. This means that * the entries in the laplacian matrix for dummy vars and other vars would be 0 - thus, sparse * matrix techniques can be used in laplacian operations. * The constraints are specified by a two lists of pairs of variable indexes and required separation. * The two lists are: * leftof: variables to which left must be to the left of, * rightof: variables to which right must be to the right of. */ 00077 class DummyVarPair { public: DummyVarPair(double desiredDist) : dist(desiredDist), lap2(1.0/(desiredDist*desiredDist)) { } CList leftof; // variables to which left dummy var must be to the left of CList rightof; // variables to which right dummy var must be to the right of double place_l; double place_r; void computeLinearTerm(double euclideanDistance) { if(euclideanDistance > 1e-30) { b = place_r - place_l; b /= euclideanDistance * dist; } else { b=0; } } double stress(double euclideanDistance) { double diff = dist - euclideanDistance; return diff*diff / (dist*dist); } private: friend class GradientProjection; /** * Setup vars and constraints for an instance of the VPSC problem. * Adds generated vars and constraints to the argument vectors. */ 00100 void setupVPSC(Variables &vars, Constraints &cs) { double weight=1; left = new vpsc::Variable(vars.size(),place_l,weight); vars.push_back(left); right = new vpsc::Variable(vars.size(),place_r,weight); vars.push_back(right); for(CList::iterator cit=leftof.begin(); cit!=leftof.end(); ++cit) { vpsc::Variable* v = vars[(*cit).first]; cs.push_back(new vpsc::Constraint(left,v,(*cit).second)); } for(CList::iterator cit=rightof.begin(); cit!=rightof.end(); ++cit) { vpsc::Variable* v = vars[(*cit).first]; cs.push_back(new vpsc::Constraint(v,right,(*cit).second)); } } /** * Extract the result of a VPSC solution to the variable positions */ 00120 void updatePosition() { place_l=left->position(); place_r=right->position(); } /** * Compute the descent vector, also copying the current position to old_place for * future reference */ 00128 void computeDescentVector() { old_place_l=place_l; old_place_r=place_r; g = 2.0 * ( b + lap2 * ( place_l - place_r ) ); } /** * move in the direction of steepest descent (based on g computed by computeGradient) * alpha is the step size. */ 00137 void steepestDescent(double alpha) { place_l -= alpha*g; place_r += alpha*g; left->desiredPosition=place_l; right->desiredPosition=place_r; } /** * add dummy vars' contribution to numerator and denominator for * beta step size calculation */ 00147 void betaCalc(double &numerator, double &denominator) { double dl = place_l-old_place_l, dr = place_r-old_place_r, r = 2.0 * lap2 * ( dr - dl ); numerator += g * ( dl - dr ); denominator += r*dl - r * dr; } /** * move by stepsize beta from old_place to place */ 00157 void feasibleDescent(double beta) { left->desiredPosition = place_l = old_place_l + beta*(place_l - old_place_l); right->desiredPosition = place_r = old_place_r + beta*(place_r - old_place_r); } double absoluteDisplacement() { return fabs(place_l - old_place_l) + fabs(place_r - old_place_r); } double dist; // ideal distance between vars double b; // linear coefficient in quad form for left (b_right = -b) vpsc::Variable* left; // Variables used in constraints vpsc::Variable* right; double lap2; // laplacian entry double g; // descent vec for quad form for left (g_right = -g) double old_place_l; // old_place is where the descent vec g was computed double old_place_r; }; typedef vector<DummyVarPair*> DummyVars; enum Dim { HORIZONTAL, VERTICAL }; 00177 class GradientProjection { public: GradientProjection( const Dim k, unsigned n, double** A, double* x, double tol, unsigned max_iterations, AlignmentConstraints* acs=NULL, bool nonOverlapConstraints=false, vpsc::Rectangle** rs=NULL, PageBoundaryConstraints *pbc = NULL, SimpleConstraints *sc = NULL) : k(k), n(n), A(A), place(x), rs(rs), nonOverlapConstraints(nonOverlapConstraints), tolerance(tol), acs(acs), max_iterations(max_iterations), g(new double[n]), d(new double[n]), old_place(new double[n]), constrained(false) { for(unsigned i=0;i<n;i++) { vars.push_back(new vpsc::Variable(i,1,1)); } if(acs) { for(AlignmentConstraints::iterator iac=acs->begin(); iac!=acs->end();++iac) { AlignmentConstraint* ac=*iac; vpsc::Variable *v=ac->variable=new vpsc::Variable(vars.size(),ac->position,0.0001); vars.push_back(v); for(OffsetList::iterator o=ac->offsets.begin(); o!=ac->offsets.end(); o++) { gcs.push_back(new vpsc::Constraint(v,vars[o->first],o->second,true)); } } } if (pbc) { pbc->createVarsAndConstraints(vars,gcs); } if (sc) { for(SimpleConstraints::iterator c=sc->begin(); c!=sc->end();++c) { gcs.push_back(new vpsc::Constraint( vars[(*c)->left],vars[(*c)->right],(*c)->gap)); } } if(!gcs.empty() || nonOverlapConstraints) { constrained=true; } } virtual ~GradientProjection() { delete [] g; delete [] d; delete [] old_place; for(Constraints::iterator i(gcs.begin()); i!=gcs.end(); i++) { delete *i; } gcs.clear(); for(unsigned i=0;i<vars.size();i++) { delete vars[i]; } } void clearDummyVars(); unsigned solve(double* b); DummyVars dummy_vars; // special vars that must be considered in Lapl. private: vpsc::IncSolver* setupVPSC(); void destroyVPSC(vpsc::IncSolver *vpsc); Dim k; unsigned n; // number of actual vars double** A; // Graph laplacian matrix double* place; Variables vars; // all variables // computations Constraints gcs; /* global constraints - persist throughout all iterations */ Constraints lcs; /* local constraints - only for current iteration */ vpsc::Rectangle** rs; bool nonOverlapConstraints; double tolerance; AlignmentConstraints* acs; unsigned max_iterations; double* g; /* gradient */ double* d; double* old_place; bool constrained; }; #endif /* _GRADIENT_PROJECTION_H */ // vim: filetype=cpp:expandtab:shiftwidth=4:tabstop=4:softtabstop=4 :