00001 #ifndef SimTK_SIMMATRIX_SMALLMATRIX_MAT_H_
00002 #define SimTK_SIMMATRIX_SMALLMATRIX_MAT_H_
00003
00004
00005
00006
00007
00008
00009
00010
00011
00012
00013
00014
00015
00016
00017
00018
00019
00020
00021
00022
00023
00024
00025
00026
00027
00028
00029
00030
00031
00032
00033
00034
00039 #include "SimTKcommon/internal/common.h"
00040
00041 namespace SimTK {
00042
00045 template <int M, int N, class ELT, int CS, int RS> class Mat {
00046 typedef ELT E;
00047 typedef typename CNT<E>::TNeg ENeg;
00048 typedef typename CNT<E>::TWithoutNegator EWithoutNegator;
00049 typedef typename CNT<E>::TReal EReal;
00050 typedef typename CNT<E>::TImag EImag;
00051 typedef typename CNT<E>::TComplex EComplex;
00052 typedef typename CNT<E>::THerm EHerm;
00053 typedef typename CNT<E>::TPosTrans EPosTrans;
00054 typedef typename CNT<E>::TSqHermT ESqHermT;
00055 typedef typename CNT<E>::TSqTHerm ESqTHerm;
00056
00057 typedef typename CNT<E>::TSqrt ESqrt;
00058 typedef typename CNT<E>::TAbs EAbs;
00059 typedef typename CNT<E>::TStandard EStandard;
00060 typedef typename CNT<E>::TInvert EInvert;
00061 typedef typename CNT<E>::TNormalize ENormalize;
00062
00063 typedef typename CNT<E>::Scalar EScalar;
00064 typedef typename CNT<E>::ULessScalar EULessScalar;
00065 typedef typename CNT<E>::Number ENumber;
00066 typedef typename CNT<E>::StdNumber EStdNumber;
00067 typedef typename CNT<E>::Precision EPrecision;
00068 typedef typename CNT<E>::ScalarNormSq EScalarNormSq;
00069
00070 public:
00071
00072 enum {
00073 NRows = M,
00074 NCols = N,
00075 MinDim = N < M ? N : M,
00076 MaxDim = N > M ? N : M,
00077 RowSpacing = RS,
00078 ColSpacing = CS,
00079 NPackedElements = M * N,
00080 NActualElements = (N-1)*CS + (M-1)*RS + 1,
00081 NActualScalars = CNT<E>::NActualScalars * NActualElements,
00082 ImagOffset = NTraits<ENumber>::ImagOffset,
00083 RealStrideFactor = 1,
00084
00085 ArgDepth = ((int)CNT<E>::ArgDepth < (int)MAX_RESOLVED_DEPTH
00086 ? CNT<E>::ArgDepth + 1
00087 : MAX_RESOLVED_DEPTH),
00088 IsScalar = 0,
00089 IsULessScalar = 0,
00090 IsNumber = 0,
00091 IsStdNumber = 0,
00092 IsPrecision = 0,
00093 SignInterpretation = CNT<E>::SignInterpretation
00094 };
00095
00096 typedef Mat<M,N,E,CS,RS> T;
00097 typedef Mat<M,N,ENeg,CS,RS> TNeg;
00098 typedef Mat<M,N,EWithoutNegator,CS,RS> TWithoutNegator;
00099
00100 typedef Mat<M,N,EReal,CS*CNT<E>::RealStrideFactor,RS*CNT<E>::RealStrideFactor>
00101 TReal;
00102 typedef Mat<M,N,EImag,CS*CNT<E>::RealStrideFactor,RS*CNT<E>::RealStrideFactor>
00103 TImag;
00104 typedef Mat<M,N,EComplex,CS,RS> TComplex;
00105 typedef Mat<N,M,EHerm,RS,CS> THerm;
00106 typedef Mat<N,M,E,RS,CS> TPosTrans;
00107 typedef E TElement;
00108 typedef Row<N,E,CS> TRow;
00109 typedef Vec<M,E,RS> TCol;
00110 typedef Vec<MinDim,E,RS+CS> TDiag;
00111
00112
00113
00114 typedef Mat<M,N,ESqrt,M,1> TSqrt;
00115 typedef Mat<M,N,EAbs,M,1> TAbs;
00116 typedef Mat<M,N,EStandard,M,1> TStandard;
00117 typedef Mat<N,M,EInvert,N,1> TInvert;
00118 typedef Mat<M,N,ENormalize,M,1> TNormalize;
00119
00120 typedef SymMat<N,ESqHermT> TSqHermT;
00121 typedef SymMat<M,ESqTHerm> TSqTHerm;
00122
00123
00124
00125 typedef Mat<M,N,E,M,1> TPacked;
00126 typedef Mat<M-1,N,E,M,1> TDropRow;
00127 typedef Mat<M,N-1,E,M,1> TDropCol;
00128 typedef Mat<M-1,N-1,E,M,1> TDropRowCol;
00129 typedef Mat<M+1,N,E,M,1> TAppendRow;
00130 typedef Mat<M,N+1,E,M,1> TAppendCol;
00131 typedef Mat<M+1,N+1,E,M,1> TAppendRowCol;
00132
00133 typedef EScalar Scalar;
00134 typedef EULessScalar ULessScalar;
00135 typedef ENumber Number;
00136 typedef EStdNumber StdNumber;
00137 typedef EPrecision Precision;
00138 typedef EScalarNormSq ScalarNormSq;
00139
00140 typedef THerm TransposeType;
00141
00142 int size() const { return M*N; }
00143 int nrow() const { return M; }
00144 int ncol() const { return N; }
00145
00146
00147 ScalarNormSq scalarNormSqr() const {
00148 ScalarNormSq sum(0);
00149 for(int j=0;j<N;++j) sum += CNT<TCol>::scalarNormSqr((*this)(j));
00150 return sum;
00151 }
00152
00153
00154
00155
00156 TSqrt sqrt() const {
00157 TSqrt msqrt;
00158 for(int j=0;j<N;++j) msqrt(j) = (*this)(j).sqrt();
00159 return msqrt;
00160 }
00161
00162
00163
00164
00165 TAbs abs() const {
00166 TAbs mabs;
00167 for(int j=0;j<N;++j) mabs(j) = (*this)(j).abs();
00168 return mabs;
00169 }
00170
00171 TStandard standardize() const {
00172 TStandard mstd;
00173 for(int j=0;j<N;++j) mstd(j) = (*this)(j).standardize();
00174 return mstd;
00175 }
00176
00177
00178
00179
00180
00181 template <class P> struct EltResult {
00182 typedef Mat<M,N, typename CNT<E>::template Result<P>::Mul, M, 1> Mul;
00183 typedef Mat<M,N, typename CNT<E>::template Result<P>::Dvd, M, 1> Dvd;
00184 typedef Mat<M,N, typename CNT<E>::template Result<P>::Add, M, 1> Add;
00185 typedef Mat<M,N, typename CNT<E>::template Result<P>::Sub, M, 1> Sub;
00186 };
00187
00188
00189
00190 template <class P> struct Result {
00191 typedef MulCNTs<M,N,ArgDepth,Mat,ColSpacing,RowSpacing,
00192 CNT<P>::NRows, CNT<P>::NCols, CNT<P>::ArgDepth,
00193 P, CNT<P>::ColSpacing, CNT<P>::RowSpacing> MulOp;
00194 typedef typename MulOp::Type Mul;
00195
00196 typedef MulCNTsNonConforming<M,N,ArgDepth,Mat,ColSpacing,RowSpacing,
00197 CNT<P>::NRows, CNT<P>::NCols, CNT<P>::ArgDepth,
00198 P, CNT<P>::ColSpacing, CNT<P>::RowSpacing> MulOpNonConforming;
00199 typedef typename MulOpNonConforming::Type MulNon;
00200
00201 typedef DvdCNTs<M,N,ArgDepth,Mat,ColSpacing,RowSpacing,
00202 CNT<P>::NRows, CNT<P>::NCols, CNT<P>::ArgDepth,
00203 P, CNT<P>::ColSpacing, CNT<P>::RowSpacing> DvdOp;
00204 typedef typename DvdOp::Type Dvd;
00205
00206 typedef AddCNTs<M,N,ArgDepth,Mat,ColSpacing,RowSpacing,
00207 CNT<P>::NRows, CNT<P>::NCols, CNT<P>::ArgDepth,
00208 P, CNT<P>::ColSpacing, CNT<P>::RowSpacing> AddOp;
00209 typedef typename AddOp::Type Add;
00210
00211 typedef SubCNTs<M,N,ArgDepth,Mat,ColSpacing,RowSpacing,
00212 CNT<P>::NRows, CNT<P>::NCols, CNT<P>::ArgDepth,
00213 P, CNT<P>::ColSpacing, CNT<P>::RowSpacing> SubOp;
00214 typedef typename SubOp::Type Sub;
00215 };
00216
00217
00218 template <class P> struct Substitute {
00219 typedef Mat<M,N,P> Type;
00220 };
00221
00222
00223
00224 Mat(){
00225 #ifndef NDEBUG
00226 setToNaN();
00227 #endif
00228 }
00229
00230
00231
00232
00233 Mat(const Mat& src) {
00234 for (int j=0; j<N; ++j)
00235 (*this)(j) = src(j);
00236 }
00237 Mat& operator=(const Mat& src) {
00238 for (int j=0; j<N; ++j)
00239 (*this)(j) = src(j);
00240 return *this;
00241 }
00242
00243
00244
00245 explicit Mat(const SymMat<M, ELT>& src) {
00246 diag() = src.diag();
00247 for (int j = 0; j < M; ++j)
00248 for (int i = j+1; i < M; ++i) {
00249 (*this)(i, j) = src.getEltLower(i, j);
00250 (*this)(j, i) = src.getEltUpper(j, i);
00251 }
00252 }
00253
00254
00255
00256 template <int CSS, int RSS> Mat(const Mat<M,N,E,CSS,RSS>& src) {
00257 for (int j=0; j<N; ++j)
00258 (*this)(j) = src(j);
00259 }
00260
00261
00262
00263 template <int CSS, int RSS> Mat(const Mat<M,N,ENeg,CSS,RSS>& src) {
00264 for (int j=0; j<N; ++j)
00265 (*this)(j) = src(j);
00266 }
00267
00268
00269
00270 template <class EE, int CSS, int RSS> explicit Mat(const Mat<M,N,EE,CSS,RSS>& mm)
00271 { for (int j=0;j<N;++j) (*this)(j) = mm(j);}
00272
00273
00274
00275 explicit Mat(const E& e)
00276 { for (int j=0;j<N;++j) (*this)(j) = E(0); diag()=e; }
00277
00278
00279 Mat(const E& e0,const E& e1)
00280 {assert(M*N==2);d[rIx(0)]=e0;d[rIx(1)]=e1;}
00281 Mat(const E& e0,const E& e1,const E& e2)
00282 {assert(M*N==3);d[rIx(0)]=e0;d[rIx(1)]=e1;d[rIx(2)]=e2;}
00283 Mat(const E& e0,const E& e1,const E& e2,const E& e3)
00284 {assert(M*N==4);d[rIx(0)]=e0;d[rIx(1)]=e1;d[rIx(2)]=e2;d[rIx(3)]=e3;}
00285 Mat(const E& e0,const E& e1,const E& e2,const E& e3,const E& e4)
00286 {assert(M*N==5);d[rIx(0)]=e0;d[rIx(1)]=e1;d[rIx(2)]=e2;d[rIx(3)]=e3;d[rIx(4)]=e4;}
00287 Mat(const E& e0,const E& e1,const E& e2,const E& e3,const E& e4,
00288 const E& e5)
00289 {assert(M*N==6);d[rIx(0)]=e0;d[rIx(1)]=e1;d[rIx(2)]=e2;d[rIx(3)]=e3;d[rIx(4)]=e4;
00290 d[rIx(5)]=e5;}
00291 Mat(const E& e0,const E& e1,const E& e2,const E& e3,const E& e4,
00292 const E& e5,const E& e6)
00293 {assert(M*N==7);d[rIx(0)]=e0;d[rIx(1)]=e1;d[rIx(2)]=e2;d[rIx(3)]=e3;d[rIx(4)]=e4;
00294 d[rIx(5)]=e5;d[rIx(6)]=e6;}
00295 Mat(const E& e0,const E& e1,const E& e2,const E& e3,const E& e4,
00296 const E& e5,const E& e6,const E& e7)
00297 {assert(M*N==8);d[rIx(0)]=e0;d[rIx(1)]=e1;d[rIx(2)]=e2;d[rIx(3)]=e3;d[rIx(4)]=e4;
00298 d[rIx(5)]=e5;d[rIx(6)]=e6;d[rIx(7)]=e7;}
00299 Mat(const E& e0,const E& e1,const E& e2,const E& e3,const E& e4,
00300 const E& e5,const E& e6,const E& e7,const E& e8)
00301 {assert(M*N==9);d[rIx(0)]=e0;d[rIx(1)]=e1;d[rIx(2)]=e2;d[rIx(3)]=e3;d[rIx(4)]=e4;
00302 d[rIx(5)]=e5;d[rIx(6)]=e6;d[rIx(7)]=e7;d[rIx(8)]=e8;}
00303 Mat(const E& e0,const E& e1,const E& e2,const E& e3,const E& e4,
00304 const E& e5,const E& e6,const E& e7,const E& e8,const E& e9)
00305 {assert(M*N==10);d[rIx(0)]=e0;d[rIx(1)]=e1;d[rIx(2)]=e2;d[rIx(3)]=e3;d[rIx(4)]=e4;
00306 d[rIx(5)]=e5;d[rIx(6)]=e6;d[rIx(7)]=e7;d[rIx(8)]=e8;d[rIx(9)]=e9;}
00307 Mat(const E& e0,const E& e1,const E& e2,const E& e3,const E& e4,
00308 const E& e5,const E& e6,const E& e7,const E& e8,const E& e9,
00309 const E& e10)
00310 {assert(M*N==11);d[rIx(0)]=e0;d[rIx(1)]=e1;d[rIx(2)]=e2;d[rIx(3)]=e3;d[rIx(4)]=e4;
00311 d[rIx(5)]=e5;d[rIx(6)]=e6;d[rIx(7)]=e7;d[rIx(8)]=e8;d[rIx(9)]=e9;d[rIx(10)]=e10;}
00312 Mat(const E& e0,const E& e1,const E& e2,const E& e3,const E& e4,
00313 const E& e5,const E& e6,const E& e7,const E& e8,const E& e9,
00314 const E& e10, const E& e11)
00315 {assert(M*N==12);d[rIx(0)]=e0;d[rIx(1)]=e1;d[rIx(2)]=e2;d[rIx(3)]=e3;d[rIx(4)]=e4;
00316 d[rIx(5)]=e5;d[rIx(6)]=e6;d[rIx(7)]=e7;d[rIx(8)]=e8;d[rIx(9)]=e9;d[rIx(10)]=e10;
00317 d[rIx(11)]=e11;}
00318 Mat(const E& e0,const E& e1,const E& e2,const E& e3,const E& e4,
00319 const E& e5,const E& e6,const E& e7,const E& e8,const E& e9,
00320 const E& e10, const E& e11, const E& e12)
00321 {assert(M*N==13);d[rIx(0)]=e0;d[rIx(1)]=e1;d[rIx(2)]=e2;d[rIx(3)]=e3;d[rIx(4)]=e4;
00322 d[rIx(5)]=e5;d[rIx(6)]=e6;d[rIx(7)]=e7;d[rIx(8)]=e8;d[rIx(9)]=e9;d[rIx(10)]=e10;
00323 d[rIx(11)]=e11;d[rIx(12)]=e12;}
00324 Mat(const E& e0,const E& e1,const E& e2,const E& e3,const E& e4,
00325 const E& e5,const E& e6,const E& e7,const E& e8,const E& e9,
00326 const E& e10, const E& e11, const E& e12, const E& e13)
00327 {assert(M*N==14);d[rIx(0)]=e0;d[rIx(1)]=e1;d[rIx(2)]=e2;d[rIx(3)]=e3;d[rIx(4)]=e4;
00328 d[rIx(5)]=e5;d[rIx(6)]=e6;d[rIx(7)]=e7;d[rIx(8)]=e8;d[rIx(9)]=e9;d[rIx(10)]=e10;
00329 d[rIx(11)]=e11;d[rIx(12)]=e12;d[rIx(13)]=e13;}
00330 Mat(const E& e0,const E& e1,const E& e2,const E& e3,const E& e4,
00331 const E& e5,const E& e6,const E& e7,const E& e8,const E& e9,
00332 const E& e10, const E& e11, const E& e12, const E& e13, const E& e14)
00333 {assert(M*N==15);d[rIx(0)]=e0;d[rIx(1)]=e1;d[rIx(2)]=e2;d[rIx(3)]=e3;d[rIx(4)]=e4;
00334 d[rIx(5)]=e5;d[rIx(6)]=e6;d[rIx(7)]=e7;d[rIx(8)]=e8;d[rIx(9)]=e9;d[rIx(10)]=e10;
00335 d[rIx(11)]=e11;d[rIx(12)]=e12;d[rIx(13)]=e13;d[rIx(14)]=e14;}
00336 Mat(const E& e0,const E& e1,const E& e2,const E& e3,const E& e4,
00337 const E& e5,const E& e6,const E& e7,const E& e8,const E& e9,
00338 const E& e10, const E& e11, const E& e12, const E& e13, const E& e14,
00339 const E& e15)
00340 {assert(M*N==16);d[rIx(0)]=e0;d[rIx(1)]=e1;d[rIx(2)]=e2;d[rIx(3)]=e3;d[rIx(4)]=e4;
00341 d[rIx(5)]=e5;d[rIx(6)]=e6;d[rIx(7)]=e7;d[rIx(8)]=e8;d[rIx(9)]=e9;d[rIx(10)]=e10;
00342 d[rIx(11)]=e11;d[rIx(12)]=e12;d[rIx(13)]=e13;d[rIx(14)]=e14;d[rIx(15)]=e15;}
00343
00344
00345 explicit Mat(const TRow& r0)
00346 { assert(M==1); (*this)[0]=r0; }
00347 Mat(const TRow& r0,const TRow& r1)
00348 { assert(M==2);(*this)[0]=r0;(*this)[1]=r1; }
00349 Mat(const TRow& r0,const TRow& r1,const TRow& r2)
00350 { assert(M==3);(*this)[0]=r0;(*this)[1]=r1;(*this)[2]=r2; }
00351 Mat(const TRow& r0,const TRow& r1,const TRow& r2,
00352 const TRow& r3)
00353 { assert(M==4);(*this)[0]=r0;(*this)[1]=r1;(*this)[2]=r2;(*this)[3]=r3; }
00354 Mat(const TRow& r0,const TRow& r1,const TRow& r2,
00355 const TRow& r3,const TRow& r4)
00356 { assert(M==5);(*this)[0]=r0;(*this)[1]=r1;(*this)[2]=r2;
00357 (*this)[3]=r3;(*this)[4]=r4; }
00358 Mat(const TRow& r0,const TRow& r1,const TRow& r2,
00359 const TRow& r3,const TRow& r4,const TRow& r5)
00360 { assert(M==6);(*this)[0]=r0;(*this)[1]=r1;(*this)[2]=r2;
00361 (*this)[3]=r3;(*this)[4]=r4;(*this)[5]=r5; }
00362
00363
00364 template <class EE, int SS> explicit Mat(const Row<N,EE,SS>& r0)
00365 { assert(M==1); (*this)[0]=r0; }
00366 template <class EE, int SS> Mat(const Row<N,EE,SS>& r0,const Row<N,EE,SS>& r1)
00367 { assert(M==2);(*this)[0]=r0;(*this)[1]=r1; }
00368 template <class EE, int SS>
00369 Mat(const Row<N,EE,SS>& r0,const Row<N,EE,SS>& r1,const Row<N,EE,SS>& r2)
00370 { assert(M==3);(*this)[0]=r0;(*this)[1]=r1;(*this)[2]=r2; }
00371 template <class EE, int SS>
00372 Mat(const Row<N,EE,SS>& r0,const Row<N,EE,SS>& r1,const Row<N,EE,SS>& r2,
00373 const Row<N,EE,SS>& r3)
00374 { assert(M==4);(*this)[0]=r0;(*this)[1]=r1;(*this)[2]=r2;(*this)[3]=r3; }
00375 template <class EE, int SS>
00376 Mat(const Row<N,EE,SS>& r0,const Row<N,EE,SS>& r1,const Row<N,EE,SS>& r2,
00377 const Row<N,EE,SS>& r3,const Row<N,EE,SS>& r4)
00378 { assert(M==5);(*this)[0]=r0;(*this)[1]=r1;(*this)[2]=r2;
00379 (*this)[3]=r3;(*this)[4]=r4; }
00380 template <class EE, int SS>
00381 Mat(const Row<N,EE,SS>& r0,const Row<N,EE,SS>& r1,const Row<N,EE,SS>& r2,
00382 const Row<N,EE,SS>& r3,const Row<N,EE,SS>& r4,const Row<N,EE,SS>& r5)
00383 { assert(M==6);(*this)[0]=r0;(*this)[1]=r1;(*this)[2]=r2;
00384 (*this)[3]=r3;(*this)[4]=r4;(*this)[5]=r5; }
00385
00386
00387
00388 explicit Mat(const TCol& r0)
00389 { assert(N==1); (*this)(0)=r0; }
00390 Mat(const TCol& r0,const TCol& r1)
00391 { assert(N==2);(*this)(0)=r0;(*this)(1)=r1; }
00392 Mat(const TCol& r0,const TCol& r1,const TCol& r2)
00393 { assert(N==3);(*this)(0)=r0;(*this)(1)=r1;(*this)(2)=r2; }
00394 Mat(const TCol& r0,const TCol& r1,const TCol& r2,
00395 const TCol& r3)
00396 { assert(N==4);(*this)(0)=r0;(*this)(1)=r1;(*this)(2)=r2;(*this)(3)=r3; }
00397 Mat(const TCol& r0,const TCol& r1,const TCol& r2,
00398 const TCol& r3,const TCol& r4)
00399 { assert(N==5);(*this)(0)=r0;(*this)(1)=r1;(*this)(2)=r2;
00400 (*this)(3)=r3;(*this)(4)=r4; }
00401 Mat(const TCol& r0,const TCol& r1,const TCol& r2,
00402 const TCol& r3,const TCol& r4,const TCol& r5)
00403 { assert(N==6);(*this)(0)=r0;(*this)(1)=r1;(*this)(2)=r2;
00404 (*this)(3)=r3;(*this)(4)=r4;(*this)(5)=r5; }
00405
00406
00407 template <class EE, int SS> explicit Mat(const Vec<M,EE,SS>& r0)
00408 { assert(N==1); (*this)(0)=r0; }
00409 template <class EE, int SS> Mat(const Vec<M,EE,SS>& r0,const Vec<M,EE,SS>& r1)
00410 { assert(N==2);(*this)(0)=r0;(*this)(1)=r1; }
00411 template <class EE, int SS>
00412 Mat(const Vec<M,EE,SS>& r0,const Vec<M,EE,SS>& r1,const Vec<M,EE,SS>& r2)
00413 { assert(N==3);(*this)(0)=r0;(*this)(1)=r1;(*this)(2)=r2; }
00414 template <class EE, int SS>
00415 Mat(const Vec<M,EE,SS>& r0,const Vec<M,EE,SS>& r1,const Vec<M,EE,SS>& r2,
00416 const Vec<M,EE,SS>& r3)
00417 { assert(N==4);(*this)(0)=r0;(*this)(1)=r1;(*this)(2)=r2;(*this)(3)=r3; }
00418 template <class EE, int SS>
00419 Mat(const Vec<M,EE,SS>& r0,const Vec<M,EE,SS>& r1,const Vec<M,EE,SS>& r2,
00420 const Vec<M,EE,SS>& r3,const Vec<M,EE,SS>& r4)
00421 { assert(N==5);(*this)(0)=r0;(*this)(1)=r1;(*this)(2)=r2;
00422 (*this)(3)=r3;(*this)(4)=r4; }
00423 template <class EE, int SS>
00424 Mat(const Vec<M,EE,SS>& r0,const Vec<M,EE,SS>& r1,const Vec<M,EE,SS>& r2,
00425 const Vec<M,EE,SS>& r3,const Vec<M,EE,SS>& r4,const Vec<M,EE,SS>& r5)
00426 { assert(N==6);(*this)(0)=r0;(*this)(1)=r1;(*this)(2)=r2;
00427 (*this)(3)=r3;(*this)(4)=r4;(*this)(5)=r5; }
00428
00429
00430
00431 template <class EE> explicit Mat(const EE* p)
00432 { assert(p); for(int i=0;i<M;++i) (*this)[i]=&p[i*N]; }
00433
00434
00435
00436
00437 template <class EE, int CSS, int RSS> Mat& operator=(const Mat<M,N,EE,CSS,RSS>& mm) {
00438 for (int j=0; j<N; ++j) (*this)(j) = mm(j);
00439 return *this;
00440 }
00441
00442 template <class EE> Mat& operator=(const EE* p) {
00443 assert(p); for(int i=0;i<M;++i) (*this)[i]=&p[i*N];
00444 return *this;
00445 }
00446
00447
00448 template <class EE, int CSS, int RSS> Mat&
00449 operator+=(const Mat<M,N,EE,CSS,RSS>& mm) {
00450 for (int j=0; j<N; ++j) (*this)(j) += mm(j);
00451 return *this;
00452 }
00453 template <class EE, int CSS, int RSS> Mat&
00454 operator+=(const Mat<M,N,negator<EE>,CSS,RSS>& mm) {
00455 for (int j=0; j<N; ++j) (*this)(j) -= -(mm(j));
00456 return *this;
00457 }
00458
00459 template <class EE, int CSS, int RSS> Mat&
00460 operator-=(const Mat<M,N,EE,CSS,RSS>& mm) {
00461 for (int j=0; j<N; ++j) (*this)(j) -= mm(j);
00462 return *this;
00463 }
00464 template <class EE, int CSS, int RSS> Mat&
00465 operator-=(const Mat<M,N,negator<EE>,CSS,RSS>& mm) {
00466 for (int j=0; j<N; ++j) (*this)(j) += -(mm(j));
00467 return *this;
00468 }
00469
00470
00471
00472 template <class EE, int CSS, int RSS> Mat&
00473 operator*=(const Mat<N,N,EE,CSS,RSS>& mm) {
00474 const Mat t(*this);
00475 for (int j=0; j<N; ++j)
00476 for (int i=0; i<M; ++i)
00477 (*this)(i,j) = t[i] * mm(j);
00478 return *this;
00479 }
00480
00481
00482
00483
00484
00485 template <class E2, int CS2, int RS2>
00486 typename Result<Mat<M,N,E2,CS2,RS2> >::Add
00487 conformingAdd(const Mat<M,N,E2,CS2,RS2>& r) const {
00488 typename Result<Mat<M,N,E2,CS2,RS2> >::Add result;
00489 for (int j=0;j<N;++j) result(j) = (*this)(j) + r(j);
00490 return result;
00491 }
00492
00493 template <class E2, int CS2, int RS2>
00494 typename Result<Mat<M,N,E2,CS2,RS2> >::Sub
00495 conformingSubtract(const Mat<M,N,E2,CS2,RS2>& r) const {
00496 typename Result<Mat<M,N,E2,CS2,RS2> >::Sub result;
00497 for (int j=0;j<N;++j) result(j) = (*this)(j) - r(j);
00498 return result;
00499 }
00500
00501 template <class E2, int CS2, int RS2>
00502 typename Mat<M,N,E2,CS2,RS2>::template Result<Mat>::Sub
00503 conformingSubtractFromLeft(const Mat<M,N,E2,CS2,RS2>& l) const {
00504 return l.conformingSubtract(*this);
00505 }
00506
00507
00508
00509 template <class E2, int RS2>
00510 typename Result<SymMat<M,E2,RS2> >::Add
00511 conformingAdd(const SymMat<M,E2,RS2>& sy) const {
00512 assert(M==N);
00513 return sy.conformingAdd(*this);
00514 }
00515
00516 template <class E2, int RS2>
00517 typename Result<SymMat<M,E2,RS2> >::Sub
00518 conformingSubtract(const SymMat<M,E2,RS2>& sy) const {
00519 assert(M==N);
00520 return sy.conformingSubtractFromLeft(*this);
00521 }
00522
00523 template <class E2, int RS2>
00524 typename SymMat<M,E2,RS2>::template Result<Mat>::Sub
00525 conformingSubtractFromLeft(const SymMat<M,E2,RS2>& sy) const {
00526 assert(M==N);
00527 return sy.conformingSubtract(*this);
00528 }
00529
00530
00531 template <int N2, class E2, int CS2, int RS2>
00532 typename Result<Mat<N,N2,E2,CS2,RS2> >::Mul
00533 conformingMultiply(const Mat<N,N2,E2,CS2,RS2>& m) const {
00534 typename Result<Mat<N,N2,E2,CS2,RS2> >::Mul result;
00535 for (int j=0;j<N2;++j)
00536 for (int i=0;i<M;++i)
00537 result(i,j) = (*this)[i].conformingMultiply(m(j));
00538 return result;
00539 }
00540
00541 template <int M2, class E2, int CS2, int RS2>
00542 typename Mat<M2,M,E2,CS2,RS2>::template Result<Mat>::Mul
00543 conformingMultiplyFromLeft(const Mat<M2,M,E2,CS2,RS2>& m) const {
00544 return m.conformingMultiply(*this);
00545 }
00546
00547
00548 template <int M2, class E2, int CS2, int RS2>
00549 typename Result<Mat<M2,N,E2,CS2,RS2> >::Dvd
00550 conformingDivide(const Mat<M2,N,E2,CS2,RS2>& m) const {
00551 return conformingMultiply(m.invert());
00552 }
00553
00554 template <int M2, class E2, int CS2, int RS2>
00555 typename Mat<M2,N,E2,CS2,RS2>::template Result<Mat>::Dvd
00556 conformingDivideFromLeft(const Mat<M2,N,E2,CS2,RS2>& m) const {
00557 return m.conformingMultiply((*this).invert());
00558 }
00559
00560 const TRow& operator[](int i) const { return row(i); }
00561 TRow& operator[](int i) { return row(i); }
00562 const TCol& operator()(int j) const { return col(j); }
00563 TCol& operator()(int j) { return col(j); }
00564
00565 const E& operator()(int i,int j) const { return d[i*RS+j*CS]; }
00566 E& operator()(int i,int j) { return d[i*RS+j*CS]; }
00567
00568
00569 ScalarNormSq normSqr() const { return scalarNormSqr(); }
00570 typename CNT<ScalarNormSq>::TSqrt
00571 norm() const { return CNT<ScalarNormSq>::sqrt(scalarNormSqr()); }
00572
00573
00574
00575
00576
00577
00578
00579
00580
00581
00582
00583
00584
00585
00586 TNormalize normalize() const {
00587 if (CNT<E>::IsScalar) {
00588 return castAwayNegatorIfAny() / (SignInterpretation*norm());
00589 } else {
00590 TNormalize elementwiseNormalized;
00591
00592 for (int j=0; j<N; ++j)
00593 elementwiseNormalized(j) = (*this)(j).normalize();
00594 return elementwiseNormalized;
00595 }
00596 }
00597
00598
00599
00600 TInvert invert() const;
00601
00602 const Mat& operator+() const { return *this; }
00603 const TNeg& operator-() const { return negate(); }
00604 TNeg& operator-() { return updNegate(); }
00605 const THerm& operator~() const { return transpose(); }
00606 THerm& operator~() { return updTranspose(); }
00607
00608 const TNeg& negate() const { return *reinterpret_cast<const TNeg*>(this); }
00609 TNeg& updNegate() { return *reinterpret_cast<TNeg*>(this); }
00610
00611 const THerm& transpose() const { return *reinterpret_cast<const THerm*>(this); }
00612 THerm& updTranspose() { return *reinterpret_cast<THerm*>(this); }
00613
00614 const TPosTrans& positionalTranspose() const
00615 { return *reinterpret_cast<const TPosTrans*>(this); }
00616 TPosTrans& updPositionalTranspose()
00617 { return *reinterpret_cast<TPosTrans*>(this); }
00618
00619
00620
00621
00622
00623 const TReal& real() const { return *reinterpret_cast<const TReal*>(this); }
00624 TReal& real() { return *reinterpret_cast< TReal*>(this); }
00625
00626
00627
00628
00629
00630
00631
00632
00633 const TImag& imag() const {
00634 const int offs = ImagOffset;
00635 const Precision* p = reinterpret_cast<const Precision*>(this);
00636 return *reinterpret_cast<const TImag*>(p+offs);
00637 }
00638 TImag& imag() {
00639 const int offs = ImagOffset;
00640 Precision* p = reinterpret_cast<Precision*>(this);
00641 return *reinterpret_cast<TImag*>(p+offs);
00642 }
00643
00644 const TWithoutNegator& castAwayNegatorIfAny() const {return *reinterpret_cast<const TWithoutNegator*>(this);}
00645 TWithoutNegator& updCastAwayNegatorIfAny() {return *reinterpret_cast<TWithoutNegator*>(this);}
00646
00647 const TRow& row(int i) const
00648 { assert(0<=i&&i<M); return *reinterpret_cast<const TRow*>(&d[i*RS]); }
00649 TRow& row(int i)
00650 { assert(0<=i&&i<M); return *reinterpret_cast< TRow*>(&d[i*RS]); }
00651
00652 const TCol& col(int j) const
00653 { assert(0<=j&&j<N); return *reinterpret_cast<const TCol*>(&d[j*CS]); }
00654 TCol& col(int j)
00655 { assert(0<=j&&j<N); return *reinterpret_cast< TCol*>(&d[j*CS]); }
00656
00657
00658 const TDiag& diag() const { return *reinterpret_cast<const TDiag*>(d); }
00659 TDiag& diag() { return *reinterpret_cast<TDiag*>(d); }
00660
00661 EStandard trace() const {return diag().sum();}
00662
00663
00664
00665
00666
00667
00668
00669
00670 template <class EE> Mat<M,N, typename CNT<E>::template Result<EE>::Mul>
00671 scalarMultiply(const EE& e) const {
00672 Mat<M,N, typename CNT<E>::template Result<EE>::Mul> result;
00673 for (int j=0; j<N; ++j) result(j) = (*this)(j).scalarMultiply(e);
00674 return result;
00675 }
00676 template <class EE> Mat<M,N, typename CNT<EE>::template Result<E>::Mul>
00677 scalarMultiplyFromLeft(const EE& e) const {
00678 Mat<M,N, typename CNT<EE>::template Result<E>::Mul> result;
00679 for (int j=0; j<N; ++j) result(j) = (*this)(j).scalarMultiplyFromLeft(e);
00680 return result;
00681 }
00682
00683
00684
00685 template <class EE> Mat<M,N, typename CNT<E>::template Result<EE>::Dvd>
00686 scalarDivide(const EE& e) const {
00687 Mat<M,N, typename CNT<E>::template Result<EE>::Dvd> result;
00688 for (int j=0; j<N; ++j) result(j) = (*this)(j).scalarDivide(e);
00689 return result;
00690 }
00691 template <class EE> Mat<M,N, typename CNT<EE>::template Result<E>::Dvd>
00692 scalarDivideFromLeft(const EE& e) const {
00693 Mat<M,N, typename CNT<EE>::template Result<E>::Dvd> result;
00694 for (int j=0; j<N; ++j) result(j) = (*this)(j).scalarDivideFromLeft(e);
00695 return result;
00696 }
00697
00698
00699 template <class EE> Mat<M,N, typename CNT<E>::template Result<EE>::Add>
00700 scalarAdd(const EE& e) const {
00701 Mat<M,N, typename CNT<E>::template Result<EE>::Add> result(*this);
00702 result.diag() += e;
00703 return result;
00704 }
00705
00706
00707 template <class EE> Mat<M,N, typename CNT<E>::template Result<EE>::Sub>
00708 scalarSubtract(const EE& e) const {
00709 Mat<M,N, typename CNT<E>::template Result<EE>::Sub> result(*this);
00710 result.diag() -= e;
00711 return result;
00712 }
00713
00714 template <class EE> Mat<M,N, typename CNT<EE>::template Result<E>::Sub>
00715 scalarSubtractFromLeft(const EE& e) const {
00716 Mat<M,N, typename CNT<EE>::template Result<E>::Sub> result(-(*this));
00717 result.diag() += e;
00718 return result;
00719 }
00720
00721
00722
00723
00724 template <class EE> Mat& operator =(const EE& e) {return scalarEq(e);}
00725 template <class EE> Mat& operator+=(const EE& e) {return scalarPlusEq(e);}
00726 template <class EE> Mat& operator-=(const EE& e) {return scalarMinusEq(e);}
00727 template <class EE> Mat& operator*=(const EE& e) {return scalarTimesEq(e);}
00728 template <class EE> Mat& operator/=(const EE& e) {return scalarDivideEq(e);}
00729
00730
00731
00732 template <class EE> Mat& scalarEq(const EE& ee)
00733 { for(int j=0; j<N; ++j) (*this)(j).scalarEq(EE(0));
00734 diag().scalarEq(ee);
00735 return *this; }
00736
00737 template <class EE> Mat& scalarPlusEq(const EE& ee)
00738 { diag().scalarPlusEq(ee); return *this; }
00739
00740 template <class EE> Mat& scalarMinusEq(const EE& ee)
00741 { diag().scalarMinusEq(ee); return *this; }
00742
00743 template <class EE> Mat& scalarMinusEqFromLeft(const EE& ee)
00744 { scalarTimesEq(E(-1)); diag().scalarAdd(ee); return *this; }
00745
00746 template <class EE> Mat& scalarTimesEq(const EE& ee)
00747 { for(int j=0; j<N; ++j) (*this)(j).scalarTimesEq(ee); return *this; }
00748 template <class EE> Mat& scalarTimesEqFromLeft(const EE& ee)
00749 { for(int j=0; j<N; ++j) (*this)(j).scalarTimesEqFromLeft(ee); return *this; }
00750
00751 template <class EE> Mat& scalarDivideEq(const EE& ee)
00752 { for(int j=0; j<N; ++j) (*this)(j).scalarDivideEq(ee); return *this; }
00753 template <class EE> Mat& scalarDivideEqFromLeft(const EE& ee)
00754 { for(int j=0; j<N; ++j) (*this)(j).scalarDivideEqFromLeft(ee); return *this; }
00755
00756 void setToNaN() {
00757 for (int j=0; j<N; ++j)
00758 (*this)(j).setToNaN();
00759 }
00760
00761
00762
00763
00764 template <int MM, int NN> struct SubMat {
00765 typedef Mat<MM,NN,ELT,CS,RS> Type;
00766 };
00767
00768 template <int MM, int NN>
00769 const typename SubMat<MM,NN>::Type& getSubMat(int i, int j) const {
00770 assert(0 <= i && i + MM <= M);
00771 assert(0 <= j && j + NN <= N);
00772 return SubMat<MM,NN>::Type::getAs(&(*this)(i,j));
00773 }
00774 template <int MM, int NN>
00775 typename SubMat<MM,NN>::Type& updSubMat(int i, int j) {
00776 assert(0 <= i && i + MM <= M);
00777 assert(0 <= j && j + NN <= N);
00778 return SubMat<MM,NN>::Type::updAs(&(*this)(i,j));
00779 }
00780
00781
00784 TDropRow dropRow(int i) const {
00785 assert(0 <= i && i < M);
00786 TDropRow out;
00787 for (int r=0, nxt=0; r<M-1; ++r, ++nxt) {
00788 if (nxt==i) ++nxt;
00789 out[r] = (*this)[nxt];
00790 }
00791 return out;
00792 }
00793
00796 TDropCol dropCol(int j) const {
00797 assert(0 <= j && j < N);
00798 TDropCol out;
00799 for (int c=0, nxt=0; c<N-1; ++c, ++nxt) {
00800 if (nxt==j) ++nxt;
00801 out(c) = (*this)(nxt);
00802 }
00803 return out;
00804 }
00805
00809 TDropRowCol dropRowCol(int i, int j) const {
00810 assert(0 <= i && i < M);
00811 assert(0 <= j && j < N);
00812 TDropRowCol out;
00813 for (int c=0, nxtc=0; c<N-1; ++c, ++nxtc) {
00814 if (nxtc==j) ++nxtc;
00815 for (int r=0, nxtr=0; r<M-1; ++r, ++nxtr) {
00816 if (nxtr==i) ++nxtr;
00817 out(r,c) = (*this)(nxtr,nxtc);
00818 }
00819 }
00820 return out;
00821 }
00822
00826 template <class EE, int SS>
00827 TAppendRow appendRow(const Row<N,EE,SS>& row) const {
00828 TAppendRow out;
00829 out.updSubMat<M,N>(0,0) = (*this);
00830 out[M] = row;
00831 return out;
00832 }
00833
00837 template <class EE, int SS>
00838 TAppendCol appendCol(const Vec<M,EE,SS>& col) const {
00839 TAppendCol out;
00840 out.updSubMat<M,N>(0,0) = (*this);
00841 out(N) = col;
00842 return out;
00843 }
00844
00851 template <class ER, int SR, class EC, int SC>
00852 TAppendRowCol appendRowCol(const Row<N+1,ER,SR>& row,
00853 const Vec<M+1,EC,SC>& col) const
00854 {
00855 TAppendRowCol out;
00856 out.updSubMat<M,N>(0,0) = (*this);
00857 out[M].updSubRow<N>(0) = row.updSubRow<N>(0);
00858 out(N) = col;
00859 return out;
00860 }
00861
00867 template <class EE, int SS>
00868 TAppendRow insertRow(int i, const Row<N,EE,SS>& row) const {
00869 assert(0 <= i && i <= M);
00870 if (i==M) return appendRow(row);
00871 TAppendRow out;
00872 for (int r=0, nxt=0; r<M; ++r, ++nxt) {
00873 if (nxt==i) out[nxt++] = row;
00874 out[nxt] = (*this)[r];
00875 }
00876 return out;
00877 }
00878
00884 template <class EE, int SS>
00885 TAppendCol insertCol(int j, const Vec<M,EE,SS>& col) const {
00886 assert(0 <= j && j <= N);
00887 if (j==N) return appendCol(col);
00888 TAppendCol out;
00889 for (int c=0, nxt=0; c<N; ++c, ++nxt) {
00890 if (nxt==j) out(nxt++) = col;
00891 out(nxt) = (*this)(c);
00892 }
00893 return out;
00894 }
00895
00903 template <class ER, int SR, class EC, int SC>
00904 TAppendRowCol insertRowCol(int i, int j, const Row<N+1,ER,SR>& row,
00905 const Vec<M+1,EC,SC>& col) const {
00906 assert(0 <= i && i <= M);
00907 assert(0 <= j && j <= N);
00908 TAppendRowCol out;
00909 for (int c=0, nxtc=0; c<N; ++c, ++nxtc) {
00910 if (nxtc==i) ++nxtc;
00911 for (int r=0, nxtr=0; r<M; ++r, ++nxtr) {
00912 if (nxtr==i) ++nxtr;
00913 out(nxtr,nxtc) = (*this)(r,c);
00914 }
00915 }
00916 out[i] = row;
00917 out(j) = col;
00918 return out;
00919 }
00920
00921
00922 static const Mat& getAs(const ELT* p) {return *reinterpret_cast<const Mat*>(p);}
00923 static Mat& updAs(ELT* p) {return *reinterpret_cast<Mat*>(p);}
00924
00925
00926 static Mat<M,N,ELT,M,1> getNaN() {
00927 Mat<M,N,ELT,M,1> m;
00928 m.setToNaN();
00929 return m;
00930 }
00931
00932 TRow sum() const {
00933 TRow temp;
00934 for (int i = 0; i < N; ++i)
00935 temp[i] = col(i).sum();
00936 return temp;
00937 }
00938
00939 private:
00940 E d[NActualElements];
00941
00942
00943
00944
00945
00946 int rIx(int k) const {
00947 const int row = k / N;
00948 const int col = k % N;
00949 return row*RS + col*CS;
00950 }
00951 };
00952
00954
00955
00957
00958 template <int M, int N, class EL, int CSL, int RSL, class ER, int CSR, int RSR> inline
00959 typename Mat<M,N,EL,CSL,RSL>::template Result<Mat<M,N,ER,CSR,RSR> >::Add
00960 operator+(const Mat<M,N,EL,CSL,RSL>& l, const Mat<M,N,ER,CSR,RSR>& r) {
00961 return Mat<M,N,EL,CSL,RSL>::template Result<Mat<M,N,ER,CSR,RSR> >
00962 ::AddOp::perform(l,r);
00963 }
00964
00965 template <int M, int N, class EL, int CSL, int RSL, class ER, int CSR, int RSR> inline
00966 typename Mat<M,N,EL,CSL,RSL>::template Result<Mat<M,N,ER,CSR,RSR> >::Sub
00967 operator-(const Mat<M,N,EL,CSL,RSL>& l, const Mat<M,N,ER,CSR,RSR>& r) {
00968 return Mat<M,N,EL,CSL,RSL>::template Result<Mat<M,N,ER,CSR,RSR> >
00969 ::SubOp::perform(l,r);
00970 }
00971
00972
00973 template <int M, int N, class EL, int CSL, int RSL, int P, class ER, int CSR, int RSR> inline
00974 typename Mat<M,N,EL,CSL,RSL>::template Result<Mat<N,P,ER,CSR,RSR> >::Mul
00975 operator*(const Mat<M,N,EL,CSL,RSL>& l, const Mat<N,P,ER,CSR,RSR>& r) {
00976 return Mat<M,N,EL,CSL,RSL>::template Result<Mat<N,P,ER,CSR,RSR> >
00977 ::MulOp::perform(l,r);
00978 }
00979
00980
00981
00982 template <int M, int N, class EL, int CSL, int RSL, int MM, int NN, class ER, int CSR, int RSR> inline
00983 typename Mat<M,N,EL,CSL,RSL>::template Result<Mat<MM,NN,ER,CSR,RSR> >::MulNon
00984 operator*(const Mat<M,N,EL,CSL,RSL>& l, const Mat<MM,NN,ER,CSR,RSR>& r) {
00985 return Mat<M,N,EL,CSL,RSL>::template Result<Mat<MM,NN,ER,CSR,RSR> >
00986 ::MulOpNonConforming::perform(l,r);
00987 }
00988
00989 template <int M, int N, class EL, int CSL, int RSL, class ER, int CSR, int RSR> inline
00990 bool operator==(const Mat<M,N,EL,CSL,RSL>& l, const Mat<M,N,ER,CSR,RSR>& r) {
00991 for (int j=0; j<N; ++j)
00992 if (l(j) != r(j)) return false;
00993 return true;
00994 }
00995 template <int M, int N, class EL, int CSL, int RSL, class ER, int CSR, int RSR> inline
00996 bool operator!=(const Mat<M,N,EL,CSL,RSL>& l, const Mat<M,N,ER,CSR,RSR>& r) {
00997 return !(l==r);
00998 }
00999
01000
01002
01004
01005
01006
01007
01008 template <int M, int N, class E, int CS, int RS> inline
01009 typename Mat<M,N,E,CS,RS>::template Result<float>::Mul
01010 operator*(const Mat<M,N,E,CS,RS>& l, const float& r)
01011 { return Mat<M,N,E,CS,RS>::template Result<float>::MulOp::perform(l,r); }
01012 template <int M, int N, class E, int CS, int RS> inline
01013 typename Mat<M,N,E,CS,RS>::template Result<float>::Mul
01014 operator*(const float& l, const Mat<M,N,E,CS,RS>& r) {return r*l;}
01015
01016 template <int M, int N, class E, int CS, int RS> inline
01017 typename Mat<M,N,E,CS,RS>::template Result<double>::Mul
01018 operator*(const Mat<M,N,E,CS,RS>& l, const double& r)
01019 { return Mat<M,N,E,CS,RS>::template Result<double>::MulOp::perform(l,r); }
01020 template <int M, int N, class E, int CS, int RS> inline
01021 typename Mat<M,N,E,CS,RS>::template Result<double>::Mul
01022 operator*(const double& l, const Mat<M,N,E,CS,RS>& r) {return r*l;}
01023
01024 template <int M, int N, class E, int CS, int RS> inline
01025 typename Mat<M,N,E,CS,RS>::template Result<long double>::Mul
01026 operator*(const Mat<M,N,E,CS,RS>& l, const long double& r)
01027 { return Mat<M,N,E,CS,RS>::template Result<long double>::MulOp::perform(l,r); }
01028 template <int M, int N, class E, int CS, int RS> inline
01029 typename Mat<M,N,E,CS,RS>::template Result<long double>::Mul
01030 operator*(const long double& l, const Mat<M,N,E,CS,RS>& r) {return r*l;}
01031
01032
01033 template <int M, int N, class E, int CS, int RS> inline
01034 typename Mat<M,N,E,CS,RS>::template Result<typename CNT<E>::Precision>::Mul
01035 operator*(const Mat<M,N,E,CS,RS>& l, int r) {return l * (typename CNT<E>::Precision)r;}
01036 template <int M, int N, class E, int CS, int RS> inline
01037 typename Mat<M,N,E,CS,RS>::template Result<typename CNT<E>::Precision>::Mul
01038 operator*(int l, const Mat<M,N,E,CS,RS>& r) {return r * (typename CNT<E>::Precision)l;}
01039
01040
01041
01042
01043 template <int M, int N, class E, int CS, int RS, class R> inline
01044 typename Mat<M,N,E,CS,RS>::template Result<std::complex<R> >::Mul
01045 operator*(const Mat<M,N,E,CS,RS>& l, const std::complex<R>& r)
01046 { return Mat<M,N,E,CS,RS>::template Result<std::complex<R> >::MulOp::perform(l,r); }
01047 template <int M, int N, class E, int CS, int RS, class R> inline
01048 typename Mat<M,N,E,CS,RS>::template Result<std::complex<R> >::Mul
01049 operator*(const std::complex<R>& l, const Mat<M,N,E,CS,RS>& r) {return r*l;}
01050
01051
01052 template <int M, int N, class E, int CS, int RS, class R> inline
01053 typename Mat<M,N,E,CS,RS>::template Result<std::complex<R> >::Mul
01054 operator*(const Mat<M,N,E,CS,RS>& l, const conjugate<R>& r) {return l*(std::complex<R>)r;}
01055 template <int M, int N, class E, int CS, int RS, class R> inline
01056 typename Mat<M,N,E,CS,RS>::template Result<std::complex<R> >::Mul
01057 operator*(const conjugate<R>& l, const Mat<M,N,E,CS,RS>& r) {return r*(std::complex<R>)l;}
01058
01059
01060 template <int M, int N, class E, int CS, int RS, class R> inline
01061 typename Mat<M,N,E,CS,RS>::template Result<typename negator<R>::StdNumber>::Mul
01062 operator*(const Mat<M,N,E,CS,RS>& l, const negator<R>& r) {return l * (typename negator<R>::StdNumber)(R)r;}
01063 template <int M, int N, class E, int CS, int RS, class R> inline
01064 typename Mat<M,N,E,CS,RS>::template Result<typename negator<R>::StdNumber>::Mul
01065 operator*(const negator<R>& l, const Mat<M,N,E,CS,RS>& r) {return r * (typename negator<R>::StdNumber)(R)l;}
01066
01067
01068
01069
01070
01071
01072
01073 template <int M, int N, class E, int CS, int RS> inline
01074 typename Mat<M,N,E,CS,RS>::template Result<float>::Dvd
01075 operator/(const Mat<M,N,E,CS,RS>& l, const float& r)
01076 { return Mat<M,N,E,CS,RS>::template Result<float>::DvdOp::perform(l,r); }
01077 template <int M, int N, class E, int CS, int RS> inline
01078 typename CNT<float>::template Result<Mat<M,N,E,CS,RS> >::Dvd
01079 operator/(const float& l, const Mat<M,N,E,CS,RS>& r)
01080 { return CNT<float>::template Result<Mat<M,N,E,CS,RS> >::DvdOp::perform(l,r); }
01081
01082 template <int M, int N, class E, int CS, int RS> inline
01083 typename Mat<M,N,E,CS,RS>::template Result<double>::Dvd
01084 operator/(const Mat<M,N,E,CS,RS>& l, const double& r)
01085 { return Mat<M,N,E,CS,RS>::template Result<double>::DvdOp::perform(l,r); }
01086 template <int M, int N, class E, int CS, int RS> inline
01087 typename CNT<double>::template Result<Mat<M,N,E,CS,RS> >::Dvd
01088 operator/(const double& l, const Mat<M,N,E,CS,RS>& r)
01089 { return CNT<double>::template Result<Mat<M,N,E,CS,RS> >::DvdOp::perform(l,r); }
01090
01091 template <int M, int N, class E, int CS, int RS> inline
01092 typename Mat<M,N,E,CS,RS>::template Result<long double>::Dvd
01093 operator/(const Mat<M,N,E,CS,RS>& l, const long double& r)
01094 { return Mat<M,N,E,CS,RS>::template Result<long double>::DvdOp::perform(l,r); }
01095 template <int M, int N, class E, int CS, int RS> inline
01096 typename CNT<long double>::template Result<Mat<M,N,E,CS,RS> >::Dvd
01097 operator/(const long double& l, const Mat<M,N,E,CS,RS>& r)
01098 { return CNT<long double>::template Result<Mat<M,N,E,CS,RS> >::DvdOp::perform(l,r); }
01099
01100
01101 template <int M, int N, class E, int CS, int RS> inline
01102 typename Mat<M,N,E,CS,RS>::template Result<typename CNT<E>::Precision>::Dvd
01103 operator/(const Mat<M,N,E,CS,RS>& l, int r) {return l / (typename CNT<E>::Precision)r;}
01104 template <int M, int N, class E, int CS, int RS> inline
01105 typename CNT<typename CNT<E>::Precision>::template Result<Mat<M,N,E,CS,RS> >::Dvd
01106 operator/(int l, const Mat<M,N,E,CS,RS>& r) {return (typename CNT<E>::Precision)l / r;}
01107
01108
01109
01110
01111
01112 template <int M, int N, class E, int CS, int RS, class R> inline
01113 typename Mat<M,N,E,CS,RS>::template Result<std::complex<R> >::Dvd
01114 operator/(const Mat<M,N,E,CS,RS>& l, const std::complex<R>& r)
01115 { return Mat<M,N,E,CS,RS>::template Result<std::complex<R> >::DvdOp::perform(l,r); }
01116 template <int M, int N, class E, int CS, int RS, class R> inline
01117 typename CNT<std::complex<R> >::template Result<Mat<M,N,E,CS,RS> >::Dvd
01118 operator/(const std::complex<R>& l, const Mat<M,N,E,CS,RS>& r)
01119 { return CNT<std::complex<R> >::template Result<Mat<M,N,E,CS,RS> >::DvdOp::perform(l,r); }
01120
01121
01122 template <int M, int N, class E, int CS, int RS, class R> inline
01123 typename Mat<M,N,E,CS,RS>::template Result<std::complex<R> >::Dvd
01124 operator/(const Mat<M,N,E,CS,RS>& l, const conjugate<R>& r) {return l/(std::complex<R>)r;}
01125 template <int M, int N, class E, int CS, int RS, class R> inline
01126 typename CNT<std::complex<R> >::template Result<Mat<M,N,E,CS,RS> >::Dvd
01127 operator/(const conjugate<R>& l, const Mat<M,N,E,CS,RS>& r) {return (std::complex<R>)l/r;}
01128
01129
01130 template <int M, int N, class E, int CS, int RS, class R> inline
01131 typename Mat<M,N,E,CS,RS>::template Result<typename negator<R>::StdNumber>::Dvd
01132 operator/(const Mat<M,N,E,CS,RS>& l, const negator<R>& r) {return l/(typename negator<R>::StdNumber)(R)r;}
01133 template <int M, int N, class E, int CS, int RS, class R> inline
01134 typename CNT<R>::template Result<Mat<M,N,E,CS,RS> >::Dvd
01135 operator/(const negator<R>& l, const Mat<M,N,E,CS,RS>& r) {return (typename negator<R>::StdNumber)(R)l/r;}
01136
01137
01138
01139
01140
01141
01142
01143
01144
01145 template <int M, int N, class E, int CS, int RS> inline
01146 typename Mat<M,N,E,CS,RS>::template Result<float>::Add
01147 operator+(const Mat<M,N,E,CS,RS>& l, const float& r)
01148 { return Mat<M,N,E,CS,RS>::template Result<float>::AddOp::perform(l,r); }
01149 template <int M, int N, class E, int CS, int RS> inline
01150 typename Mat<M,N,E,CS,RS>::template Result<float>::Add
01151 operator+(const float& l, const Mat<M,N,E,CS,RS>& r) {return r+l;}
01152
01153 template <int M, int N, class E, int CS, int RS> inline
01154 typename Mat<M,N,E,CS,RS>::template Result<double>::Add
01155 operator+(const Mat<M,N,E,CS,RS>& l, const double& r)
01156 { return Mat<M,N,E,CS,RS>::template Result<double>::AddOp::perform(l,r); }
01157 template <int M, int N, class E, int CS, int RS> inline
01158 typename Mat<M,N,E,CS,RS>::template Result<double>::Add
01159 operator+(const double& l, const Mat<M,N,E,CS,RS>& r) {return r+l;}
01160
01161 template <int M, int N, class E, int CS, int RS> inline
01162 typename Mat<M,N,E,CS,RS>::template Result<long double>::Add
01163 operator+(const Mat<M,N,E,CS,RS>& l, const long double& r)
01164 { return Mat<M,N,E,CS,RS>::template Result<long double>::AddOp::perform(l,r); }
01165 template <int M, int N, class E, int CS, int RS> inline
01166 typename Mat<M,N,E,CS,RS>::template Result<long double>::Add
01167 operator+(const long double& l, const Mat<M,N,E,CS,RS>& r) {return r+l;}
01168
01169
01170 template <int M, int N, class E, int CS, int RS> inline
01171 typename Mat<M,N,E,CS,RS>::template Result<typename CNT<E>::Precision>::Add
01172 operator+(const Mat<M,N,E,CS,RS>& l, int r) {return l + (typename CNT<E>::Precision)r;}
01173 template <int M, int N, class E, int CS, int RS> inline
01174 typename Mat<M,N,E,CS,RS>::template Result<typename CNT<E>::Precision>::Add
01175 operator+(int l, const Mat<M,N,E,CS,RS>& r) {return r + (typename CNT<E>::Precision)l;}
01176
01177
01178
01179
01180 template <int M, int N, class E, int CS, int RS, class R> inline
01181 typename Mat<M,N,E,CS,RS>::template Result<std::complex<R> >::Add
01182 operator+(const Mat<M,N,E,CS,RS>& l, const std::complex<R>& r)
01183 { return Mat<M,N,E,CS,RS>::template Result<std::complex<R> >::AddOp::perform(l,r); }
01184 template <int M, int N, class E, int CS, int RS, class R> inline
01185 typename Mat<M,N,E,CS,RS>::template Result<std::complex<R> >::Add
01186 operator+(const std::complex<R>& l, const Mat<M,N,E,CS,RS>& r) {return r+l;}
01187
01188
01189 template <int M, int N, class E, int CS, int RS, class R> inline
01190 typename Mat<M,N,E,CS,RS>::template Result<std::complex<R> >::Add
01191 operator+(const Mat<M,N,E,CS,RS>& l, const conjugate<R>& r) {return l+(std::complex<R>)r;}
01192 template <int M, int N, class E, int CS, int RS, class R> inline
01193 typename Mat<M,N,E,CS,RS>::template Result<std::complex<R> >::Add
01194 operator+(const conjugate<R>& l, const Mat<M,N,E,CS,RS>& r) {return r+(std::complex<R>)l;}
01195
01196
01197 template <int M, int N, class E, int CS, int RS, class R> inline
01198 typename Mat<M,N,E,CS,RS>::template Result<typename negator<R>::StdNumber>::Add
01199 operator+(const Mat<M,N,E,CS,RS>& l, const negator<R>& r) {return l + (typename negator<R>::StdNumber)(R)r;}
01200 template <int M, int N, class E, int CS, int RS, class R> inline
01201 typename Mat<M,N,E,CS,RS>::template Result<typename negator<R>::StdNumber>::Add
01202 operator+(const negator<R>& l, const Mat<M,N,E,CS,RS>& r) {return r + (typename negator<R>::StdNumber)(R)l;}
01203
01204
01205
01206
01207 template <int M, int N, class E, int CS, int RS> inline
01208 typename Mat<M,N,E,CS,RS>::template Result<float>::Sub
01209 operator-(const Mat<M,N,E,CS,RS>& l, const float& r)
01210 { return Mat<M,N,E,CS,RS>::template Result<float>::SubOp::perform(l,r); }
01211 template <int M, int N, class E, int CS, int RS> inline
01212 typename CNT<float>::template Result<Mat<M,N,E,CS,RS> >::Sub
01213 operator-(const float& l, const Mat<M,N,E,CS,RS>& r)
01214 { return CNT<float>::template Result<Mat<M,N,E,CS,RS> >::SubOp::perform(l,r); }
01215
01216 template <int M, int N, class E, int CS, int RS> inline
01217 typename Mat<M,N,E,CS,RS>::template Result<double>::Sub
01218 operator-(const Mat<M,N,E,CS,RS>& l, const double& r)
01219 { return Mat<M,N,E,CS,RS>::template Result<double>::SubOp::perform(l,r); }
01220 template <int M, int N, class E, int CS, int RS> inline
01221 typename CNT<double>::template Result<Mat<M,N,E,CS,RS> >::Sub
01222 operator-(const double& l, const Mat<M,N,E,CS,RS>& r)
01223 { return CNT<double>::template Result<Mat<M,N,E,CS,RS> >::SubOp::perform(l,r); }
01224
01225 template <int M, int N, class E, int CS, int RS> inline
01226 typename Mat<M,N,E,CS,RS>::template Result<long double>::Sub
01227 operator-(const Mat<M,N,E,CS,RS>& l, const long double& r)
01228 { return Mat<M,N,E,CS,RS>::template Result<long double>::SubOp::perform(l,r); }
01229 template <int M, int N, class E, int CS, int RS> inline
01230 typename CNT<long double>::template Result<Mat<M,N,E,CS,RS> >::Sub
01231 operator-(const long double& l, const Mat<M,N,E,CS,RS>& r)
01232 { return CNT<long double>::template Result<Mat<M,N,E,CS,RS> >::SubOp::perform(l,r); }
01233
01234
01235 template <int M, int N, class E, int CS, int RS> inline
01236 typename Mat<M,N,E,CS,RS>::template Result<typename CNT<E>::Precision>::Sub
01237 operator-(const Mat<M,N,E,CS,RS>& l, int r) {return l - (typename CNT<E>::Precision)r;}
01238 template <int M, int N, class E, int CS, int RS> inline
01239 typename CNT<typename CNT<E>::Precision>::template Result<Mat<M,N,E,CS,RS> >::Sub
01240 operator-(int l, const Mat<M,N,E,CS,RS>& r) {return (typename CNT<E>::Precision)l - r;}
01241
01242
01243
01244
01245
01246 template <int M, int N, class E, int CS, int RS, class R> inline
01247 typename Mat<M,N,E,CS,RS>::template Result<std::complex<R> >::Sub
01248 operator-(const Mat<M,N,E,CS,RS>& l, const std::complex<R>& r)
01249 { return Mat<M,N,E,CS,RS>::template Result<std::complex<R> >::SubOp::perform(l,r); }
01250 template <int M, int N, class E, int CS, int RS, class R> inline
01251 typename CNT<std::complex<R> >::template Result<Mat<M,N,E,CS,RS> >::Sub
01252 operator-(const std::complex<R>& l, const Mat<M,N,E,CS,RS>& r)
01253 { return CNT<std::complex<R> >::template Result<Mat<M,N,E,CS,RS> >::SubOp::perform(l,r); }
01254
01255
01256 template <int M, int N, class E, int CS, int RS, class R> inline
01257 typename Mat<M,N,E,CS,RS>::template Result<std::complex<R> >::Sub
01258 operator-(const Mat<M,N,E,CS,RS>& l, const conjugate<R>& r) {return l-(std::complex<R>)r;}
01259 template <int M, int N, class E, int CS, int RS, class R> inline
01260 typename CNT<std::complex<R> >::template Result<Mat<M,N,E,CS,RS> >::Sub
01261 operator-(const conjugate<R>& l, const Mat<M,N,E,CS,RS>& r) {return (std::complex<R>)l-r;}
01262
01263
01264 template <int M, int N, class E, int CS, int RS, class R> inline
01265 typename Mat<M,N,E,CS,RS>::template Result<typename negator<R>::StdNumber>::Sub
01266 operator-(const Mat<M,N,E,CS,RS>& l, const negator<R>& r) {return l-(typename negator<R>::StdNumber)(R)r;}
01267 template <int M, int N, class E, int CS, int RS, class R> inline
01268 typename CNT<R>::template Result<Mat<M,N,E,CS,RS> >::Sub
01269 operator-(const negator<R>& l, const Mat<M,N,E,CS,RS>& r) {return (typename negator<R>::StdNumber)(R)l-r;}
01270
01271
01272
01273 template <int M, int N, class E, int CS, int RS, class CHAR, class TRAITS> inline
01274 std::basic_ostream<CHAR,TRAITS>&
01275 operator<<(std::basic_ostream<CHAR,TRAITS>& o, const Mat<M,N,E,CS,RS>& m) {
01276 for (int i=0;i<M;++i) {
01277 o << std::endl << "[";
01278 for (int j=0;j<N;++j)
01279 o << (j>0?",":"") << m(i,j);
01280 o << "]";
01281 }
01282 if (M) o << std::endl;
01283 return o;
01284 }
01285
01286 template <int M, int N, class E, int CS, int RS, class CHAR, class TRAITS> inline
01287 std::basic_istream<CHAR,TRAITS>&
01288 operator>>(std::basic_istream<CHAR,TRAITS>& is, Mat<M,N,E,CS,RS>& m) {
01289
01290 assert(false);
01291 return is;
01292 }
01293
01294 }
01295
01296
01297 #endif //SimTK_SIMMATRIX_SMALLMATRIX_MAT_H_