tmtd.cpp
Go to the documentation of this file.
1 
6 
7 
8 //#define WANT_STREAM
9 
10 #include "include.h"
11 #include "newmatap.h"
12 
13 #include "tmt.h"
14 
15 #ifdef use_namespace
16 using namespace NEWMAT;
17 #endif
18 
20 {
21  Matrix Y = X.i();
22  Y.Release();
23  return Y.ForReturn();
24 }
25 
26 // this version forces a copy
28 {
29  Matrix Y = X.i();
30  Y.Release();
31  return Y.ForReturn();
32 }
33 
35 {
36  Matrix Y = X.i();
37  Y.Release();
38  return Y.ForReturn();
39 }
40 
41 // this version forces a copy
43 {
44  Matrix Y = X.i();
45  Y.Release();
46  return Y.ForReturn();
47 }
48 
50 {
51  Tracer et1("LU1 - Crout");
52  CroutMatrix X = A;
53  return X.for_return();
54 }
55 
57 {
58  Tracer et1("LU2 - Crout");
59  CroutMatrix X = A; X.release();
60  return X.for_return();
61 }
62 
64 {
65  Tracer et1("LU3 - Crout");
67  X->release_and_delete();
68  return X->for_return();
69 }
70 
72 {
73  Tracer et1("LU1 - BandLU");
74  BandLUMatrix X = A;
75  return X.for_return();
76 }
77 
79 {
80  Tracer et1("LU2 - BandLU");
81  BandLUMatrix X = A; X.release();
82  return X.for_return();
83 }
84 
86 {
87  Tracer et1("LU3 - BandLU");
89  X->release_and_delete();
90  return X->for_return();
91 }
92 
93 
94 
95 void CircularShift(const Matrix& X1, int first, int last)
96 {
97  Matrix X; UpperTriangularMatrix U1, U2;
98  int n = X1.Ncols();
99 
100  // Try right circular shift of columns
101  X = X1; QRZ(X, U1);
102  RightCircularUpdateCholesky(U1, first, last);
103  X = X1.Columns(1,first-1) | X1.Column(last)
104  | X1.Columns(first,last-1) | X1.Columns(last+1,n);
105  QRZ(X, U2);
106  X = U1 - U2; Clean(X, 0.000000001); Print(X);
107 
108  // Try left circular shift of columns
109  X = X1; QRZ(X, U1);
110  LeftCircularUpdateCholesky(U1, first, last);
111  X = X1.Columns(1,first-1) | X1.Columns(first+1,last)
112  | X1.Column(first) | X1.Columns(last+1,n);
113  QRZ(X, U2);
114  X = U1 - U2; Clean(X, 0.000000001); Print(X);
115 }
116 
118 {
119  int m,n1,n2,n3;
120  Matrix X1, X2, X3;
121  MultWithCarry mwc; // Uniform random number generator
122 public:
123  void Reset();
124  TestUpdateQRZ(int mx, int n1x, int n2x=0, int n3x=0)
125  : m(mx), n1(n1x), n2(n2x), n3(n3x) { Reset(); }
126  void DoTest();
127  void ClearRow(int i) { X1.Row(i) = 0.0; }
128  void SetRow(int i, int j) { X1.Row(i) = X1.Row(j); }
129 };
130 
131 
132 
133 void trymatd()
134 {
135  Tracer et("Thirteenth test of Matrix package");
137  Matrix X(5,20);
138  int i,j;
139  for (j=1;j<=20;j++) X(1,j) = j+1;
140  for (i=2;i<=5;i++) for (j=1;j<=20; j++) X(i,j) = (long)X(i-1,j) * j % 1001;
141  SymmetricMatrix S; S << X * X.t();
142  Matrix SM = X * X.t() - S;
143  Print(SM);
145  Matrix Diff = L*L.t()-S; Clean(Diff, 0.000000001);
146  Print(Diff);
147  {
148  Tracer et1("Stage 1");
149  LowerTriangularMatrix L1(5);
150  Matrix Xt = X.t(); Matrix Xt2 = Xt;
151  QRZT(X,L1);
152  Diff = L - L1; Clean(Diff,0.000000001); Print(Diff);
153  UpperTriangularMatrix Ut(5);
154  QRZ(Xt,Ut);
155  Diff = L - Ut.t(); Clean(Diff,0.000000001); Print(Diff);
156  Matrix Y(3,20);
157  for (j=1;j<=20;j++) Y(1,j) = 22-j;
158  for (i=2;i<=3;i++) for (j=1;j<=20; j++)
159  Y(i,j) = (long)Y(i-1,j) * j % 101;
160  Matrix Yt = Y.t(); Matrix M,Mt; Matrix Y2=Y;
161  QRZT(X,Y,M); QRZ(Xt,Yt,Mt);
162  Diff = Xt - X.t(); Clean(Diff,0.000000001); Print(Diff);
163  Diff = Yt - Y.t(); Clean(Diff,0.000000001); Print(Diff);
164  Diff = Mt - M.t(); Clean(Diff,0.000000001); Print(Diff);
165  Diff = Y2 * Xt2 * S.i() - M * L.i();
166  Clean(Diff,0.000000001); Print(Diff);
167  }
168 
169  ColumnVector C1(5);
170  {
171  Tracer et1("Stage 2");
172  X.ReSize(5,5);
173  for (j=1;j<=5;j++) X(1,j) = j+1;
174  for (i=2;i<=5;i++) for (j=1;j<=5; j++)
175  X(i,j) = (long)X(i-1,j) * j % 1001;
176  for (i=1;i<=5;i++) C1(i) = i*i;
177  CroutMatrix A = X;
178  ColumnVector C2 = A.i() * C1; C1 = X.i() * C1;
179  X = C1 - C2; Clean(X,0.000000001); Print(X);
180  }
181 
182  {
183  Tracer et1("Stage 3");
184  X.ReSize(7,7);
185  for (j=1;j<=7;j++) X(1,j) = j+1;
186  for (i=2;i<=7;i++) for (j=1;j<=7; j++)
187  X(i,j) = (long)X(i-1,j) * j % 1001;
188  C1.ReSize(7);
189  for (i=1;i<=7;i++) C1(i) = i*i;
190  RowVector R1 = C1.t();
191  Diff = R1 * X.i() - ( X.t().i() * R1.t() ).t(); Clean(Diff,0.000000001);
192  Print(Diff);
193  }
194 
195  {
196  Tracer et1("Stage 4");
197  X.ReSize(5,5);
198  for (j=1;j<=5;j++) X(1,j) = j+1;
199  for (i=2;i<=5;i++) for (j=1;j<=5; j++)
200  X(i,j) = (long)X(i-1,j) * j % 1001;
201  C1.ReSize(5);
202  for (i=1;i<=5;i++) C1(i) = i*i;
203  CroutMatrix A1 = X*X;
204  ColumnVector C2 = A1.i() * C1; C1 = X.i() * C1; C1 = X.i() * C1;
205  X = C1 - C2; Clean(X,0.000000001); Print(X);
206  }
207 
208 
209  {
210  Tracer et1("Stage 5");
211  int n = 40;
212  SymmetricBandMatrix B(n,2); B = 0.0;
213  for (i=1; i<=n; i++)
214  {
215  B(i,i) = 6;
216  if (i<=n-1) B(i,i+1) = -4;
217  if (i<=n-2) B(i,i+2) = 1;
218  }
219  B(1,1) = 5; B(n,n) = 5;
220  SymmetricMatrix A = B;
221  ColumnVector X(n);
222  X(1) = 429;
223  for (i=2;i<=n;i++) X(i) = (long)X(i-1) * 31 % 1001;
224  X = X / 100000L;
225  // the matrix B is rather ill-conditioned so the difficulty is getting
226  // good agreement (we have chosen X very small) may not be surprising;
227  // maximum element size in B.i() is around 1400
228  ColumnVector Y1 = A.i() * X;
230  ColumnVector Y2 = C1.t().i() * (C1.i() * X) - Y1;
231  Clean(Y2, 0.000000001); Print(Y2);
232  UpperTriangularMatrix CU = C1.t().i();
233  LowerTriangularMatrix CL = C1.i();
234  Y2 = CU * (CL * X) - Y1;
235  Clean(Y2, 0.000000001); Print(Y2);
236  Y2 = B.i() * X - Y1; Clean(Y2, 0.000000001); Print(Y2);
237 
238  LowerBandMatrix C2 = Cholesky(B);
239  Matrix M = C2 - C1; Clean(M, 0.000000001); Print(M);
240  ColumnVector Y3 = C2.t().i() * (C2.i() * X) - Y1;
241  Clean(Y3, 0.000000001); Print(Y3);
242  CU = C1.t().i();
243  CL = C1.i();
244  Y3 = CU * (CL * X) - Y1;
245  Clean(Y3, 0.000000001); Print(Y3);
246 
247  Y3 = B.i() * X - Y1; Clean(Y3, 0.000000001); Print(Y3);
248 
249  SymmetricMatrix AI = A.i();
250  Y2 = AI*X - Y1; Clean(Y2, 0.000000001); Print(Y2);
251  SymmetricMatrix BI = B.i();
252  BandMatrix C = B; Matrix CI = C.i();
253  M = A.i() - CI; Clean(M, 0.000000001); Print(M);
254  M = B.i() - CI; Clean(M, 0.000000001); Print(M);
255  M = AI-BI; Clean(M, 0.000000001); Print(M);
256  M = AI-CI; Clean(M, 0.000000001); Print(M);
257 
258  M = A; AI << M; M = AI-A; Clean(M, 0.000000001); Print(M);
259  C = B; BI << C; M = BI-B; Clean(M, 0.000000001); Print(M);
260  }
261 
262  {
263  Tracer et1("Stage 5");
264  SymmetricMatrix A(4), B(4);
265  A << 5
266  << 1 << 4
267  << 2 << 1 << 6
268  << 1 << 0 << 1 << 7;
269  B << 8
270  << 1 << 5
271  << 1 << 0 << 9
272  << 2 << 1 << 0 << 6;
274  Matrix M = Cholesky(A) * B * Cholesky(A).t() - AB*AB.t();
275  Clean(M, 0.000000001); Print(M);
276  M = A * Cholesky(B); M = M * M.t() - A * B * A;
277  Clean(M, 0.000000001); Print(M);
278  }
279 
280  {
281  Tracer et1("Stage 6");
282  int N=49;
283  int i;
284  SymmetricBandMatrix S(N,1);
285  Matrix B(N,N+1); B=0;
286  for (i=1;i<=N;i++) { S(i,i)=1; B(i,i)=1; B(i,i+1)=-1; }
287  for (i=1;i<N; i++) S(i,i+1)=-.5;
288  DiagonalMatrix D(N+1); D = 1;
289  B = B.t()*S.i()*B - (D-1.0/(N+1))*2.0;
290  Clean(B, 0.000000001); Print(B);
291  }
292 
293  {
294  Tracer et1("Stage 7");
295  // Copying and moving CroutMatrix
296  Matrix A(7,7);
297  A.Row(1) << 3 << 2 << -1 << 4 << -3 << 5 << 9;
298  A.Row(2) << -8 << 7 << 2 << 0 << 7 << 0 << -1;
299  A.Row(3) << 2 << -2 << 3 << 1 << 9 << 0 << 3;
300  A.Row(4) << -1 << 5 << 2 << 2 << 5 << -1 << 2;
301  A.Row(5) << 4 << -4 << 1 << 9 << -8 << 7 << 5;
302  A.Row(6) << 1 << -2 << 5 << -1 << -2 << 5 << 1;
303  A.Row(7) << -6 << 3 << -1 << 8 << -1 << 2 << 2;
304  RowVector D(30); D = 0;
305  Real x = determinant(A);
306  CroutMatrix B = A;
307  D(1) = determinant(B) / x - 1;
308  Matrix C = A * Inverter1(B) - IdentityMatrix(7);
309  Clean(C, 0.000000001); Print(C);
310  // Test copy constructor (in Inverter2 and ordinary copy)
311  CroutMatrix B1; B1 = B;
312  D(2) = determinant(B1) / x - 1;
313  C = A * Inverter2(B1) - IdentityMatrix(7);
314  Clean(C, 0.000000001); Print(C);
315  // Do it again with release
316  B.release(); B1 = B;
317  D(2) = B.nrows(); D(3) = B.ncols(); D(4) = B.size();
318  D(5) = determinant(B1) / x - 1;
319  B1.release();
320  C = A * Inverter2(B1) - IdentityMatrix(7);
321  D(6) = B1.nrows(); D(7) = B1.ncols(); D(8) = B1.size();
322  Clean(C, 0.000000001); Print(C);
323  // see if we get an implicit invert
324  B1 = -A;
325  D(9) = determinant(B1) / x + 1; // odd number of rows - sign will change
326  C = -A * Inverter2(B1) - IdentityMatrix(7);
327  Clean(C, 0.000000001); Print(C);
328  // check for_return
329  B = LU1(A); B1 = LU2(A); CroutMatrix B2 = LU3(A);
330  C = A * B.i() - IdentityMatrix(7); Clean(C, 0.000000001); Print(C);
331  D(10) = (B == B1 ? 0 : 1) + (B == B2 ? 0 : 1);
332  // check lengths
333  D(13) = B.size()-49;
334  // check release(2)
335  B1.release(2);
336  B2 = B1; D(15) = B == B2 ? 0 : 1;
337  CroutMatrix B3 = B1; D(16) = B == B3 ? 0 : 1;
338  D(17) = B1.size();
339  // some oddments
340  B1 = B; B1 = B1.i(); C = A - B1.i(); Clean(C, 0.000000001); Print(C);
341  B1 = B; B1.release(); B1 = B1; B2 = B1;
342  D(19) = B == B1 ? 0 : 1; D(20) = B == B2 ? 0 : 1;
343  B1.cleanup(); B2 = B1; D(21) = B1.size(); D(22) = B2.size();
344  GenericMatrix GM = B; C = A.i() - GM.i(); Clean(C, 0.000000001); Print(C);
345  B1 = GM; D(23) = B == B1 ? 0 : 1;
346  B1 = A * 0; B2 = B1; D(24) = B2.is_singular() ? 0 : 1;
347  // check release again - see if memory moves
348  const Real* d = B.const_data();
349  const int* i = B.const_data_indx();
350  B1 = B;
351  const Real* d1 = B1.const_data();
352  const int* i1 = B1.const_data_indx();
353  B1.release(); B2 = B1;
354  const Real* d2 = B2.const_data();
355  const int* i2 = B2.const_data_indx();
356  D(25) = (d != d1 ? 0 : 1) + (d1 == d2 ? 0 : 1)
357  + (i != i1 ? 0 : 1) + (i1 == i2 ? 0 : 1);
358 
359  Clean(D, 0.000000001); Print(D);
360  }
361 
362  {
363  Tracer et1("Stage 8");
364  // Same for BandLUMatrix
365  BandMatrix A(7,3,2);
366  A.Row(1) << 3 << 2 << -1;
367  A.Row(2) << -8 << 7 << 2 << 0;
368  A.Row(3) << 2 << -2 << 3 << 1 << 9;
369  A.Row(4) << -1 << 5 << 2 << 2 << 5 << -1;
370  A.Row(5) << -4 << 1 << 9 << -8 << 7 << 5;
371  A.Row(6) << 5 << -1 << -2 << 5 << 1;
372  A.Row(7) << 8 << -1 << 2 << 2;
373  RowVector D(30); D = 0;
374  Real x = determinant(A);
375  BandLUMatrix B = A;
376  D(1) = determinant(B) / x - 1;
377  Matrix C = A * Inverter1(B) - IdentityMatrix(7);
378  Clean(C, 0.000000001); Print(C);
379  // Test copy constructor (in Inverter2 and ordinary copy)
380  BandLUMatrix B1; B1 = B;
381  D(2) = determinant(B1) / x - 1;
382  C = A * Inverter2(B1) - IdentityMatrix(7);
383  Clean(C, 0.000000001); Print(C);
384  // Do it again with release
385  B.release(); B1 = B;
386  D(2) = B.nrows(); D(3) = B.ncols(); D(4) = B.size();
387  D(5) = determinant(B1) / x - 1;
388  B1.release();
389  C = A * Inverter2(B1) - IdentityMatrix(7);
390  D(6) = B1.nrows(); D(7) = B1.ncols(); D(8) = B1.size();
391  Clean(C, 0.000000001); Print(C);
392  // see if we get an implicit invert
393  B1 = -A;
394  D(9) = determinant(B1) / x + 1; // odd number of rows - sign will change
395  C = -A * Inverter2(B1) - IdentityMatrix(7);
396  Clean(C, 0.000000001); Print(C);
397  // check for_return
398  B = LU1(A); B1 = LU2(A); BandLUMatrix B2 = LU3(A);
399  C = A * B.i() - IdentityMatrix(7); Clean(C, 0.000000001); Print(C);
400  D(10) = (B == B1 ? 0 : 1) + (B == B2 ? 0 : 1);
401  // check lengths
402  D(11) = B.bandwidth().lower()-3;
403  D(12) = B.bandwidth().upper()-2;
404  D(13) = B.size()-42;
405  D(14) = B.size2()-21;
406  // check release(2)
407  B1.release(2);
408  B2 = B1; D(15) = B == B2 ? 0 : 1;
409  BandLUMatrix B3 = B1; D(16) = B == B3 ? 0 : 1;
410  D(17) = B1.size();
411  // Compare with CroutMatrix
412  CroutMatrix CM = A;
413  C = CM.i() - B.i(); Clean(C, 0.000000001); Print(C);
414  D(18) = determinant(CM) / x - 1;
415  // some oddments
416  B1 = B; CM = B1.i(); C = A - CM.i(); Clean(C, 0.000000001); Print(C);
417  B1 = B; B1.release(); B1 = B1; B2 = B1;
418  D(19) = B == B1 ? 0 : 1; D(20) = B == B2 ? 0 : 1;
419  B1.cleanup(); B2 = B1; D(21) = B1.size(); D(22) = B2.size();
420  GenericMatrix GM = B; C = A.i() - GM.i(); Clean(C, 0.000000001); Print(C);
421  B1 = GM; D(23) = B == B1 ? 0 : 1;
422  B1 = A * 0; B2 = B1; D(24) = B2.is_singular() ? 0 : 1;
423  // check release again - see if memory moves
424  const Real* d = B.const_data(); const Real* dd = B.const_data();
425  const int* i = B.const_data_indx();
426  B1 = B;
427  const Real* d1 = B1.const_data(); const Real* dd1 = B1.const_data();
428  const int* i1 = B1.const_data_indx();
429  B1.release(); B2 = B1;
430  const Real* d2 = B2.const_data(); const Real* dd2 = B2.const_data();
431  const int* i2 = B2.const_data_indx();
432  D(25) = (d != d1 ? 0 : 1) + (d1 == d2 ? 0 : 1)
433  + (dd != dd1 ? 0 : 1) + (dd1 == dd2 ? 0 : 1)
434  + (i != i1 ? 0 : 1) + (i1 == i2 ? 0 : 1);
435 
436  Clean(D, 0.000000001); Print(D);
437  }
438 
439  {
440  Tracer et1("Stage 9");
441  // Modification of Cholesky decomposition
442 
443  int i, j;
444 
445  // Build test matrix
446  Matrix X(100, 10);
447  MultWithCarry mwc; // Uniform random number generator
448  for (i = 1; i <= 100; ++i) for (j = 1; j <= 10; ++j)
449  X(i, j) = 2.0 * (mwc.Next() - 0.5);
450  Matrix X1 = X; // save copy
451 
452  // Form sums of squares and products matrix and Cholesky decompose
453  SymmetricMatrix A; A << X.t() * X;
454  UpperTriangularMatrix U1 = Cholesky(A).t();
455 
456  // Do QR decomposition of X and check we get same triangular matrix
458  QRZ(X, U2);
459  Matrix Diff = U1 - U2; Clean(Diff, 0.000000001); Print(Diff);
460 
461  // Try adding new row to X and updating triangular matrix
462  RowVector NewRow(10);
463  for (j = 1; j <= 10; ++j) NewRow(j) = 2.0 * (mwc.Next() - 0.5);
464  UpdateCholesky(U2, NewRow);
465  X = X1 & NewRow; QRZ(X, U1);
466  Diff = U1 - U2; Clean(Diff, 0.000000001); Print(Diff);
467 
468  // Try removing two rows and updating triangular matrix
469  DowndateCholesky(U2, X1.Row(20));
470  DowndateCholesky(U2, X1.Row(35));
471  X = X1.Rows(1,19) & X1.Rows(21,34) & X1.Rows(36,100) & NewRow; QRZ(X, U1);
472  Diff = U1 - U2; Clean(Diff, 0.000000001); Print(Diff);
473 
474  // Circular shifts
475 
476  CircularShift(X, 3,6);
477  CircularShift(X, 5,5);
478  CircularShift(X, 4,5);
479  CircularShift(X, 1,6);
480  CircularShift(X, 6,10);
481  }
482 
483  {
484  Tracer et1("Stage 10");
485  // Try updating QRZ, QRZT decomposition
486  TestUpdateQRZ tuqrz1(10, 100, 50, 25); tuqrz1.DoTest();
487  tuqrz1.Reset(); tuqrz1.ClearRow(1); tuqrz1.DoTest();
488  tuqrz1.Reset(); tuqrz1.ClearRow(1); tuqrz1.ClearRow(2); tuqrz1.DoTest();
489  tuqrz1.Reset(); tuqrz1.ClearRow(5); tuqrz1.ClearRow(6); tuqrz1.DoTest();
490  tuqrz1.Reset(); tuqrz1.ClearRow(10); tuqrz1.DoTest();
491  TestUpdateQRZ tuqrz2(15, 100, 0, 0); tuqrz2.DoTest();
492  tuqrz2.Reset(); tuqrz2.ClearRow(1); tuqrz2.DoTest();
493  tuqrz2.Reset(); tuqrz2.ClearRow(1); tuqrz2.ClearRow(2); tuqrz2.DoTest();
494  tuqrz2.Reset(); tuqrz2.ClearRow(5); tuqrz2.ClearRow(6); tuqrz2.DoTest();
495  tuqrz2.Reset(); tuqrz2.ClearRow(15); tuqrz2.DoTest();
496  TestUpdateQRZ tuqrz3(5, 0, 10, 0); tuqrz3.DoTest();
497 
498  }
499 
500 // cout << "\nEnd of Thirteenth test\n";
501 }
502 
504 {
505  int i, j;
506  // set up Matrices
507  X1.ReSize(m, n1); X2.ReSize(m, n2); X3.ReSize(m, n3);
508  for (i = 1; i <= m; ++i) for (j = 1; j <= n1; ++j)
509  X1(i, j) = 2.0 * (mwc.Next() - 0.5);
510  for (i = 1; i <= m; ++i) for (j = 1; j <= n2; ++j)
511  X2(i, j) = 2.0 * (mwc.Next() - 0.5);
512  for (i = 1; i <= m; ++i) for (j = 1; j <= n3; ++j)
513  X3(i, j) = 2.0 * (mwc.Next() - 0.5);
514 }
515 
517 {
518  Matrix XT1 = X1.t(), XT2 = X2.t(), XT3 = X3.t();
519  Matrix X = X1 | X2 | X3;
520  SymmetricMatrix SM; SM << X * X.t();
521  LowerTriangularMatrix L, LX, LY, L0;
522  QRZT(X, L);
523  X = X1 | X2 | X3;
524  LY.ReSize(m); LY = 0.0;
525  UpdateQRZT(X, LY);
526  QRZT(X1, LX); UpdateQRZT(X2, LX); UpdateQRZT(X3, LX);
527  Matrix Diff = L * L.t() - SM; Clean(Diff, 0.000000001); Print(Diff);
528  Diff = SM - LX * LX.t(); Clean(Diff, 0.000000001); Print(Diff);
529  Diff = SM - LY * LY.t(); Clean(Diff, 0.000000001); Print(Diff);
530  UpperTriangularMatrix U, UX, UY;
531  X = XT1 & XT2 & XT3;
532  QRZ(X, U);
533  Diff = U.t() * U - SM; Clean(Diff, 0.000000001); Print(Diff);
534  UY.ReSize(m); UY = 0.0;
535  X = XT1 & XT2 & XT3;
536  UpdateQRZ(X, UY);
537  Diff = SM - UY.t() * UY; Clean(Diff, 0.000000001); Print(Diff);
538  QRZ(XT1, UX); UpdateQRZ(XT2, UX); UpdateQRZ(XT3, UX);
539  Diff = SM - UX.t() * UX; Clean(Diff, 0.000000001); Print(Diff);
540 }
541 
void Release()
Definition: newmat.h:514
void QRZT(Matrix &X, LowerTriangularMatrix &L)
Definition: hholder.cpp:33
bool is_singular() const
Definition: newmat.h:1336
void UpdateCholesky(UpperTriangularMatrix &chol, const RowVector &x)
Definition: newmatap.h:61
TestUpdateQRZ(int mx, int n1x, int n2x=0, int n3x=0)
Definition: tmtd.cpp:124
void ReSize(int m)
Definition: newmat.h:1034
bool is_singular() const
Definition: newmat.h:1085
void UpdateQRZT(Matrix &X, LowerTriangularMatrix &L)
Definition: newmatap.h:40
ReturnMatrix LU2(const Matrix &A)
Definition: tmtd.cpp:56
ReturnMatrix Inverter2(CroutMatrix X)
Definition: tmtd.cpp:27
void DowndateCholesky(UpperTriangularMatrix &chol, const RowVector &x)
Definition: newmatap.h:67
ReturnMatrix for_return() const
Definition: newmat4.cpp:249
GetSubMatrix Columns(int f, int l) const
Definition: newmat.h:2153
void trymatd()
Definition: tmtd.cpp:133
ReturnMatrix Inverter1(const CroutMatrix &X)
Definition: tmtd.cpp:19
void release_and_delete()
Definition: newmat.h:519
const int * const_data_indx() const
Definition: newmat.h:1087
virtual void ReSize(int m, int n)
Definition: newmat.h:662
int size2() const
Definition: newmat.h:1339
const int * const_data_indx() const
Definition: newmat.h:1340
double Real
Definition: include.h:307
void LeftCircularUpdateCholesky(UpperTriangularMatrix &chol, int k, int l)
Definition: newmatap.h:81
ReturnMatrix LU1(const Matrix &A)
Definition: tmtd.cpp:49
int lower() const
Definition: newmat.h:219
MatrixBandWidth bandwidth() const
Definition: newmat4.cpp:684
int upper() const
Definition: newmat.h:217
void SetRow(int i, int j)
Definition: tmtd.cpp:128
Matrix X3
Definition: tmtd.cpp:120
Upper triangular matrix.
Definition: newmat.h:799
void RightCircularUpdateCholesky(UpperTriangularMatrix &chol, int k, int l)
Definition: newmatap.h:74
Real Next()
Definition: tmt.cpp:187
ReturnMatrix LU3(const Matrix &A)
Definition: tmtd.cpp:63
void Clean(Matrix &A, Real c)
Definition: tmt.cpp:128
GetSubMatrix Column(int f) const
Definition: newmat.h:2152
FloatVector * d
void Reset()
Definition: tmtd.cpp:503
int size() const
Definition: newmat.h:501
Band matrix.
Definition: newmat.h:1096
TransposedMatrix t() const
Definition: newmat6.cpp:320
void ReSize(int m)
Definition: newmat.h:881
const Real * const_data() const
Definition: newmat.h:504
static void PrintTrace()
Definition: myexcept.cpp:109
The usual rectangular matrix.
Definition: newmat.h:625
InvertedMatrix i() const
Definition: newmat6.cpp:329
ReturnMatrix Cholesky(const SymmetricMatrix &S)
Definition: cholesky.cpp:36
void DoTest()
Definition: tmtd.cpp:516
int Ncols() const
Definition: newmat.h:495
void QRZ(Matrix &X, UpperTriangularMatrix &U)
Definition: hholder.cpp:117
LU decomposition of a band matrix.
Definition: newmat.h:1306
void MatrixErrorNoSpace(const void *)
test for allocation fails
Definition: newmatex.cpp:301
Diagonal matrix.
Definition: newmat.h:896
Lower triangular matrix.
Definition: newmat.h:848
void cleanup()
Definition: newmat4.cpp:1128
GetSubMatrix Row(int f) const
Definition: newmat.h:2150
int ncols() const
Definition: newmat.h:500
Symmetric band matrix.
Definition: newmat.h:1245
Lower triangular band matrix.
Definition: newmat.h:1206
void cleanup()
Definition: newmat4.cpp:1142
ReturnMatrix ForReturn() const
Definition: newmat.h:2157
void UpdateQRZ(Matrix &X, UpperTriangularMatrix &U)
Definition: newmatap.h:43
MultWithCarry mwc
Definition: tmtd.cpp:121
int nrows() const
Definition: newmat.h:499
Real determinant(const BaseMatrix &B)
Definition: newmat.h:2090
void CircularShift(const Matrix &X1, int first, int last)
Definition: tmtd.cpp:95
Row vector.
Definition: newmat.h:953
void Print(const Matrix &X)
Definition: tmt.cpp:42
Column vector.
Definition: newmat.h:1008
void ClearRow(int i)
Definition: tmtd.cpp:127
void ReSize(int m)
Definition: newmat.h:833
A matrix which can be of any GeneralMatrix type.
Definition: newmat.h:1397
GetSubMatrix Rows(int f, int l) const
Definition: newmat.h:2151
void release()
Definition: newmat.h:517
Identity matrix.
Definition: newmat.h:1350
Symmetric matrix.
Definition: newmat.h:753


kni
Author(s): Martin Günther
autogenerated on Fri Jun 7 2019 22:06:45