Main Page | Class Hierarchy | Class List | File List | Class Members | Related Pages

matrix.h

00001 // 00002 // matrix.h 00003 // 00004 // Copyright (C) 1996 Limit Point Systems, Inc. 00005 // 00006 // Author: Curtis Janssen <cljanss@limitpt.com> 00007 // Maintainer: LPS 00008 // 00009 // This file is part of the SC Toolkit. 00010 // 00011 // The SC Toolkit is free software; you can redistribute it and/or modify 00012 // it under the terms of the GNU Library General Public License as published by 00013 // the Free Software Foundation; either version 2, or (at your option) 00014 // any later version. 00015 // 00016 // The SC Toolkit is distributed in the hope that it will be useful, 00017 // but WITHOUT ANY WARRANTY; without even the implied warranty of 00018 // MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the 00019 // GNU Library General Public License for more details. 00020 // 00021 // You should have received a copy of the GNU Library General Public License 00022 // along with the SC Toolkit; see the file COPYING.LIB. If not, write to 00023 // the Free Software Foundation, 675 Mass Ave, Cambridge, MA 02139, USA. 00024 // 00025 // The U.S. Government is granted a limited license as per AL 91-7. 00026 // 00027 00028 #ifndef _math_scmat_matrix_h 00029 #define _math_scmat_matrix_h 00030 #ifdef __GNUC__ 00031 #pragma interface 00032 #endif 00033 00034 #include <iostream> 00035 00036 #include <math/scmat/abstract.h> 00037 00038 namespace sc { 00039 00040 class SCVectordouble; 00041 class SCMatrixdouble; 00042 class SymmSCMatrixdouble; 00043 class DiagSCMatrixdouble; 00044 00045 class SCMatrixBlockIter; 00046 class SCMatrixRectBlock; 00047 class SCMatrixLTriBlock; 00048 class SCMatrixDiagBlock; 00049 class SCVectorSimpleBlock; 00050 00051 class RefSCMatrix; 00052 class RefSymmSCMatrix; 00055 class RefSCVector: public Ref<SCVector> { 00056 // standard overrides 00057 public: 00060 RefSCVector(); 00062 RefSCVector(const RefSCVector& v); 00064 RefSCVector(SCVector *v); 00065 // don't allow automatic conversion from any reference to a 00066 // described class 00067 ~RefSCVector(); 00069 RefSCVector& operator=(SCVector* v); 00071 RefSCVector& operator=(const RefSCVector& v); 00072 00073 // vector specific members 00074 public: 00077 RefSCVector(const RefSCDimension& dim,const Ref<SCMatrixKit>&); 00078 00080 SCVectordouble operator()(int) const; 00082 SCVectordouble operator[](int) const; 00084 RefSCVector operator+(const RefSCVector&a) const; 00086 RefSCVector operator-(const RefSCVector&a) const; 00088 RefSCVector operator*(double) const; 00090 RefSCMatrix outer_product(const RefSCVector& v) const; 00092 RefSymmSCMatrix symmetric_outer_product() const; 00093 00094 void set_element(int i,double val) const; 00095 void accumulate_element(int i,double val) const; 00096 double get_element(int) const; 00097 int n() const; 00098 RefSCDimension dim() const; 00099 Ref<SCMatrixKit> kit() const; 00100 RefSCVector clone() const; 00101 RefSCVector copy() const; 00102 double maxabs() const; 00103 double scalar_product(const RefSCVector&) const; 00104 double dot(const RefSCVector&) const; 00105 void normalize() const; 00106 void randomize() const; 00107 void assign(const RefSCVector& v) const; 00108 void assign(double val) const; 00109 void assign(const double* v) const; 00110 void convert(double*) const; 00111 void scale(double val) const; 00112 void accumulate(const RefSCVector& v) const; 00113 void accumulate_product(const RefSymmSCMatrix&, const RefSCVector&); 00114 void accumulate_product(const RefSCMatrix&, const RefSCVector&); 00115 void element_op(const Ref<SCElementOp>& op) const; 00116 void element_op(const Ref<SCElementOp2>&, 00117 const RefSCVector&) const; 00118 void element_op(const Ref<SCElementOp3>&, 00119 const RefSCVector&, 00120 const RefSCVector&) const; 00121 void print(std::ostream&out) const; 00122 void print(const char*title=0, 00123 std::ostream&out=ExEnv::out0(), int precision=10) const; 00124 void save(StateOut&); 00125 void restore(StateIn&); 00126 }; 00127 RefSCVector operator*(double,const RefSCVector&); 00128 00129 class RefSymmSCMatrix; 00130 class RefDiagSCMatrix; 00134 class RefSCMatrix: public Ref<SCMatrix> { 00135 // standard overrides 00136 public: 00139 RefSCMatrix(); 00141 RefSCMatrix(const RefSCMatrix& m); 00143 RefSCMatrix(SCMatrix* m); 00144 ~RefSCMatrix(); 00146 RefSCMatrix& operator=(SCMatrix* m); 00148 RefSCMatrix& operator=(const RefSCMatrix& m); 00149 00150 // matrix specific members 00151 public: 00154 RefSCMatrix(const RefSCDimension& d1,const RefSCDimension& d2, 00155 const Ref<SCMatrixKit>&); 00156 00158 RefSCVector operator*(const RefSCVector&) const; 00159 00161 RefSCMatrix operator*(const RefSCMatrix&) const; 00162 RefSCMatrix operator*(const RefSymmSCMatrix&) const; 00163 RefSCMatrix operator*(const RefDiagSCMatrix&) const; 00164 00166 RefSCMatrix operator*(double) const; 00167 00169 RefSCMatrix operator+(const RefSCMatrix&) const; 00171 RefSCMatrix operator-(const RefSCMatrix&) const; 00172 00174 RefSCMatrix t() const; 00176 RefSCMatrix i() const; 00178 RefSCMatrix gi() const; 00179 00182 RefSCMatrix clone() const; 00183 RefSCMatrix copy() const; 00184 00185 RefSCMatrix get_subblock(int br, int er, int bc, int ec); 00186 void assign_subblock(const RefSCMatrix&, int br, int er, int bc, int ec, 00187 int source_br = 0, int source_bc = 0); 00188 void accumulate_subblock(const RefSCMatrix&, int, int, int, int, 00189 int source_br = 0, int source_bc = 0); 00190 RefSCVector get_row(int) const; 00191 RefSCVector get_column(int) const; 00192 void assign_row(const RefSCVector&, int) const; 00193 void assign_column(const RefSCVector&, int) const; 00194 void accumulate_row(const RefSCVector&, int) const; 00195 void accumulate_column(const RefSCVector&, int) const; 00196 00197 void accumulate_outer_product(const RefSCVector&,const RefSCVector&) const; 00198 void accumulate_product(const RefSCMatrix&,const RefSCMatrix&) const; 00199 void assign(const RefSCMatrix&) const; 00200 void scale(double) const; 00201 void randomize() const; 00202 void assign(double) const; 00203 void assign(const double*) const; 00204 void assign(const double**) const; 00205 void convert(double*) const; 00206 void convert(double**) const; 00207 void accumulate(const RefSCMatrix&) const; 00208 void accumulate(const RefSymmSCMatrix&) const; 00209 void accumulate(const RefDiagSCMatrix&) const; 00210 void element_op(const Ref<SCElementOp>&) const; 00211 void element_op(const Ref<SCElementOp2>&, 00212 const RefSCMatrix&) const; 00213 void element_op(const Ref<SCElementOp3>&, 00214 const RefSCMatrix&, 00215 const RefSCMatrix&) const; 00216 int nrow() const; 00217 int ncol() const; 00218 RefSCDimension rowdim() const; 00219 RefSCDimension coldim() const; 00220 Ref<SCMatrixKit> kit() const; 00221 void set_element(int,int,double) const; 00222 void accumulate_element(int,int,double) const; 00223 double get_element(int,int) const; 00224 void print(std::ostream&) const; 00225 void print(const char*title=0, 00226 std::ostream&out=ExEnv::out0(), int =10) const; 00227 double trace() const; 00228 void save(StateOut&); 00229 void restore(StateIn&); 00230 00235 void svd(const RefSCMatrix &U, 00236 const RefDiagSCMatrix &sigma, 00237 const RefSCMatrix &V); 00239 double solve_lin(const RefSCVector& v) const; 00241 double determ() const; 00243 SCMatrixdouble operator()(int i,int j) const; 00244 00248 int nblock() const; 00252 RefSCMatrix block(int i) const; 00253 }; 00255 RefSCMatrix operator*(double,const RefSCMatrix&); 00256 00259 class RefSymmSCMatrix: public Ref<SymmSCMatrix> { 00260 // standard overrides 00261 public: 00264 RefSymmSCMatrix(); 00266 RefSymmSCMatrix(const RefSymmSCMatrix& m); 00268 RefSymmSCMatrix(SymmSCMatrix *m); 00269 ~RefSymmSCMatrix(); 00271 RefSymmSCMatrix& operator=(SymmSCMatrix* m); 00273 RefSymmSCMatrix& operator=(const RefSymmSCMatrix& m); 00274 00275 // matrix specific members 00276 public: 00279 RefSymmSCMatrix(const RefSCDimension& d,const Ref<SCMatrixKit>&); 00281 RefSCMatrix operator*(const RefSCMatrix&) const; 00282 RefSCMatrix operator*(const RefSymmSCMatrix&) const; 00284 RefSCVector operator*(const RefSCVector&a) const; 00285 RefSymmSCMatrix operator*(double) const; 00287 RefSymmSCMatrix operator+(const RefSymmSCMatrix&) const; 00288 RefSymmSCMatrix operator-(const RefSymmSCMatrix&) const; 00290 RefSymmSCMatrix i() const; 00292 RefSymmSCMatrix gi() const; 00295 RefSymmSCMatrix clone() const; 00296 RefSymmSCMatrix copy() const; 00297 void set_element(int,int,double) const; 00298 void accumulate_element(int,int,double) const; 00299 double get_element(int,int) const; 00300 00301 RefSCMatrix get_subblock(int br, int er, int bc, int ec); 00302 RefSymmSCMatrix get_subblock(int br, int er); 00303 void assign_subblock(const RefSCMatrix&, int br, int er, int bc, int ec); 00304 void assign_subblock(const RefSymmSCMatrix&, int br, int er); 00305 void accumulate_subblock(const RefSCMatrix&, int, int, int, int); 00306 void accumulate_subblock(const RefSymmSCMatrix&, int, int); 00307 RefSCVector get_row(int); 00308 void assign_row(const RefSCVector&, int); 00309 void accumulate_row(const RefSCVector&, int); 00310 00311 void accumulate_symmetric_outer_product(const RefSCVector&) const; 00312 double scalar_product(const RefSCVector&) const; 00313 void accumulate_symmetric_product(const RefSCMatrix&) const; 00314 void accumulate_symmetric_sum(const RefSCMatrix&) const; 00316 void accumulate_transform(const RefSCMatrix&a,const RefSymmSCMatrix&b, 00317 SCMatrix::Transform = SCMatrix::NormalTransform) const; 00318 void accumulate_transform(const RefSCMatrix&a,const RefDiagSCMatrix&b, 00319 SCMatrix::Transform = SCMatrix::NormalTransform) const; 00320 void accumulate_transform(const RefSymmSCMatrix&a, 00321 const RefSymmSCMatrix&b) const; 00322 00323 void randomize() const; 00324 void assign(const RefSymmSCMatrix&) const; 00325 void scale(double) const; 00326 void assign(double) const; 00327 void assign(const double*) const; 00328 void assign(const double**) const; 00329 void convert(double*) const; 00330 void convert(double**) const; 00331 void accumulate(const RefSymmSCMatrix&) const; 00332 void element_op(const Ref<SCElementOp>&) const; 00333 void element_op(const Ref<SCElementOp2>&, 00334 const RefSymmSCMatrix&) const; 00335 void element_op(const Ref<SCElementOp3>&, 00336 const RefSymmSCMatrix&, 00337 const RefSymmSCMatrix&) const; 00338 double trace() const; 00339 int n() const; 00340 RefSCDimension dim() const; 00341 Ref<SCMatrixKit> kit() const; 00342 void print(std::ostream&) const; 00343 void print(const char*title=0, 00344 std::ostream&out=ExEnv::out0(), int =10) const; 00345 void save(StateOut&); 00346 void restore(StateIn&); 00347 00349 double solve_lin(const RefSCVector&) const; 00351 double determ() const; 00353 RefDiagSCMatrix eigvals() const; 00355 RefSCMatrix eigvecs() const; 00358 void diagonalize(const RefDiagSCMatrix& eigvals, 00359 const RefSCMatrix& eigvecs) const; 00361 SymmSCMatrixdouble operator()(int i,int j) const; 00365 int nblock() const; 00369 RefSymmSCMatrix block(int i) const; 00370 }; 00372 RefSymmSCMatrix operator*(double,const RefSymmSCMatrix&); 00373 00376 class RefDiagSCMatrix: public Ref<DiagSCMatrix> { 00377 // standard overrides 00378 public: 00381 RefDiagSCMatrix(); 00383 RefDiagSCMatrix(const RefDiagSCMatrix& m); 00385 RefDiagSCMatrix(DiagSCMatrix *m); 00386 ~RefDiagSCMatrix(); 00388 RefDiagSCMatrix& operator=(DiagSCMatrix* m); 00390 RefDiagSCMatrix& operator=(const RefDiagSCMatrix & m); 00391 00392 // matrix specific members 00393 public: 00396 RefDiagSCMatrix(const RefSCDimension&,const Ref<SCMatrixKit>&); 00398 RefSCMatrix operator*(const RefSCMatrix&) const; 00399 RefDiagSCMatrix operator*(double) const; 00401 RefDiagSCMatrix operator+(const RefDiagSCMatrix&) const; 00402 RefDiagSCMatrix operator-(const RefDiagSCMatrix&) const; 00404 RefDiagSCMatrix i() const; 00406 RefDiagSCMatrix gi() const; 00409 RefDiagSCMatrix clone() const; 00410 RefDiagSCMatrix copy() const; 00411 void set_element(int,double) const; 00412 void accumulate_element(int,double) const; 00413 double get_element(int) const; 00414 void randomize() const; 00415 void assign(const RefDiagSCMatrix&) const; 00416 void scale(double) const; 00417 void assign(double) const; 00418 void assign(const double*) const; 00419 void convert(double*) const; 00420 void accumulate(const RefDiagSCMatrix&) const; 00421 void element_op(const Ref<SCElementOp>&) const; 00422 void element_op(const Ref<SCElementOp2>&, 00423 const RefDiagSCMatrix&) const; 00424 void element_op(const Ref<SCElementOp3>&, 00425 const RefDiagSCMatrix&, 00426 const RefDiagSCMatrix&) const; 00427 int n() const; 00428 RefSCDimension dim() const; 00429 Ref<SCMatrixKit> kit() const; 00430 double trace() const; 00431 void print(std::ostream&) const; 00432 void print(const char*title=0, 00433 std::ostream&out=ExEnv::out0(), int =10) const; 00434 void save(StateOut&); 00435 void restore(StateIn&); 00437 double determ() const; 00439 DiagSCMatrixdouble operator()(int i) const; 00443 int nblock() const; 00447 RefDiagSCMatrix block(int i) const; 00448 }; 00450 RefDiagSCMatrix operator*(double,const RefDiagSCMatrix&); 00451 00452 class SCVectordouble { 00453 friend class RefSCVector; 00454 private: 00455 RefSCVector vector; 00456 int i; 00457 00458 SCVectordouble(SCVector*,int); 00459 public: 00460 SCVectordouble(const SCVectordouble&); 00461 ~SCVectordouble(); 00462 double operator=(double a); 00463 double operator=(const SCVectordouble&); 00464 operator double(); 00465 double val() const; 00466 }; 00467 00468 class SCMatrixdouble { 00469 friend class RefSCMatrix; 00470 private: 00471 RefSCMatrix matrix; 00472 int i; 00473 int j; 00474 00475 SCMatrixdouble(SCMatrix*,int,int); 00476 public: 00477 SCMatrixdouble(const SCMatrixdouble&); 00478 ~SCMatrixdouble(); 00479 double operator=(double a); 00480 double operator=(const SCMatrixdouble&); 00481 operator double(); 00482 double val() const; 00483 }; 00484 00485 class SymmSCMatrixdouble { 00486 friend class RefSymmSCMatrix; 00487 private: 00488 RefSymmSCMatrix matrix; 00489 int i; 00490 int j; 00491 00492 SymmSCMatrixdouble(SymmSCMatrix*,int,int); 00493 public: 00494 SymmSCMatrixdouble(const SCMatrixdouble&); 00495 ~SymmSCMatrixdouble(); 00496 double operator=(double a); 00497 double operator=(const SymmSCMatrixdouble&); 00498 operator double(); 00499 double val() const; 00500 }; 00501 00502 class DiagSCMatrixdouble { 00503 friend class RefDiagSCMatrix; 00504 private: 00505 RefDiagSCMatrix matrix; 00506 int i; 00507 int j; 00508 00509 DiagSCMatrixdouble(DiagSCMatrix*,int,int); 00510 public: 00511 DiagSCMatrixdouble(const SCMatrixdouble&); 00512 ~DiagSCMatrixdouble(); 00513 double operator=(double a); 00514 double operator=(const DiagSCMatrixdouble&); 00515 operator double(); 00516 double val() const; 00517 }; 00518 00519 } 00520 00521 #ifdef INLINE_FUNCTIONS 00522 #include <math/scmat/matrix_i.h> 00523 #endif 00524 00525 #endif 00526 00527 // Local Variables: 00528 // mode: c++ 00529 // c-file-style: "CLJ" 00530 // End:

Generated at Sat Aug 7 00:04:21 2004 for MPQC 2.2.2 using the documentation package Doxygen 1.3.8.