Reputation: 34
I am looking through GSL functions to calculate Z*Z^T, where Z is n*1 column vector, but I could not find any fit function, every help is much appreciated.
Upvotes: 1
Views: 594
Reputation: 347
GSL supports BLAS (basic linear algebra subprograms), see [http://www.gnu.org/software/gsl/manual/html_node/GSL-BLAS-Interface.html][1].
The functions are classified by the complexity of the operation:
Most functions come in different versions for float, double and complex numbers. Your operation is basically an outer product of the vector Z with itself. You can initialize the vector as a column vector (here double precision numbers):
gsl_matrix * Z = gsl_matrix_calloc (n,1);
and then use the BLAS function gsl_blas_dgemm to compute Z * Z^T. The first arguments of this function determine, whether or not the input matrices should be transposed before the matrix multiplication:
gsl_blas_dgemm (CblasNoTrans, CblasTrans, 1.0, Z, Z, 0.0, C);
Here's a working test program (you may need to link it against gsl and blas):
#include <gsl/gsl_matrix.h>
#include <gsl/gsl_blas.h>
int main(int argc, char ** argv)
{
size_t n = 4;
gsl_matrix * Z = gsl_matrix_calloc (n,1);
gsl_matrix * C = gsl_matrix_calloc (n,n);
gsl_matrix_set(Z,0,0,1);
gsl_matrix_set(Z,1,0,2);
gsl_matrix_set(Z,2,0,0);
gsl_matrix_set(Z,3,0,1);
gsl_blas_dgemm (CblasNoTrans,
CblasTrans, 1.0, Z, Z, 0.0, C);
int i,j;
for (i = 0; i < n; i++)
{
for (j = 0; j < n; j++)
{
printf ("%g\t", gsl_matrix_get (C, i, j));
}
printf("\n");
}
gsl_matrix_free(Z);
gsl_matrix_free(C);
return 0;
}
Upvotes: 1