$search
00001 00002 00003 00006 00007 00008 //#define WANT_STREAM 00009 00010 #include "include.h" 00011 00012 #include "newmat.h" 00013 00014 #include "tmt.h" 00015 00016 #ifdef use_namespace 00017 using namespace NEWMAT; 00018 #endif 00019 00020 00021 // **************************** test program ****************************** 00022 00023 00024 00025 ReturnMatrix Returner0(const GenericMatrix& GM) 00026 { Matrix M = GM; M.Release(); return M; } 00027 00028 ReturnMatrix Returner1(const GenericMatrix& GM) 00029 { Matrix M = GM+1; M.Release(); return M; } 00030 00031 ReturnMatrix Returner2(const GenericMatrix& GM) 00032 { UpperBandMatrix M = GM*2; M.Release(); return M; } 00033 00034 ReturnMatrix Returner3(const GenericMatrix& GM) 00035 { LowerBandMatrix M = GM*3; M.Release(); return M; } 00036 00037 ReturnMatrix Returner4(const GenericMatrix& GM) 00038 { SymmetricMatrix M = GM+4; M.Release(); return M; } 00039 00040 ReturnMatrix Returner5(const GenericMatrix& GM) 00041 { SymmetricBandMatrix M = GM*5; M.Release(); return M; } 00042 00043 ReturnMatrix Returner6(const GenericMatrix& GM) 00044 { BandMatrix M = GM*6; M.Release(); return M; } 00045 00046 ReturnMatrix Returner7(const GenericMatrix& GM) 00047 { DiagonalMatrix M = GM*7; M.Release(); return M; } 00048 00049 void trymat5() 00050 { 00051 // cout << "\nFifth test of Matrix package\n"; 00052 Tracer et("Fifth test of Matrix package"); 00053 Tracer::PrintTrace(); 00054 00055 int i,j; 00056 00057 Matrix A(5,6); 00058 for (i=1;i<=5;i++) for (j=1;j<=6;j++) A(i,j)=1+i*j+i*i+j*j; 00059 ColumnVector CV(6); 00060 for (i=1;i<=6;i++) CV(i)=i*i+3; 00061 ColumnVector CV2(5); for (i=1;i<=5;i++) CV2(i)=1.0; 00062 ColumnVector CV1=CV; 00063 00064 { 00065 CV=A*CV; 00066 RowVector RV=CV.t(); // RowVector RV; RV=CV.t(); 00067 RV=RV-1.0; 00068 CV=(RV*A).t()+A.t()*CV2; CV1=(A.t()*A)*CV1 - CV; 00069 Print(CV1); 00070 } 00071 00072 CV1.ReSize(6); 00073 CV2.ReSize(6); 00074 CV.ReSize(6); 00075 for (i=1;i<=6;i++) { CV1(i)=i*3+1; CV2(i)=10-i; CV(i)=11+i*2; } 00076 ColumnVector CX=CV2-CV; { CX=CX+CV1; Print(CX); } 00077 Print(ColumnVector(CV1+CV2-CV)); 00078 RowVector RV=CV.t(); RowVector RV1=CV1.t(); 00079 RowVector R=RV-RV1; Print(RowVector(R-CV2.t())); 00080 00081 // test loading of list 00082 00083 RV.ReSize(10); 00084 for (i=1;i<=10;i++) RV(i) = i*i; 00085 RV1.ReSize(10); 00086 RV1 << 1 << 4 << 9 << 16 << 25 << 36 << 49 << 64 << 81 << 100; // << 121; 00087 Print(RowVector(RV-RV1)); 00088 00089 et.ReName("Fifth test of Matrix package - almost at end"); 00090 00091 Matrix X(2,3); 00092 X << 11 << 12 << 13 00093 << 21 << 22 << 23; 00094 00095 Matrix Y = X.t(); // check simple transpose 00096 00097 X(1,1) -= 11; X(1,2) -= 12; X(1,3) -= 13; 00098 X(2,1) -= 21; X(2,2) -= 22; X(2,3) -= 23; 00099 Print(X); 00100 00101 Y(1,1) -= 11; Y(2,1) -= 12; Y(3,1) -= 13; 00102 Y(1,2) -= 21; Y(2,2) -= 22; Y(3,2) -= 23; 00103 Print(Y); 00104 00105 et.ReName("Fifth test of Matrix package - at end"); 00106 00107 RV = Returner1(RV)-1; Print(RowVector(RV-RV1)); 00108 CV1 = Returner1(RV.t())-1; Print(ColumnVector(RV.t()-CV1)); 00109 #ifndef DONT_DO_NRIC 00110 nricMatrix AA = A; 00111 X = Returner1(AA)-A-1; Print(X); 00112 #endif 00113 UpperTriangularMatrix UT(31); 00114 for (i=1; i<=31; i++) for (j=i; j<=31; j++) UT(i,j) = i+j+(i-j)*(i-2*j); 00115 UpperBandMatrix UB(31,5); UB.Inject(UT); 00116 LowerTriangularMatrix LT = UT.t(); 00117 LowerBandMatrix LB(31,5); LB.Inject(LT); 00118 A = Returner0(UB)-LB.t(); Print(A); 00119 A = Returner2(UB).t()-LB*2; Print(A); 00120 A = Returner3(LB).t()-UB*3; Print(A); 00121 SymmetricMatrix SM; SM << (UT+LT); 00122 A = Returner4(SM)-UT-LT-4; Print(A); 00123 SymmetricBandMatrix SB(31,5); SB.Inject(SM); 00124 A = Returner5(SB)/5-UB-LB; Print(A); 00125 BandMatrix B = UB+LB*LB; A = LB; 00126 A = Returner6(B)/6 - UB - A*A; Print(A); 00127 DiagonalMatrix D; D << UT; 00128 D << (Returner7(D)/7 - UT); Print(D); 00129 00130 // cout << "\nEnd of fifth test\n"; 00131 } 00132 00133 00134