Hi Mark,
Bug 4751232 - GMatrix SVD() computation is incorrect
is submitted for investigation.
Thanks for your bug report.
- Kelvin
----------------------
Java 3D Team
Sun Microsystems Inc.
Mark Whitehorn wrote:
> Kelvin and Josh,
>
> I believe there is a problem with the SVD implementation in javax.vecmath.
> The following code demonstrates that the vecmath decomposition U W V is not
> the correct decomposition of the original matrix. (after having problems
> with the GMatrix implementation, I just switched to using the JNL package)
>
> Mark
>
> ----------------------------------------------------------------------------
> -
>
> package svdTest;
>
> import javax.vecmath.GMatrix;
>
> import VisualNumerics.math.DoubleSVD;
>
> public class SVDtest {
>
> //Main method
> public static void main(String[] args) {
> final int N=3;
> // generate NxN array of random values
> double[][] rv2 = new double[N][N];
> for (int i=0; i<N; i++)
> for (int j=0; j<N; j++)
> rv2[i][j] = Math.random();
>
> // construct a GMatrix from the NxN array
> GMatrix Q = new GMatrix(N, N);
> for (int i=0; i<N; i++) Q.setRow(i, rv2[i]);
>
> // invoke the JNL SVD routine
> final int M = 10000;
> DoubleSVD[] qsvda = new DoubleSVD[M];
> long startTime = System.currentTimeMillis();
> for (int i=0; i<M; i++) qsvda[i] = new DoubleSVD(rv2);
> long duration = System.currentTimeMillis() - startTime;
> double interval = (double)duration / M;
> System.out.println("JNL SVD computation time: " + interval + " msec");
> DoubleSVD qsvd = qsvda[0];
> try {
> System.out.println("JNL rank: " + qsvd.rank());
> } catch (Exception e) {e.printStackTrace();};
>
> // construct result matrices for SVD
> GMatrix U = new GMatrix(N, N);
> GMatrix S = new GMatrix(N, N);
> GMatrix V = new GMatrix(N, N);
>
> // test results by computing X = U*S*transpose(V)
> for (int i=0; i<N; i++) U.setRow(i, qsvd.U()[i]);
> for (int i=0; i<N; i++) S.setElement(i, i, qsvd.S()[i]);
> for (int i=0; i<N; i++) V.setRow(i, qsvd.V()[i]);
>
> GMatrix X = new GMatrix(U);
> X.mul(S);
> V.transpose();
> X.mul(V);
>
> // test X = Q
> X.sub(Q);
> double ssq = 0;
> for (int i=0; i<N; i++) {
> for (int j=0; j<N; j++) {
> ssq += Math.pow(X.getElement(i,j), 2);
> }
> }
> System.out.println("JNL error \n" + X);
> System.out.println("sum of squared errors: " + ssq);
>
>
> // repeat test using GMatrix SVD routine
> int rank = 0;
> startTime = System.currentTimeMillis();
> for (int i=0; i<M; i++) rank = Q.SVD(U, S, V);
> duration = System.currentTimeMillis() - startTime;
> interval = (double)duration / M;
> System.out.println("GMatrix rank: " + rank);
> System.out.println("vecmath SVD computation time: " + interval + "
> msec");
>
> X = new GMatrix(U);
> X.mul(S);
> V.transpose();
> X.mul(V);
> X.sub(Q);
> ssq = 0;
> for (int i=0; i<N; i++) {
> for (int j=0; j<N; j++) {
> ssq += Math.pow(X.getElement(i,j), 2);
> }
> }
> System.out.println("GMatrix error \n" + X);
> System.out.println("sum of squared errors: " + ssq);
> }
> }
>
> ----------------------------------------------------------------------------
> -
>
>
> C:\j2sdk1.4.0-rc\bin\javaw -classpath
> "C:\local\markw\markwJava\jbprojects\svdTest\classes;C:\j2sdk1.4.0-rc\demo\j
> fc\Java2D\Java2Demo.jar;C:\j2sdk1.4.0-rc\demo\plugin\jfc\Java2D\Java2Demo.ja
> r;C:\j2sdk1.4.0-rc\jre\lib\charsets.jar;C:\j2sdk1.4.0-rc\jre\lib\ext\dnsns.j
> ar;C:\j2sdk1.4.0-rc\jre\lib\ext\j3daudio.jar;C:\j2sdk1.4.0-rc\jre\lib\ext\j3
> dcore.jar;C:\j2sdk1.4.0-rc\jre\lib\ext\j3dtree.jar;C:\j2sdk1.4.0-rc\jre\lib\
> ext\j3dutils.jar;C:\j2sdk1.4.0-rc\jre\lib\ext\jai_codec.jar;C:\j2sdk1.4.0-rc
> \jre\lib\ext\jai_core.jar;C:\j2sdk1.4.0-rc\jre\lib\ext\jmf.jar;C:\j2sdk1.4.0
> -rc\jre\lib\ext\ldapsec.jar;C:\j2sdk1.4.0-rc\jre\lib\ext\localedata.jar;C:\j
> 2sdk1.4.0-rc\jre\lib\ext\mlibwrapper_jai.jar;C:\j2sdk1.4.0-rc\jre\lib\ext\po
> rtfolio\portfolio.jar;C:\j2sdk1.4.0-rc\jre\lib\ext\portfolio.jar;C:\j2sdk1.4
> .0-rc\jre\lib\ext\sunjce_provider.jar;C:\j2sdk1.4.0-rc\jre\lib\ext\vecmath.j
> ar;C:\j2sdk1.4.0-rc\jre\lib\ext;C:\j2sdk1.4.0-rc\jre\lib\jaws.jar;C:\j2sdk1.
> 4.0-rc\jre\lib\jce.jar;C:\j2sdk1.4.0-rc\jre\lib\jsse.jar;C:\j2sdk1.4.0-rc\jr
> e\lib\rt.jar;C:\j2sdk1.4.0-rc\jre\lib\sunrsasign.jar;C:\j2sdk1.4.0-rc\lib\dt
> .jar;C:\j2sdk1.4.0-rc\lib\htmlconverter.jar;C:\j2sdk1.4.0-rc\lib\tools.jar"
> svdTest.SVDtest
> JNL SVD computation time: 0.0844 msec
>
> JNL rank: 3
>
> JNL error
> 2.0816681711721685E-16 9.992007221626409E-16 2.220446049250313E-16
> 3.3306690738754696E-16 1.8041124150158794E-16 3.3306690738754696E-16
> 1.6653345369377348E-16 5.551115123125783E-16 -2.914335439641036E-16
>
>
> sum of squared errors: 1.766270351997534E-30
>
> GMatrix rank: 3
>
> vecmath SVD computation time: 0.0313 msec
>
> GMatrix error
> -0.004779392884258549 -0.5042163801665813 0.5760133790720257
> -0.7285576864679053 -0.6700460809138684 -1.0894544843687317
> -0.2314124499011076 0.08774924790691563 -0.7340488647018298
>
>
> sum of squared errors: 3.3527969283127526
>
>
> ----------------------------------------------------------------------------
> -
>
> ----- Original Message -----
> From: "Kelvin Chung" <[EMAIL PROTECTED]>
> To: <[EMAIL PROTECTED]>
> Sent: Thursday, September 19, 2002 12:18 PM
> Subject: Re: [JAVA3D] Matrix4d bug!
>
>
>
>>Josh Richmond wrote:
>>
>>>Hi everyone,
>>>
>>>Either I misunderstand the definition of Matrix4d.get(Matrix3d) or I've
>>
> discovered a surprising bug.
>
>>>Run the attached program to see, but it looks like (at least using my
>>
> initial values), the 3x3 matrix returned is not the 3x3 rotation matrix of
> the 4x4 matrix, even with no affine transformation.
>
>>>E:\temp>java -cp . MatrixTest
>>>4x4:
>>>0.996823, 0.036243, 0.070927, -0.245034
>>>-0.036412, 0.997481, 0.001096, 0.034035
>>>-0.07084, -0.003675, 0.999336, 0.036529
>>>0.0, 0.0, 0.0, 1.0
>>>
>>>
>>>3x3:
>>>0.03627668530489642, 0.9968262981253496, 0.07086136795886099
>>>0.9993350269113018, -0.03644589019989963, 0.0010959360038394175
>>>-0.0036750634657327997, -0.0707744901306144, 0.997485571552428
>>>
>>>As you can see it is vastly re-arranged. Am I misinterpreting the
>>
> definition of the method? Is the problem with my values? They're practically
> identity.
>
>>>josh
>>>
>>
>>
>>I don't see a problem here, if you pass the m3 matrix to Transform3D and
>>get its type :
>>
>> Transform3D t = new Transform3D();
>> t.set(m3);
>> int type = t.getType();
>>
>> if ((type & Transform3D.ZERO) > 0 )
>>System.out.print(" ZERO");
>> if ((type & Transform3D.IDENTITY) > 0 )
>
> System.out.print("
>
>>IDENTITY");
>> if ((type & Transform3D.SCALE) > 0 )
>
> System.out.print("
>
>>SCALE");
>> if ((type & Transform3D.TRANSLATION) > 0 )
>
> System.out.print("
>
>>TRANSLATION");
>> if ((type & Transform3D.ORTHOGONAL) > 0 )
>
> System.out.print("
>
>>ORTHOGONAL");
>> if ((type & Transform3D.RIGID) > 0 )
>
> System.out.print("
>
>>RIGID");
>> if ((type & Transform3D.CONGRUENT) > 0 )
>
> System.out.print("
>
>>CONGRUENT");
>> if ((type & Transform3D.AFFINE) > 0 )
>
> System.out.print("
>
>>AFFINE");
>> if ((type & Transform3D.NEGATIVE_DETERMINANT) > 0 )
>
> System.out.print("
>
>>NEGATIVE_DETERMINANT");
>>
>>
>>
>>The output is :
>>
>>ORTHOGONAL RIGID CONGRUENT AFFINE NEGATIVE_DETERMINANT4x4
>>
>>
>>- Kelvin
>>
>>
>
> ===========================================================================
>
>>To unsubscribe, send email to [EMAIL PROTECTED] and include in the
>
> body
>
>>of the message "signoff JAVA3D-INTEREST". For general help, send email to
>>[EMAIL PROTECTED] and include in the body of the message "help".
>>
>
>
>
===========================================================================
To unsubscribe, send email to [EMAIL PROTECTED] and include in the body
of the message "signoff JAVA3D-INTEREST". For general help, send email to
[EMAIL PROTECTED] and include in the body of the message "help".
