Skip to content
Snippets Groups Projects
Commit 14ad2448 authored by Christophe Geuzaine's avatar Christophe Geuzaine
Browse files

fix compile without gsl

parent bb8e1bce
No related branches found
No related tags found
No related merge requests found
# $Id: Makefile,v 1.46 2008-02-17 08:47:56 geuzaine Exp $ # $Id: Makefile,v 1.47 2008-02-21 12:42:40 geuzaine Exp $
# #
# Copyright (C) 1997-2008 C. Geuzaine, J.-F. Remacle # Copyright (C) 1997-2008 C. Geuzaine, J.-F. Remacle
# #
...@@ -71,8 +71,9 @@ Box.o: Box.cpp ../Geo/GModel.h ../Geo/GVertex.h ../Geo/GEntity.h \ ...@@ -71,8 +71,9 @@ Box.o: Box.cpp ../Geo/GModel.h ../Geo/GVertex.h ../Geo/GEntity.h \
../Post/PView.h ../Post/PViewData.h ../Post/PViewOptions.h \ ../Post/PView.h ../Post/PViewData.h ../Post/PViewOptions.h \
../Post/ColorTable.h ../Post/PViewDataList.h ../Post/PViewData.h \ ../Post/ColorTable.h ../Post/PViewDataList.h ../Post/PViewData.h \
../Post/AdaptiveViews.h ../Common/GmshMatrix.h ../Mesh/Field.h \ ../Post/AdaptiveViews.h ../Common/GmshMatrix.h ../Mesh/Field.h \
../Geo/Geo.h ../Common/GmshDefines.h ../Geo/gmshSurface.h ../Geo/Pair.h \ ../contrib/ANN/include/ANN/ANN.h ../Geo/Geo.h ../Common/GmshDefines.h \
../Geo/Range.h ../Geo/SPoint2.h ../Geo/SPoint3.h ../Geo/SVector3.h \ ../Geo/gmshSurface.h ../Geo/Pair.h ../Geo/Range.h ../Geo/SPoint2.h \
../Geo/SBoundingBox3d.h ../Geo/SPoint2.h ../Geo/ExtrudeParams.h \ ../Geo/SPoint3.h ../Geo/SVector3.h ../Geo/SBoundingBox3d.h \
../Common/SmoothData.h ../Post/OctreePost.h ../Common/Octree.h \ ../Geo/SPoint2.h ../Geo/ExtrudeParams.h ../Common/SmoothData.h \
../Common/OctreeInternals.h ../Mesh/BackgroundMesh.h ../Post/OctreePost.h ../Common/Octree.h ../Common/OctreeInternals.h \
../Mesh/BackgroundMesh.h
...@@ -30,17 +30,15 @@ public: ...@@ -30,17 +30,15 @@ public:
inline int size() const { return r; } inline int size() const { return r; }
SCALAR *data; SCALAR *data;
~Gmsh_Vector() { delete [] data; } ~Gmsh_Vector() { delete [] data; }
Gmsh_Vector(int R) Gmsh_Vector(int R) : r(R)
: r(R)
{ {
data = new SCALAR[r]; data = new SCALAR[r];
scal(0); scale(0);
} }
Gmsh_Vector(const Gmsh_Vector<SCALAR> &other) Gmsh_Vector(const Gmsh_Vector<SCALAR> &other) : r(other.r)
: r(other.r)
{ {
data = new double[r]; data = new double[r];
copy(other.data); for (int i = 0; i < r; ++i) data[i] = other.data[i];
} }
inline SCALAR operator () (int i) const inline SCALAR operator () (int i) const
{ {
...@@ -50,22 +48,16 @@ public: ...@@ -50,22 +48,16 @@ public:
{ {
return data[i]; return data[i];
} }
inline SCALAR operator *(const Gmsh_Vector<SCALAR> &other) inline double norm()
{ {
throw; double n = 0.;
for(int i = 0; i < r; ++i) n += data[i] * data[i];
return sqrt(n);
} }
inline void scal(const SCALAR s) inline void scale(const SCALAR s)
{ {
for (int i = 0; i < r; ++i) data[i] *= s; for (int i = 0; i < r; ++i) data[i] *= s;
} }
inline void copy(const SCALAR **other)
{
for (int i = 0; i < r; ++i) data[i] = other.data[i];
}
inline void lu_solve (const Gmsh_Vector<SCALAR> &rhs, Gmsh_Vector<SCALAR> &result)
{
throw;
}
}; };
template <class SCALAR> template <class SCALAR>
...@@ -78,17 +70,20 @@ public: ...@@ -78,17 +70,20 @@ public:
inline int size2() const { return c; } inline int size2() const { return c; }
SCALAR *data; SCALAR *data;
~Gmsh_Matrix() { delete [] data; } ~Gmsh_Matrix() { delete [] data; }
Gmsh_Matrix(int R,int C) Gmsh_Matrix(int R,int C) : r(R), c(C)
: r(R), c(C)
{ {
data = new SCALAR[r * c]; data = new SCALAR[r * c];
scal(0); scale(0);
} }
Gmsh_Matrix(const Gmsh_Matrix<SCALAR> &other) Gmsh_Matrix(const Gmsh_Matrix<SCALAR> &other) : r(other.r), c(other.c)
: r(other.r), c(other.c)
{ {
data = new double[r * c]; data = new double[r * c];
copy(other.data); memcpy(other);
}
Gmsh_Matrix() : r(0), c(0), data(0) {}
void memcpy(const Gmsh_Matrix &other)
{
for (int i = 0; i < r * c; ++i) data[i] = other.data[i];
} }
inline SCALAR operator () (int i, int j) const inline SCALAR operator () (int i, int j) const
{ {
...@@ -98,47 +93,50 @@ public: ...@@ -98,47 +93,50 @@ public:
{ {
return data[i + r * j]; return data[i + r * j];
} }
inline Gmsh_Matrix operator *(const Gmsh_Matrix<SCALAR> &other) inline void mult(const Gmsh_Matrix<SCALAR> &x, const Gmsh_Matrix<SCALAR> &b)
{ {
throw; throw;
} }
inline void scal(const SCALAR s) inline void mult(const Gmsh_Vector<SCALAR> &x, Gmsh_Vector<SCALAR> &b)
{ {
for (int i = 0; i < r * c; ++i) data[i] *= s; throw;
} }
inline void copy(const SCALAR **other) inline void set_all(const double &m)
{ {
for (int i = 0; i < r * c; ++i) data[i] = other.data[i]; throw;
} }
inline void mult(const Gmsh_Matrix<SCALAR> &x, const Gmsh_Matrix<SCALAR> &b) inline void lu_solve(const Gmsh_Vector<SCALAR> &rhs, Gmsh_Vector<SCALAR> &result)
{ {
throw; throw;
} }
inline void mult (const Gmsh_Vector<SCALAR> &x, Gmsh_Vector<SCALAR> &b) Gmsh_Matrix cofactor(int i, int j) const
{ {
throw; throw;
} }
inline void least_squares(const Gmsh_Vector<SCALAR> &rhs, Gmsh_Vector<SCALAR> &result) inline void invert()
{ {
throw; throw;
} }
inline void lu_solve (const Gmsh_Vector<SCALAR> &rhs, Gmsh_Vector<SCALAR> &result) double determinant() const
{ {
throw; throw;
} }
Gmsh_Matrix cofactor(int i, int j) const inline Gmsh_Matrix touchSubmatrix(int i0, int ni, int j0, int nj)
{ {
throw; throw;
}
inline void scale(const SCALAR s)
{
for (int i = 0; i < r * c; ++i) data[i] *= s;
} }
inline void invert() inline void add(const double &a)
{ {
throw; throw;
} }
double determinant() const inline void add(const Gmsh_Matrix &m)
{ {
throw; throw;
} }
}; };
#ifdef HAVE_GSL #ifdef HAVE_GSL
...@@ -173,7 +171,8 @@ public: ...@@ -173,7 +171,8 @@ public:
{ {
return *gsl_vector_ptr(data, i); return *gsl_vector_ptr(data, i);
} }
inline double norm(){ inline double norm()
{
return gsl_blas_dnrm2(data); return gsl_blas_dnrm2(data);
} }
inline void scale(const double &y) inline void scale(const double &y)
...@@ -187,22 +186,20 @@ class GSL_Matrix ...@@ -187,22 +186,20 @@ class GSL_Matrix
{ {
private: private:
gsl_matrix_view view; gsl_matrix_view view;
int r, c;
public: public:
inline size_t size1() const { return data->size1; } inline size_t size1() const { return data->size1; }
inline size_t size2() const { return data->size2; } inline size_t size2() const { return data->size2; }
gsl_matrix *data; gsl_matrix *data;
GSL_Matrix(gsl_matrix_view _data) : view(_data),data(&view.matrix) {} GSL_Matrix(gsl_matrix_view _data) : view(_data), data(&view.matrix) {}
GSL_Matrix(size_t R, size_t C) { data = gsl_matrix_calloc (R, C); } GSL_Matrix(size_t R, size_t C) { data = gsl_matrix_calloc(R, C); }
GSL_Matrix() : data(0) {} GSL_Matrix() : r(0), c(0), data(0) {}
GSL_Matrix(const GSL_Matrix&other) : data(0) GSL_Matrix(const GSL_Matrix &other) : data(0)
{ {
if(data) gsl_matrix_free(data); if(data) gsl_matrix_free(data);
data = gsl_matrix_calloc(other.data->size1, other.data->size2); data = gsl_matrix_calloc(other.data->size1, other.data->size2);
gsl_matrix_memcpy(data, other.data); gsl_matrix_memcpy(data, other.data);
} }
virtual ~GSL_Matrix() { if(data && data->owner == 1) gsl_matrix_free(data); } virtual ~GSL_Matrix() { if(data && data->owner == 1) gsl_matrix_free(data); }
GSL_Matrix & operator = (const GSL_Matrix&other) GSL_Matrix & operator = (const GSL_Matrix&other)
{ {
if (&other != this){ if (&other != this){
...@@ -232,19 +229,6 @@ public: ...@@ -232,19 +229,6 @@ public:
{ {
gsl_matrix_set_all(data, m); gsl_matrix_set_all(data, m);
} }
inline void least_squares (const GSL_Vector &rhs, GSL_Vector &result)
{
assert(r > c);
assert(rhs.size() == r);
assert(result.size() == c);
GSL_Matrix *ls = new GSL_Matrix(c, c);
GSL_Vector *ls_rhs = new GSL_Vector(c);
gsl_blas_dgemm(CblasTrans, CblasNoTrans, 1.0, data, data, 1.0, ls->data);
gsl_blas_dgemv(CblasTrans, 1.0, data, rhs.data, 1.0, ls_rhs->data);
ls->lu_solve (*ls_rhs,result);
delete ls;
delete ls_rhs;
}
inline void lu_solve(const GSL_Vector &rhs, GSL_Vector &result) inline void lu_solve(const GSL_Vector &rhs, GSL_Vector &result)
{ {
int s; int s;
...@@ -262,7 +246,7 @@ public: ...@@ -262,7 +246,7 @@ public:
gsl_linalg_LU_invert(data, p, data_inv) ; gsl_linalg_LU_invert(data, p, data_inv) ;
gsl_matrix_memcpy(data, data_inv); gsl_matrix_memcpy(data, data_inv);
gsl_matrix_free(data_inv); gsl_matrix_free(data_inv);
gsl_permutation_free (p); gsl_permutation_free(p);
} }
inline bool invertSecure(double &det) inline bool invertSecure(double &det)
{ {
...@@ -307,7 +291,7 @@ public: ...@@ -307,7 +291,7 @@ public:
} }
return cof; return cof;
} }
inline void mult(const GSL_Vector &x, GSL_Vector &b) inline void mult(const GSL_Vector &x, GSL_Vector &b)
{ {
gsl_blas_dgemv(CblasNoTrans, 1.0, data, x.data, 1.0, b.data); gsl_blas_dgemv(CblasNoTrans, 1.0, data, x.data, 1.0, b.data);
} }
...@@ -333,12 +317,18 @@ public: ...@@ -333,12 +317,18 @@ public:
gsl_matrix_add (data, m.data); gsl_matrix_add (data, m.data);
} }
}; };
typedef GSL_Matrix Double_Matrix; typedef GSL_Matrix Double_Matrix;
typedef GSL_Vector Double_Vector; typedef GSL_Vector Double_Vector;
#else #else
typedef Gmsh_Matrix<double> Double_Matrix; typedef Gmsh_Matrix<double> Double_Matrix;
typedef Gmsh_Vector<double> Double_Vector; typedef Gmsh_Vector<double> Double_Vector;
#endif #endif
typedef Gmsh_Matrix<int> Int_Matrix; typedef Gmsh_Matrix<int> Int_Matrix;
typedef Gmsh_Vector<int> Int_Vector; typedef Gmsh_Vector<int> Int_Vector;
#endif #endif
lc = 0.01;
Point(1) = {0, 0, 0, lc/5};
Point(2) = {.1, 0, 0, lc/5};
Point(3) = {.1, .3, 0, lc};
Point(4) = {0, .3, 0, lc};
Line(1) = {1,2} ;
Line(2) = {3,2} ;
Line(3) = {3,4} ;
Line(4) = {4,1} ;
Line Loop(5) = {4,1,-2,3} ;
Plane Surface(6) = {5} ;
Recombine Surface{6};
0% Loading or .
You are about to add 0 people to the discussion. Proceed with caution.
Please register or to comment