synchronization with 2geom library

(bzr r5723)

3 files changed, 481 insertions, 480 deletions

diff --git a/src/2geom/numeric/linear_system.h b/src/2geom/numeric/linear_system.h
index 2b4963998..9cb521eb2 100644
--- a/src/2geom/numeric/linear_system.h
+++ b/src/2geom/numeric/linear_system.h
@@ -1,89 +1,89 @@
-#ifndef _NL_LINEAR_SYSTEM_H_

-#define _NL_LINEAR_SYSTEM_H_

-

-

-#include <cassert>

-

-#include <gsl/gsl_linalg.h>

-

-#include "numeric/matrix.h"

-#include "numeric/vector.h"

-

-

-namespace Geom { namespace NL {

-

-

-class LinearSystem

-{

-public:

-	LinearSystem(Matrix & _matrix, Vector & _vector)

-		: m_matrix(_matrix), m_vector(_vector), m_solution(_matrix.columns())

-	{

-	}

-	

-	const Vector & LU_solve()

-	{

-		assert( matrix().rows() == matrix().columns() 

-				&& matrix().rows() == vector().size() );

-		int s;

-		gsl_permutation * p = gsl_permutation_alloc(matrix().rows());

-		gsl_linalg_LU_decomp (matrix().get_gsl_matrix(), p, &s);

-		gsl_linalg_LU_solve( matrix().get_gsl_matrix(), 

-							 p, 

-				             vector().get_gsl_vector(), 

-				             m_solution.get_gsl_vector()

-				           );

-		gsl_permutation_free(p);

-		return solution();

-	}

-	

-	const Vector & SV_solve()

-	{

-		assert( matrix().rows() >= matrix().columns()

-				&& matrix().rows() == vector().size() );

-		

-		gsl_matrix* U = matrix().get_gsl_matrix();

-		gsl_matrix* V = gsl_matrix_alloc(matrix().columns(), matrix().columns());

-		gsl_vector* S = gsl_vector_alloc(matrix().columns());

-		gsl_vector* work = gsl_vector_alloc(matrix().columns());

-		

-		gsl_linalg_SV_decomp( U, V, S, work );

-		

-		gsl_vector* b = vector().get_gsl_vector();

-		gsl_vector* x = m_solution.get_gsl_vector();

-		

-		gsl_linalg_SV_solve( U, V, S, b, x);

-		

-		gsl_matrix_free(V);

-		gsl_vector_free(S);

-		gsl_vector_free(work);

-		

-		return solution();			  

-	}

-	

-	Matrix & matrix()

-	{

-		return m_matrix;

-	}

-	

-	Vector & vector()

-	{

-		return m_vector;

-	}

-	

-	const Vector & solution() const

-	{

-		return m_solution;

-	}

-	

-private:

-	Matrix & m_matrix;

-	Vector & m_vector;

-	Vector m_solution;

-};

-

-

-} } // end namespaces

-

-

-#endif /*_NL_LINEAR_SYSTEM_H_*/

+#ifndef _NL_LINEAR_SYSTEM_H_
+#define _NL_LINEAR_SYSTEM_H_
+
+
+#include <cassert>
+
+#include <gsl/gsl_linalg.h>
+
+#include "numeric/matrix.h"
+#include "numeric/vector.h"
+
+
+namespace Geom { namespace NL {
+
+
+class LinearSystem
+{
+public:
+	LinearSystem(Matrix & _matrix, Vector & _vector)
+		: m_matrix(_matrix), m_vector(_vector), m_solution(_matrix.columns())
+	{
+	}
+	
+	const Vector & LU_solve()
+	{
+		assert( matrix().rows() == matrix().columns() 
+				&& matrix().rows() == vector().size() );
+		int s;
+		gsl_permutation * p = gsl_permutation_alloc(matrix().rows());
+		gsl_linalg_LU_decomp (matrix().get_gsl_matrix(), p, &s);
+		gsl_linalg_LU_solve( matrix().get_gsl_matrix(), 
+							 p, 
+				             vector().get_gsl_vector(), 
+				             m_solution.get_gsl_vector()
+				           );
+		gsl_permutation_free(p);
+		return solution();
+	}
+	
+	const Vector & SV_solve()
+	{
+		assert( matrix().rows() >= matrix().columns()
+				&& matrix().rows() == vector().size() );
+		
+		gsl_matrix* U = matrix().get_gsl_matrix();
+		gsl_matrix* V = gsl_matrix_alloc(matrix().columns(), matrix().columns());
+		gsl_vector* S = gsl_vector_alloc(matrix().columns());
+		gsl_vector* work = gsl_vector_alloc(matrix().columns());
+		
+		gsl_linalg_SV_decomp( U, V, S, work );
+		
+		gsl_vector* b = vector().get_gsl_vector();
+		gsl_vector* x = m_solution.get_gsl_vector();
+		
+		gsl_linalg_SV_solve( U, V, S, b, x);
+		
+		gsl_matrix_free(V);
+		gsl_vector_free(S);
+		gsl_vector_free(work);
+		
+		return solution();			  
+	}
+	
+	Matrix & matrix()
+	{
+		return m_matrix;
+	}
+	
+	Vector & vector()
+	{
+		return m_vector;
+	}
+	
+	const Vector & solution() const
+	{
+		return m_solution;
+	}
+	
+private:
+	Matrix & m_matrix;
+	Vector & m_vector;
+	Vector m_solution;
+};
+
+
+} } // end namespaces
+
+
+#endif /*_NL_LINEAR_SYSTEM_H_*/
diff --git a/src/2geom/numeric/matrix.h b/src/2geom/numeric/matrix.h
index bdfab75ef..bdd647e53 100644
--- a/src/2geom/numeric/matrix.h
+++ b/src/2geom/numeric/matrix.h
@@ -1,207 +1,207 @@
-#ifndef _NL_MATRIX_H_

-#define _NL_MATRIX_H_

-

-

-#include <cassert>

-#include <utility>

-#include <algorithm>

-

-#include <gsl/gsl_matrix.h>

-

-

-

-namespace Geom { namespace NL {

-

-class Matrix;

-void swap(Matrix & m1, Matrix & m2);

-

-

-class Matrix

-{

-public:

-	// the matrix is not inizialized

-	Matrix(size_t n1, size_t n2)

-		: m_rows(n1), m_columns(n2)

-	{

-		m_matrix = gsl_matrix_alloc(n1, n2);

-	}

-	

-	Matrix(size_t n1, size_t n2, double x)

-		:m_rows(n1), m_columns(n2)

-	{

-		m_matrix = gsl_matrix_alloc(n1, n2);

-		gsl_matrix_set_all(m_matrix, x);

-	}

-	

-	Matrix( Matrix const& _matrix )

-		: m_rows(_matrix.rows()), m_columns(_matrix.columns())

-	{

-		m_matrix = gsl_matrix_alloc(rows(), columns());

-		gsl_matrix_memcpy(m_matrix, _matrix.m_matrix);

-	}

-	

-//	explicit

-//	Matrix( gsl_matrix* m, size_t n1, size_t n2)

-//		: m_rows(n1), m_columns(n2), m_matrix(m)

-//	{	

-//	}

-	

-	Matrix & operator=(Matrix const& _matrix)

-	{

-		assert( rows() == _matrix.rows() && columns() ==  _matrix.columns() );

-		gsl_matrix_memcpy(m_matrix, _matrix.m_matrix);

-		

-		return *this;

-	}

-	

-	virtual ~Matrix()

-	{

-		gsl_matrix_free(m_matrix);

-	}

-	

-	void set_all( double x = 0)

-	{

-		gsl_matrix_set_all(m_matrix, x);

-	}

-	

-	void set_identity()

-	{

-		gsl_matrix_set_identity(m_matrix);

-	}

-	

-	const double & operator() (size_t i, size_t j) const

-	{

-		return *gsl_matrix_const_ptr(m_matrix, i, j);

-	}

-	

-	double & operator() (size_t i, size_t j)

-	{

-		return *gsl_matrix_ptr(m_matrix, i, j);

-	}

-	

-	gsl_matrix* get_gsl_matrix()

-	{

-		return m_matrix;

-	}

-	

-	void swap_rows(size_t i, size_t j)

-	{

-		 gsl_matrix_swap_rows(m_matrix, i, j);

-	}

-	

-	void swap_columns(size_t i, size_t j)

-	{

-		gsl_matrix_swap_columns(m_matrix, i, j);

-	}

-	

-	Matrix & transpose()

-	{

-		gsl_matrix_transpose(m_matrix);

-		return (*this);

-	}

-	

-	Matrix & scale(double x)

-	{

-		gsl_matrix_scale(m_matrix, x);

-		return (*this);

-	}

-	

-	Matrix & translate(double x)

-	{

-		 gsl_matrix_add_constant(m_matrix, x);

-		 return (*this);

-	}

-	

-	Matrix & operator+=(Matrix const& _matrix)

-	{

-		gsl_matrix_add(m_matrix, _matrix.m_matrix);

-		return (*this);

-	}

-	

-	Matrix & operator-=(Matrix const& _matrix)

-	{

-		gsl_matrix_sub(m_matrix, _matrix.m_matrix);

-		return (*this);

-	}

-	

-	bool is_zero() const

-	{

-		return gsl_matrix_isnull(m_matrix);

-	}

-	

-	bool is_positive() const

-	{

-		return gsl_matrix_ispos(m_matrix);

-	}

-	

-	bool is_negative() const

-	{

-		return gsl_matrix_isneg(m_matrix);

-	}

-	

-	bool is_non_negative() const

-	{

-		for ( unsigned int i = 0; i < rows(); ++i )

-		{

-			for ( unsigned int j = 0; j < columns(); ++j )

-			{

-				if ( (*this)(i,j) < 0 ) return false;

-			}

-		}

-		return true;

-	}

-	

-	double max() const

-	{

-		return gsl_matrix_max(m_matrix);

-	}

-	

-	double min() const

-	{

-		return gsl_matrix_min(m_matrix);

-	}

-	

-	std::pair<size_t, size_t>

-	max_index() const

-	{

-		std::pair<size_t, size_t> indices;

-		gsl_matrix_max_index(m_matrix, &(indices.first), &(indices.second));

-		return indices;

-	}

-	

-	std::pair<size_t, size_t>

-	min_index() const

-	{

-		std::pair<size_t, size_t> indices;

-		gsl_matrix_min_index(m_matrix, &(indices.first), &(indices.second));

-		return indices;

-	}

-	

-	friend 

-	void swap(Matrix & m1, Matrix & m2);

-	

-	size_t rows() const

-	{

-		return m_rows;

-	}

-	

-	size_t columns() const

-	{

-		return m_columns;

-	}

-	

-private:

-	size_t m_rows, m_columns;

-	gsl_matrix* m_matrix;

-};

-

-void swap(Matrix & m1, Matrix & m2)

-{

-	assert( m1.rows() == m2.rows() && m1.columns() ==  m2.columns() );

-	std::swap(m1.m_matrix, m2.m_matrix);

-}

-

-

-} } // end namespaces

-

-#endif /*_NL_MATRIX_H_*/

+#ifndef _NL_MATRIX_H_
+#define _NL_MATRIX_H_
+
+
+#include <cassert>
+#include <utility>
+#include <algorithm>
+
+#include <gsl/gsl_matrix.h>
+
+
+
+namespace Geom { namespace NL {
+
+class Matrix;
+void swap(Matrix & m1, Matrix & m2);
+
+
+class Matrix
+{
+public:
+	// the matrix is not inizialized
+	Matrix(size_t n1, size_t n2)
+		: m_rows(n1), m_columns(n2)
+	{
+		m_matrix = gsl_matrix_alloc(n1, n2);
+	}
+	
+	Matrix(size_t n1, size_t n2, double x)
+		:m_rows(n1), m_columns(n2)
+	{
+		m_matrix = gsl_matrix_alloc(n1, n2);
+		gsl_matrix_set_all(m_matrix, x);
+	}
+	
+	Matrix( Matrix const& _matrix )
+		: m_rows(_matrix.rows()), m_columns(_matrix.columns())
+	{
+		m_matrix = gsl_matrix_alloc(rows(), columns());
+		gsl_matrix_memcpy(m_matrix, _matrix.m_matrix);
+	}
+	
+//	explicit
+//	Matrix( gsl_matrix* m, size_t n1, size_t n2)
+//		: m_rows(n1), m_columns(n2), m_matrix(m)
+//	{	
+//	}
+	
+	Matrix & operator=(Matrix const& _matrix)
+	{
+		assert( rows() == _matrix.rows() && columns() ==  _matrix.columns() );
+		gsl_matrix_memcpy(m_matrix, _matrix.m_matrix);
+		
+		return *this;
+	}
+	
+	virtual ~Matrix()
+	{
+		gsl_matrix_free(m_matrix);
+	}
+	
+	void set_all( double x = 0)
+	{
+		gsl_matrix_set_all(m_matrix, x);
+	}
+	
+	void set_identity()
+	{
+		gsl_matrix_set_identity(m_matrix);
+	}
+	
+	const double & operator() (size_t i, size_t j) const
+	{
+		return *gsl_matrix_const_ptr(m_matrix, i, j);
+	}
+	
+	double & operator() (size_t i, size_t j)
+	{
+		return *gsl_matrix_ptr(m_matrix, i, j);
+	}
+	
+	gsl_matrix* get_gsl_matrix()
+	{
+		return m_matrix;
+	}
+	
+	void swap_rows(size_t i, size_t j)
+	{
+		 gsl_matrix_swap_rows(m_matrix, i, j);
+	}
+	
+	void swap_columns(size_t i, size_t j)
+	{
+		gsl_matrix_swap_columns(m_matrix, i, j);
+	}
+	
+	Matrix & transpose()
+	{
+		gsl_matrix_transpose(m_matrix);
+		return (*this);
+	}
+	
+	Matrix & scale(double x)
+	{
+		gsl_matrix_scale(m_matrix, x);
+		return (*this);
+	}
+	
+	Matrix & translate(double x)
+	{
+		 gsl_matrix_add_constant(m_matrix, x);
+		 return (*this);
+	}
+	
+	Matrix & operator+=(Matrix const& _matrix)
+	{
+		gsl_matrix_add(m_matrix, _matrix.m_matrix);
+		return (*this);
+	}
+	
+	Matrix & operator-=(Matrix const& _matrix)
+	{
+		gsl_matrix_sub(m_matrix, _matrix.m_matrix);
+		return (*this);
+	}
+	
+	bool is_zero() const
+	{
+		return gsl_matrix_isnull(m_matrix);
+	}
+	
+	bool is_positive() const
+	{
+		return gsl_matrix_ispos(m_matrix);
+	}
+	
+	bool is_negative() const
+	{
+		return gsl_matrix_isneg(m_matrix);
+	}
+	
+	bool is_non_negative() const
+	{
+		for ( unsigned int i = 0; i < rows(); ++i )
+		{
+			for ( unsigned int j = 0; j < columns(); ++j )
+			{
+				if ( (*this)(i,j) < 0 ) return false;
+			}
+		}
+		return true;
+	}
+	
+	double max() const
+	{
+		return gsl_matrix_max(m_matrix);
+	}
+	
+	double min() const
+	{
+		return gsl_matrix_min(m_matrix);
+	}
+	
+	std::pair<size_t, size_t>
+	max_index() const
+	{
+		std::pair<size_t, size_t> indices;
+		gsl_matrix_max_index(m_matrix, &(indices.first), &(indices.second));
+		return indices;
+	}
+	
+	std::pair<size_t, size_t>
+	min_index() const
+	{
+		std::pair<size_t, size_t> indices;
+		gsl_matrix_min_index(m_matrix, &(indices.first), &(indices.second));
+		return indices;
+	}
+	
+	friend 
+	void swap(Matrix & m1, Matrix & m2);
+	
+	size_t rows() const
+	{
+		return m_rows;
+	}
+	
+	size_t columns() const
+	{
+		return m_columns;
+	}
+	
+private:
+	size_t m_rows, m_columns;
+	gsl_matrix* m_matrix;
+};
+
+void swap(Matrix & m1, Matrix & m2)
+{
+	assert( m1.rows() == m2.rows() && m1.columns() ==  m2.columns() );
+	std::swap(m1.m_matrix, m2.m_matrix);
+}
+
+
+} } // end namespaces
+
+#endif /*_NL_MATRIX_H_*/
diff --git a/src/2geom/numeric/vector.h b/src/2geom/numeric/vector.h
index 136b2cc3f..b25861e22 100644
--- a/src/2geom/numeric/vector.h
+++ b/src/2geom/numeric/vector.h
@@ -1,184 +1,185 @@
-#ifndef _NL_VECTOR_H_

-#define _NL_VECTOR_H_

-

-#include <cassert>

-#include <utility>

-

-#include <gsl/gsl_vector.h>

-

-

-namespace Geom { namespace NL {

-

-class Vector;

-void swap(Vector & v1, Vector & v2);

-

-

-class Vector

-{

-public:

-	Vector(size_t n)

-		: m_size(n)

-	{

-		m_vector = gsl_vector_alloc(n);

-	}

-	

-	Vector(size_t n, double x)

-		: m_size(n)

-	{

-		m_vector = gsl_vector_alloc(n);

-		gsl_vector_set_all(m_vector, x);

-	}

-	

-	// create a vector with n elements all set to zero 

-	// but the i-th that is set to 1

-	Vector(size_t n, size_t i)

-		: m_size(n)

-	{

-		m_vector = gsl_vector_alloc(n);

-		gsl_vector_set_basis(m_vector, i);

-	}

-	

-	Vector(Vector const& _vector)

-		: m_size(_vector.size())

-	{

-		m_vector = gsl_vector_alloc(size());

-		gsl_vector_memcpy(m_vector, _vector.m_vector);

-	}

-	

-	virtual ~Vector()

-	{

-		gsl_vector_free(m_vector);

-	}

-	

-	Vector & operator=(Vector const& _vector)

-	{

-		assert( size() == _vector.size() );

-		gsl_vector_memcpy(m_vector, _vector.m_vector);

-	}

-	

-	void set_all(double x = 0)

-	{

-		gsl_vector_set_all(m_vector, x);

-	}

-

-	void set_basis(size_t i)

-	{

-		gsl_vector_set_basis(m_vector, i);

-	}

-	

-	double const& operator[](size_t i) const

-	{

-		return *gsl_vector_const_ptr(m_vector, i);

-	}

-	

-	double & operator[](size_t i)

-	{

-		return *gsl_vector_ptr(m_vector, i);

-	}

-	

-	gsl_vector* get_gsl_vector()

-	{

-		return m_vector;

-	}

-	

-	void swap_elements(size_t i, size_t j)

-	{

-		gsl_vector_swap_elements(m_vector, i, j);

-	}

-	

-	void reverse()

-	{

-		gsl_vector_reverse(m_vector);

-	}

-	

-	Vector & scale(double x)

-	{

-		gsl_vector_scale(m_vector, x);

-		return (*this);

-	}

-	

-	Vector & translate(double x)

-	{

-		gsl_vector_add_constant(m_vector, x);

-		return (*this);

-	}

-	

-	Vector & operator+=(Vector const& _vector)

-	{

-		gsl_vector_add(m_vector, _vector.m_vector);

-		return (*this);

-	}

-	

-	Vector & operator-=(Vector const& _vector)

-	{

-		gsl_vector_sub(m_vector, _vector.m_vector);

-		return (*this);

-	}

-	

-	bool is_zero() const

-	{

-		return gsl_vector_isnull(m_vector);

-	}

-	

-	bool is_positive() const

-	{

-		return gsl_vector_ispos(m_vector);

-	}

-	

-	bool is_negative() const

-	{

-		return gsl_vector_isneg(m_vector);

-	}

-	

-	bool is_non_negative() const

-	{

-		for ( size_t i = 0; i < size(); ++i )

-		{

-			if ( (*this)[i] < 0 ) return false;

-		}

-		return true;

-	}

-	

-	double max() const

-	{

-		return gsl_vector_max(m_vector);

-	}

-	

-	double min() const

-	{

-		return gsl_vector_min(m_vector);

-	}

-	

-	size_t max_index() const

-	{

-		return gsl_vector_max_index(m_vector);

-	}

-	

-	size_t min_index() const

-	{

-		return gsl_vector_min_index(m_vector);

-	}

-	

-	friend

-	void swap(Vector & v1, Vector & v2);

-	

-	size_t size() const

-	{

-		return m_size;

-	}

-	

-private:

-	size_t m_size;

-	gsl_vector* m_vector;

-};

-

-void swap(Vector & v1, Vector & v2)

-{

-	assert( v1.size() == v2.size() );

-	std::swap(v1.m_vector, v2.m_vector);

-}

-

-} } // end namespaces

-

-

-#endif /*_NL_VECTOR_H_*/

+#ifndef _NL_VECTOR_H_
+#define _NL_VECTOR_H_
+
+#include <cassert>
+#include <utility>
+
+#include <gsl/gsl_vector.h>
+
+
+namespace Geom { namespace NL {
+
+class Vector;
+void swap(Vector & v1, Vector & v2);
+
+
+class Vector
+{
+public:
+	Vector(size_t n)
+		: m_size(n)
+	{
+		m_vector = gsl_vector_alloc(n);
+	}
+	
+	Vector(size_t n, double x)
+		: m_size(n)
+	{
+		m_vector = gsl_vector_alloc(n);
+		gsl_vector_set_all(m_vector, x);
+	}
+	
+	// create a vector with n elements all set to zero 
+	// but the i-th that is set to 1
+	Vector(size_t n, size_t i)
+		: m_size(n)
+	{
+		m_vector = gsl_vector_alloc(n);
+		gsl_vector_set_basis(m_vector, i);
+	}
+	
+	Vector(Vector const& _vector)
+		: m_size(_vector.size())
+	{
+		m_vector = gsl_vector_alloc(size());
+		gsl_vector_memcpy(m_vector, _vector.m_vector);
+	}
+	
+	virtual ~Vector()
+	{
+		gsl_vector_free(m_vector);
+	}
+	
+	Vector & operator=(Vector const& _vector)
+	{
+		assert( size() == _vector.size() );
+		gsl_vector_memcpy(m_vector, _vector.m_vector);
+		return (*this);
+	}
+	
+	void set_all(double x = 0)
+	{
+		gsl_vector_set_all(m_vector, x);
+	}
+
+	void set_basis(size_t i)
+	{
+		gsl_vector_set_basis(m_vector, i);
+	}
+	
+	double const& operator[](size_t i) const
+	{
+		return *gsl_vector_const_ptr(m_vector, i);
+	}
+	
+	double & operator[](size_t i)
+	{
+		return *gsl_vector_ptr(m_vector, i);
+	}
+	
+	gsl_vector* get_gsl_vector()
+	{
+		return m_vector;
+	}
+	
+	void swap_elements(size_t i, size_t j)
+	{
+		gsl_vector_swap_elements(m_vector, i, j);
+	}
+	
+	void reverse()
+	{
+		gsl_vector_reverse(m_vector);
+	}
+	
+	Vector & scale(double x)
+	{
+		gsl_vector_scale(m_vector, x);
+		return (*this);
+	}
+	
+	Vector & translate(double x)
+	{
+		gsl_vector_add_constant(m_vector, x);
+		return (*this);
+	}
+	
+	Vector & operator+=(Vector const& _vector)
+	{
+		gsl_vector_add(m_vector, _vector.m_vector);
+		return (*this);
+	}
+	
+	Vector & operator-=(Vector const& _vector)
+	{
+		gsl_vector_sub(m_vector, _vector.m_vector);
+		return (*this);
+	}
+	
+	bool is_zero() const
+	{
+		return gsl_vector_isnull(m_vector);
+	}
+	
+	bool is_positive() const
+	{
+		return gsl_vector_ispos(m_vector);
+	}
+	
+	bool is_negative() const
+	{
+		return gsl_vector_isneg(m_vector);
+	}
+	
+	bool is_non_negative() const
+	{
+		for ( size_t i = 0; i < size(); ++i )
+		{
+			if ( (*this)[i] < 0 ) return false;
+		}
+		return true;
+	}
+	
+	double max() const
+	{
+		return gsl_vector_max(m_vector);
+	}
+	
+	double min() const
+	{
+		return gsl_vector_min(m_vector);
+	}
+	
+	size_t max_index() const
+	{
+		return gsl_vector_max_index(m_vector);
+	}
+	
+	size_t min_index() const
+	{
+		return gsl_vector_min_index(m_vector);
+	}
+	
+	friend
+	void swap(Vector & v1, Vector & v2);
+	
+	size_t size() const
+	{
+		return m_size;
+	}
+	
+private:
+	size_t m_size;
+	gsl_vector* m_vector;
+};
+
+void swap(Vector & v1, Vector & v2)
+{
+	assert( v1.size() == v2.size() );
+	std::swap(v1.m_vector, v2.m_vector);
+}
+
+} } // end namespaces
+
+
+#endif /*_NL_VECTOR_H_*/