~centralelyon2010/inkscape/imagelinks2

Viewing changes to src/libcola/conjugate_gradient.cpp

Committer: tgdwyer
Date: 2006-07-12 00:55:58 UTC
Revision ID: tgdwyer@users.sourceforge.net-20060712005558-4pqys3ou7f5er3dm

Previously graph layout was done using the Kamada-Kawai layout algorithm
implemented in Boost. I am replacing this with a custom implementation of
a constrained stress-majorization algorithm.

The stress-majorization algorithm is more robust and has better convergence
characteristics than Kamada-Kawai, and also simple constraints can be placed
on node position (for example, to enforce downward-pointing edges, non-overlap constraints, or cluster constraints).

Another big advantage is that we no longer need Boost.

I've tested the basic functionality, but I have yet to properly handle
disconnected graphs or to properly scale the resulting layout.

This commit also includes significant refactoring... the quadratic program solver - libvpsc (Variable Placement with Separation Constraints) has been moved to src/libvpsc and the actual graph layout algorithm is in libcola.

files added:
src/libcola

src/libcola/Makefile_insert

src/libcola/cola.cpp

src/libcola/cola.h

src/libcola/conjugate_gradient.cpp

src/libcola/conjugate_gradient.h

src/libcola/cycle_detector.cpp

src/libcola/cycle_detector.h

src/libcola/defs.h

src/libcola/gradient_projection.cpp

src/libcola/gradient_projection.h

src/libcola/shortest_paths.cpp

src/libcola/shortest_paths.h

src/libcola/straightener.cpp

src/libcola/straightener.h

src/libvpsc

src/libvpsc/COPYING

src/libvpsc/Makefile_insert

src/libvpsc/block.cpp

src/libvpsc/block.h

src/libvpsc/blocks.cpp

src/libvpsc/blocks.h

src/libvpsc/constraint.cpp

src/libvpsc/constraint.h

src/libvpsc/csolve_VPSC.cpp

src/libvpsc/csolve_VPSC.h

src/libvpsc/generate-constraints.cpp

src/libvpsc/generate-constraints.h

src/libvpsc/isnan.h

src/libvpsc/pairingheap

src/libvpsc/pairingheap/.dirstamp

src/libvpsc/pairingheap/PairingHeap.cpp

src/libvpsc/pairingheap/PairingHeap.h

src/libvpsc/pairingheap/dsexceptions.h

src/libvpsc/placement_SolveVPSC.h

src/libvpsc/remove_rectangle_overlap.cpp

src/libvpsc/remove_rectangle_overlap.h

src/libvpsc/solve_VPSC.cpp

src/libvpsc/solve_VPSC.h

src/libvpsc/variable.cpp

src/libvpsc/variable.h

files removed:
src/removeoverlap/block.cpp

src/removeoverlap/block.h

src/removeoverlap/blocks.cpp

src/removeoverlap/blocks.h

src/removeoverlap/constraint.cpp

src/removeoverlap/constraint.h

src/removeoverlap/generate-constraints.cpp

src/removeoverlap/generate-constraints.h

src/removeoverlap/pairingheap

src/removeoverlap/pairingheap/.cvsignore

src/removeoverlap/pairingheap/PairingHeap.cpp

src/removeoverlap/pairingheap/PairingHeap.h

src/removeoverlap/pairingheap/dsexceptions.h

src/removeoverlap/placement_SolveVPSC.cpp

src/removeoverlap/placement_SolveVPSC.h

src/removeoverlap/remove_rectangle_overlap-test.cpp

src/removeoverlap/remove_rectangle_overlap.cpp

src/removeoverlap/remove_rectangle_overlap.h

src/removeoverlap/solve_VPSC.cpp

src/removeoverlap/solve_VPSC.h

src/removeoverlap/variable.cpp

src/removeoverlap/variable.h

files modified:
configure.ac

mkinstalldirs

src/Makefile.am

src/Makefile_insert

src/graphlayout/graphlayout.cpp

src/removeoverlap/Makefile_insert

src/removeoverlap/removeoverlap.cpp

Show diffs side-by-side

added added

removed removed

src/libcola/conjugate_gradient.cpp

#include <math.h>

#include <stdlib.h>

#include <valarray>

#include <cassert>

#include "conjugate_gradient.h"

* Authors:

* Nathan Hurst <njh@njhurst.com>

* Tim Dwyer <tgdwyer@gmail.com>

* Released under GNU LGPL.

/* lifted wholely from wikipedia. Well, apart from the bug in the wikipedia version. */

using std::valarray;

static void

matrix_times_vector(valarray<double> const &matrix, /* m * n */

valarray<double> const &vec, /* n */

valarray<double> &result) /* m */

{

unsigned n = vec.size();

unsigned m = result.size();

assert(m*n == matrix.size());

const double* mp = &matrix[0];

for (unsigned i = 0; i < m; i++) {

double res = 0;

for (unsigned j = 0; j < n; j++)

res += *mp++ * vec[j];

result[i] = res;

}

static double Linfty(valarray<double> const &vec) {

return std::max(vec.max(), -vec.min());

}

double

inner(valarray<double> const &x,

valarray<double> const &y) {

double total = 0;

for(unsigned i = 0; i < x.size(); i++)

total += x[i]*y[i];

return total;// (x*y).sum(); <- this is more concise, but ineff

}

void

conjugate_gradient(double **A,

double *x,

double *b,

unsigned n,

double tol,

unsigned max_iterations) {

valarray<double> vA(n*n);

valarray<double> vx(n);

valarray<double> vb(n);

for(unsigned i=0;i<n;i++) {

vx[i]=x[i];

vb[i]=b[i];

for(unsigned j=0;j<n;j++) {

vA[i*n+j]=A[i][j];

}

conjugate_gradient(vA,vx,vb,n,tol,max_iterations);

for(unsigned i=0;i<n;i++) {

x[i]=vx[i];

}

void

conjugate_gradient(valarray<double> const &A,

valarray<double> &x,

valarray<double> const &b,

unsigned n, double tol,

unsigned max_iterations) {

valarray<double> Ap(n), p(n), r(n);

matrix_times_vector(A,x,Ap);

r=b-Ap;

double r_r = inner(r,r);

unsigned k = 0;

tol *= tol;

while(k < max_iterations && r_r > tol) {

k++;

double r_r_new = r_r;

if(k == 1)

p = r;

else {

r_r_new = inner(r,r);

p = r + (r_r_new/r_r)*p;

}

matrix_times_vector(A, p, Ap);

double alpha_k = r_r_new / inner(p, Ap);

x += alpha_k*p;

r -= alpha_k*Ap;

r_r = r_r_new;

}

100

printf("njh: %d iters, Linfty = %g L2 = %g\n", k,

101

std::max(-r.min(), r.max()), sqrt(r_r));

102

// x is solution

103

}

104

105

Local Variables:

106

mode:c++

107

c-file-style:"stroustrup"

108

c-file-offsets:((innamespace . 0)(inline-open . 0)(case-label . +))

109

indent-tabs-mode:nil

110

fill-column:99

111

End:

112

113

// vim: filetype=cpp:expandtab:shiftwidth=4:tabstop=4:softtabstop=4

Older »