$search
00001 00002 00003 00006 00007 // Copyright (C) 1991,2,3,4: R B Davies 00008 00009 #include "include.h" 00010 00011 #include "newmat.h" 00012 #include "newmatrc.h" 00013 00014 #ifdef use_namespace 00015 namespace NEWMAT { 00016 #endif 00017 00018 #ifdef DO_REPORT 00019 #define REPORT { static ExeCounter ExeCount(__LINE__,11); ++ExeCount; } 00020 #else 00021 #define REPORT {} 00022 #endif 00023 00024 00025 /****************************** submatrices *********************************/ 00026 00027 GetSubMatrix BaseMatrix::submatrix(int first_row, int last_row, int first_col, 00028 int last_col) const 00029 { 00030 REPORT 00031 Tracer tr("submatrix"); 00032 int a = first_row - 1; int b = last_row - first_row + 1; 00033 int c = first_col - 1; int d = last_col - first_col + 1; 00034 if (a<0 || b<0 || c<0 || d<0) Throw(SubMatrixDimensionException()); 00035 // allow zero rows or columns 00036 return GetSubMatrix(this, a, b, c, d, false); 00037 } 00038 00039 GetSubMatrix BaseMatrix::sym_submatrix(int first_row, int last_row) const 00040 { 00041 REPORT 00042 Tracer tr("sym_submatrix"); 00043 int a = first_row - 1; int b = last_row - first_row + 1; 00044 if (a<0 || b<0) Throw(SubMatrixDimensionException()); 00045 // allow zero rows or columns 00046 return GetSubMatrix( this, a, b, a, b, true); 00047 } 00048 00049 GetSubMatrix BaseMatrix::row(int first_row) const 00050 { 00051 REPORT 00052 Tracer tr("SubMatrix(row)"); 00053 int a = first_row - 1; 00054 if (a<0) Throw(SubMatrixDimensionException()); 00055 return GetSubMatrix(this, a, 1, 0, -1, false); 00056 } 00057 00058 GetSubMatrix BaseMatrix::rows(int first_row, int last_row) const 00059 { 00060 REPORT 00061 Tracer tr("SubMatrix(rows)"); 00062 int a = first_row - 1; int b = last_row - first_row + 1; 00063 if (a<0 || b<0) Throw(SubMatrixDimensionException()); 00064 // allow zero rows or columns 00065 return GetSubMatrix(this, a, b, 0, -1, false); 00066 } 00067 00068 GetSubMatrix BaseMatrix::column(int first_col) const 00069 { 00070 REPORT 00071 Tracer tr("SubMatrix(column)"); 00072 int c = first_col - 1; 00073 if (c<0) Throw(SubMatrixDimensionException()); 00074 return GetSubMatrix(this, 0, -1, c, 1, false); 00075 } 00076 00077 GetSubMatrix BaseMatrix::columns(int first_col, int last_col) const 00078 { 00079 REPORT 00080 Tracer tr("SubMatrix(columns)"); 00081 int c = first_col - 1; int d = last_col - first_col + 1; 00082 if (c<0 || d<0) Throw(SubMatrixDimensionException()); 00083 // allow zero rows or columns 00084 return GetSubMatrix(this, 0, -1, c, d, false); 00085 } 00086 00087 void GetSubMatrix::SetUpLHS() 00088 { 00089 REPORT 00090 Tracer tr("SubMatrix(LHS)"); 00091 const BaseMatrix* bm1 = bm; 00092 GeneralMatrix* gm1 = ((BaseMatrix*&)bm)->Evaluate(); 00093 if ((BaseMatrix*)gm1!=bm1) 00094 Throw(ProgramException("Invalid LHS")); 00095 if (row_number < 0) row_number = gm1->Nrows(); 00096 if (col_number < 0) col_number = gm1->Ncols(); 00097 if (row_skip+row_number > gm1->Nrows() 00098 || col_skip+col_number > gm1->Ncols()) 00099 Throw(SubMatrixDimensionException()); 00100 } 00101 00102 void GetSubMatrix::operator<<(const BaseMatrix& bmx) 00103 { 00104 REPORT 00105 Tracer tr("SubMatrix(<<)"); GeneralMatrix* gmx = 0; 00106 Try 00107 { 00108 SetUpLHS(); gmx = ((BaseMatrix&)bmx).Evaluate(); 00109 if (row_number != gmx->Nrows() || col_number != gmx->Ncols()) 00110 Throw(IncompatibleDimensionsException()); 00111 MatrixRow mrx(gmx, LoadOnEntry); 00112 MatrixRow mr(gm, LoadOnEntry+StoreOnExit+DirectPart, row_skip); 00113 // do need LoadOnEntry 00114 MatrixRowCol sub; int i = row_number; 00115 while (i--) 00116 { 00117 mr.SubRowCol(sub, col_skip, col_number); // put values in sub 00118 sub.Copy(mrx); mr.Next(); mrx.Next(); 00119 } 00120 gmx->tDelete(); 00121 } 00122 00123 CatchAll 00124 { 00125 if (gmx) gmx->tDelete(); 00126 ReThrow; 00127 } 00128 } 00129 00130 void GetSubMatrix::operator=(const BaseMatrix& bmx) 00131 { 00132 REPORT 00133 Tracer tr("SubMatrix(=)"); GeneralMatrix* gmx = 0; 00134 // MatrixConversionCheck mcc; // Check for loss of info 00135 Try 00136 { 00137 SetUpLHS(); gmx = ((BaseMatrix&)bmx).Evaluate(); 00138 if (row_number != gmx->Nrows() || col_number != gmx->Ncols()) 00139 Throw(IncompatibleDimensionsException()); 00140 LoadAndStoreFlag lasf = 00141 ( row_skip == col_skip 00142 && gm->type().is_symmetric() 00143 && gmx->type().is_symmetric() ) 00144 ? LoadOnEntry+DirectPart 00145 : LoadOnEntry; 00146 MatrixRow mrx(gmx, lasf); 00147 MatrixRow mr(gm, LoadOnEntry+StoreOnExit+DirectPart, row_skip); 00148 // do need LoadOnEntry 00149 MatrixRowCol sub; int i = row_number; 00150 while (i--) 00151 { 00152 mr.SubRowCol(sub, col_skip, col_number); // put values in sub 00153 sub.CopyCheck(mrx); mr.Next(); mrx.Next(); 00154 } 00155 gmx->tDelete(); 00156 } 00157 00158 CatchAll 00159 { 00160 if (gmx) gmx->tDelete(); 00161 ReThrow; 00162 } 00163 } 00164 00165 void GetSubMatrix::operator<<(const double* r) 00166 { 00167 REPORT 00168 Tracer tr("SubMatrix(<<double*)"); 00169 SetUpLHS(); 00170 if (row_skip+row_number > gm->Nrows() || col_skip+col_number > gm->Ncols()) 00171 Throw(SubMatrixDimensionException()); 00172 MatrixRow mr(gm, LoadOnEntry+StoreOnExit+DirectPart, row_skip); 00173 // do need LoadOnEntry 00174 MatrixRowCol sub; int i = row_number; 00175 while (i--) 00176 { 00177 mr.SubRowCol(sub, col_skip, col_number); // put values in sub 00178 sub.Copy(r); mr.Next(); 00179 } 00180 } 00181 00182 void GetSubMatrix::operator<<(const float* r) 00183 { 00184 REPORT 00185 Tracer tr("SubMatrix(<<float*)"); 00186 SetUpLHS(); 00187 if (row_skip+row_number > gm->Nrows() || col_skip+col_number > gm->Ncols()) 00188 Throw(SubMatrixDimensionException()); 00189 MatrixRow mr(gm, LoadOnEntry+StoreOnExit+DirectPart, row_skip); 00190 // do need LoadOnEntry 00191 MatrixRowCol sub; int i = row_number; 00192 while (i--) 00193 { 00194 mr.SubRowCol(sub, col_skip, col_number); // put values in sub 00195 sub.Copy(r); mr.Next(); 00196 } 00197 } 00198 00199 void GetSubMatrix::operator<<(const int* r) 00200 { 00201 REPORT 00202 Tracer tr("SubMatrix(<<int*)"); 00203 SetUpLHS(); 00204 if (row_skip+row_number > gm->Nrows() || col_skip+col_number > gm->Ncols()) 00205 Throw(SubMatrixDimensionException()); 00206 MatrixRow mr(gm, LoadOnEntry+StoreOnExit+DirectPart, row_skip); 00207 // do need LoadOnEntry 00208 MatrixRowCol sub; int i = row_number; 00209 while (i--) 00210 { 00211 mr.SubRowCol(sub, col_skip, col_number); // put values in sub 00212 sub.Copy(r); mr.Next(); 00213 } 00214 } 00215 00216 void GetSubMatrix::operator=(Real r) 00217 { 00218 REPORT 00219 Tracer tr("SubMatrix(=Real)"); 00220 SetUpLHS(); 00221 MatrixRow mr(gm, LoadOnEntry+StoreOnExit+DirectPart, row_skip); 00222 // do need LoadOnEntry 00223 MatrixRowCol sub; int i = row_number; 00224 while (i--) 00225 { 00226 mr.SubRowCol(sub, col_skip, col_number); // put values in sub 00227 sub.Copy(r); mr.Next(); 00228 } 00229 } 00230 00231 void GetSubMatrix::inject(const GeneralMatrix& gmx) 00232 { 00233 REPORT 00234 Tracer tr("SubMatrix(inject)"); 00235 SetUpLHS(); 00236 if (row_number != gmx.Nrows() || col_number != gmx.Ncols()) 00237 Throw(IncompatibleDimensionsException()); 00238 MatrixRow mrx((GeneralMatrix*)(&gmx), LoadOnEntry); 00239 MatrixRow mr(gm, LoadOnEntry+StoreOnExit+DirectPart, row_skip); 00240 // do need LoadOnEntry 00241 MatrixRowCol sub; int i = row_number; 00242 while (i--) 00243 { 00244 mr.SubRowCol(sub, col_skip, col_number); // put values in sub 00245 sub.Inject(mrx); mr.Next(); mrx.Next(); 00246 } 00247 } 00248 00249 void GetSubMatrix::operator+=(const BaseMatrix& bmx) 00250 { 00251 REPORT 00252 Tracer tr("SubMatrix(+=)"); GeneralMatrix* gmx = 0; 00253 // MatrixConversionCheck mcc; // Check for loss of info 00254 Try 00255 { 00256 SetUpLHS(); gmx = ((BaseMatrix&)bmx).Evaluate(); 00257 if (row_number != gmx->Nrows() || col_number != gmx->Ncols()) 00258 Throw(IncompatibleDimensionsException()); 00259 MatrixRow mrx(gmx, LoadOnEntry); 00260 MatrixRow mr(gm, LoadOnEntry+StoreOnExit+DirectPart, row_skip); 00261 // do need LoadOnEntry 00262 MatrixRowCol sub; int i = row_number; 00263 while (i--) 00264 { 00265 mr.SubRowCol(sub, col_skip, col_number); // put values in sub 00266 sub.Check(mrx); // check for loss of info 00267 sub.Add(mrx); mr.Next(); mrx.Next(); 00268 } 00269 gmx->tDelete(); 00270 } 00271 00272 CatchAll 00273 { 00274 if (gmx) gmx->tDelete(); 00275 ReThrow; 00276 } 00277 } 00278 00279 void GetSubMatrix::operator-=(const BaseMatrix& bmx) 00280 { 00281 REPORT 00282 Tracer tr("SubMatrix(-=)"); GeneralMatrix* gmx = 0; 00283 // MatrixConversionCheck mcc; // Check for loss of info 00284 Try 00285 { 00286 SetUpLHS(); gmx = ((BaseMatrix&)bmx).Evaluate(); 00287 if (row_number != gmx->Nrows() || col_number != gmx->Ncols()) 00288 Throw(IncompatibleDimensionsException()); 00289 MatrixRow mrx(gmx, LoadOnEntry); 00290 MatrixRow mr(gm, LoadOnEntry+StoreOnExit+DirectPart, row_skip); 00291 // do need LoadOnEntry 00292 MatrixRowCol sub; int i = row_number; 00293 while (i--) 00294 { 00295 mr.SubRowCol(sub, col_skip, col_number); // put values in sub 00296 sub.Check(mrx); // check for loss of info 00297 sub.Sub(mrx); mr.Next(); mrx.Next(); 00298 } 00299 gmx->tDelete(); 00300 } 00301 00302 CatchAll 00303 { 00304 if (gmx) gmx->tDelete(); 00305 ReThrow; 00306 } 00307 } 00308 00309 void GetSubMatrix::operator+=(Real r) 00310 { 00311 REPORT 00312 Tracer tr("SubMatrix(+= or -= Real)"); 00313 // MatrixConversionCheck mcc; // Check for loss of info 00314 Try 00315 { 00316 SetUpLHS(); 00317 MatrixRow mr(gm, LoadOnEntry+StoreOnExit+DirectPart, row_skip); 00318 // do need LoadOnEntry 00319 MatrixRowCol sub; int i = row_number; 00320 while (i--) 00321 { 00322 mr.SubRowCol(sub, col_skip, col_number); // put values in sub 00323 sub.Check(); // check for loss of info 00324 sub.Add(r); mr.Next(); 00325 } 00326 } 00327 00328 CatchAll 00329 { 00330 ReThrow; 00331 } 00332 } 00333 00334 void GetSubMatrix::operator*=(Real r) 00335 { 00336 REPORT 00337 Tracer tr("SubMatrix(*= or /= Real)"); 00338 // MatrixConversionCheck mcc; // Check for loss of info 00339 Try 00340 { 00341 SetUpLHS(); 00342 MatrixRow mr(gm, LoadOnEntry+StoreOnExit+DirectPart, row_skip); 00343 // do need LoadOnEntry 00344 MatrixRowCol sub; int i = row_number; 00345 while (i--) 00346 { 00347 mr.SubRowCol(sub, col_skip, col_number); // put values in sub 00348 sub.Multiply(r); mr.Next(); 00349 } 00350 } 00351 00352 CatchAll 00353 { 00354 ReThrow; 00355 } 00356 } 00357 00358 #ifdef use_namespace 00359 } 00360 #endif 00361 00363