#11248: SVD matrix decomposition may return a non-invertible "unitary" matrix
------------------------------+---------------------------------------------
Reporter: rbeezer | Owner: jason, was
Type: defect | Status: new
Priority: major | Milestone: sage-4.7
Component: linear algebra | Keywords:
Work_issues: | Upstream: N/A
Reviewer: | Author:
Merged: | Dependencies:
------------------------------+---------------------------------------------
Comment(by jason):
I guess we need to make one more change to the program to print out the
last row of V**H. Here is my patch to the program:
{{{
--- test_zgesdd.c.orig 2011-04-26 11:12:41.000000000 -0500
+++ test_zgesdd.c 2011-04-26 11:13:26.000000000 -0500
@@ -105,14 +105,15 @@
};
/* Executable statements */
printf( " ZGESDD Example Program Results\n" );
+ print_matrix("Original Matrix", m, n, a, LDA);
/* Query and allocate the optimal workspace */
lwork = -1;
- zgesdd( "Singular vectors", &m, &n, a, &lda, s, u, &ldu, vt,
&ldvt, &wkopt,
+ zgesdd( "A", &m, &n, a, &lda, s, u, &ldu, vt, &ldvt, &wkopt,
&lwork, rwork, iwork, &info );
lwork = (int)wkopt.re;
work = (dcomplex*)malloc( lwork*sizeof(dcomplex) );
/* Compute SVD */
- zgesdd( "Singular vectors", &m, &n, a, &lda, s, u, &ldu, vt,
&ldvt, work,
+ zgesdd( "A", &m, &n, a, &lda, s, u, &ldu, vt, &ldvt, work,
&lwork, rwork, iwork, &info );
/* Check for convergence */
if( info > 0 ) {
@@ -124,7 +125,7 @@
/* Print left singular vectors */
print_matrix( "Left singular vectors (stored columnwise)", m, m,
u, ldu );
/* Print right singular vectors */
- print_matrix( "Right singular vectors (stored rowwise)", m, n,
vt, ldvt );
+ print_matrix( "Right singular vectors (stored rowwise)", n, n,
vt, ldvt );
/* Free workspace */
free( (void*)work );
exit( 0 );
}}}
Notice that I changed the `print_matrix` call to print out n rows, instead
of m rows (that's the first argument after the string).
So: for someone that had a problem with numpy, could you run the example
program:
1. Copy
http://software.intel.com/sites/products/documentation/hpc/mkl/lapack/mkl_lapack_examples/zgesdd_ex.c.htm
to a file `test_zgesdd.c`
2. Apply the patch above
3. Compile with `gcc test_zgesdd.c -o test_zgesdd -framework Accelerate`
4. Run: `./test_zgesdd`
Then compare the resulting matrices with your output for the following
input to `sage -python`
{{{
I=1j
import numpy
a=numpy.asarray([( -5.40+ I* 7.40), ( 6.00+ I* 6.38), ( 9.91+ I*
0.16), ( -5.28+ I* -4.16), ( 1.09+ I* 1.55), ( 2.60+ I* 0.07), (
3.98+ I* -5.26), ( 2.03+ I* 1.11), ( 9.88+ I* 1.91), ( 4.92+ I*
6.31), ( -2.11+ I* 7.39), ( -9.81+ I* -8.98)],dtype=complex).reshape(3,4)
numpy.linalg.svd(a)
}}}
--
Ticket URL: <http://trac.sagemath.org/sage_trac/ticket/11248#comment:24>
Sage <http://www.sagemath.org>
Sage: Creating a Viable Open Source Alternative to Magma, Maple, Mathematica,
and MATLAB
--
You received this message because you are subscribed to the Google Groups
"sage-trac" group.
To post to this group, send email to [email protected].
To unsubscribe from this group, send email to
[email protected].
For more options, visit this group at
http://groups.google.com/group/sage-trac?hl=en.
