Singular value decomposition: Difference between revisions
Added C
(Added Go) |
(Added C) |
||
Line 49:
2. The algorithm should be applicable for general case(<math>m\times n</math>).
=={{header|C}}==
{{libheader|GNU Scientific Library}}
The gsl_linalg_SV_decomp function can decompose any m x n matrix though the example here is for a 2 x 2 matrix.
Requires a C99 or later compiler.
<syntaxhighlight lang="c">#include <stdio.h>
#include <gsl/gsl_linalg.h>
/* Custom function for printing a gsl_matrix in matrix form. */
void gsl_matrix_print(const gsl_matrix *M) {
int rows = M->size1;
int cols = M->size2;
for (int i = 0; i < rows; i++) {
printf("|");
for (int j = 0; j < cols; j++) {
printf("% 12.10f ", gsl_matrix_get(M, i, j));
}
printf("\b|\n");
}
printf("\n");
}
int main(){
double a[] = {3, 0, 4, 5};
gsl_matrix_view A = gsl_matrix_view_array(a, 2, 2);
gsl_matrix *V = gsl_matrix_alloc(2, 2);
gsl_vector *S = gsl_vector_alloc(2);
gsl_vector *work = gsl_vector_alloc(2);
/* V is returned here in untransposed form. */
gsl_linalg_SV_decomp(&A.matrix, V, S, work);
gsl_matrix_transpose(V);
double s[] = {S->data[0], 0, 0, S->data[1]};
gsl_matrix_view SM = gsl_matrix_view_array(s, 2, 2);
printf("U:\n");
gsl_matrix_print(&A.matrix);
printf("S:\n");
gsl_matrix_print(&SM.matrix);
printf("V:\n");
gsl_matrix_print(V);
gsl_matrix_free(V);
gsl_vector_free(S);
gsl_vector_free(work);
return 0;
}</syntaxhighlight>
{{out}}
<pre>
U:
|-0.3162277660 -0.9486832981|
|-0.9486832981 0.3162277660|
S:
| 6.7082039325 0.0000000000|
| 0.0000000000 2.2360679775|
V:
|-0.7071067812 -0.7071067812|
|-0.7071067812 0.7071067812|
</pre>
=={{header|Go}}==
| |||