26 #if defined HAVE_CONFIG_H 30 #include <gsl/gsl_linalg.h> 31 #include <gsl/gsl_eigen.h> 32 #include <gsl/gsl_permute_vector.h> 33 #include <gsl/gsl_sort_vector.h> 34 #include <gsl/gsl_sort.h> 35 #include <gsl/gsl_permutation.h> 42 #include <blitz/tinyvec-et.h> 82 ptrdiff_t
modulo(ptrdiff_t a, ptrdiff_t b);
102 blitz::TinyVector<ptrdiff_t,Dim>
modulo(
103 blitz::TinyVector<ptrdiff_t,Dim>
const &a,
104 blitz::TinyVector<ptrdiff_t,Dim>
const &b);
115 template<
typename DataT,
int Dim>
116 blitz::TinyVector<DataT,Dim+1>
128 template<
typename DataT,
int Dim>
129 blitz::TinyVector<DataT,Dim-1>
149 blitz::TinyVector<double,3>
151 blitz::Array<double,2>
const &trafo);
165 blitz::TinyVector<double,3>
186 blitz::TinyVector<double,3>
188 blitz::Array<double,2>
const &trafo);
202 blitz::TinyVector<double,3>
217 template<
typename DataT>
218 blitz::TinyVector<DataT,3>
219 rotate(blitz::TinyVector<DataT,3>
const &vec,
220 blitz::TinyVector<DataT,3>
const &axis,
234 template<
typename DataT>
235 double dot(blitz::Array<DataT,1>
const &vec1,
236 blitz::Array<DataT,1>
const &vec2);
248 template<
typename DataT,
int Dim1,
int Dim2>
249 blitz::TinyMatrix<DataT,Dim1,Dim2>
251 blitz::TinyVector<DataT,Dim2>
const &vec2);
262 template<
typename DataT,
int Dim1,
int Dim2>
263 blitz::TinyMatrix<DataT,Dim2,Dim1>
264 transpose(blitz::TinyMatrix<DataT,Dim1,Dim2>
const &in);
277 template<
typename DataT>
278 blitz::Array<DataT,2>
transpose(blitz::Array<DataT,2>
const &in);
292 template<
typename MatrixT,
typename BaseT,
int Dim>
293 blitz::TinyVector<BaseT,Dim>
294 mvMult(MatrixT
const &m, blitz::TinyVector<BaseT,Dim>
const &v);
306 template<
typename BaseT,
int Dim>
307 blitz::TinyVector<BaseT,Dim>
308 mvMult(blitz::Array<BaseT,2>
const &m, blitz::TinyVector<BaseT,Dim>
const &v);
322 template<
typename BaseT>
323 blitz::Array<BaseT,1>
324 mvMult(blitz::Array<BaseT,2>
const &m, blitz::Array<BaseT,1>
const &v);
338 template<
typename DataT>
339 blitz::Array<DataT,2>
340 mmMult(blitz::Array<DataT,2>
const &A, blitz::Array<DataT,2>
const &B);
352 template<
typename DataT,
int Dim1,
int Dim2,
int Dim3>
353 blitz::TinyMatrix<DataT,Dim1,Dim3>
354 mmMult(blitz::TinyMatrix<DataT,Dim1,Dim2>
const &A,
355 blitz::TinyMatrix<DataT,Dim2,Dim3>
const &B);
365 void invert(blitz::Array<double,2>
const &A, blitz::Array<double,2> &Ainv);
377 blitz::Array<double,2>
378 invert(blitz::Array<double,2>
const &A);
389 blitz::TinyMatrix<double,Dim,Dim>
390 invert(blitz::TinyMatrix<double,Dim,Dim>
const &A);
402 double determinant(blitz::TinyMatrix<double,Dim,Dim>
const &A);
417 blitz::TinyMatrix<double,Dim,Dim>
const &A,
418 blitz::TinyVector<double,Dim> &lambda,
434 blitz::TinyMatrix<double,Dim,Dim>
const &A,
435 blitz::TinyMatrix<double,Dim,Dim> &U,
436 blitz::TinyVector<double,Dim>
const &lambda);
453 blitz::TinyMatrix<double,Dim,Dim>
const &A,
454 blitz::TinyMatrix<double,Dim,Dim> &U,
455 blitz::TinyVector<double,Dim> &lambda,
467 template<
typename DataT>
484 blitz::Array<double,2>
const &A,
485 blitz::Array<double,2> &U, blitz::Array<double,1> &lambda,
490 #include "ATBLinAlg.icc"
ptrdiff_t modulo(ptrdiff_t a, ptrdiff_t b)
Strictly positive modulo, to compute valid indices for Array access modulo Array size.
blitz::TinyVector< DataT, 3 > rotate(blitz::TinyVector< DataT, 3 > const &vec, blitz::TinyVector< DataT, 3 > const &axis, DataT angle)
Rotate the given 3-D vector by a specified angle angle around the given axis.
void eigenvalueDecompositionRealSymmetric(blitz::TinyMatrix< double, Dim, Dim > const &A, blitz::TinyVector< double, Dim > &lambda, SortingType sort=Descending)
Computation of the eigen value decomposition of real symmetric matrices.
blitz::TinyVector< double, 3 > euclideanToSpherical(blitz::TinyVector< double, 3 > const &eucl, blitz::Array< double, 2 > const &trafo)
Get the spherical (theta, phi, r) coordinates of the given euclidean point.
double determinant(blitz::TinyMatrix< double, Dim, Dim > const &A)
Compute the determinant of matrix A.
Query specific information about different data types.
ATBGSLWrapper.hh provides wrappers and setup routines to use GSL functions with ArrayToolbox data str...
blitz::TinyVector< double, 3 > sphericalToEuclidean(blitz::TinyVector< double, 3 > const &spherical, blitz::Array< double, 2 > const &trafo)
Get the euclidean (z, y, x) coordinates of the given point on the 2-sphere.
void invert(blitz::Array< double, 2 > const &A, blitz::Array< double, 2 > &Ainv)
Compute the inverse of matrix A.
blitz::TinyMatrix< DataT, Dim1, Dim2 > outerProduct(blitz::TinyVector< DataT, Dim1 > const &vec1, blitz::TinyVector< DataT, Dim2 > const &vec2)
Calculation of the outer product of two vectors.
blitz::TinyVector< BaseT, Dim > mvMult(MatrixT const &m, blitz::TinyVector< BaseT, Dim > const &v)
Quadratic Matrix * Vector product (generic).
blitz::TinyMatrix< DataT, Dim2, Dim1 > transpose(blitz::TinyMatrix< DataT, Dim1, Dim2 > const &in)
Compute the transpose of the given Matrix.
double factorialFraction(int n, int k)
Compute the value of expression .
blitz::Array< DataT, 2 > mmMult(blitz::Array< DataT, 2 > const &A, blitz::Array< DataT, 2 > const &B)
Naive matrix multiplication A * B = C.
blitz::TinyVector< DataT, Dim+1 > euclideanToHomogeneous(blitz::TinyVector< DataT, Dim > const &pos)
Convert the given vector to homogeneous coordinates.
blitz::TinyVector< DataT, Dim-1 > homogeneousToEuclidean(blitz::TinyVector< DataT, Dim > const &pos)
Convert the given homogeneous vector to euclidean coordinates.
double dot(blitz::Array< DataT, 1 > const &vec1, blitz::Array< DataT, 1 > const &vec2)
Calculation of the inner product (dot product) of the given vectors.
void computeEigenvectors(blitz::TinyMatrix< double, Dim, Dim > const &A, blitz::TinyMatrix< double, Dim, Dim > &U, blitz::TinyVector< double, Dim > const &lambda)
Computates the eigenvectors of the real symmetric matrix A given the eigenvalues. ...
void sortElements(blitz::TinyVector< DataT, 3 > &v, SortingType sort)
Sort the values of the given TinyVector according to the given SortingType.