00001
00002
00003
00004 #include "include.h"
00005
00006 #include "newmatap.h"
00007
00008
00009 #include "tmt.h"
00010
00011 #ifdef use_namespace
00012 using namespace NEWMAT;
00013 #endif
00014
00015
00016 void trymatj()
00017 {
00018 Tracer et("Nineteenth test of Matrix package");
00019 Tracer::PrintTrace();
00020
00021
00022 {
00023 Tracer et1("Stage 1");
00024 Matrix A(13,7), B(13,7), C(13,7);
00025 int i,j;
00026 for (i=1;i<=13;i++) for (j=1; j<=7; j++)
00027 {
00028 Real a = (i+j*j)/2, b = (i*j-i/4);
00029 A(i,j)=a; B(i,j)=b; C(i,j)=a*b;
00030 }
00031
00032 Matrix X = SP(A,B)-C; Print(X);
00033 X = SP(A,B+1.0)-A-C; Print(X);
00034 X = SP(A-1,B)+B-C; Print(X);
00035 X = SP(A-1,B+1)+B-A-C+1; Print(X);
00036
00037 A = A.Rows(7,13); B = B.Rows(7,13); C = C.Rows(7,13);
00038 LowerTriangularMatrix LTA; LTA << A;
00039 UpperTriangularMatrix UTB; UTB << B;
00040 DiagonalMatrix DC; DC << C;
00041 X = SP(LTA,UTB) - DC; Print(X);
00042 X = SP(LTA*2,UTB) - DC*2; Print(X);
00043 X = SP(LTA, UTB /2) - DC/2; Print(X);
00044 X = SP(LTA/2, UTB*2) - DC; Print(X);
00045 DiagonalMatrix DX;
00046 DX << SP(A,B); DX << (DX-C); Print(DX);
00047 DX << SP(A*4,B); DX << (DX-C*4); Print(DX);
00048 DX << SP(A,B*2); DX << (DX-C*2); Print(DX);
00049 DX << SP(A/4,B/4); DX << (DX-C/16); Print(DX);
00050 LowerTriangularMatrix LX;
00051 LX = SP(LTA,B); LX << (LX-C); Print(LX);
00052 LX = SP(LTA*3,B); LX << (LX-C*3); Print(LX);
00053 LX = SP(LTA,B*5); LX << (LX-C*5); Print(LX);
00054 LX = SP(-LTA,-B); LX << (LX-C); Print(LX);
00055 }
00056 {
00057
00058 Tracer et1("Stage 2");
00059 SymmetricMatrix A(25), B(25), C(25);
00060 int i,j;
00061 for (i=1;i<=25;i++) for (j=i;j<=25;j++)
00062 {
00063 Real a = i*j +i - j + 3;
00064 Real b = i * i + j;
00065 A(i,j)=a; B(i,j)=b; C(i,j)=a*b;
00066 }
00067 UpperTriangularMatrix UT;
00068 UT << SP(A,B); UT << (UT - C); Print(UT);
00069 Matrix MA = A, X;
00070 X = SP(MA,B)-C; Print(X);
00071 X = SP(A,B)-C; Print(X);
00072 SymmetricBandMatrix BA(25,5), BB(25,5), BC(25,5);
00073 BA.Inject(A); BB.Inject(B); BC.Inject(C);
00074 X = SP(BA,BB)-BC; Print(X);
00075 X = SP(BA*7,BB)-BC*7; Print(X);
00076 X = SP(BA,BB/8)-BC/8; Print(X);
00077 X = SP(BA*16,BB/16)-BC; Print(X);
00078 X = SP(BA,BB); X=X-BC; Print(X);
00079 X = SP(BA*2, BB/2)-BC; Print(X);
00080 X = SP(BA, BB/2)-BC/2; Print(X);
00081 X = SP(BA*2, BB)-BC*2; Print(X);
00082 }
00083 {
00084
00085 Tracer et1("Stage 3");
00086 Matrix A(19,19), B(19,19), C(19,19);
00087 int i,j;
00088 for (i=1;i<=19;i++) for (j=1;j<=19;j++)
00089 {
00090 Real a = i*j +i - 1.5*j + 3;
00091 Real b = i * i + j;
00092 A(i,j)=a; B(i,j)=b; C(i,j)=a*b;
00093 }
00094 BandMatrix BA(19,10,7), BB(19,8,15), BC(19,8,7);
00095 BA.Inject(A); BB.Inject(B); BC.Inject(C);
00096 Matrix X; BandMatrix BX; ColumnVector BW(2);
00097 X = SP(BA,BB); X=X-BC; Print(X);
00098 X = SP(BA/8,BB); X=X-BC/8; Print(X);
00099 X = SP(BA,BB*17); X=X-BC*17; Print(X);
00100 X = SP(BA/4,BB*7); X=X-BC*7/4; Print(X);
00101 X = SP(BA,BB)-BC; Print(X);
00102 X = SP(BA/8,BB)-BC/8; Print(X);
00103 X = SP(BA,BB*17)-BC*17; Print(X);
00104 X = SP(BA/4,BB*7)-BC*7/4; Print(X);
00105 BX = SP(BA,BB);
00106 BW(1)=BX.upper-7; BW(2)=BX.lower-8; Print(BW);
00107
00108 BA.ReSize(19,7,10); BB.ReSize(19,15,8);
00109 BC.ReSize(19,7,8);
00110 BA.Inject(A); BB.Inject(B); BC.Inject(C);
00111
00112 X = SP(BA,BB); X=X-BC; Print(X);
00113 X = SP(BA/8,BB); X=X-BC/8; Print(X);
00114 X = SP(BA,BB*17); X=X-BC*17; Print(X);
00115 X = SP(BA/4,BB*7); X=X-BC*7/4; Print(X);
00116 X = SP(BA,BB)-BC; Print(X);
00117 X = SP(BA/8,BB)-BC/8; Print(X);
00118 X = SP(BA,BB*17)-BC*17; Print(X);
00119 X = SP(BA/4,BB*7)-BC*7/4; Print(X);
00120 BX = SP(BA,BB);
00121 BW(1)=BX.upper-8; BW(2)=BX.lower-7; Print(BW);
00122 }
00123 {
00124
00125 Tracer et1("Stage 4");
00126 Matrix A(7,7), B(7,7);
00127 int i,j;
00128 for (i=1;i<=7;i++) for (j=1;j<=7;j++)
00129 {
00130 Real a = i*j +i - 1.5*j + 3;
00131 Real b = i * i + j;
00132 A(i,j)=a; B(i,j)=b;
00133 }
00134 BandMatrix BA(7,2,4), BB(7,3,1), BC(7,2,1);
00135 BA.Inject(A);
00136 SymmetricBandMatrix SB(7,3);
00137 SymmetricMatrix S; S << (B+B.t());
00138 SB.Inject(S); A = BA; S = SB;
00139 Matrix X;
00140 X = SP(BA,SB); X=X-SP(A,S); Print(X);
00141 X = SP(BA*2,SB); X=X-SP(A,S*2); Print(X);
00142 X = SP(BA,SB/4); X=X-SP(A/4,S); Print(X);
00143 X = SP(BA*4,SB/4); X=X-SP(A,S); Print(X);
00144 X = SP(BA,SB)-SP(A,S); Print(X);
00145 X = SP(BA*2,SB)-SP(A,S*2); Print(X);
00146 X = SP(BA,SB/4)-SP(A/4,S); Print(X);
00147 X = SP(BA*4,SB/4)-SP(A,S); Print(X);
00148 }
00149
00150 }
00151
00152