From: "Jim Love" <Jim.Love@asml.com>
To: <gsl-discuss@sources.redhat.com>
Subject: Problem with Singular Value Decomposition Algorithm
Date: Wed, 19 Dec 2001 13:20:00 -0000 [thread overview]
Message-ID: <sb9f3956.087@wiltonhub.svg.com> (raw)
I have downloaded the latest beta release and the SVD algorithm produces the wrong answers. It appears that columns are swapped in the output and possibly a sign problem.
Here was my test array:
1 1 0.975
1 -1 0.975
-1 -1 -0.925
-1 1 -1.025
The correct answer for the S vector is: 2.7940 2.0000 0.0358
The gls output was: 2.0000 2.7940 0.0358
The Correct Q matrix is:
-0.7155 0.0256 -0.6981
0.0183 0.9997 0.0179
-0.6983 -0.0000 0.7158
The gls output was:
-0.025633 -0.715538 -0.698103
-0.999671 0.018347 0.017900
-0.000000 -0.698332 0.715774
A similar problem was seen in the U matrix.
Any ideas? Is this caused by my implementation or is it a real bug?
Here is the source code to the program I made to test the LIB:
#include "stdafx.h"
#include <math.h>
#include <stdio.h>
#include <C:\Program Files\gsl-0-9-2\gsl\gsl_matrix.h>
#include <C:\Program Files\gsl-0-9-2\gsl\gsl_blas.h>
#include <C:\Program Files\gsl-0-9-2\gsl\gsl_linalg.h>
int main(void)
{
gsl_matrix *a = gsl_matrix_alloc (4, 3);
gsl_matrix *u = gsl_matrix_alloc (4, 4);
gsl_matrix *v = gsl_matrix_alloc (3, 3);
gsl_vector *s = gsl_vector_alloc (3);
gsl_vector *work = gsl_vector_alloc (3);
gsl_matrix_set_zero(a);
gsl_matrix_set (a, 0, 0, 1);
gsl_matrix_set (a, 0, 1, 1);
gsl_matrix_set (a, 0, 2, 0.975);
gsl_matrix_set (a, 1, 0, 1);
gsl_matrix_set (a, 1, 1, -1);
gsl_matrix_set (a, 1, 2, 0.975);
gsl_matrix_set (a, 2, 0, -1);
gsl_matrix_set (a, 2, 1, -1);
gsl_matrix_set (a, 2, 2, -0.925);
gsl_matrix_set (a, 3, 0, -1);
gsl_matrix_set (a, 3, 1, 1);
gsl_matrix_set (a, 3, 2, -1.025);
printf ("%f %f %f \n",gsl_matrix_get (a, 0, 0), gsl_matrix_get (a, 0, 1), gsl_matrix_get (a, 0, 2));
printf ("%f %f %f \n",gsl_matrix_get (a, 1, 0), gsl_matrix_get (a, 1, 1), gsl_matrix_get (a, 1, 2));
printf ("%f %f %f \n",gsl_matrix_get (a, 2, 0), gsl_matrix_get (a, 2, 1), gsl_matrix_get (a, 2, 2));
printf ("%f %f %f \n",gsl_matrix_get (a, 3, 0), gsl_matrix_get (a, 3, 1), gsl_matrix_get (a, 3, 2));
gsl_linalg_SV_decomp (a, v, s, work);
printf ("\n");
printf ("%f %f %f \n",gsl_matrix_get (a, 0, 0), gsl_matrix_get (a, 0, 1), gsl_matrix_get (a, 0, 2));
printf ("%f %f %f \n",gsl_matrix_get (a, 1, 0), gsl_matrix_get (a, 1, 1), gsl_matrix_get (a, 1, 2));
printf ("%f %f %f \n",gsl_matrix_get (a, 2, 0), gsl_matrix_get (a, 2, 1), gsl_matrix_get (a, 2, 2));
printf ("%f %f %f \n",gsl_matrix_get (a, 3, 0), gsl_matrix_get (a, 3, 1), gsl_matrix_get (a, 3, 2));
printf ("\n");
printf ("%f %f %f \n",gsl_vector_get (s, 0), gsl_vector_get (s, 1), gsl_vector_get (s, 2));
printf ("\n");
printf ("%f %f %f \n",gsl_matrix_get (v, 0, 0), gsl_matrix_get (v, 0, 1), gsl_matrix_get (v, 0, 2));
printf ("%f %f %f \n",gsl_matrix_get (v, 1, 0), gsl_matrix_get (v, 1, 1), gsl_matrix_get (v, 1, 2));
printf ("%f %f %f \n",gsl_matrix_get (v, 2, 0), gsl_matrix_get (v, 2, 1), gsl_matrix_get (v, 2, 2));
gsl_matrix_free (a);
gsl_matrix_free (u);
gsl_matrix_free (v);
gsl_vector_free (work);
gsl_vector_free (s);
return 0;
}
James A. Love
X4477
Pager 1-800-286-1188 Pin# 400659
next reply other threads:[~2001-12-19 13:20 UTC|newest]
Thread overview: 15+ messages / expand[flat|nested] mbox.gz Atom feed top
2001-12-19 13:20 Jim Love [this message]
2001-12-19 13:20 ` James Theiler
-- strict thread matches above, loose matches on Subject: below --
2001-12-19 13:20 Jim Love
2001-12-19 13:20 ` James Theiler
2001-12-19 13:20 ` Brian Gough
2001-12-19 13:20 ` Peter Schmitteckert
2001-12-19 13:20 ` Brian Gough
2001-12-19 13:20 Jim Love
2001-12-19 13:20 Jim Love
2001-12-19 13:20 ` Brian Gough
2001-12-19 13:20 ` Brian Gough
2001-12-19 13:20 Jim Love
2001-12-19 13:20 ` Brian Gough
[not found] <sba09ba1.005@wiltonhub.svg.com>
2001-12-19 13:20 ` Brian Gough
2001-12-19 13:20 Jim Love
Reply instructions:
You may reply publicly to this message via plain-text email
using any one of the following methods:
* Save the following mbox file, import it into your mail client,
and reply-to-all from there: mbox
Avoid top-posting and favor interleaved quoting:
https://en.wikipedia.org/wiki/Posting_style#Interleaved_style
* Reply using the --to, --cc, and --in-reply-to
switches of git-send-email(1):
git send-email \
--in-reply-to=sb9f3956.087@wiltonhub.svg.com \
--to=jim.love@asml.com \
--cc=gsl-discuss@sources.redhat.com \
/path/to/YOUR_REPLY
https://kernel.org/pub/software/scm/git/docs/git-send-email.html
* If your mail client supports setting the In-Reply-To header
via mailto: links, try the mailto: link
Be sure your reply has a Subject: header at the top and a blank line
before the message body.
This is a public inbox, see mirroring instructions
for how to clone and mirror all data and code used for this inbox;
as well as URLs for read-only IMAP folder(s) and NNTP newsgroup(s).