src/libcola/cola.cpp

   1 #include "cola.h"
   2 #include "conjugate_gradient.h"
   3 #include "straightener.h"
   4
   5 namespace cola {
   6
   7 /**
   8  * Find the euclidean distance between a pair of dummy variables
   9  */
  10 inline double dummy_var_euclidean_dist(GradientProjection* gpx, GradientProjection* gpy, unsigned i) {
  11     double dx = gpx->dummy_vars[i]->place_r - gpx->dummy_vars[i]->place_l,
  12         dy = gpy->dummy_vars[i]->place_r - gpy->dummy_vars[i]->place_l;
  13     return sqrt(dx*dx + dy*dy);
  14 }
  15
  16 void
  17 ConstrainedMajorizationLayout
  18 ::setupDummyVars() {
  19     if(clusters==NULL) return;
  20     double* coords[2]={X,Y};
  21     GradientProjection* gp[2]={gpX,gpY};
  22     for(unsigned k=0;k<2;k++) {
  23         gp[k]->clearDummyVars();
  24         if(clusters) {
  25             for(Clusters::iterator cit=clusters->begin();
  26                     cit!=clusters->end(); ++cit) {
  27                 Cluster *c = *cit;
  28                 DummyVarPair* p = new DummyVarPair(edge_length);
  29                 gp[k]->dummy_vars.push_back(p);
  30                 double minPos=DBL_MAX, maxPos=-DBL_MAX;
  31                 for(Cluster::iterator vit=c->begin();
  32                         vit!=c->end(); ++vit) {
  33                     double pos = coords[k][*vit];
  34                     minPos=min(pos,minPos);
  35                     maxPos=max(pos,maxPos);
  36                     p->leftof.push_back(make_pair(*vit,0));
  37                     p->rightof.push_back(make_pair(*vit,0));
  38                 }
  39                 p->place_l = minPos;
  40                 p->place_r = maxPos;
  41             }
  42         }
  43     }
  44     for(unsigned k=0;k<2;k++) {
  45         unsigned n_d = gp[k]->dummy_vars.size();
  46         if(n_d > 0) {
  47             for(unsigned i=0; i<n_d; i++) {
  48                 gp[k]->dummy_vars[i]->computeLinearTerm(dummy_var_euclidean_dist(gpX, gpY, i));
  49             }
  50         }
  51     }
  52 }
  53 void ConstrainedMajorizationLayout::majlayout(
  54         double** Dij, GradientProjection* gp, double* coords)
  55 {
  56     double b[n];
  57     fill(b,b+n,0);
  58     majlayout(Dij,gp,coords,b);
  59 }
  60 void ConstrainedMajorizationLayout::majlayout(
  61         double** Dij, GradientProjection* gp, double* coords, double* b)
  62 {
  63     double L_ij,dist_ij,degree;
  64     /* compute the vector b */
  65     /* multiply on-the-fly with distance-based laplacian */
  66     for (unsigned i = 0; i < n; i++) {
  67         degree = 0;
  68         if(i<lapSize) {
  69             for (unsigned j = 0; j < lapSize; j++) {
  70                 if (j == i) continue;
  71                 dist_ij = euclidean_distance(i, j);
  72                 if (dist_ij > 1e-30 && Dij[i][j] > 1e-30) {     /* skip zero distances */
  73                     /* calculate L_ij := w_{ij}*d_{ij}/dist_{ij} */
  74                     L_ij = 1.0 / (dist_ij * Dij[i][j]);
  75                     degree -= L_ij;
  76                     b[i] += L_ij * coords[j];
  77                 }
  78             }
  79             b[i] += degree * coords[i];
  80         }
  81         assert(!isnan(b[i]));
  82     }
  83     if(constrainedLayout) {
  84         setupDummyVars();
  85         gp->solve(b);
  86     } else {
  87         conjugate_gradient(lap2, coords, b, n, tol, n);
  88     }
  89     moveBoundingBoxes();
  90 }
  91 inline double ConstrainedMajorizationLayout
  92 ::compute_stress(double **Dij) {
  93     double sum = 0, d, diff;
  94     for (unsigned i = 1; i < lapSize; i++) {
  95         for (unsigned j = 0; j < i; j++) {
  96             d = Dij[i][j];
  97             diff = d - euclidean_distance(i,j);
  98             sum += diff*diff / (d*d);
  99         }
 100     }
 101     if(clusters!=NULL) {
 102         for(unsigned i=0; i<gpX->dummy_vars.size(); i++) {
 103             sum += gpX->dummy_vars[i]->stress(dummy_var_euclidean_dist(gpX, gpY, i));
 104         }
 105     }
 106     return sum;
 107 }
 108 /*
 109 void ConstrainedMajorizationLayout
 110 ::addLinearConstraints(LinearConstraints* linearConstraints) {
 111     n=lapSize+linearConstraints->size();
 112     Q=new double*[n];
 113     X=new double[n];
 114     Y=new double[n];
 115     for(unsigned i = 0; i<n; i++) {
 116         X[i]=rs[i]->getCentreX();
 117         Y[i]=rs[i]->getCentreY();
 118         Q[i]=new double[n];
 119         for(unsigned j=0; j<n; j++) {
 120             if(i<lapSize&&j<lapSize) {
 121                 Q[i][j]=lap2[i][j];
 122             } else {
 123                 Q[i][j]=0;
 124             }
 125         }
 126     }
 127     for(LinearConstraints::iterator i=linearConstraints->begin();
 128            i!= linearConstraints->end();i++) {
 129         LinearConstraint* c=*i;
 130         Q[c->u][c->u]+=c->w*c->duu;
 131         Q[c->u][c->v]+=c->w*c->duv;
 132         Q[c->u][c->b]+=c->w*c->dub;
 133         Q[c->v][c->u]+=c->w*c->duv;
 134         Q[c->v][c->v]+=c->w*c->dvv;
 135         Q[c->v][c->b]+=c->w*c->dvb;
 136         Q[c->b][c->b]+=c->w*c->dbb;
 137         Q[c->b][c->u]+=c->w*c->dub;
 138         Q[c->b][c->v]+=c->w*c->dvb;
 139     }
 140 }
 141 */
 142
 143 bool ConstrainedMajorizationLayout::run() {
 144     /*
 145     for(unsigned i=0;i<n;i++) {
 146         for(unsigned j=0;j<n;j++) {
 147             cout << lap2[i][j] << " ";
 148         }
 149         cout << endl;
 150     }
 151     */
 152     do {
 153         /* Axis-by-axis optimization: */
 154         if(straightenEdges) {
 155             straighten(*straightenEdges,HORIZONTAL);
 156             straighten(*straightenEdges,VERTICAL);
 157         } else {
 158             majlayout(Dij,gpX,X);
 159             majlayout(Dij,gpY,Y);
 160         }
 161     } while(!done(compute_stress(Dij),X,Y));
 162     return true;
 163 }
 164 static bool straightenToProjection=true;
 165 void ConstrainedMajorizationLayout::straighten(vector<straightener::Edge*>& sedges, Dim dim) {
 166         vector<straightener::Node*> snodes;
 167         for (unsigned i=0;i<lapSize;i++) {
 168                 snodes.push_back(new straightener::Node(i,boundingBoxes[i]));
 169         }
 170     SimpleConstraints cs;
 171     straightener::generateConstraints(snodes,sedges,cs,dim);
 172     n=snodes.size();
 173     Q=new double*[n];
 174     delete [] X;
 175     delete [] Y;
 176     X=new double[n];
 177     Y=new double[n];
 178     for(unsigned i = 0; i<n; i++) {
 179         X[i]=snodes[i]->x;
 180         Y[i]=snodes[i]->y;
 181         Q[i]=new double[n];
 182         for(unsigned j=0; j<n; j++) {
 183             if(i<lapSize&&j<lapSize) {
 184                 Q[i][j]=lap2[i][j];
 185             } else {
 186                 Q[i][j]=0;
 187             }
 188         }
 189     }
 190     LinearConstraints linearConstraints;
 191     for(unsigned i=0;i<sedges.size();i++) {
 192         sedges[i]->nodePath(snodes);
 193         vector<unsigned>& path=sedges[i]->path;
 194         // take u and v as the ends of the line
 195         //unsigned u=path[0];
 196         //unsigned v=path[path.size()-1];
 197         double total_length=0;
 198         for(unsigned j=1;j<path.size();j++) {
 199             unsigned u=path[j-1], v=path[j];
 200             total_length+=euclidean_distance(u,v);
 201         }
 202         for(unsigned j=1;j<path.size()-1;j++) {
 203             // find new u and v for each line segment
 204             unsigned u=path[j-1];
 205             unsigned b=path[j];
 206             unsigned v=path[j+1];
 207             double weight=-0.01;
 208             double wub=euclidean_distance(u,b)/total_length;
 209             double wbv=euclidean_distance(b,v)/total_length;
 210             linearConstraints.push_back(new cola::LinearConstraint(u,v,b,weight,wub,wbv,X,Y));
 211         }
 212     }
 213     //cout << "Generated "<<linearConstraints.size()<< " linear constraints"<<endl;
 214     assert(snodes.size()==lapSize+linearConstraints.size());
 215     double b[n],*coords=dim==HORIZONTAL?X:Y,dist_ub,dist_bv;
 216     fill(b,b+n,0);
 217     for(LinearConstraints::iterator i=linearConstraints.begin();
 218            i!= linearConstraints.end();i++) {
 219         LinearConstraint* c=*i;
 220         if(straightenToProjection) {
 221             Q[c->u][c->u]+=c->w*c->duu;
 222             Q[c->u][c->v]+=c->w*c->duv;
 223             Q[c->u][c->b]+=c->w*c->dub;
 224             Q[c->v][c->u]+=c->w*c->duv;
 225             Q[c->v][c->v]+=c->w*c->dvv;
 226             Q[c->v][c->b]+=c->w*c->dvb;
 227             Q[c->b][c->b]+=c->w*c->dbb;
 228             Q[c->b][c->u]+=c->w*c->dub;
 229             Q[c->b][c->v]+=c->w*c->dvb;
 230         } else {
 231             double wub=edge_length*c->frac_ub;
 232             double wbv=edge_length*c->frac_bv;
 233             dist_ub=euclidean_distance(c->u,c->b)*wub;
 234             dist_bv=euclidean_distance(c->b,c->v)*wbv;
 235             wub=max(wub,0.00001);
 236             wbv=max(wbv,0.00001);
 237             dist_ub=max(dist_ub,0.00001);
 238             dist_bv=max(dist_bv,0.00001);
 239             wub=1/(wub*wub);
 240             wbv=1/(wbv*wbv);
 241             Q[c->u][c->u]-=wub;
 242             Q[c->u][c->b]+=wub;
 243             Q[c->v][c->v]-=wbv;
 244             Q[c->v][c->b]+=wbv;
 245             Q[c->b][c->b]-=wbv+wub;
 246             Q[c->b][c->u]+=wub;
 247             Q[c->b][c->v]+=wbv;
 248
 249             b[c->u]+=(coords[c->b]-coords[c->u]) / dist_ub;
 250             b[c->v]+=(coords[c->b]-coords[c->v]) / dist_bv;
 251             b[c->b]+=coords[c->u] / dist_ub + coords[c->v] / dist_bv
 252                    - coords[c->b] / dist_ub - coords[c->b] / dist_bv;
 253         }
 254     }
 255         GradientProjection gp(dim,n,Q,coords,tol,100,
 256             (AlignmentConstraints*)NULL,false,(Rectangle**)NULL,(PageBoundaryConstraints*)NULL,&cs);
 257     constrainedLayout = true;
 258     majlayout(Dij,&gp,coords,b);
 259     for(unsigned i=0;i<sedges.size();i++) {
 260         sedges[i]->createRouteFromPath(X,Y);
 261         sedges[i]->dummyNodes.clear();
 262         sedges[i]->path.clear();
 263     }
 264     for(unsigned i=0;i<cs.size();i++) {
 265         delete cs[i];
 266     }
 267     for(unsigned i=0;i<linearConstraints.size();i++) {
 268         delete linearConstraints[i];
 269     }
 270     for(unsigned i=0;i<snodes.size();i++) {
 271         delete snodes[i];
 272     }
 273     for(unsigned i = 0; i<n; i++) {
 274         delete [] Q[i];
 275     }
 276     delete [] Q;
 277     snodes.resize(lapSize);
 278 }
 279
 280 void ConstrainedMajorizationLayout::setupConstraints(
 281         AlignmentConstraints* acsx, AlignmentConstraints* acsy,
 282         bool avoidOverlaps,
 283         PageBoundaryConstraints* pbcx, PageBoundaryConstraints* pbcy,
 284         SimpleConstraints* scx, SimpleConstraints* scy,
 285         Clusters* cs,
 286         vector<straightener::Edge*>* straightenEdges) {
 287     constrainedLayout = true;
 288     this->avoidOverlaps = avoidOverlaps;
 289     if(cs) {
 290         clusters=cs;
 291     }
 292         gpX=new GradientProjection(
 293             HORIZONTAL,n,Q,X,tol,100,acsx,avoidOverlaps,boundingBoxes,pbcx,scx);
 294         gpY=new GradientProjection(
 295             VERTICAL,n,Q,Y,tol,100,acsy,avoidOverlaps,boundingBoxes,pbcy,scy);
 296     this->straightenEdges = straightenEdges;
 297 }
 298 } // namespace cola
 299 // vim: filetype=cpp:expandtab:shiftwidth=4:tabstop=4:softtabstop=4