From d02e66cd8528bc0f98d1fe2439ff5b2ed25075c2 Mon Sep 17 00:00:00 2001 From: fred Date: Tue, 26 Mar 2002 22:48:17 +0000 Subject: [PATCH] lilypond-1.2.10 --- flower/choleski.cc | 209 -------------- flower/full-storage.cc | 168 ----------- flower/include/choleski.hh | 52 ---- flower/include/full-storage.hh | 171 ----------- flower/include/matrix-storage.hh | 153 ---------- flower/include/matrix.hh | 159 ---------- flower/include/vector.hh | 117 -------- flower/matrix-debug.cc | 64 ----- flower/matrix.cc | 299 ------------------- flower/vector.cc | 35 --- input/bugs/auto-melismas.ly | 17 -- input/bugs/hara.ly | 88 ------ input/bugs/l.ly | 13 - input/bugs/utremi_old.ly | 262 ----------------- lily/col-info.cc | 70 ----- lily/include/ineq-constrained-qp.hh | 91 ------ lily/include/qlp.hh | 44 --- lily/include/qlpsolve.hh | 95 ------ lily/include/spring-spacer.hh | 99 ------- lily/ineq-constrained-qp.cc | 294 ------------------- lily/mixed-qp.cc | 82 ------ lily/p-col.cc | 144 ---------- lily/qlpsolve.cc | 169 ----------- lily/spring-spacer.cc | 432 ---------------------------- 24 files changed, 3327 deletions(-) delete mode 100644 flower/choleski.cc delete mode 100644 flower/full-storage.cc delete mode 100644 flower/include/choleski.hh delete mode 100644 flower/include/full-storage.hh delete mode 100644 flower/include/matrix-storage.hh delete mode 100644 flower/include/matrix.hh delete mode 100644 flower/include/vector.hh delete mode 100644 flower/matrix-debug.cc delete mode 100644 flower/matrix.cc delete mode 100644 flower/vector.cc delete mode 100644 input/bugs/auto-melismas.ly delete mode 100644 input/bugs/hara.ly delete mode 100644 input/bugs/l.ly delete mode 100644 input/bugs/utremi_old.ly delete mode 100644 lily/col-info.cc delete mode 100644 lily/include/ineq-constrained-qp.hh delete mode 100644 lily/include/qlp.hh delete mode 100644 lily/include/qlpsolve.hh delete mode 100644 lily/include/spring-spacer.hh delete mode 100644 lily/ineq-constrained-qp.cc delete mode 100644 lily/mixed-qp.cc delete mode 100644 lily/p-col.cc delete mode 100644 lily/qlpsolve.cc delete mode 100644 lily/spring-spacer.cc diff --git a/flower/choleski.cc b/flower/choleski.cc deleted file mode 100644 index ae40c45613..0000000000 --- a/flower/choleski.cc +++ /dev/null @@ -1,209 +0,0 @@ -/* - choleski.cc -- implement Choleski_decomposition - - source file of the Flower Library - - (c) 1997--1999 Han-Wen Nienhuys -*/ - -#include "choleski.hh" -const Real EPS = 1e-7; // so sue me. Hard coded - -// for testing new Matrix_storage. -//#define PARANOID - -void -Choleski_decomposition::full_matrix_solve (Vector &out, Vector const &rhs) const -{ - int n= rhs.dim(); - assert (n == L.dim()); - Vector y; - y.set_dim (n); - out.set_dim (n); - - // forward substitution - for (int i=0; i < n; i++) - { - Real sum (0.0); - for (int j=0; j < i; j++) - sum += y (j) * L(i,j); - y (i) = (rhs (i) - sum)/L(i,i); - } - - for (int i=0; i < n; i++) - y (i) /= D(i); - - // backward subst - Vector &x (out); // using input as return val. - for (int i=n-1; i >= 0; i--) - { - Real sum (0.0); - for (int j=i+1; j < n; j++) - sum += L(j,i)*x (j); - x (i) = (y (i) - sum)/L(i,i); - } -} - -void -Choleski_decomposition::band_matrix_solve (Vector &out, Vector const &rhs) const -{ - int n= rhs.dim(); - int b = L.band_i(); - assert (n == L.dim()); - - out.set_dim (n); - - Vector y; - y.set_dim (n); - - // forward substitution - for (int i=0; i < n; i++) - { - Real sum (0.0); - for (int j= 0 >? i - b; j < i; j++) - sum += y (j) * L(i,j); - y (i) = (rhs (i) - sum)/L(i,i); - } - for (int i=0; i < n; i++) - y (i) /= D(i); - - // backward subst - Vector &x (out); // using input as return val. - for (int i=n-1; i >= 0; i--) - { - Real sum (0.0); - for (int j=i+1; j <= i + b&&j < n ; j++) - sum += L(j,i)*x (j); - x (i) = (y (i) - sum)/L(i,i); - } -} - -void -Choleski_decomposition::solve (Vector &x, Vector const &rhs) const -{ - if (band_b_) - { - band_matrix_solve (x,rhs); - } - else - full_matrix_solve (x,rhs); -} - -Vector -Choleski_decomposition::solve (Vector rhs) const -{ - Vector r; - solve (r, rhs); - return r; -} - -void -Choleski_decomposition::full_matrix_decompose (Matrix const & P) -{ - - int n = P.dim(); - L.unit(); - for (int k= 0; k < n; k++) - { - for (int j = 0; j < k; j++) - { - Real sum (0.0); - for (int l=0; l < j; l++) - sum += L(k,l)*L(j,l)*D(l); - L(k,j) = (P(k,j) - sum)/D(j); - } - Real sum=0.0; - - for (int l=0; l < k; l++) - sum += sqr (L(k,l))*D(l); - Real d = P(k,k) - sum; - D(k) = d; - } - -} - -void -Choleski_decomposition::band_matrix_decompose (Matrix const &P) -{ - int n = P.dim(); - int b = P.band_i(); - L.unit(); - - for (int i= 0; i < n; i++) - { - for (int j = 0 >? i - b; j < i; j++) - { - Real sum (0.0); - for (int l=0 >? i - b; l < j; l++) - sum += L(i,l)*L(j,l)*D(l); - L(i,j) = (P(i,j) - sum)/D(j); - } - Real sum=0.0; - - for (int l=0 >? i - b; l < i; l++) - sum += sqr (L(i,l))*D(l); - Real d = P(i,i) - sum; - D(i) = d; - } - L.set_band(); - band_b_ = true; -} - - - - -/* - Standard matrix algorithm. - */ - -Choleski_decomposition::Choleski_decomposition (Matrix const & P) - : L(P.dim()), D(P.dim ()) -{ -#ifdef PARANOID - assert ((P-P.transposed()).norm ()/P.norm () < EPS); -#endif - band_b_ = false; - - int b = P.calc_band_i (); - - if (b <= P.dim ()/2) - band_matrix_decompose (P); - else - full_matrix_decompose (P); - -#ifdef PARANOID - assert ((original()-P).norm () / P.norm () < EPS); -#endif -} - -Matrix -Choleski_decomposition::original() const -{ - Matrix T(L.dim()); - T.set_diag (D); - return L*T*L.transposed(); -} - -Matrix -Choleski_decomposition::inverse() const -{ - int n=L.dim(); - Matrix invm (n); - Vector e_i (n); - Vector inv (n); - for (int i = 0; i < n; i++) - { - e_i.set_unit (i); - solve (inv, e_i); - for (int j = 0 ; j -*/ - -#include "full-storage.hh" - - -void -Full_storage::operator=(Full_storage const &fs) -{ - resize (fs.height_i_, fs.width_i_); - OK(); - fs.OK(); - for (int i=0; i= height_i_ && max_width_i_ >= width_i_); - assert (height_i_ >= 0 && width_i_ >= 0); - assert (els_p_p_||!max_height_i_); -#endif -} - - - - -Full_storage::~Full_storage() -{ - for (int i=0; i < max_height_i_; i++) - delete [] els_p_p_[i]; - delete[] els_p_p_; -} - -void -Full_storage::resize (int rows, int cols) -{ - OK(); - resize_cols (rows); - resize_rows (cols); - band_i_ = rows >? cols; -} - -void -Full_storage::delete_column (int k) -{ - assert (0 <= k &&k k ; i--) - for (int j=0; j -Full_storage::row (int n) const -{ - Array r; - for (int j = 0; j < cols(); j++) - r.push (elem (n,j)); - return r; -} - -Array -Full_storage::column (int n) const -{ - Array r; - for (int i = 0; i < rows(); i++) - r.push (elem (i,n)); - return r; -} - -void -Full_storage::set_size (int rows, int cols) -{ - resize (rows,cols); -} - -void -Full_storage::set_size (int rows) -{ - - resize (rows); -} - - - -void -Full_storage::resize_cols (int newh) -{ - if (newh <= max_height_i_) - { - height_i_=newh; - return; - } - - Real ** newa=new Real*[newh]; - int j=0; - for (; j < height_i_; j++) - newa[j] = els_p_p_[j]; - for (; j < newh; j++) - newa[j] = new Real[max_width_i_]; - delete[] els_p_p_; - els_p_p_=newa; - - height_i_ = max_height_i_ = newh; -} - - - -void -Full_storage::resize_rows (int neww) -{ - if (neww <= max_width_i_) - { - width_i_=neww; - return; - } - for (int i=0; i < max_height_i_ ; i++) - { - Real* newa = new Real[neww]; - for (int k=0; k < width_i_; k++) - newa[k] = els_p_p_[i][k]; - - delete[] els_p_p_[i]; - els_p_p_[i] = newa; - } - width_i_ = max_width_i_ = neww; -} - -#ifdef INLINE -#undef INLINE -#endif -#define INLINE - - -#include "full-storage.icc" diff --git a/flower/include/choleski.hh b/flower/include/choleski.hh deleted file mode 100644 index c6db243ba8..0000000000 --- a/flower/include/choleski.hh +++ /dev/null @@ -1,52 +0,0 @@ -#ifndef CHOLESKI_HH -#define CHOLESKI_HH - -#include "matrix.hh" - -/** - Choleski decomposition of a matrix - structure for using the LU decomposition of a positive definite matrix. - - #P# is split into - - LD transpose (L) - */ -struct Choleski_decomposition { - - /// lower triangle of Choleski decomposition - Matrix L; - - bool band_b_; - - /// diagonal - Vector D; - - /** Create decomposition of P. - PRE - P needs to be symmetric positive definite - */ - - Choleski_decomposition (Matrix const &P); - - /** - solve Px = rhs - */ - Vector solve (Vector rhs) const; - void solve (Vector &dest, Vector const &rhs) const; - Vector operator * (Vector rhs) const { return solve (rhs); } - /** - return the inverse of the matrix P. - */ - Matrix inverse() const; - /** - return P, calc'ed from L and D - */ - Matrix original() const; -private: - void full_matrix_solve (Vector &,Vector const&) const; - void band_matrix_solve (Vector &, Vector const&) const; - void full_matrix_decompose (Matrix const & P); - void band_matrix_decompose (Matrix const &P); - -}; -#endif diff --git a/flower/include/full-storage.hh b/flower/include/full-storage.hh deleted file mode 100644 index d666953d6f..0000000000 --- a/flower/include/full-storage.hh +++ /dev/null @@ -1,171 +0,0 @@ -/* - full-storage.hh -- declare Full_storage - - source file of the Flower Library - - (c) 1997--1999 Han-Wen Nienhuys -*/ - - -#ifndef FULL_STORAGE_HH -#define FULL_STORAGE_HH - -#include "array.hh" -#include "real.hh" - -#ifndef INLINE -#define INLINE inline -#endif - -/// simplest matrix storage. refer to its baseclass for the doco. -class Full_storage -{ - /// height, width - int height_i_; - int width_i_; - /// maxima. - int max_height_i_; - int max_width_i_; - - /// the storage - Real** els_p_p_; - - INLINE void init() ; - - INLINE bool valid (int i, int j) const ; - - void resize_rows (int); - void resize_cols (int); - -public: - int band_i_; // ugh - - /// check invariants - void OK() const; - - /// height of matrix - INLINE int rows() const; - - /// width of matrix - INLINE int cols() const; - - /// size if square - INLINE int dim() const; - - /** set the size. contents lost. - PRE - i >=0, j>=0 - */ - void set_size (int rows, int cols) ; - - - /**set the size to square dimen. contents lost - PRE - i>=0 - */ - void set_size (int i) ; - /**set the size to i. - - keep contents. If enlarged contents unspecified - - PRE - i>=0, j>=0 - - */ - void resize (int rows, int cols); - - /** - set the size to square dimen. contents kept - Keep contents. If enlarged contents are unspecified - - PRE - i>=0 - */ - void resize (int i); - - - /** - access an element. - - Generate an errormessage, if this happens - in the 0-part of a sparse matrix. - */ - - INLINE Real& elem (int i,int j); - - /// access a element, no modify - INLINE Real elem (int i, int j) const; - - Array row (int i) const ; - Array column (int j) const; - - - /** - add a row to the matrix before row k. Contents - of added row are unspecified - - 0 <= k <= rows() - */ - void insert_row (int k); - - - /** - delete a row from this matrix. - - PRE - 0 <= k < rows(); - */ - void delete_row (int k); - void delete_column (int k); - - - - /** - at end of matrix?. when doing loop - - for (i=0; i -*/ - - -#ifndef MATRIX_STORAGE_HH -#define MATRIX_STORAGE_HH - -#include "array.hh" -#include "real.hh" - -/** - - base class for interface with matrix storageclasses. There are no - iterators for matrixclasses, since matrices are (like arrays) - explicitly int-indexed. - - Iteration is provided by *_next, *_ok, which update and check both - index variables simultaneously. - - TODO - determine type of product matrix. - -*/ -class Matrix_storage { - - -public: - /// check invariants - void OK() const; - - /// height of matrix - int rows() const; - - /// width of matrix - int cols() const; - - /// size if square - int dim() const; - - /** set the size. contents lost. - PRE - i >=0, j>=0 - */ - void set_size (int rows, int cols) ; - - - /**set the size to square dimen. contents lost - PRE - i>=0 - */ - void set_size (int i) ; - /**set the size to i. - - keep contents. If enlarged contents unspecified - - PRE - i>=0, j>=0 - - */ - void resize (int rows, int cols); - - /** - set the size to square dimen. contents kept - Keep contents. If enlarged contents are unspecified - - PRE - i>=0 - */ - void resize (int i); - - - /** - access an element. - - Generate an errormessage, if this happens - in the 0-part of a sparse matrix. - */ - - Real& elem (int i,int j); - - /// access a element, no modify - Real elem (int i, int j) const; - - Array row (int i) const ; - Array column (int j) const; - - - /** - add a row to the matrix before row k. Contents - of added row are unspecified - - 0 <= k <= rows() - */ - void insert_row (int k); - - - /** - delete a row from this matrix. - - PRE - 0 <= k < rows(); - */ - void delete_row (int k); - void delete_column (int k); - ~Matrix_storage() { } - - - /** - at end of matrix?. when doing loop - - for (i=0; i -*/ - - -#ifndef MATRIX_HH -#define MATRIX_HH - -#include "full-storage.hh" -#include "vector.hh" -#include "offset.hh" - -/** a Real matrix. This is a class for a nonsquare block of #Real#s. The - implementation of sparse matrices is done in the appropriate #smat# - class. Matrix only does the mathematical actions (adding, - multiplying, etc.) - - - TODO - implement ref counting? */ - - -class Matrix { - friend Matrix operator *(Matrix const &m1, Matrix const &m2); - -protected: - Full_storage *dat_; - void set (Matrix_storage*); -public: - void OK() const { dat_->OK(); } - void set_band (); - int calc_band_i () const; - int cols() const { return dat_->cols (); } - int rows() const { return dat_->rows (); } - - /** return the size of a matrix. - PRE - the matrix needs to be square. - */ - int dim() const; - - /** - the band size of the matrix. - @ret - - 0 <= band_i() <= dim - */ - int band_i() const; - - /// set entries to r - void fill (Real r); - - /// set diagonal to d - void set_diag (Real d); - - void set_diag (Vector d); - /// set unit matrix - void unit() { set_diag (1.0); } - - void operator+=(Matrix const &m); - void operator-=(Matrix const &m); - void operator*=(Real a); - void operator/=(Real a) { (*this) *= 1/a; } - - /** add a row. - add a row to the matrix before row k - - PRE - v.dim() == cols () - 0 <= k <= rows() - */ - void insert_row (Vector v,int k); - /** . - delete a row from this matrix. - - PRE - 0 <= k < rows(); - */ - void delete_row (int k) { dat_->delete_row (k); } - void delete_column (int k) { dat_->delete_column (k); } - - /** - square n matrix, initialised to null - */ - Matrix (int n); - - /** - n x m matrix, init to 0 - */ - Matrix (int n, int m); - Matrix (Matrix const &m); - - /// dyadic product: v * w.transpose - Matrix (Vector v, Vector w); - void operator=(Matrix const &m); - - /// access an element - Real operator()(int i,int j) const { return dat_->elem (i,j); } - - /// access an element - Real &operator()(int i, int j) { return dat_->elem (i,j); } - - /// Matrix multiply with vec (from right) - Vector operator *(Vector const &v) const; - - /// set this to m1*m2. - void set_product (Matrix const &m1, Matrix const &m2); - - Vector left_multiply (Vector const &) const; - - Matrix operator-() const; - - /// transpose this. - void transpose(); - - /// return a transposed copy. - Matrix transposed() const ; - - Real norm() const; - /** swap. - PRE - 0 <= c1,c2 < cols() - */ - void swap_columns (int c1, int c2); - - /** swap. - PRE - 0 <= c1,c2 < rows() - */ - void swap_rows (int c1, int c2); - - - Vector row (int) const; - Vector col (int) const; - - String str () const; - void print() const; - ~Matrix (); -}; - -inline Vector -operator *(Vector &v, Matrix const & m) { return m.left_multiply (v); } -Matrix operator *(Matrix const & m1,Matrix const &m2); -Matrix operator /(Matrix const &m1,Real a); -inline Matrix operator -(Matrix m1,const Matrix m2) -{ - m1 -= m2; - return m1; -} -inline Matrix operator +(Matrix m1,const Matrix m2) -{ - m1 += m2; - return m1; -} -#endif diff --git a/flower/include/vector.hh b/flower/include/vector.hh deleted file mode 100644 index b4eeb97597..0000000000 --- a/flower/include/vector.hh +++ /dev/null @@ -1,117 +0,0 @@ -#ifndef VECTOR_HH -#define VECTOR_HH - -#include -#include "real.hh" -#include "array.hh" - -class Dstream; -class String; - -/** a row of numbers. - a vector. Storage is handled in Array, Vector only does the mathematics. - */ -class Vector { - Array dat; -public: - void OK() const { dat.OK();} - int dim() const { return dat.size (); } - void fill (Real r) { - for (int i=0; i < dim(); i++) - dat[i] =r; - } - - Vector() { } - Vector (Array d); - Vector (Vector const &n); - - Real &operator()(int i) { return dat[i]; } - Real operator()(int i) const { return dat[i]; } - Real elem (int i) { return dat[i]; } - Vector (int n) { - dat.set_size (n); - fill (0); - } - void set_dim (int i) - { - dat.set_size (i); - } - - void insert (Real v, int i) { - dat.insert (v,i); - } - void del (int i) { dat.del (i); } - - String str () const; - - void operator +=(Vector v) { - assert (v.dim() == dim ()); - for (int i=0; i < dim(); i++) - dat[i] += v.dat[i]; - } - - void operator -=(Vector v) { - assert (v.dim() == dim ()); - for (int i=0; i < dim(); i++) - dat[i] -= v (i); - } - - Real operator *(Vector v) const { - Real ip=0; - assert (v.dim() == dim ()); - for (int i=0; i < dim(); i++) - ip += dat[i] *v (i); - return ip; - } - void operator *=(Real a) { - for (int i=0; i < dim(); i++) - dat[i] *= a; - } - void operator /=(Real a) { - (*this) *= 1/a; - } - Vector operator-() const; - Real norm_sq() { - return ((*this) * (*this)); - } - Real norm() { - return sqrt (norm_sq()); - } - operator Array() { return dat; } - Array const &to_array()const { return dat; } - void print() const; - /// set to j-th element of unit-base - void set_unit (int j) ; -}; - -inline Vector -operator+(Vector a, Vector const &b) { - a += b; - return a; -} - -inline Vector -operator-(Vector a, Vector const &b) { - a -= b; - return a; -} - -inline Vector -operator*(Vector v, Real a) { - v *= a; - return v; -} - -inline Vector -operator*(Real a,Vector v) { - v *= a; - return v; -} - -inline Vector -operator/(Vector v,Real a) { - v *= 1/a; - return v; -} - -#endif diff --git a/flower/matrix-debug.cc b/flower/matrix-debug.cc deleted file mode 100644 index 323962e2c2..0000000000 --- a/flower/matrix-debug.cc +++ /dev/null @@ -1,64 +0,0 @@ -/* - matrix-debug.cc -- implement Matrix print routines - - source file of the Flower Library - - (c) 1997--1999 Han-Wen Nienhuys -*/ - - -#include "flower-debug.hh" -#include "matrix.hh" - -String -Matrix::str () const -{ - String s; -#ifndef NPRINT - Full_storage const * stor_c_l = dat_; - s = String ("matrix {"); - for (int i=0; i< rows(); i++) - { - for (int j = 0; j < cols(); j++) - { - s+= to_str (stor_c_l->elem (i,j), "%6f "); - } - s+="\n"; - } - s+="}\n"; -#endif - return s; -} - - -void -Matrix::print () const -{ -#ifndef NPRINT - fdebug << *this; -#endif -} - -String -Vector::str () const -{ - String s; -#ifndef NPRINT - s = String ("vector (") + to_str (dim ()) + ") ["; - for (int i=0; i < dim(); i++) - { - s += to_str (dat[i], "%6f") + to_str (' '); - } - s += "]\n"; -#endif - return s; -} - - -void -Vector::print() const -{ -#ifndef NDEBUG - fdebug << *this << '\n'; -#endif -} diff --git a/flower/matrix.cc b/flower/matrix.cc deleted file mode 100644 index b671238ac5..0000000000 --- a/flower/matrix.cc +++ /dev/null @@ -1,299 +0,0 @@ -/* - matrix.cc -- implement Matrix - - source file of the Flower Library - - (c) 1997--1999 Han-Wen Nienhuys -*/ - -#include "matrix.hh" -#include "full-storage.hh" - -Real -Matrix::norm() const -{ - Real r =0.0; - for (int i=0, j=0; dat_->mult_ok (i,j); dat_->mult_next (i,j)) - r += sqr (dat_->elem (i,j)); - return sqrt (r); -} - -void -Matrix::fill (Real r) -{ - for (int i=0, j=0; dat_->mult_ok (i,j); dat_->mult_next (i,j)) - dat_->elem (i,j)=r; -} - -void -Matrix::set_diag (Real r) -{ - for (int i=0, j=0; dat_->mult_ok (i,j); dat_->mult_next (i,j)) - dat_->elem (i,j)=(i==j) ? r: 0.0; -} - -void -Matrix::set_diag (Vector d) -{ - for (int i=0, j=0; dat_->mult_ok (i,j); dat_->mult_next (i,j)) - dat_->elem (i,j)=(i==j) ? d (i): 0.0; -} - -void -Matrix::operator+=(Matrix const &m) -{ - assert (m.cols() == cols ()); - assert (m.rows() == rows ()); - for (int i=0, j=0; dat_->mult_ok (i,j); dat_->mult_next (i,j)) - dat_->elem (i,j) += m (i,j); -} - -void -Matrix::operator-=(Matrix const &m) -{ - assert (m.cols() == cols ()); - assert (m.rows() == rows ()); - for (int i=0, j=0; dat_->mult_ok (i,j); dat_->mult_next (i,j)) - dat_->elem (i,j) -= m (i,j); -} - - -void -Matrix::operator*=(Real a) -{ - for (int i=0, j=0; dat_->mult_ok (i,j); dat_->mult_next (i,j)) - dat_->elem (i,j) *= a; -} - -void -Matrix::operator=(Matrix const &m) -{ - if (&m == this) - return ; - delete dat_; - dat_ = new Full_storage (*m.dat_); -} - - -int -Matrix::band_i () const -{ - return dat_->band_i_; -} - -void -Matrix::set_band () -{ - dat_->band_i_ = calc_band_i (); -} - -int -Matrix::calc_band_i() const -{ - int starty = dim(); - while (starty >= 0) - { - for (int i = starty, j = 0; i < dim(); i++, j++) - if (dat_->elem (i,j)) - goto gotcha; - for (int i=0, j = starty; j < dim(); i++,j++) - if (dat_->elem (i,j)) - goto gotcha; - starty --; - } - gotcha: - return starty; -} - -Matrix::Matrix (Matrix const &m) -{ - dat_ = new Full_storage (*m.dat_); -} - - -Matrix::Matrix (int n, int m) -{ - dat_ = new Full_storage(n,m); - fill (0); -} - - -Matrix::Matrix (int n) -{ - dat_= new Full_storage (n,n); - fill (0); -} - -Matrix::Matrix (Vector v, Vector w) -{ - int n = v.dim(); - int m = w.dim (); - dat_= new Full_storage (n,m); - for (int i=0; i < n; i++) - for (int j=0; j < m ; j++) - dat_->elem (i,j)=v (i)*w (j); -} - - -Vector -Matrix::row (int k) const -{ - int n=cols(); - - - Vector v (n); - for (int i=0; i < n; i++) - v (i)=dat_->elem (k,i); - - return v; -} - -Vector -Matrix::col (int k) const -{ - int n=rows(); - Vector v (n); - for (int i=0; i < n; i++) - v (i)=dat_->elem (i,k); - return v; -} - -Vector -Matrix::left_multiply (Vector const & v) const -{ - Vector dest (v.dim()); - assert (dat_->cols()==v.dim ()); - for (int i=0, j=0; dat_->mult_ok (i,j); dat_->mult_next (i,j)) - dest (i)+= dat_->elem (j,i)*v (j); - return dest; -} - -Vector -Matrix::operator *(Vector const & v) const -{ - Vector dest (rows()); - assert (dat_->cols()==v.dim ()); - for (int i=0, j=0; dat_->mult_ok (i,j); dat_->mult_next (i,j)) - dest (i)+= dat_->elem (i,j)*v (j); - return dest; -} - -Matrix -operator /(Matrix const& m1,Real a) -{ - Matrix m (m1); - m /= a; - return m; -} - -/* - ugh. Only works for square matrices. - */ -void -Matrix::transpose() -{ - for (int i=0, j=0; dat_->mult_ok (i,j); dat_->mult_next (i,j)) - { - if (i >= j) - continue; - Real r=dat_->elem (i,j); - dat_->elem (i,j) = dat_->elem (j,i); - dat_->elem (j,i)=r; - } -} - -Matrix -Matrix::operator-() const -{ - OK(); - Matrix m (*this); - m*=-1.0; - return m; -} - -Matrix -Matrix::transposed() const -{ - Matrix m (*this); - m.transpose(); - return m; -} - -Matrix -operator *(Matrix const &m1, Matrix const &m2) -{ - Matrix result (m1.rows (), m2.cols ()); - result.set_product (m1,m2); - return result; -} - - - -void -Matrix::set_product (Matrix const &m1, Matrix const &m2) -{ - assert (m1.cols()==m2.rows ()); - assert (cols()==m2.cols () && rows ()==m1.rows ()); - - for (int i=0; i < rows (); i++) - for (int j=0; j < cols (); j++) - { - Real r=0.0; - for (int k = 0 >? i - m1.band_i () >? j - m2.band_i (); - k < m1.cols(); k++) - { - r += m1(i,k)*m2(k,j); - } - dat_->elem (i,j)=r; - } -} - -void -Matrix::insert_row (Vector v, int k) -{ - int c = cols(); - assert (v.dim()==cols ()); - dat_->insert_row (k); - for (int j=0; j < c; j++) - dat_->elem (k,j)=v (j); -} - - -void -Matrix::swap_columns (int c1, int c2) -{ - assert (c1>=0&& c1 < cols()&&c2 < cols () && c2 >=0); - int r = rows(); - for (int i=0; i< r; i++) - { - Real r=dat_->elem (i,c1); - dat_->elem (i,c1) = dat_->elem (i,c2); - dat_->elem (i,c2)=r; - } -} - -void -Matrix::swap_rows (int c1, int c2) -{ - assert (c1>=0&& c1 < rows()&&c2 < rows () && c2 >=0); - int c = cols(); - for (int i=0; i< c; i++) - { - Real r=dat_->elem (c1,i); - dat_->elem (c1,i) = dat_->elem (c2,i); - dat_->elem (c2,i)=r; - } -} - - -int -Matrix::dim() const -{ - assert (cols() == rows ()); - return rows(); -} - -Matrix::~Matrix () -{ - delete dat_; -} diff --git a/flower/vector.cc b/flower/vector.cc deleted file mode 100644 index e4eb654d69..0000000000 --- a/flower/vector.cc +++ /dev/null @@ -1,35 +0,0 @@ -/* - vector.cc -- implement Vector - - source file of the Flower Library - - (c) 1996-98 Han-Wen Nienhuys - - */ - -#include "vector.hh" - -Vector::Vector (Array d) - : dat (d) -{ - -} -Vector::Vector (Vector const &n) - : dat (n.dat) -{ -} - -Vector -Vector::operator-() const -{ - Vector v (*this); - v*=-1; - return v; -} - -void -Vector::set_unit (int j) -{ - fill (0.0); - dat[j] = 1.0; -} diff --git a/input/bugs/auto-melismas.ly b/input/bugs/auto-melismas.ly deleted file mode 100644 index bafa169559..0000000000 --- a/input/bugs/auto-melismas.ly +++ /dev/null @@ -1,17 +0,0 @@ -\score{ - \addlyrics - \context Staff { - \notes\relative c''{ c () d r e f } - \notes\relative c''{ c ~ c r c c } - } - { - \context Lyrics \lyrics { foo __ bar baz } - \context Lyrics \lyrics { foo __ bar baz } - } - \paper{ - \translator{ - \VoiceContext - automaticMelismas=1; - } - } -} diff --git a/input/bugs/hara.ly b/input/bugs/hara.ly deleted file mode 100644 index 85ff14adf8..0000000000 --- a/input/bugs/hara.ly +++ /dev/null @@ -1,88 +0,0 @@ - - -The following example shows 2 bugs with HaraKiriStaffContext. - -The first bug is that when the second staff disappears, the instr -string in the margin "2" is still printed, and surprisingly appears -above staff 1. - -The second and more serious bug is that when staff 2 reappears, -the bar separators are no longer being printed for that staff. - -Regards, -Roy Rankin - -\header { - composer = "Music: Trad."; - crossRefNumber = "4"; - title = "HaraKiri Bugs"; -} -voice1 = \notes { -\property Staff.instrument = "Part1" -\property Staff.instr = "1" -\property Staff.timeSignatureStyle="C" - \time 4/4; \key c; - e''4. e''8 e''4 r4 | e''4 d''4 - c''2 | - e''4 d''4 c''2 ( | - e''2 c''2 ~ | - ) c''1 \bar "|."; \break - e''4. e''8 e''4 r4 | e''4 d''4 - c''2 | - e''4 d''4 c''2 ( | - e''2 c''2 ~ | - ) c''1 \bar "|."; \break - e''4. e''8 e''4 r4 | e''4 d''4 - c''2 | - e''4 d''4 c''2 ( | - e''2 c''2 ~ | - ) c''1 \bar "|."; \break - e''4. e''8 e''4 r4 | e''4 d''4 - c''2 | - e''4 d''4 c''2 ( | - e''2 c''2 ~ | - ) c''1 \bar "|."; -} -voice2 = \notes { -\property Staff.instrument = "Part2" -\property Staff.instr = "2" -\property Staff.timeSignatureStyle="C" - \time 4/4; \key c; - e''4. e''8 e''4 r4 | e''4 d''4 - c''2 | - e''4 d''4 c''2 ( | - e''2 c''2 ~ | - ) c''1 \bar "|."; -s1 * 10 - e''4. e''8 e''4 r4 | e''4 d''4 - c''2 | - e''4 d''4 c''2 ( | - e''2 c''2 ~ | - ) c''1 \bar "|."; -} -\score{ - \notes < - - \context Staff="1" - { - \$voice1 - } - \context Staff="2" - { - \$voice2 - } - - > - \paper - { - \translator - { - \HaraKiriStaffContext - \consists Staff_margin_engraver; - } - } - \midi {} -} - - - diff --git a/input/bugs/l.ly b/input/bugs/l.ly deleted file mode 100644 index f042e73068..0000000000 --- a/input/bugs/l.ly +++ /dev/null @@ -1,13 +0,0 @@ -\score { - \notes { - \repeat volta 5 { - c' c' c' c' - } - \alternative { - {d' d' d' d'} - {e' e' e' e'} - } - f' - } - } - diff --git a/input/bugs/utremi_old.ly b/input/bugs/utremi_old.ly deleted file mode 100644 index 60a62403fa..0000000000 --- a/input/bugs/utremi_old.ly +++ /dev/null @@ -1,262 +0,0 @@ -\header { -filename = "utremi.ly"; -enteredby = "Christian Mondrup"; -%composer = "Thomas Ravenscroft"; -opus = "Thomas Ravenscroft"; -%opus = "Pammelia, 1609, no. 31"; -arranger = "Pammelia, 1609, no. 31"; -title = "Vt, re, me, fa, sol, la, sol, fa, re, vt"; -tagline = "Typeset with GNU LilyPond by Christian Mondrup. Non-commercial copying welcome."; -} - -papersize = "letter" - -\include "paper16.ly" -%\include "a4.ly" - -lyricsOne = \lyrics { - ""\breve ""8 - "Vt,"\breve re, me, fa, sol, la, "" - sol, la, fa, me, re, vt. -} - -lyricsTwo = \lyrics { - "" "" - Hey downe downe hey downe downe downe downe hey down - hey down "" down down a. - My heart of gold as true as steele - as I me leant "" vn -- to the bowres. - but if my La -- dy loue me well, - Lord so Ro -- bin bowres, -} - -lyricsThree = \lyrics { - "" "" "" "" "" "" "" "" "" "" "" "" - heaue and hoe Rum -- be -- lo, hey tro -- lo tro -- ly lo, - hey tro -- ly trol -- ly hey tro -- ly trol -- ly - hey tro -- ly trol -- ly hey trol -- "" ly trol -- ly lo. - My La -- dies gone to Can -- ter -- bu -- ry - S. Tho -- mas be her boote. - Shee met with Kate of Malms -- "" bu -- ry, - why weepst thou ma -- ple -} - -lyricsFour = \lyrics { - "" - roote: - O sleepst thou or wakst thon Ies -- se -- ry, Cooke, - the rost it burnes, turne round turne round "" a -- bout "" - turne round "" a -- bout, "" - turne round "" a -- bout, "" turne round. "" - O Fri -- er how fares thy ban -- de -- low ban -- de -- low - Fri -- er, how fares thy San -- de -- low, San -- de -- low. -} - -global = \notes { - \key c; - %\time 4/2; - %\time 17/8; - \time 33/16; - \clef "tenor"; -} - -dummyBeat = \notes { - \property Staff.defaultBarType = "||" - \property Staff.barAlways = "1" - s16 - \property Staff.barAlways = "0" - \property Staff.defaultBarType = "|" -} - -incipOne = \notes\relative c' { - g\breve - \dummyBeat - \bar ""; -} - -incipTwo = \notes\relative c'' { - \property Voice.noteHeadStyle = "diamond" - g1. - \property Voice.noteHeadStyle = "" - f2 - \dummyBeat - \bar ""; -} - -incipThree = \notes\relative c' { - \property Voice.noteHeadStyle = "harmonic" - \property Voice.tupletVisibility = 0 - \times 2/3 { d4. e8 d4 } \times 2/3 { b4. a8 g4 } - \times 2/3 { b4. c8 b4 } \times 2/3 { b4. a8 g4 } - \dummyBeat - \property Voice.tupletVisibility = 3 - \bar ""; -} - -incipFour = \notes\relative c' { - \property Voice.noteHeadStyle = "diamond" - %\property Voice.restStyle = "mensural" - b1 r1 - \dummyBeat - \bar ""; -} - -partOne = \notes\relative c' { - %\property Score.currentBarNumber = "1" - \property Voice.noteHeadStyle = "" - %\property Staff.barNumberScriptPadding = 1 - g1 ~ | g | a ~ | a | bes ~ | bes | c ~ | c | d ~ | d | e ~ | e | - r | r | - e ~ | e | d ~ | d | c ~ | c | bes ~ | bes | a ~ | a | g ~ | g - \bar "|."; -} - -partTwo = \notes\relative c'' { - \property Voice.noteHeadStyle = "" - g2. f4 | e2 d | c4. b8 a4 g | f2 e | - d g ~ | g4 a g2 | c,1 | r2 c | - g'2. a4 | b2 g | c c | c c | - b4 c d2 ~ | d4 d b g | c1 | r2 c | - g g | d2. e4 | f2 e | c1 | - g'2. g4 | g2 g | c,1 ~ | c | r | r - \bar "|."; -} - -partThree = \notes\relative c' { - \property Voice.noteHeadStyle = "" - \times 2/3 { d4. e8 d4 } \times 2/3 {b4. a8 g4 } | - \times 2/3 { b4. c8 b4 } \times 2/3 { b4. a8 g4 } - c2 c4 c | c2 c | - d2 d4 d | d2 d | e2. f4 | e2 d4 c | - b2 g ~ | g4 f e d | c1 | r2 c | - g'2. g4 | g2 g | c,2. c4 | c c c'2 | - bes g | bes bes | a1 | r2 e' | - d d | g g | a2. g4 | fis e d2 | - g g | d2. c4 - \bar "|."; -} - -partFour = \notes\relative c' { - \property Voice.noteHeadStyle = "" - b1^\fermata | r2 g' | e e4 e | a2. g4 | - fis2 e4 d | g1 | g2 g | g g | - g g | g g ~ | g4 f e d | c2 d ~ | - d4 c b a | g2 g' ~ | g4 f e d | c2 g' ~ | - g r | \times 2/3 { r1 d2 } \times 2/3 { a'2 a g } | - \times 2/3 { e1 f2 } \times 2/3 { g2. g4 g2 } - \times 2/3 { d2. d4 e2 } \times 2/3 { f2. f4 e2 } - \times 2/3 { d1 c2 } \times 2/3 { b2. a4 g2 } - \times 2/3 { d'2. d4 d2 } - \bar "|."; -} - -partOneStaff = < - \context Staff = vocal1 < - \property Staff.clefStyle = "fullSizeChanges" - \property Staff.instrument = "\large{1.}" - %\property Staff.instr = "" - \notes { - \global - \property Staff.timeSignatureStyle = "old4/4" - \incipOne - \clef "G_8"; - \property Staff.timeSignatureStyle = "C2/2" - \time 4/4; - \partOne - } - {\context Lyrics = lyrOne \lyricsOne } - > -> - -partTwoStaff = < - \context Staff = vocal2 < - \property Staff.clefStyle = "fullSizeChanges" - \property Staff.instrument = "\large{2.}" - %\property Staff.instr = "" - %\property Voice.automaticMelismata = "1" - \addlyrics - \notes { - \global - \property Staff.timeSignatureStyle = "old4/4" - \incipTwo - \clef "G_8"; - \property Staff.timeSignatureStyle = "C2/2" - \time 4/4; - \partTwo - } - {\context Lyrics = lyrTwo \lyricsTwo } - > -> - - -partThreeStaff = < - \context Staff = vocal3 < - \property Staff.clefStyle = "fullSizeChanges" - \property Staff.instrument = "\large{3.}" - %\property Staff.instr = "" - %\property Voice.automaticMelismata = "1" - \addlyrics - \notes { - \context Voice = vthree - \property Voice.tupletDirection = \up - \global - \property Staff.timeSignatureStyle = "old4/4" - \incipThree - \clef "G_8"; - \property Staff.timeSignatureStyle = "C2/2" - \time 4/4; - \partThree - } - {\context Lyrics = lyrThree \lyricsThree } - > -> - -partFourStaff = < - \context Staff = vocal4 < - \property Staff.clefStyle = "fullSizeChanges" - \property Staff.timeSignatureStyle = "C2/2" - \property Staff.barNumberDirection = \up - \property Staff.instrument = "\large{4.}" - \property Staff.instr = "" - %\property Voice.automaticMelismata = "1" - \addlyrics - \notes { - \context Voice = vfour - \property Voice.tupletDirection = \up - \global - %\property Staff.instrument = "\large{4.}" - %\property Staff.instr = "" - \property Staff.timeSignatureStyle = "old4/4" - \incipFour - \clef "G_8"; - \property Staff.timeSignatureStyle = "C2/2" - \time 4/4; - \partFour - } - {\context Lyrics = lyrFour \lyricsFour} - > -> - -\score { - \context StaffGroup < - \partFourStaff - \partThreeStaff - \partTwoStaff - \partOneStaff - > - %\paper {\translator {\BarNumberingStaffContext } } - \paper { - indent = 0.\mm; - linewidth = 17.0\cm; - textheight = 27.0\cm; - gourlay_maxmeasures=6.0; - \translator { \StaffContext - \consists "Staff_margin_engraver"; } - %\translator { \ScoreContext - % minVerticalAlign = 1.5*\staffheight; } - } - %\midi { - % output = "utremi.mid"; - % \tempo 2 = 80; - %} -} diff --git a/lily/col-info.cc b/lily/col-info.cc deleted file mode 100644 index 2ebfd15d4f..0000000000 --- a/lily/col-info.cc +++ /dev/null @@ -1,70 +0,0 @@ -/* - col-info.cc -- implement Column_info - - source file of the GNU LilyPond music typesetter - - (c) 1997--1999 Han-Wen Nienhuys -*/ - -#include "paper-column.hh" -#include "col-info.hh" -#include "debug.hh" - -void -Column_info::print() const -{ -#ifndef NPRINT - DOUT << "column { "; - if (fixed_b()) - DOUT << "fixed at " << fixed_position() << ", "; - assert (pcol_l_); - DOUT << width_.str(); - Direction d = LEFT; - do { - for (int i=0; i < rods_[d].size (); i++) - rods_[d][i].print (); - } while (flip (&d) != LEFT); - - DOUT <<"}\n"; -#endif -} - -Column_info::Column_info (Paper_column *col_l, Real const *fixed_C) -{ - if (fixed_C) - fixpos_p_.set_l (fixed_C); - pcol_l_ = col_l; - width_ = pcol_l_->extent(X_AXIS); - if (width_.empty_b()) - width_ = Interval(0,0); -} - - -Column_info::Column_info() -{ - pcol_l_ =0; -} - -bool -Column_info::fixed_b () const -{ - return fixpos_p_.get_C(); -} - -Real -Column_info::fixed_position () const -{ - return *fixpos_p_; -} - -int -Column_info::rank_i () const -{ - return pcol_l_->rank_i (); -} - -void -Spacer_rod::print ()const -{ - DOUT << "Other " << other_idx_ << "dist = " << distance_f_ << '\n'; -} diff --git a/lily/include/ineq-constrained-qp.hh b/lily/include/ineq-constrained-qp.hh deleted file mode 100644 index 247f2a7688..0000000000 --- a/lily/include/ineq-constrained-qp.hh +++ /dev/null @@ -1,91 +0,0 @@ -/* - ineq-constrained-qp.hh -- declare Ineq_constrained_qp - - source file of the GNU LilyPond music typesetter - - (c) 1997--1999 Han-Wen Nienhuys -*/ - - -#ifndef INEQ_CONSTRAINED_QP_HH -#define INEQ_CONSTRAINED_QP_HH - - - -#include "matrix.hh" - -/** inequality constrained quadratic program - - It takes the form of - - optimise for x : x*quad_ *x + lin_* x + const_term_ - - subject to for all i: cons_[i] * x >= consrhs_[i] - - - @usage: - instantiate Ineq_constrained_qp. - - Modify quad_, lin_ and const_term_ directly. Use - add_inequality_cons () to add constraints. Call solve () with a - feasible solution - - - */ -class Ineq_constrained_qp { - friend class Active_constraints; - - Array cons_; - Array consrhs_; -public: - Matrix quad_; - Vector lin_; - Real const_term_; - - - /** - use a KKT method to assert optimality of sol - */ - void assert_solution (Vector sol) const; - /// solve the problem using a projected gradient method - Vector constraint_solve (Vector) const; - /** - Solve it. First try it the easy way. - */ - Vector solve (Vector start) const; - - /** - @return the number of variables in the problem - */ - int dim() const; - - /** - add a constraint - - - c*vars >= r - - PRE - c.dim() == dim (); - - */ - void add_inequality_cons (Vector c, double r); - - /** set up matrices to go with the problem. */ - Ineq_constrained_qp (int novars); - - /** - evaluate the quadratic function for input #v# - */ - Real eval (Vector v); - - void eliminate_var (int idx, Real value); - void OK() const; - void print() const; - -}; - -// ugh -const Real EPS=1e-7; - -#endif // INEQ_CONSTRAINED_QP_HH diff --git a/lily/include/qlp.hh b/lily/include/qlp.hh deleted file mode 100644 index af7513e5d0..0000000000 --- a/lily/include/qlp.hh +++ /dev/null @@ -1,44 +0,0 @@ -/* - qlp.hh -- declare Ineq_constrained_qp, Mixed_qp - - source file of the GNU LilyPond music typesetter - - (c) 1997--1999 Han-Wen Nienhuys -*/ - -#ifndef QLP_HH -#define QLP_HH - -#include "ineq-constrained-qp.hh" - -/** - Quadratic programming with mixed linear constraints. - problem definition of a quadratic optimisation problem with linear - inequality and equality constraints - - - x^T QUAD x /2 + b^T x -*/ -class Mixed_qp :public Ineq_constrained_qp { - Array eq_cons; - Array eq_consrhs; -public: - Mixed_qp (int n); - void OK() const; - void print() const; - - Vector solve (Vector start) const; - void add_fixed_var (int i , Real value); - - - /** - add a constraint, - - c*vars == r - - PRE - c.dim()==dim (); - */ - void add_equality_cons (Vector c, double r); -}; -#endif diff --git a/lily/include/qlpsolve.hh b/lily/include/qlpsolve.hh deleted file mode 100644 index 9a342527ee..0000000000 --- a/lily/include/qlpsolve.hh +++ /dev/null @@ -1,95 +0,0 @@ -/* - qlpsolve.hh -- declare Active_constraints, Inactive_iter - - source file of the GNU LilyPond music typesetter - - (c) 1997--1999 Han-Wen Nienhuys -*/ - - -#ifndef QLPSOLVE_HH -#define QLPSOLVE_HH - -#include "matrix.hh" - - -/** - This class represents the set of active (binding) constraints - which can be active while the QLP algorithm is in a feasible - point. The active constraints are numbered. - If the constraints are of the form - - A^T*x >= b - - then the binding constraints are those where the >= is equality. - - */ - -class Active_constraints { - friend class Inactive_iter; - - - Matrix A,H; - Array active; - Array inactive; // actually this is a set, not an array. - Ineq_constrained_qp const *opt; - -public: - - /** This counts the number of errors the algorithms makes. The - optimum search should be abandoned if it becomes too high. */ - int degenerate_count_i_; - String status() const; - - Vector vec (int k) const { return opt->cons_[k]; } - Real rhs (int k) const { return opt->consrhs_[k]; } - - - /** drop constraint. drop constraint k from the active set. k is the index of the - constraint in #active# - - */ - void drop_constraint (int k); - - - /** add constraint j. - add constraint j to the active set j is the index of the - constraint in #inactive# - */ - void add_constraint (int j); - - /// exchange in and out. - void exchange (int in, int out) { add_constraint (in); drop_constraint (out); } - - - Vector find_active_optimum (Vector g); - - /// get lagrange multipliers. - Vector get_lagrange (Vector v); - - Active_constraints (Ineq_constrained_qp const *op); - /** construct: no constraints active, n vars. Put the equalities - into the constraints. */ - - /// check invariants - void OK(); -}; - - -/** - loop through the inactive constraints. - */ -class Inactive_iter { - int j; - Active_constraints const* ac; -public: - Inactive_iter (Active_constraints const &c) { ac=&c; j=0; } - int idx() const { return j; } - void operator ++(int) { j++; } - int constraint_id() const { return ac->inactive[j]; } - Vector vec() const { return ac->vec (constraint_id ()); } - Real rhs() const { return ac->rhs (constraint_id ()); } - bool ok() const { return j < ac->inactive.size (); } -}; - -#endif // QLPSOLVE_HH diff --git a/lily/include/spring-spacer.hh b/lily/include/spring-spacer.hh deleted file mode 100644 index 7359f6a1a6..0000000000 --- a/lily/include/spring-spacer.hh +++ /dev/null @@ -1,99 +0,0 @@ -/* - spring-spacer.hh -- declare Spring_spacer - - source file of the GNU LilyPond music typesetter - - (c) 1997--1999 Han-Wen Nienhuys -*/ - - -#ifndef SPRING_SPACER_HH -#define SPRING_SPACER_HH - -#include "line-spacer.hh" -#include "cons.hh" -#include "col-info.hh" -#include "column-x-positions.hh" -#include "moment.hh" - - - -/** - Determine positions of columns connected by springs and held apart by rods. - - - Generate a spacing which - \begin{itemize} - \item - Satisfies spacing constraints (notes can't be printed through each other) - \item - Looks good, ie tries to conform to an ideal spacing as much as possible. - \end{itemize} - - This is converted by regarding idealspacing as "springs" attached - to columns. The equilibrium of one spring is the ideal - distance. The columns have a size, this imposes "hard" constraints - on the distances. This transforms the problem into a quadratic - programming problem with linear constraints. - - The quality is given by the total potential energy in the - springs. The lower the energy, the better the configuration. - - - TODO: too complicated. Revise. - Use force iso. energy. Also optimise for uniform density -*/ - -class Spring_spacer : public Line_spacer { - friend class Durations_iter; -private: - // can't copy me. - Spring_spacer (Spring_spacer const&s); - Cons *ideal_p_list_; - Array cols_; - Real indent_f_; - - /// the index of #c# in #cols# - int col_id (Paper_column const *c) const; - - /// generate an (nonoptimal) solution - Vector find_initial_solution() const; - - /// does #this# contain the column #w#? - bool contains_b (Paper_column const *w); - - /// make the energy function - void make_matrices (Matrix &quad, Vector &lin,Real&) const; - - - /// generate the LP constraints - void make_constraints (Mixed_qp& lp) const; - - - void handle_loose_cols(); - bool try_initial_solution_and_tell (Vector&)const; - Vector try_initial_solution() const; - - void set_fixed_cols (Mixed_qp&) const; - - Score_column* scol_l (int); - void connect (int i,int j, Real,Real); - Real calculate_energy_f (Vector) const; -public: - static Line_spacer *constructor(); - Real energy_normalisation_f_; - Spring_spacer (); - virtual ~Spring_spacer (); - virtual void solve (Column_x_positions*) const; - virtual void lower_bound_solution (Column_x_positions*) const; - virtual void add_columns (Link_array); - void add_column (Paper_column *, bool, Real); - - virtual Vector default_solution() const; - virtual bool check_constraints (Vector v) const; - virtual void OK() const; - virtual void print() const; - virtual void prepare(); -}; - -#endif // SPRING_SPACER_HH diff --git a/lily/ineq-constrained-qp.cc b/lily/ineq-constrained-qp.cc deleted file mode 100644 index 563199c938..0000000000 --- a/lily/ineq-constrained-qp.cc +++ /dev/null @@ -1,294 +0,0 @@ -/* - ineq-constrained-qp.cc -- implement Ineq_constrained_qp - - source file of the GNU LilyPond music typesetter - - (c) 1996, 1997--1999 Han-Wen Nienhuys -*/ -#include "ineq-constrained-qp.hh" -#include "qlpsolve.hh" -#include "debug.hh" -#include "choleski.hh" -#include "main.hh" - -/* - May be this also should go into a library - */ - -const int MAXITER=100; // qlpsolve.hh - -const int MAXDEGEN=5; - -/* - assume x (idx) == value, and adjust constraints, lin and quad accordingly - - TODO: add const_term - */ -void -Ineq_constrained_qp::eliminate_var (int idx, Real value) -{ - Vector row (quad_.row (idx)); - row*= value; - - quad_.delete_row (idx); - - quad_.delete_column (idx); - - lin_.del (idx); - row.del (idx); - lin_ +=row ; - - for (int i=0; i < cons_.size(); i++) - { - consrhs_[i] -= cons_[i](idx) *value; - cons_[i].del (idx); - } -} - -void -Ineq_constrained_qp::add_inequality_cons (Vector c, double r) -{ - cons_.push (c); - consrhs_.push (r); -} - -Ineq_constrained_qp::Ineq_constrained_qp (int novars): - quad_ (novars), - lin_ (novars), - const_term_ (0.0) -{ -} - -void -Ineq_constrained_qp::OK() const -{ -#if !defined (NDEBUG) && defined (PARANOID) - assert (cons_.size() == consrhs_.size ()); - Matrix Qdif= quad_ - quad_.transposed(); - assert (Qdif.norm()/quad_.norm () < EPS); -#endif -} - - -Real -Ineq_constrained_qp::eval (Vector v) -{ - return v * quad_ * v + lin_ * v + const_term_; -} - - -int -min_elt_index (Vector v) -{ - Real m=infinity_f; - int idx=-1; - for (int i = 0; i < v.dim(); i++) - { - if (v (i) < m) - { - idx = i; - m = v (i); - } - assert (v (i) <= infinity_f); - } - return idx; -} - - -/**the numerical solving. Mordecai Avriel, Nonlinear Programming: analysis and methods (1976) - Prentice Hall. - - Section 13.3 - - This is a "projected gradient" algorithm. Starting from a point x - the next point is found in a direction determined by projecting - the gradient onto the active constraints. (well, not really the - gradient. The optimal solution obeying the active constraints is - tried. This is why H = Q^-1 in initialisation)) - - - */ -Vector -Ineq_constrained_qp::constraint_solve (Vector start) const -{ - if (!dim()) - return Vector (0); - - Active_constraints act (this); - act.OK(); - - - Vector x (start); - Vector gradient=quad_*x+lin_; - - - // Real fvalue = x*quad_*x/2 + lin*x + const_term;// it's no use. - Vector last_gradient (gradient); - int iterations=0; - - while (iterations++ < MAXITER && act.degenerate_count_i_ < MAXDEGEN) - { - if (experimental_features_global_b) - assert_solution (x); - - Vector direction= - act.find_active_optimum (gradient); - - DOUT << "gradient "<< gradient<< "\ndirection " << direction<< '\n'; - - if (direction.norm() > EPS) - { - DOUT << act.status() << '\n'; - - Real unbounded_alfa = 1.0; - Real minalf = unbounded_alfa; - - Inactive_iter minidx (act); - - - /* - we know the optimum on this "hyperplane". Check if we - bump into the edges of the simplex - */ - - for (Inactive_iter ia (act); ia.ok(); ia++) - { - Real dot = ia.vec() * direction; - Real mindot = (experimental_features_global_b) - ? -EPS - : 0; - - if (dot >= mindot) - continue; - - - Real numerator = ia.rhs () - ia.vec()*x; - if (numerator >= 0) - { - if (numerator > EPS) - { - warning (_f ("Ineq_constrained_qp::solve (): Constraint off by %f", numerator)); - act.degenerate_count_i_ ++; - } - minalf = -numerator; - minidx = ia; - break; - } - - Real alfa = numerator / dot; - - - if (minalf > alfa) - { - minidx = ia; - minalf = alfa; - } - } - - Real optimal_step = minalf; - - Vector deltax = direction * optimal_step; - x += deltax; - gradient += optimal_step * (quad_ * deltax); - - DOUT << "step = " << optimal_step << " (|dx| = " << - to_str (deltax.norm()) << ")\n"; - - if (minalf < unbounded_alfa) - { - /* bumped into an edge. try again, in smaller space. */ - act.add_constraint (minidx.idx()); - DOUT << "adding cons "<< minidx.idx () << '\n'; - continue; - } - /*ASSERT: we are at the optimal solution for this "plane"*/ - } - - Vector lagrange_mult=act.get_lagrange (gradient); - int m= min_elt_index (lagrange_mult); - - if (m>=0 && lagrange_mult (m) > 0) - { - break; // optimal sol. - } - else if (m<0) - { - Real n =gradient.norm(); - if (n >= EPS) - { - programming_error ("Huh? Gradient should be zero "); - act.degenerate_count_i_ ++; - } - else if (n <0) - { - programming_error ("Huh? Norm should be positive"); - act.degenerate_count_i_++; - } - break; - } - - DOUT << "dropping cons " << m << '\n'; - act.drop_constraint (m); - } - if (iterations >= MAXITER) - WARN << _ ("Didn't converge!") << '\n'; - if (act.degenerate_count_i_ >= MAXDEGEN) - WARN << _ ("too much degeneracy") << '\n'; - DOUT << ": found " << x << " in " << iterations <<" iterations\n"; - assert_solution (x); - return x; -} - - -Vector -Ineq_constrained_qp::solve (Vector start) const -{ - /* no hassle if no constraints*/ - if (! cons_.size()) - { - Choleski_decomposition chol (quad_); - return - chol.solve (lin_); - } - else - { - return constraint_solve (start); - } -} - -int -Ineq_constrained_qp::dim () const -{ - return lin_.dim(); -} - - - - -void -Ineq_constrained_qp::assert_solution (Vector sol) const -{ - bool sol_b=true; - - for (int i=0; sol_b && i < cons_.size(); i++) - { - Real R=cons_[i] * sol- consrhs_[i]; - if (R> -EPS) - sol_b = false; - } -} - -void -Ineq_constrained_qp::print() const -{ -#ifndef NPRINT - DOUT << "Ineq_constrained_qp { " << '\n'; - DOUT << "Quad " << quad_; - DOUT << "lin " << lin_ << '\n' - << "const " << const_term_<< '\n'; - for (int i=0; i < cons_.size(); i++) - { - DOUT << "constraint["<= " << consrhs_[i]; - DOUT << '\n'; - } - DOUT << "}\n"; -#endif -} diff --git a/lily/mixed-qp.cc b/lily/mixed-qp.cc deleted file mode 100644 index 25946d444f..0000000000 --- a/lily/mixed-qp.cc +++ /dev/null @@ -1,82 +0,0 @@ -/* - qlp.cc -- implement Mixed_qp - - source file of the GNU LilyPond music typesetter - - (c) 1997--1999 Han-Wen Nienhuys -*/ - -#include "debug.hh" -#include "qlp.hh" - - -void -Mixed_qp::add_equality_cons (Vector , double) -{ - assert (false); -} - -void -Mixed_qp::add_fixed_var (int i, Real r) -{ - eq_cons.push (i); - eq_consrhs.push (r); -} - - -/** - eliminate appropriate variables, until we have a Ineq_constrained_qp - then solve that. - - PRE - cons should be ascending - */ -Vector -Mixed_qp::solve (Vector start) const -{ - if (!dim()) - return Vector (0); - - print(); - Ineq_constrained_qp pure (*this); - - for (int i= eq_cons.size()-1; i>=0; i--) - { - pure.eliminate_var (eq_cons[i], eq_consrhs[i]); - start.del (eq_cons[i]); - } - Vector sol = pure.solve (start); - for (int i= 0; i < eq_cons.size(); i++) - { - sol.insert (eq_consrhs[i],eq_cons[i]); - } - return sol; -} - - -Mixed_qp::Mixed_qp (int n) - : Ineq_constrained_qp (n) -{ -} - -void -Mixed_qp::OK() const -{ -#ifndef NDEBUG - Ineq_constrained_qp::OK(); - assert (eq_consrhs.size() == eq_cons.size ()); -#endif -} - -void -Mixed_qp::print() const -{ -#ifndef NPRINT - Ineq_constrained_qp::print(); - for (int i=0; i < eq_cons.size(); i++) - { - DOUT << "eq cons "< -*/ - -#include "paper-column.hh" -#include "paper-score.hh" -#include "debug.hh" - -void -Paper_column::add_rod (Paper_column * p, Real d) -{ - Direction dir = Direction (sign (p->rank_i () - rank_i ())); - - if (!dir) - { - programming_error ("Must set minimum distance between differing columns."); - return; - } - - for (int i=0; i < minimal_dists_arr_drul_[dir].size (); i++) - { - Column_rod &rod = minimal_dists_arr_drul_[dir][i]; - if (rod.other_l_ == p) - { - rod.distance_f_ = rod.distance_f_ >? d; - return ; - } - } - - Column_rod cr; - cr.distance_f_ = d; - cr.other_l_ = p; - - minimal_dists_arr_drul_[dir].push (cr); -} - -void -Paper_column::add_spring (Paper_column * p, Real d, Real s) -{ - Direction dir = Direction (sign (p->rank_i () - rank_i ())); - - if (!dir) - { - warning (_ ("Must set spring between differing columns")); - return; - } - - for (int i=0; i < spring_arr_drul_[dir].size (); i++) - { - Column_spring &spring = spring_arr_drul_[dir][i]; - if (spring.other_l_ == p) - { - spring.distance_f_ = spring.distance_f_ >? d; - return ; - } - } - - Column_spring cr; - cr.distance_f_ = d; - cr.strength_f_ = s; - cr.other_l_ = p; - - spring_arr_drul_[dir].push (cr); -} - -int -Paper_column::rank_i() const -{ - return rank_i_; -} - -void -Paper_column::set_rank (int i) -{ - rank_i_ = i; -} - -void -Paper_column::do_print() const -{ -#ifndef NPRINT - DOUT << "rank: " << rank_i_ << '\n'; - Direction d = LEFT; - do - { - for (int i=0; i < minimal_dists_arr_drul_[d].size (); i++) - { - minimal_dists_arr_drul_[d][i].print (); - } - for (int i=0; i < spring_arr_drul_[d].size (); i++) - { - spring_arr_drul_[d][i].print (); - } - - } - while ((flip (&d))!=LEFT); - Item::do_print (); -#endif -} - -bool -Paper_column::breakpoint_b() const -{ - return !line_l_; -} - -Paper_column::Paper_column() -{ - set_axes (X_AXIS, X_AXIS); - - line_l_=0; - rank_i_ = -1; -} - -Line_of_score* -Paper_column::line_l() const -{ - return line_l_; -} - - - - -Paper_column* -Paper_column::column_l () const -{ - return (Paper_column*)(this); -} - -/* - ugh. - */ -void -Paper_column::preprocess () -{ - minimal_dists_arr_drul_[LEFT].sort (Column_rod::compare); - minimal_dists_arr_drul_[RIGHT].sort (Column_rod::compare); - spring_arr_drul_[LEFT].sort (Column_spring::compare); - spring_arr_drul_[RIGHT].sort (Column_spring::compare); -} diff --git a/lily/qlpsolve.cc b/lily/qlpsolve.cc deleted file mode 100644 index 3d64b37792..0000000000 --- a/lily/qlpsolve.cc +++ /dev/null @@ -1,169 +0,0 @@ -/* - qlpsolve.cc -- implement Active_constraints, Inactive_iter - - source file of the GNU LilyPond music typesetter - - (c) 1996, 1997--1999 Han-Wen Nienhuys - - TODO: - try fixed point arithmetic, to speed up lily. - */ - -#include "ineq-constrained-qp.hh" -#include "qlpsolve.hh" -#include "debug.hh" -#include "choleski.hh" - -const Real TOL=1e-1; // roughly 1/30 mm - -String -Active_constraints::status() const -{ - String s ("Active|Inactive ["); - for (int i=0; i< active.size(); i++) - { - s += to_str (active[i]) + " "; - } - - s+="| "; - for (int i=0; i< inactive.size(); i++) - { - s += to_str (inactive[i]) + " "; - } - s+="]"; - - return s; -} - -void -Active_constraints::OK() -{ -#ifndef NDEBUG - H.OK(); - A.OK(); - assert (active.size() +inactive.size () == opt->cons_.size ()); - assert (H.dim() == opt->dim ()); - assert (active.size() == A.rows ()); - Array allcons; - - for (int i=0; i < opt->cons_.size(); i++) - allcons.push (0); - for (int i=0; i < active.size(); i++) - { - int j = active[i]; - allcons[j]++; - } - for (int i=0; i < inactive.size(); i++) - { - int j = inactive[i]; - allcons[j]++; - } - for (int i=0; i < allcons.size(); i++) - assert (allcons[i] == 1); -#endif -} - -Vector -Active_constraints::get_lagrange (Vector gradient) -{ - return (A*gradient); -} - -void -Active_constraints::add_constraint (int k) -{ - // add indices - int cidx=inactive[k]; - - Vector a (opt->cons_[cidx]); - // update of matrices - Vector Ha = H*a; - Real aHa = a*Ha; - Vector addrow (Ha.dim()); - bool degenerate = (abs (aHa) < EPS); - - if (degenerate) - { - warning (String ("Active_constraints::add ():") - + _("degenerate constraints")); - DOUT << "Ha = " << Ha.str () << '\n'; - /* - a != 0, so if Ha = O(EPS), then - Ha * aH / aHa = O(EPS^2/EPS) - - if H*a == 0, the constraints are dependent. - */ - degenerate_count_i_ ++; - } - if (!degenerate) - { - active.push (cidx); - inactive.swap (k,inactive.size()-1); - inactive.pop(); - - H -= Matrix (Ha/aHa , Ha); - - addrow=Ha; - addrow /= aHa; - A -= Matrix (A*a, addrow); - A.insert_row (addrow,A.rows()); - } -} - -void -Active_constraints::drop_constraint (int k) -{ - int q=active.size()-1; - - - Vector a (A.row (q)); - if (a.norm() > EPS) - { - // drop indices - inactive.push (active[k]); - active.swap (k,q); - A.swap_rows (k,q); - active.pop(); - /* - - */ - Real aqa = a*opt->quad_*a; - Matrix aaq (a,a/aqa); - H += aaq; - A -= A*opt->quad_*aaq; - - A.delete_row (q); - }else { - degenerate_count_i_ ++; - warning (String ("Active_constraints::drop ():") - + _("degenerate constraints")); - } -} - - -Active_constraints::Active_constraints (Ineq_constrained_qp const *op) - : A(0,op->dim()), - H(op->dim()), - opt (op) -{ - for (int i=0; i < op->cons_.size(); i++) - inactive.push (i); - Choleski_decomposition chol (op->quad_); - - /* - ugh. - */ - H=chol.inverse(); - OK(); - - degenerate_count_i_ = 0; -} - -/** Find the optimum which is in the planes generated by the active - constraints. - */ -Vector -Active_constraints::find_active_optimum (Vector g) -{ - return H*g; -} diff --git a/lily/spring-spacer.cc b/lily/spring-spacer.cc deleted file mode 100644 index 45dcb76a8e..0000000000 --- a/lily/spring-spacer.cc +++ /dev/null @@ -1,432 +0,0 @@ -/* - spring-spacer.cc -- implement Spring_spacer - - source file of the GNU LilyPond music typesetter - - (c) 1996--1999 Han-Wen Nienhuys -*/ - - -#include -#include -#include "killing-cons.tcc" -#include "spring-spacer.hh" -#include "paper-column.hh" -#include "debug.hh" -#include "dimensions.hh" -#include "qlp.hh" -#include "unionfind.hh" -#include "idealspacing.hh" -#include "pointer.tcc" -#include "score-column.hh" -#include "paper-def.hh" -#include "column-x-positions.hh" -#include "main.hh" - -Vector -Spring_spacer::default_solution() const -{ - return try_initial_solution(); -} - -Score_column* -Spring_spacer::scol_l (int i) -{ - return dynamic_cast(cols_[i].pcol_l_); -} - -const Real COLFUDGE=1e-3; -template class P; // ugh. - -bool -Spring_spacer::contains_b (Paper_column const *w) -{ - for (int i=0; i< cols_.size(); i++) - if (cols_[i].pcol_l_ == w) - return true; - return false; -} - - -void -Spring_spacer::OK() const -{ -#ifndef NDEBUG - for (int i = 1; i < cols_.size(); i++) - assert (cols_[i].rank_i_ > cols_[i-1].rank_i_); -#endif -} - -/** - Make sure no unconnected columns happen. - */ -void -Spring_spacer::handle_loose_cols() -{ - Union_find connected (cols_.size()); - Array fixed; - - for (Cons *i = ideal_p_list_; i; i = i->next_) - { - connected.connect (i->car_->cols_drul_[LEFT],i->car_->cols_drul_[RIGHT]); - } - for (int i = 0; i < cols_.size(); i++) - if (cols_[i].fixed_b()) - fixed.push (i); - for (int i=1; i < fixed.size(); i++) - connected.connect (fixed[i-1], fixed[i]); - - /* - If columns do not have spacing information set, we need to supply our own. - */ - Real d = default_space_f_; - for (int i = cols_.size(); i--;) - { - if (! connected.equiv (fixed[0], i)) - { - connected.connect (i-1, i); - connect (i-1, i, d, 1.0); - } - } -} - -bool -Spring_spacer::check_constraints (Vector v) const -{ - int dim=v.dim(); - assert (dim == cols_.size()); - DOUT << "checking " << v; - for (int i=0; i < dim; i++) - { - if (cols_[i].fixed_b() && - abs (cols_[i].fixed_position() - v (i)) > COLFUDGE) - { - DOUT << "Fixpos broken\n"; - return false; - } - Array const &rods (cols_[i].rods_[RIGHT]); - for (int j =0; j < rods.size (); j++) - { - int other =rods[j].other_idx_; - Real diff =v (other) - v (i) ; - if (COLFUDGE +diff < rods[j].distance_f_) - { - DOUT << "i, other_i: " << i << " " << other << '\n'; - DOUT << "dist, minimal = " << diff << " " - << rods[j].distance_f_ << '\n'; - return false; - } - } - - } - return true; -} - -/** try to generate a solution which obeys the min - distances and fixed positions - */ -Vector -Spring_spacer::try_initial_solution() const -{ - Vector v; - if (!try_initial_solution_and_tell (v)) - { - warning (_ ("I'm too fat; call Oprah")); - } - return v; - -} - -bool -Spring_spacer::try_initial_solution_and_tell (Vector &v) const -{ - int dim=cols_.size(); - bool succeeded = true; - Vector initsol (dim); - - assert (cols_[0].fixed_b ()); - DOUT << "fixpos 0 " << cols_[0].fixed_position (); - for (int i=0; i < dim; i++) - { - Real min_x = i ? initsol (i-1) : cols_[0].fixed_position (); - Array const &sr_arr(cols_[i].rods_[LEFT]); - for (int j=0; j < sr_arr.size (); j++) - { - min_x = min_x >? (initsol (sr_arr[j].other_idx_) + sr_arr[j].distance_f_); - } - initsol (i) = min_x; - - if (cols_[i].fixed_b()) - { - initsol (i)=cols_[i].fixed_position(); - if (initsol (i) < min_x ) - { - DOUT << "failing: init, min : " << initsol (i) << " " << min_x << '\n'; - initsol (i) = min_x; - succeeded = false; - } - } - } - v = initsol; - - DOUT << "tried and told solution: " << v; - if (!succeeded) - DOUT << "(failed)\n"; - return succeeded; -} - - - -// generate the matrices -void -Spring_spacer::make_matrices (Matrix &quad, Vector &lin, Real &c) const -{ - quad.fill (0); - lin.fill (0); - c = 0; - - for (Cons *p =ideal_p_list_; p; p = p->next_) - { - Idealspacing *i = p->car_; - int l = i->cols_drul_[LEFT]; - int r = i->cols_drul_[RIGHT]; - - quad (r,r) += i->hooke_f_; - quad (r,l) -= i->hooke_f_; - quad (l,r) -= i->hooke_f_; - quad (l,l) += i->hooke_f_; - - lin (r) -= i->space_f_*i->hooke_f_; - lin (l) += i->space_f_*i->hooke_f_; - - c += sqr (i->space_f_); - } - - if (quad.dim() > 10) - quad.set_band(); - - -} - -void -Spring_spacer::set_fixed_cols (Mixed_qp &qp) const -{ - for (int j=0; j < cols_.size(); j++) - if (cols_[j].fixed_b()) - qp.add_fixed_var (j,cols_[j].fixed_position()); -} - -// put the constraints into the LP problem -void -Spring_spacer::make_constraints (Mixed_qp& lp) const -{ - int dim=cols_.size(); - - for (int j=0; j < dim -1; j++) - { - Array const&rod_arr (cols_[j].rods_[RIGHT]); - for (int i = 0; i < rod_arr.size (); i++) - { - Vector c1(dim); - c1(rod_arr[i].other_idx_)=1.0 ; - c1(j)=-1.0 ; - - lp.add_inequality_cons (c1, rod_arr[i].distance_f_); - } - } -} - - -Real -Spring_spacer::calculate_energy_f (Vector solution) const -{ - Real e = 0.0; - for (Cons*p =ideal_p_list_; p; p = p->next_) - { - Idealspacing * i = p->car_; - e += i->energy_f(solution(i->cols_drul_[RIGHT]) - solution(i->cols_drul_[LEFT])); - } - - return e; -} -void -Spring_spacer::lower_bound_solution (Column_x_positions*positions) const -{ - Mixed_qp lp (cols_.size()); - make_matrices (lp.quad_,lp.lin_, lp.const_term_); - set_fixed_cols (lp); - - Vector start (cols_.size()); - start.fill (0.0); - Vector solution_vec (lp.solve (start)); - for (int i=0; i < solution_vec.dim (); i++) - solution_vec(i) += indent_f_; - - DOUT << "Lower bound sol: " << solution_vec; - positions->energy_f_ = calculate_energy_f (solution_vec); - positions->config_ = solution_vec; - positions->satisfies_constraints_b_ = check_constraints (solution_vec); -} - -Spring_spacer::Spring_spacer () -{ - ideal_p_list_ =0; - energy_normalisation_f_ = 1.0; -} - -void -Spring_spacer::solve (Column_x_positions*positions) const -{ - Vector solution_try; - - bool constraint_satisfaction = try_initial_solution_and_tell (solution_try); - if (constraint_satisfaction) - { - Mixed_qp lp (cols_.size()); - make_matrices (lp.quad_,lp.lin_, lp.const_term_); - make_constraints (lp); - set_fixed_cols (lp); - - - Vector solution_vec (lp.solve (solution_try)); - for (int i=0; i < solution_vec.dim (); i++) - solution_vec(i) += indent_f_; - - - positions->satisfies_constraints_b_ = check_constraints (solution_vec); - if (!positions->satisfies_constraints_b_) - { - warning (_("Solution doesn't satisfy constraints")); - } - - positions->energy_f_ = calculate_energy_f (solution_vec); - positions->config_ = solution_vec; - } - else - { - positions->set_stupid_solution (solution_try); - } - -} - -/** - add one column to the problem. - - TODO: ugh merge with add_columns. -*/ -void -Spring_spacer::add_column (Paper_column *col, bool fixed, Real fixpos) -{ - Column_info c (col,(fixed)? &fixpos : 0); - int this_rank = cols_.size(); - c.rank_i_ = this_rank; - - for (int i=0; i < col->minimal_dists_arr_drul_[LEFT].size (); i++) - { - Column_rod &cr = col->minimal_dists_arr_drul_[LEFT][i]; - int left_idx = cr.other_l_->rank_i () - cols_[0].pcol_l_->rank_i (); - if (left_idx < 0) - continue; - - if (cols_[left_idx].pcol_l_ != cr.other_l_) - continue; - - Spacer_rod l_rod; - l_rod.distance_f_ = cr.distance_f_; - l_rod.other_idx_ = left_idx; - c.rods_[LEFT].push (l_rod); - - Spacer_rod r_rod; - r_rod.distance_f_ = cr.distance_f_; - r_rod.other_idx_ = this_rank; - cols_[left_idx].rods_[RIGHT].push (r_rod); - } - - for (int i=0; i < col->spring_arr_drul_[LEFT].size (); i++) - { - Column_spring &cr = col->spring_arr_drul_[LEFT][i]; - int idx = cr.other_l_->rank_i () - cols_[0].pcol_l_->rank_i (); - if (idx < 0) - continue; - - if (cols_[idx].pcol_l_ != cr.other_l_) - continue; - connect (idx, this_rank, cr.distance_f_, cr.strength_f_); - } - - cols_.push (c); -} - - -void -Spring_spacer::add_columns (Link_array cols) -{ - energy_normalisation_f_ = sqr (line_len_f_); - add_column (cols[0], true, 0.0); - for (int i=1; i< cols.size ()-1; i++) - add_column (cols[i],false,0.0); - - if (line_len_f_ > 0) - add_column (cols.top (), true, line_len_f_); - else - add_column (cols.top (), false, 0); -} - - - -void -Spring_spacer::print() const -{ -#ifndef NPRINT - for (int i=0; i < cols_.size(); i++) - { - DOUT << "col " << i << " "; - cols_[i].print(); - } - - for (Cons *p =ideal_p_list_; p; p = p->next_) - { - p->car_->print(); - } -#endif -} - - -void -Spring_spacer::connect (int i, int j, Real d, Real h) -{ - if (d > 100 CM) - { - programming_error( _f ("Improbable distance: %f point, setting to 10 mm", d)); - d = 1 CM; - } - if(d < 0) - { - programming_error (_ ("Negative distance, setting to 10 mm")); - d = 10 MM; - } - - assert(h >=0); - - Idealspacing * s = new Idealspacing; - - s->cols_drul_[LEFT] = i ; - s->cols_drul_[RIGHT] = j; - s->space_f_ = d; - s->hooke_f_ = h; - - ideal_p_list_ = new Killing_cons (s, ideal_p_list_); -} - -void -Spring_spacer::prepare() -{ - handle_loose_cols(); - print(); -} - - -Spring_spacer::~Spring_spacer() -{ - delete ideal_p_list_; -} -- 2.39.5