tmt4.cpp
Go to the documentation of this file.
1 
6 
7 
8 //#define WANT_STREAM
9 
10 
11 #include "include.h"
12 
13 #include "newmat.h"
14 
15 #include "tmt.h"
16 
17 #ifdef use_namespace
18 using namespace NEWMAT;
19 #endif
20 
21 
22 /**************************** test program ******************************/
23 
24 
25 void trymat4()
26 {
27 // cout << "\nFourth test of Matrix package\n";
28  Tracer et("Fourth test of Matrix package");
30 
31 
32  {
33  Tracer et1("Stage 1");
34  int i, j;
35  Matrix M(10,10);
37  for (i=1;i<=10;i++) for (j=1;j<=10;j++) M(i,j) = 100*i+j;
38  U << -M;
39  Matrix X1 = M.Rows(2,4);
40  Matrix Y1 = U.t().Rows(2,4);
41  Matrix X = U; { Print(Matrix(X.Columns(2,4).t()-Y1)); }
42  RowVector RV = M.Row(5);
43  {
44  X.ReSize(3,10);
45  X.Row(1) << M.Row(2); X.Row(2) << M.Row(3); X.Row(3) << M.Row(4);
46  Print(Matrix(X-X1));
47  }
48  {
50  Matrix MV = U.SubMatrix(3,5,3,5); { Print(Matrix(MV-V)); }
51  Matrix X2 = M.t().Columns(2,4); { Print(Matrix(X2-X1.t())); }
52  Matrix Y2 = U.Columns(2,4); { Print(Matrix(Y2-Y1.t())); }
53  ColumnVector CV = M.t().Column(5); { Print(ColumnVector(CV-RV.t())); }
54  X.ReSize(10,3); M = M.t();
55  X.Column(1) << M.Column(2); X.Column(2) << M.Column(3);
56  X.Column(3) << M.Column(4);
57  Print(Matrix(X-X2));
58  }
59  }
60 
61  {
62  Tracer et1("Stage 2");
63  int i, j;
64  Matrix M; Matrix X; M.ReSize(5,8);
65  for (i=1;i<=5;i++) for (j=1;j<=8;j++) M(i,j) = 100*i+j;
66  {
67  X = M.Columns(5,8); M.Columns(5,8) << M.Columns(1,4);
68  M.Columns(1,4) << X;
69  X = M.Columns(3,4); M.Columns(3,4) << M.Columns(1,2);
70  M.Columns(1,2) << X;
71  X = M.Columns(7,8); M.Columns(7,8) << M.Columns(5,6);
72  M.Columns(5,6) << X;
73  }
74  {
75  X = M.Column(2); M.Column(2) = M.Column(1); M.Column(1) = X;
76  X = M.Column(4); M.Column(4) = M.Column(3); M.Column(3) = X;
77  X = M.Column(6); M.Column(6) = M.Column(5); M.Column(5) = X;
78  X = M.Column(8); M.Column(8) = M.Column(7); M.Column(7) = X;
79  X.ReSize(5,8);
80  }
81  for (i=1;i<=5;i++) for (j=1;j<=8;j++) X(i,9-j) = 100*i+j;
82  Print(Matrix(X-M));
83  }
84  {
85  Tracer et1("Stage 3");
86  // try submatrices of zero dimension
87  int i, j;
88  Matrix A(4,5); Matrix B, C;
89  for (i=1; i<=4; i++) for (j=1; j<=5; j++)
90  A(i,j) = 100+i*10+j;
91  B = A + 100;
92  C = A | B.Columns(4,3); Print(Matrix(A - C));
93  C = A | B.Columns(1,0); Print(Matrix(A - C));
94  C = A | B.Columns(6,5); Print(Matrix(A - C));
95  C = A & B.Rows(2,1); Print(Matrix(A - C));
96  }
97  {
98  Tracer et1("Stage 4");
99  BandMatrix BM(5,3,2);
100  BM(1,1) = 1; BM(1,2) = 2; BM(1,3) = 3;
101  BM(2,1) = 4; BM(2,2) = 5; BM(2,3) = 6; BM(2,4) = 7;
102  BM(3,1) = 8; BM(3,2) = 9; BM(3,3) =10; BM(3,4) =11; BM(3,5) =12;
103  BM(4,1) =13; BM(4,2) =14; BM(4,3) =15; BM(4,4) =16; BM(4,5) =17;
104  BM(5,2) =18; BM(5,3) =19; BM(5,4) =20; BM(5,5) =21;
105  SymmetricBandMatrix SM(5,3);
106  SM.Inject(BandMatrix(BM + BM.t()));
107  Matrix A = BM + 1;
108  Matrix M = A + A.t() - 2;
109  Matrix C = A.i() * BM;
110  C = A * C - BM; Clean(C, 0.000000001); Print(C);
111  C = A.i() * SM;
112  C = A * C - M; Clean(C, 0.000000001); Print(C);
113 
114  // check row-wise load
115  BandMatrix BM1(5,3,2);
116  BM1.Row(1) << 1 << 2 << 3;
117  BM1.Row(2) << 4 << 5 << 6 << 7;
118  BM1.Row(3) << 8 << 9 << 10 << 11 << 12;
119  BM1.Row(4) << 13 << 14 << 15 << 16 << 17;
120  BM1.Row(5) << 18 << 19 << 20 << 21;
121  Matrix M1 = BM1 - BM; Print(M1);
122  }
123  {
124  Tracer et1("Stage 5");
125  Matrix X(4,4);
126  X << 1 << 2 << 3 << 4
127  << 5 << 6 << 7 << 8
128  << 9 <<10 <<11 <<12
129  <<13 <<14 <<15 <<16;
130  Matrix Y(4,0);
131  Y = X | Y;
132  X -= Y; Print(X);
133 
134  DiagonalMatrix D(1);
135  D << 23; // matrix input with just one value
136  D(1) -= 23; Print(D);
137 
138  }
139  {
140  Tracer et1("Stage 6");
141  Matrix h (2,2);
142  h << 1.0 << 2.0 << 0.0 << 1.0 ;
143  RowVector c(2);
144  c << 0.0 << 1.0;
145  h -= c & c;
146  h -= c.t().Reverse() | c.Reverse().t();
147  Print(h);
148  }
149  {
150  Tracer et1("Stage 7");
151  // Check row-wise input for diagonal matrix
152  DiagonalMatrix D(4);
153  D << 18 << 23 << 31 << 17;
154  DiagonalMatrix D1(4);
155  D1.Row(1) << 18; D1.Row(2) << 23; D1.Row(3) << 31; D1.Row(4) << 17;
156  D1 -= D; Print(D1);
157  D1(1) = 18; D1(2) = 23; D1(3) = 31; D1(4) = 17;
158  D1 -= D; Print(D1);
159  }
160 
161  {
162  Tracer et1("Stage 8");
163  // test swap functions
164  MultWithCarry MWC;
165  Matrix A(3,4); Matrix B(5,6);
166  FillWithValues(MWC, A); FillWithValues(MWC, B);
167  Matrix A1 = A; Matrix B1 = B; A.Release(); B.Release(2);
168  swap(A, B);
169  int a = A.size() - B1.size(), b = B.size() - A1.size();
170  Matrix D = A - B1; Print(D);
171  D = B - A1; Print(D);
172  Print(B); // should be zero because of release
173  D = A - B1; Print(D);
174  Print(A); // now should be zero
175  D.ReSize(1,2); D(1,1) = a; D(1,2) = b; Print(D);
176 
177  A.ReSize(20,20); FillWithValues(MWC, A);
178 
179  UpperTriangularMatrix UA; UA << A; UpperTriangularMatrix UA1 = UA;
181  swap(UA, UB); Print(UA); UB -= UA1; Print(UB);
182 
183  LowerTriangularMatrix LA; LA << A; LowerTriangularMatrix LA1 = LA;
185  swap(LB, LA); Print(LA); LB -= LA1; Print(LB);
186 
187  SymmetricMatrix SA; SA << A; SymmetricMatrix SA1 = SA;
188  SymmetricMatrix SB;
189  swap(SA, SB); Print(SA); SB -= SA1; Print(SB);
190 
191  DiagonalMatrix DA; DA << A; DiagonalMatrix DA1 = DA;
192  DiagonalMatrix DB;
193  swap(DB, DA); Print(DA); DB -= DA1; Print(DB);
194 
195  RowVector RVA = A.Row(1); RowVector RVA1 = RVA;
196  RowVector RVB;
197  swap(RVB, RVA); Print(RVA); RVB -= RVA1; Print(RVB);
198 
199  ColumnVector CVA = A.Column(1); ColumnVector CVA1 = CVA;
200  ColumnVector CVB;
201  swap(CVA, CVB); Print(CVA); CVB -= CVA1; Print(CVB);
202 
203  BandMatrix BA(20, 7, 4); BA.Inject(A); BandMatrix BA1 = BA;
204  BandMatrix BB;
205  swap(BA, BB); D = BA; Print(D); BB -= BA1; D = BB; Print(D);
206 
207  LowerBandMatrix LBA(20, 6); LBA.Inject(A); LowerBandMatrix LBA1 = LBA;
208  LowerBandMatrix LBB;
209  swap(LBB, LBA); D = LBA; Print(D); LBB -= LBA1; D = LBB; Print(D);
210 
211  UpperBandMatrix UBA(20, 9); UBA.Inject(A); UpperBandMatrix UBA1 = UBA;
212  UpperBandMatrix UBB;
213  swap(UBA, UBB); D = UBA; Print(D); UBB -= UBA1; D = UBB; Print(D);
214 
215  SymmetricBandMatrix SBA(20, 4); SBA.Inject(A);
216  SymmetricBandMatrix SBA1 = SBA;
218 
219  swap(SBB, SBA); D = SBA; Print(D);
220  SBB -= SBA1; D = SBB; Print(D);
221 
222  B.ReSize(10,10); FillWithValues(MWC, B);
223 
224  CroutMatrix CA = A; IdentityMatrix IA(20);
225  CroutMatrix CB = B; IdentityMatrix IB(10);
226  swap(CA, CB); swap(IA, IB);
227  D = CA.i() * B - IA; Clean(D,0.00000001); Print(D);
228  D = CB.i() * A - IB; Clean(D,0.00000001); Print(D);
229 
230  BA.ReSize(20, 5, 7); BA.Inject(A); BandLUMatrix BLUA = BA;
231  BB.ReSize(10, 3, 4); BB.Inject(B); BandLUMatrix BLUB = BB;
232  swap(BLUA, BLUB);
233  D = BLUA.i() * BB - IA; Clean(D,0.00000001); Print(D);
234  D = BLUB.i() * BA - IB; Clean(D,0.00000001); Print(D);
235 
236 
237  SBA.ReSize(20, 5); SBA.Inject(A); BandLUMatrix SBLUA = SBA;
238  SBB.ReSize(10, 3); SBB.Inject(B); BandLUMatrix SBLUB = SBB;
239  swap(SBLUA, SBLUB);
240  D = SBLUA.i() * SBB - IA; Clean(D,0.00000001); Print(D);
241  D = SBLUB.i() * SBA - IB; Clean(D,0.00000001); Print(D);
242 
243  UA << A;
244  GenericMatrix GUA = UA; GenericMatrix GB = B; swap(GUA, GB);
245  D = GB - UA; Print(D); D = B - GUA; Print(D);
246 
247  }
248 
249 // cout << "\nEnd of fourth test\n";
250 }
251 
252 
void trymat4()
Definition: tmt4.cpp:25
void Release()
Definition: newmat.h:514
GetSubMatrix Columns(int f, int l) const
Definition: newmat.h:2153
void swap(Matrix &A, Matrix &B)
Definition: newmat.h:2159
virtual void ReSize(int m, int n)
Definition: newmat.h:662
void Inject(const GeneralMatrix &GM)
Definition: newmat.h:526
virtual void ReSize(int m, int n, int b)
Definition: newmat.h:1145
ReversedMatrix Reverse() const
Definition: newmat.h:2140
Upper triangular band matrix.
Definition: newmat.h:1167
Upper triangular matrix.
Definition: newmat.h:799
void Clean(Matrix &A, Real c)
Definition: tmt.cpp:128
GetSubMatrix Column(int f) const
Definition: newmat.h:2152
int size() const
Definition: newmat.h:501
Band matrix.
Definition: newmat.h:1096
TransposedMatrix t() const
Definition: newmat6.cpp:320
void ReSize(int m, int b)
Definition: newmat.h:1289
static void PrintTrace()
Definition: myexcept.cpp:109
The usual rectangular matrix.
Definition: newmat.h:625
InvertedMatrix i() const
Definition: newmat6.cpp:329
GetSubMatrix SymSubMatrix(int f, int l) const
Definition: newmat.h:2148
FloatVector FloatVector * a
LU decomposition of a band matrix.
Definition: newmat.h:1306
Diagonal matrix.
Definition: newmat.h:896
Lower triangular matrix.
Definition: newmat.h:848
GetSubMatrix Row(int f) const
Definition: newmat.h:2150
GetSubMatrix SubMatrix(int fr, int lr, int fc, int lc) const
Definition: newmat.h:2146
Symmetric band matrix.
Definition: newmat.h:1245
Lower triangular band matrix.
Definition: newmat.h:1206
Row vector.
Definition: newmat.h:953
void Print(const Matrix &X)
Definition: tmt.cpp:42
Column vector.
Definition: newmat.h:1008
A matrix which can be of any GeneralMatrix type.
Definition: newmat.h:1397
GetSubMatrix Rows(int f, int l) const
Definition: newmat.h:2151
Identity matrix.
Definition: newmat.h:1350
void FillWithValues(MultWithCarry &MWC, Matrix &M)
Definition: tmt.cpp:190
Symmetric matrix.
Definition: newmat.h:753


kni
Author(s): Martin Günther
autogenerated on Fri Jan 3 2020 04:01:17