00001
00002
00003
00006
00007
00008
00009 #define WANT_MATH
00010
00011 #include "include.h"
00012
00013 #include "newmatap.h"
00014
00015
00016 #include "tmt.h"
00017
00018 #ifdef use_namespace
00019 using namespace NEWMAT;
00020 #endif
00021
00022
00023 void CheckIsSorted(const DiagonalMatrix& D, bool ascending = false)
00024 {
00025 DiagonalMatrix D1 = D;
00026 if (ascending) SortAscending(D1); else SortDescending(D1);
00027 D1 -= D; Print(D1);
00028 }
00029
00030
00031
00032 void trymate()
00033 {
00034
00035
00036 Tracer et("Fourteenth test of Matrix package");
00037 Tracer::PrintTrace();
00038
00039 {
00040 Tracer et1("Stage 1");
00041 Matrix A(8,5);
00042 {
00043 Real a[] = { 22, 10, 2, 3, 7,
00044 14, 7, 10, 0, 8,
00045 -1, 13, -1,-11, 3,
00046 -3, -2, 13, -2, 4,
00047 9, 8, 1, -2, 4,
00048 9, 1, -7, 5, -1,
00049 2, -6, 6, 5, 1,
00050 4, 5, 0, -2, 2 };
00051 int ai[] = { 22, 10, 2, 3, 7,
00052 14, 7, 10, 0, 8,
00053 -1, 13, -1,-11, 3,
00054 -3, -2, 13, -2, 4,
00055 9, 8, 1, -2, 4,
00056 9, 1, -7, 5, -1,
00057 2, -6, 6, 5, 1,
00058 4, 5, 0, -2, 2 };
00059 A << a;
00060
00061 Matrix AI(8,5);
00062 AI << ai; AI -= A; Print(AI);
00063 int b[] = { 13, -1,-11,
00064 -2, 13, -2,
00065 8, 1, -2,
00066 1, -7, 5 };
00067 Matrix B(8, 5); B = 23;
00068 B.SubMatrix(3,6,2,4) << b;
00069 AI = A;
00070 AI.Rows(1,2) = 23; AI.Rows(7,8) = 23;
00071 AI.Column(1) = 23; AI.Column(5) = 23;
00072 AI -= B; Print(AI);
00073 }
00074 DiagonalMatrix D; Matrix U; Matrix V;
00075 int anc = A.Ncols(); IdentityMatrix I(anc);
00076 SymmetricMatrix S1; S1 << A.t() * A;
00077 SymmetricMatrix S2; S2 << A * A.t();
00078 Real zero = 0.0; SVD(A+zero,D,U,V); CheckIsSorted(D);
00079 DiagonalMatrix D1; SVD(A,D1); CheckIsSorted(D1);
00080 D1 -= D; Clean(D1,0.000000001);Print(D1);
00081 Matrix W;
00082 SVD(A, D1, W, W, true, false); D1 -= D; W -= U;
00083 Clean(W,0.000000001); Print(W); Clean(D1,0.000000001); Print(D1);
00084 Matrix WX;
00085 SVD(A, D1, WX, W, false, true); D1 -= D; W -= V;
00086 Clean(W,0.000000001); Print(W); Clean(D1,0.000000001); Print(D1);
00087 Matrix SU = U.t() * U - I; Clean(SU,0.000000001); Print(SU);
00088 Matrix SV = V.t() * V - I; Clean(SV,0.000000001); Print(SV);
00089 Matrix B = U * D * V.t() - A; Clean(B,0.000000001); Print(B);
00090 D1=0.0; SVD(A,D1,A); CheckIsSorted(D1);
00091 A -= U; Clean(A,0.000000001); Print(A);
00092 D(1) -= sqrt(1248.0); D(2) -= 20; D(3) -= sqrt(384.0);
00093 Clean(D,0.000000001); Print(D);
00094
00095 Jacobi(S1, D, V); CheckIsSorted(D, true);
00096 V = S1 - V * D * V.t(); Clean(V,0.000000001); Print(V);
00097 D = D.Reverse(); D(1)-=1248; D(2)-=400; D(3)-=384;
00098 Clean(D,0.000000001); Print(D);
00099
00100 Jacobi(S1, D); CheckIsSorted(D, true);
00101 D = D.Reverse(); D(1)-=1248; D(2)-=400; D(3)-=384;
00102 Clean(D,0.000000001); Print(D);
00103
00104 SymmetricMatrix JW(5);
00105 Jacobi(S1, D, JW); CheckIsSorted(D, true);
00106 D = D.Reverse(); D(1)-=1248; D(2)-=400; D(3)-=384;
00107 Clean(D,0.000000001); Print(D);
00108
00109 Jacobi(S2, D, V); CheckIsSorted(D, true);
00110 V = S2 - V * D * V.t(); Clean(V,0.000000001); Print(V);
00111 D = D.Reverse(); D(1)-=1248; D(2)-=400; D(3)-=384;
00112 Clean(D,0.000000001); Print(D);
00113
00114 EigenValues(S1, D, V); CheckIsSorted(D, true);
00115 V = S1 - V * D * V.t(); Clean(V,0.000000001); Print(V);
00116 D(5)-=1248; D(4)-=400; D(3)-=384;
00117 Clean(D,0.000000001); Print(D);
00118
00119 EigenValues(S2, D, V); CheckIsSorted(D, true);
00120 V = S2 - V * D * V.t(); Clean(V,0.000000001); Print(V);
00121 D(8)-=1248; D(7)-=400; D(6)-=384;
00122 Clean(D,0.000000001); Print(D);
00123
00124 EigenValues(S1, D); CheckIsSorted(D, true);
00125 D(5)-=1248; D(4)-=400; D(3)-=384;
00126 Clean(D,0.000000001); Print(D);
00127
00128 SymmetricMatrix EW(S2);
00129 EigenValues(S2, D, EW); CheckIsSorted(D, true);
00130 D(8)-=1248; D(7)-=400; D(6)-=384;
00131 Clean(D,0.000000001); Print(D);
00132
00133 }
00134
00135 {
00136 Tracer et1("Stage 2");
00137 Matrix A(20,21);
00138 int i,j;
00139 for (i=1; i<=20; i++) for (j=1; j<=21; j++)
00140 { if (i>j) A(i,j) = 0; else if (i==j) A(i,j) = 21-i; else A(i,j) = -1; }
00141 A = A.t();
00142 SymmetricMatrix S1; S1 << A.t() * A;
00143 SymmetricMatrix S2; S2 << A * A.t();
00144 DiagonalMatrix D; Matrix U; Matrix V;
00145 DiagonalMatrix I(A.Ncols());
00146 I=1.0;
00147 SVD(A,D,U,V); CheckIsSorted(D);
00148 Matrix SU = U.t() * U - I; Clean(SU,0.000000001); Print(SU);
00149 Matrix SV = V.t() * V - I; Clean(SV,0.000000001); Print(SV);
00150 Matrix B = U * D * V.t() - A; Clean(B,0.000000001); Print(B);
00151 for (i=1; i<=20; i++) D(i) -= sqrt((22.0-i)*(21.0-i));
00152 Clean(D,0.000000001); Print(D);
00153 Jacobi(S1, D, V); CheckIsSorted(D, true);
00154 V = S1 - V * D * V.t(); Clean(V,0.000000001); Print(V);
00155 D = D.Reverse();
00156 for (i=1; i<=20; i++) D(i) -= (22-i)*(21-i);
00157 Clean(D,0.000000001); Print(D);
00158 Jacobi(S2, D, V); CheckIsSorted(D, true);
00159 V = S2 - V * D * V.t(); Clean(V,0.000000001); Print(V);
00160 D = D.Reverse();
00161 for (i=1; i<=20; i++) D(i) -= (22-i)*(21-i);
00162 Clean(D,0.000000001); Print(D);
00163
00164 EigenValues(S1, D, V); CheckIsSorted(D, true);
00165 V = S1 - V * D * V.t(); Clean(V,0.000000001); Print(V);
00166 for (i=1; i<=20; i++) D(i) -= (i+1)*i;
00167 Clean(D,0.000000001); Print(D);
00168 EigenValues(S2, D, V); CheckIsSorted(D, true);
00169 V = S2 - V * D * V.t(); Clean(V,0.000000001); Print(V);
00170 for (i=2; i<=21; i++) D(i) -= (i-1)*i;
00171 Clean(D,0.000000001); Print(D);
00172
00173 EigenValues(S1, D); CheckIsSorted(D, true);
00174 for (i=1; i<=20; i++) D(i) -= (i+1)*i;
00175 Clean(D,0.000000001); Print(D);
00176 EigenValues(S2, D); CheckIsSorted(D, true);
00177 for (i=2; i<=21; i++) D(i) -= (i-1)*i;
00178 Clean(D,0.000000001); Print(D);
00179 }
00180
00181 {
00182 Tracer et1("Stage 3");
00183 Matrix A(30,30);
00184 int i,j;
00185 for (i=1; i<=30; i++) for (j=1; j<=30; j++)
00186 { if (i>j) A(i,j) = 0; else if (i==j) A(i,j) = 1; else A(i,j) = -1; }
00187 Real d1 = A.LogDeterminant().Value();
00188 DiagonalMatrix D; Matrix U; Matrix V;
00189 DiagonalMatrix I(A.Ncols());
00190 I=1.0;
00191 SVD(A,D,U,V); CheckIsSorted(D);
00192 Matrix SU = U.t() * U - I; Clean(SU,0.000000001); Print(SU);
00193 Matrix SV = V.t() * V - I; Clean(SV,0.000000001); Print(SV);
00194 Real d2 = D.LogDeterminant().Value();
00195 Matrix B = U * D * V.t() - A; Clean(B,0.000000001); Print(B);
00196 Real d3 = D.LogDeterminant().Value();
00197 ColumnVector Test(3);
00198 Test(1) = d1 - 1; Test(2) = d2 - 1; Test(3) = d3 - 1;
00199 Clean(Test,0.00000001); Print(Test);
00200 A.ReSize(2,2);
00201 Real a = 1.5; Real b = 2; Real c = 2 * (a*a + b*b);
00202 A << a << b << a << b;
00203 I.ReSize(2); I=1;
00204 SVD(A,D,U,V); CheckIsSorted(D);
00205 SU = U.t() * U - I; Clean(SU,0.000000001); Print(SU);
00206 SV = V.t() * V - I; Clean(SV,0.000000001); Print(SV);
00207 B = U * D * V.t() - A; Clean(B,0.000000001); Print(B);
00208 D = D*D; SortDescending(D);
00209 DiagonalMatrix D50(2); D50 << c << 0; D = D - D50;
00210 Clean(D,0.000000001);
00211 Print(D);
00212 A << a << a << b << b;
00213 SVD(A,D,U,V); CheckIsSorted(D);
00214 SU = U.t() * U - I; Clean(SU,0.000000001); Print(SU);
00215 SV = V.t() * V - I; Clean(SV,0.000000001); Print(SV);
00216 B = U * D * V.t() - A; Clean(B,0.000000001); Print(B);
00217 D = D*D; SortDescending(D);
00218 D = D - D50;
00219 Clean(D,0.000000001);
00220 Print(D);
00221 }
00222
00223 {
00224 Tracer et1("Stage 4");
00225
00226
00227
00228
00229
00230 Matrix A(22,20);
00231
00232 A = 0;
00233
00234 int a=1;
00235
00236 A(a+0,a+2) = 1; A(a+0,a+18) = -1;
00237 A(a+1,a+9) = 1; A(a+1,a+12) = -1;
00238 A(a+2,a+11) = 1; A(a+2,a+12) = -1;
00239 A(a+3,a+10) = 1; A(a+3,a+19) = -1;
00240 A(a+4,a+16) = 1; A(a+4,a+19) = -1;
00241 A(a+5,a+17) = 1; A(a+5,a+18) = -1;
00242 A(a+6,a+10) = 1; A(a+6,a+4) = -1;
00243 A(a+7,a+3) = 1; A(a+7,a+2) = -1;
00244 A(a+8,a+14) = 1; A(a+8,a+15) = -1;
00245 A(a+9,a+13) = 1; A(a+9,a+16) = -1;
00246 A(a+10,a+8) = 1; A(a+10,a+9) = -1;
00247 A(a+11,a+1) = 1; A(a+11,a+15) = -1;
00248 A(a+12,a+16) = 1; A(a+12,a+4) = -1;
00249 A(a+13,a+6) = 1; A(a+13,a+9) = -1;
00250 A(a+14,a+5) = 1; A(a+14,a+4) = -1;
00251 A(a+15,a+0) = 1; A(a+15,a+1) = -1;
00252 A(a+16,a+14) = 1; A(a+16,a+0) = -1;
00253 A(a+17,a+7) = 1; A(a+17,a+6) = -1;
00254 A(a+18,a+13) = 1; A(a+18,a+5) = -1;
00255 A(a+19,a+7) = 1; A(a+19,a+8) = -1;
00256 A(a+20,a+17) = 1; A(a+20,a+3) = -1;
00257 A(a+21,a+6) = 1; A(a+21,a+11) = -1;
00258
00259
00260 Matrix U, V; DiagonalMatrix S;
00261
00262 SVD(A, S, U, V, true, true); CheckIsSorted(S);
00263
00264 DiagonalMatrix D(20); D = 1;
00265
00266 Matrix tmp = U.t() * U - D;
00267 Clean(tmp,0.000000001); Print(tmp);
00268
00269 tmp = V.t() * V - D;
00270 Clean(tmp,0.000000001); Print(tmp);
00271
00272 tmp = U * S * V.t() - A ;
00273 Clean(tmp,0.000000001); Print(tmp);
00274
00275 }
00276
00277 {
00278 Tracer et1("Stage 5");
00279 Matrix A(10,10);
00280
00281 A.Row(1) << 1.00 << 0.07 << 0.05 << 0.00 << 0.06
00282 << 0.09 << 0.03 << 0.02 << 0.02 << -0.03;
00283 A.Row(2) << 0.07 << 1.00 << 0.05 << 0.05 << -0.03
00284 << 0.07 << 0.00 << 0.07 << 0.00 << 0.02;
00285 A.Row(3) << 0.05 << 0.05 << 1.00 << 0.05 << 0.02
00286 << 0.01 << -0.05 << 0.04 << 0.05 << -0.03;
00287 A.Row(4) << 0.00 << 0.05 << 0.05 << 1.00 << -0.05
00288 << 0.04 << 0.01 << 0.02 << -0.05 << 0.00;
00289 A.Row(5) << 0.06 << -0.03 << 0.02 << -0.05 << 1.00
00290 << -0.03 << 0.02 << -0.02 << 0.04 << 0.00;
00291 A.Row(6) << 0.09 << 0.07 << 0.01 << 0.04 << -0.03
00292 << 1.00 << -0.06 << 0.08 << -0.02 << -0.10;
00293 A.Row(7) << 0.03 << 0.00 << -0.05 << 0.01 << 0.02
00294 << -0.06 << 1.00 << 0.09 << 0.12 << -0.03;
00295 A.Row(8) << 0.02 << 0.07 << 0.04 << 0.02 << -0.02
00296 << 0.08 << 0.09 << 1.00 << 0.00 << -0.02;
00297 A.Row(9) << 0.02 << 0.00 << 0.05 << -0.05 << 0.04
00298 << -0.02 << 0.12 << 0.00 << 1.00 << 0.02;
00299 A.Row(10) << -0.03 << 0.02 << -0.03 << 0.00 << 0.00
00300 << -0.10 << -0.03 << -0.02 << 0.02 << 1.00;
00301
00302 SymmetricMatrix AS; AS << A;
00303 Matrix V; DiagonalMatrix D, D1;
00304 ColumnVector Check(6);
00305 EigenValues(AS,D,V); CheckIsSorted(D, true);
00306 Check(1) = MaximumAbsoluteValue(A - V * D * V.t());
00307 DiagonalMatrix I(10); I = 1;
00308 Check(2) = MaximumAbsoluteValue(V * V.t() - I);
00309 Check(3) = MaximumAbsoluteValue(V.t() * V - I);
00310
00311 EigenValues(AS, D1); CheckIsSorted(D1, true);
00312 D -= D1;
00313 Clean(D,0.000000001); Print(D);
00314
00315 Jacobi(AS,D,V);
00316 Check(4) = MaximumAbsoluteValue(A - V * D * V.t());
00317 Check(5) = MaximumAbsoluteValue(V * V.t() - I);
00318 Check(6) = MaximumAbsoluteValue(V.t() * V - I);
00319
00320 SortAscending(D);
00321 D -= D1;
00322 Clean(D,0.000000001); Print(D);
00323
00324 Clean(Check,0.000000001); Print(Check);
00325
00326
00327
00328 SymmetricMatrix B(10);
00329
00330 B.Row(1) << 1.00;
00331 B.Row(2) << 0.07 << 1.00;
00332 B.Row(3) << 0.05 << 0.05 << 1.00;
00333 B.Row(4) << 0.00 << 0.05 << 0.05 << 1.00;
00334 B.Row(5) << 0.06 << -0.03 << 0.02 << -0.05 << 1.00;
00335 B.Row(6) << 0.09 << 0.07 << 0.01 << 0.04 << -0.03
00336 << 1.00;
00337 B.Row(7) << 0.03 << 0.00 << -0.05 << 0.01 << 0.02
00338 << -0.06 << 1.00;
00339 B.Row(8) << 0.02 << 0.07 << 0.04 << 0.02 << -0.02
00340 << 0.08 << 0.09 << 1.00;
00341 B.Row(9) << 0.02 << 0.00 << 0.05 << -0.05 << 0.04
00342 << -0.02 << 0.12 << 0.00 << 1.00;
00343 B.Row(10) << -0.03 << 0.02 << -0.03 << 0.00 << 0.00
00344 << -0.10 << -0.03 << -0.02 << 0.02 << 1.00;
00345
00346 B -= AS; Print(B);
00347
00348 }
00349
00350 {
00351 Tracer et1("Stage 6");
00352
00353 Matrix A(9,9);
00354
00355 A.Row(1) << 1.13324e+012 << 3.68788e+011 << 3.35163e+009
00356 << 3.50193e+011 << 1.25335e+011 << 1.02212e+009
00357 << 3.16602e+009 << 1.02418e+009 << 9.42959e+006;
00358 A.Row(2) << 3.68788e+011 << 1.67128e+011 << 1.27449e+009
00359 << 1.25335e+011 << 6.05413e+010 << 4.34573e+008
00360 << 1.02418e+009 << 4.69192e+008 << 3.61098e+006;
00361 A.Row(3) << 3.35163e+009 << 1.27449e+009 << 1.25571e+007
00362 << 1.02212e+009 << 4.34573e+008 << 3.69769e+006
00363 << 9.42959e+006 << 3.61098e+006 << 3.59450e+004;
00364 A.Row(4) << 3.50193e+011 << 1.25335e+011 << 1.02212e+009
00365 << 1.43514e+011 << 5.42310e+010 << 4.15822e+008
00366 << 1.23068e+009 << 4.31545e+008 << 3.58714e+006;
00367 A.Row(5) << 1.25335e+011 << 6.05413e+010 << 4.34573e+008
00368 << 5.42310e+010 << 2.76601e+010 << 1.89102e+008
00369 << 4.31545e+008 << 2.09778e+008 << 1.51083e+006;
00370 A.Row(6) << 1.02212e+009 << 4.34573e+008 << 3.69769e+006
00371 << 4.15822e+008 << 1.89102e+008 << 1.47143e+006
00372 << 3.58714e+006 << 1.51083e+006 << 1.30165e+004;
00373 A.Row(7) << 3.16602e+009 << 1.02418e+009 << 9.42959e+006
00374 << 1.23068e+009 << 4.31545e+008 << 3.58714e+006
00375 << 1.12335e+007 << 3.54778e+006 << 3.34311e+004;
00376 A.Row(8) << 1.02418e+009 << 4.69192e+008 << 3.61098e+006
00377 << 4.31545e+008 << 2.09778e+008 << 1.51083e+006
00378 << 3.54778e+006 << 1.62552e+006 << 1.25885e+004;
00379 A.Row(9) << 9.42959e+006 << 3.61098e+006 << 3.59450e+004
00380 << 3.58714e+006 << 1.51083e+006 << 1.30165e+004
00381 << 3.34311e+004 << 1.25885e+004 << 1.28000e+002;
00382
00383
00384 SymmetricMatrix AS; AS << A;
00385 Matrix V; DiagonalMatrix D, D1;
00386 ColumnVector Check(6);
00387 EigenValues(AS,D,V); CheckIsSorted(D, true);
00388 Check(1) = MaximumAbsoluteValue(A - V * D * V.t()) / 100000;
00389 DiagonalMatrix I(9); I = 1;
00390 Check(2) = MaximumAbsoluteValue(V * V.t() - I);
00391 Check(3) = MaximumAbsoluteValue(V.t() * V - I);
00392
00393 EigenValues(AS, D1);
00394 D -= D1;
00395 Clean(D,0.001); Print(D);
00396
00397 Jacobi(AS,D,V);
00398 Check(4) = MaximumAbsoluteValue(A - V * D * V.t()) / 100000;
00399 Check(5) = MaximumAbsoluteValue(V * V.t() - I);
00400 Check(6) = MaximumAbsoluteValue(V.t() * V - I);
00401
00402 SortAscending(D);
00403 D -= D1;
00404 Clean(D,0.001); Print(D);
00405
00406 Clean(Check,0.0000001); Print(Check);
00407 }
00408
00409 {
00410 Tracer et1("Stage 7");
00411
00412 Matrix A(8,8);
00413 A.Row(1)<<-0.4343<<-0.0445<<-0.4582<<-0.1612<<-0.3191<<-0.6784<<0.1068<<0;
00414 A.Row(2)<<0.5791<<0.5517<<0.2575<<-0.1055<<-0.0437<<-0.5282<<0.0442<<0;
00415 A.Row(3)<<0.5709<<-0.5179<<-0.3275<<0.2598<<-0.196<<-0.1451<<-0.4143<<0;
00416 A.Row(4)<<0.2785<<-0.5258<<0.1251<<-0.4382<<0.0514<<-0.0446<<0.6586<<0;
00417 A.Row(5)<<0.2654<<0.3736<<-0.7436<<-0.0122<<0.0376<<0.3465<<0.3397<<0;
00418 A.Row(6)<<0.0173<<-0.0056<<-0.1903<<-0.7027<<0.4863<<-0.0199<<-0.4825<<0;
00419 A.Row(7)<<0.0434<<0.0966<<0.1083<<-0.4576<<-0.7857<<0.3425<<-0.1818<<0;
00420 A.Row(8)<<0.0<<0.0<<0.0<<0.0<<0.0<<0.0<<0.0<<-1.0;
00421 Matrix U,V; DiagonalMatrix D;
00422 SVD(A,D,U,V); CheckIsSorted(D);
00423 Matrix B = U * D * V.t() - A; Clean(B,0.000000001); Print(B);
00424 DiagonalMatrix I(8); I = 1; D -= I; Clean(D,0.0001); Print(D);
00425 U *= U.t(); U -= I; Clean(U,0.000000001); Print(U);
00426 V *= V.t(); V -= I; Clean(V,0.000000001); Print(V);
00427
00428 }
00429
00430 {
00431 Tracer et1("Stage 8");
00432
00433
00434 Matrix A(15, 10);
00435 int i, j;
00436 for (i = 1; i <= 15; ++i) for (j = 1; j <= 10; ++j)
00437 A(i, j) = i + j / 1000.0;
00438 DiagonalMatrix D(10);
00439 D << 0.2 << 0.5 << 0.1 << 0.7 << 0.8 << 0.3 << 0.4 << 0.7 << 0.9 << 0.6;
00440 Matrix U = A; Matrix V = 10 - 2 * A;
00441 Matrix Prod = U * D * V.t();
00442
00443 DiagonalMatrix D2 = D; SortDescending(D2);
00444 DiagonalMatrix D1 = D; SortSV(D1, U, V); Matrix X = D1 - D2; Print(X);
00445 X = Prod - U * D1 * V.t(); Clean(X,0.000000001); Print(X);
00446 U = A; V = 10 - 2 * A;
00447 D1 = D; SortSV(D1, U); X = D1 - D2; Print(X);
00448 D1 = D; SortSV(D1, V); X = D1 - D2; Print(X);
00449 X = Prod - U * D1 * V.t(); Clean(X,0.000000001); Print(X);
00450
00451 D2 = D; SortAscending(D2);
00452 U = A; V = 10 - 2 * A;
00453 D1 = D; SortSV(D1, U, V, true); X = D1 - D2; Print(X);
00454 X = Prod - U * D1 * V.t(); Clean(X,0.000000001); Print(X);
00455 U = A; V = 10 - 2 * A;
00456 D1 = D; SortSV(D1, U, true); X = D1 - D2; Print(X);
00457 D1 = D; SortSV(D1, V, true); X = D1 - D2; Print(X);
00458 X = Prod - U * D1 * V.t(); Clean(X,0.000000001); Print(X);
00459 }
00460
00461 {
00462 Tracer et1("Stage 9");
00463
00464 Matrix A(10,10);
00465 Matrix U;
00466 Matrix V;
00467 DiagonalMatrix Sigma;
00468 Real myVals[] =
00469 {
00470 1, 1, 1, 1, 1, 1, 1, 1, 1, 1,
00471 1, 1, 1, 1, 1, 1, 1, 1, 1, 1,
00472 1, 1, 1, 1, 1, 1, 1, 1, 1, 1,
00473 1, 1, 1, 1, 1, 1, 1, 1, 1, 1,
00474 1, 1, 1, 1, 1, 1, 1, 1, 1, 1,
00475 1, 1, 1, 1, 1, 1, 1, 1, 1, 0,
00476 1, 1, 1, 1, 1, 1, 1, 1, 1, 0,
00477 1, 1, 1, 1, 1, 1, 1, 1, 0, 0,
00478 1, 1, 1, 1, 1, 1, 1, 0, 0, 0,
00479 1, 1, 1, 1, 1, 0, 0, 0, 0, 0,
00480 };
00481
00482 A << myVals;
00483 SVD(A, Sigma, U, V); CheckIsSorted(Sigma);
00484 A -= U * Sigma * V.t();
00485 Clean(A, 0.000000001); Print(A);
00486 }
00487
00488
00489
00490 }
00491
00492