00001 #ifndef SimTK_SIMMATRIX_SMALLMATRIX_SYMMAT_H_
00002 #define SimTK_SIMMATRIX_SMALLMATRIX_SYMMAT_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
00074 #include "SimTKcommon/internal/common.h"
00075
00076 namespace SimTK {
00077
00079 template <int M, class ELT, int RS> class SymMat {
00080 typedef ELT E;
00081 typedef typename CNT<E>::TNeg ENeg;
00082 typedef typename CNT<E>::TWithoutNegator EWithoutNegator;
00083 typedef typename CNT<E>::TReal EReal;
00084 typedef typename CNT<E>::TImag EImag;
00085 typedef typename CNT<E>::TComplex EComplex;
00086 typedef typename CNT<E>::THerm EHerm;
00087 typedef typename CNT<E>::TPosTrans EPosTrans;
00088 typedef typename CNT<E>::TSqHermT ESqHermT;
00089 typedef typename CNT<E>::TSqTHerm ESqTHerm;
00090
00091 typedef typename CNT<E>::TSqrt ESqrt;
00092 typedef typename CNT<E>::TAbs EAbs;
00093 typedef typename CNT<E>::TStandard EStandard;
00094 typedef typename CNT<E>::TInvert EInvert;
00095 typedef typename CNT<E>::TNormalize ENormalize;
00096
00097 typedef typename CNT<E>::Scalar EScalar;
00098 typedef typename CNT<E>::ULessScalar EULessScalar;
00099 typedef typename CNT<E>::Number ENumber;
00100 typedef typename CNT<E>::StdNumber EStdNumber;
00101 typedef typename CNT<E>::Precision EPrecision;
00102 typedef typename CNT<E>::ScalarNormSq EScalarNormSq;
00103
00104 public:
00105
00106 enum {
00107 NRows = M,
00108 NCols = M,
00109 NPackedElements = (M*(M+1))/2,
00110 NActualElements = RS * NPackedElements,
00111 NActualScalars = CNT<E>::NActualScalars * NActualElements,
00112 RowSpacing = RS,
00113 ColSpacing = NActualElements,
00114 ImagOffset = NTraits<ENumber>::ImagOffset,
00115 RealStrideFactor = 1,
00116
00117 ArgDepth = ((int)CNT<E>::ArgDepth < (int)MAX_RESOLVED_DEPTH
00118 ? CNT<E>::ArgDepth + 1
00119 : MAX_RESOLVED_DEPTH),
00120 IsScalar = 0,
00121 IsULessScalar = 0,
00122 IsNumber = 0,
00123 IsStdNumber = 0,
00124 IsPrecision = 0,
00125 SignInterpretation = CNT<E>::SignInterpretation
00126 };
00127
00128 typedef SymMat<M,E,RS> T;
00129 typedef SymMat<M,ENeg,RS> TNeg;
00130 typedef SymMat<M,EWithoutNegator,RS> TWithoutNegator;
00131
00132 typedef SymMat<M,EReal,RS*CNT<E>::RealStrideFactor>
00133 TReal;
00134 typedef SymMat<M,EImag,RS*CNT<E>::RealStrideFactor>
00135 TImag;
00136 typedef SymMat<M,EComplex,RS> TComplex;
00137 typedef T THerm;
00138 typedef SymMat<M,EHerm,RS> TPosTrans;
00139 typedef E TElement;
00140 typedef Vec<M,E,RS> TDiag;
00141 typedef Vec<(M*(M-1))/2,E,RS> TLower;
00142 typedef Vec<(M*(M-1))/2,EHerm,RS> TUpper;
00143 typedef Vec<(M*(M+1))/2,E,RS> TAsVec;
00144
00145
00146
00147 typedef Row<M,E,1> TRow;
00148 typedef Vec<M,E,1> TCol;
00149 typedef SymMat<M,ESqrt,1> TSqrt;
00150 typedef SymMat<M,EAbs,1> TAbs;
00151 typedef SymMat<M,EStandard,1> TStandard;
00152 typedef SymMat<M,EInvert,1> TInvert;
00153 typedef SymMat<M,ENormalize,1> TNormalize;
00154 typedef SymMat<M,ESqHermT,1> TSqHermT;
00155 typedef SymMat<M,ESqTHerm,1> TSqTHerm;
00156 typedef SymMat<M,E,1> TPacked;
00157
00158 typedef EScalar Scalar;
00159 typedef EULessScalar ULessScalar;
00160 typedef ENumber Number;
00161 typedef EStdNumber StdNumber;
00162 typedef EPrecision Precision;
00163 typedef EScalarNormSq ScalarNormSq;
00164
00165 int size() const { return (M*(M+1))/2; }
00166 int nrow() const { return M; }
00167 int ncol() const { return M; }
00168
00169
00170
00171 ScalarNormSq scalarNormSqr() const {
00172 return getDiag().scalarNormSqr() + 2.*getLower().scalarNormSqr();
00173 }
00174
00175
00176
00177
00178 TSqrt sqrt() const {
00179 return TSqrt(getAsVec().sqrt());
00180 }
00181
00182
00183
00184
00185 TAbs abs() const {
00186 return TAbs(getAsVec().abs());
00187 }
00188
00189 TStandard standardize() const {
00190 return TStandard(getAsVec().standardize());
00191 }
00192
00193 EStandard trace() const {return getDiag().sum();}
00194
00195
00196
00197
00198 template <class P> struct EltResult {
00199 typedef SymMat<M, typename CNT<E>::template Result<P>::Mul, 1> Mul;
00200 typedef SymMat<M, typename CNT<E>::template Result<P>::Dvd, 1> Dvd;
00201 typedef SymMat<M, typename CNT<E>::template Result<P>::Add, 1> Add;
00202 typedef SymMat<M, typename CNT<E>::template Result<P>::Sub, 1> Sub;
00203 };
00204
00205
00206
00207 template <class P> struct Result {
00208 typedef MulCNTs<M,M,ArgDepth,SymMat,ColSpacing,RowSpacing,
00209 CNT<P>::NRows, CNT<P>::NCols, CNT<P>::ArgDepth,
00210 P, CNT<P>::ColSpacing, CNT<P>::RowSpacing> MulOp;
00211 typedef typename MulOp::Type Mul;
00212
00213 typedef MulCNTsNonConforming<M,M,ArgDepth,SymMat,ColSpacing,RowSpacing,
00214 CNT<P>::NRows, CNT<P>::NCols, CNT<P>::ArgDepth,
00215 P, CNT<P>::ColSpacing, CNT<P>::RowSpacing> MulOpNonConforming;
00216 typedef typename MulOpNonConforming::Type MulNon;
00217
00218
00219 typedef DvdCNTs<M,M,ArgDepth,SymMat,ColSpacing,RowSpacing,
00220 CNT<P>::NRows, CNT<P>::NCols, CNT<P>::ArgDepth,
00221 P, CNT<P>::ColSpacing, CNT<P>::RowSpacing> DvdOp;
00222 typedef typename DvdOp::Type Dvd;
00223
00224 typedef AddCNTs<M,M,ArgDepth,SymMat,ColSpacing,RowSpacing,
00225 CNT<P>::NRows, CNT<P>::NCols, CNT<P>::ArgDepth,
00226 P, CNT<P>::ColSpacing, CNT<P>::RowSpacing> AddOp;
00227 typedef typename AddOp::Type Add;
00228
00229 typedef SubCNTs<M,M,ArgDepth,SymMat,ColSpacing,RowSpacing,
00230 CNT<P>::NRows, CNT<P>::NCols, CNT<P>::ArgDepth,
00231 P, CNT<P>::ColSpacing, CNT<P>::RowSpacing> SubOp;
00232 typedef typename SubOp::Type Sub;
00233 };
00234
00235
00236 template <class P> struct Substitute {
00237 typedef SymMat<M,P> Type;
00238 };
00239
00240
00241
00242 SymMat(){
00243 #ifndef NDEBUG
00244 setToNaN();
00245 #endif
00246 }
00247
00248 SymMat(const SymMat& src) {
00249 updAsVec() = src.getAsVec();
00250 }
00251
00252 SymMat& operator=(const SymMat& src) {
00253 updAsVec() = src.getAsVec();
00254 return *this;
00255 }
00256
00257
00258
00259 template <class EE, int CSS, int RSS>
00260 explicit SymMat(const Mat<M,M,EE,CSS,RSS>& m) {
00261 updDiag() = m.diag().real();
00262 for (int j=0; j<M; ++j)
00263 for (int i=j+1; i<M; ++i)
00264 updEltLower(i,j) = m(i,j);
00265 }
00266
00267
00268
00269 template <int RSS> SymMat(const SymMat<M,E,RSS>& src)
00270 { updAsVec() = src.getAsVec(); }
00271
00272
00273
00274 template <int RSS> SymMat(const SymMat<M,ENeg,RSS>& src)
00275 { updAsVec() = src.getAsVec(); }
00276
00277
00278
00279 template <class EE, int RSS> explicit SymMat(const SymMat<M,EE,RSS>& src)
00280 { updAsVec() = src.getAsVec(); }
00281
00282
00283
00284
00285 explicit SymMat(const E& e)
00286 { updDiag() = e; updLower() = E(0); }
00287
00301
00302 SymMat(const E& e0,
00303 const E& e1,const E& e2)
00304 { assert(M==2); TDiag& d=updDiag(); TLower& l=updLower();
00305 d[0]=CNT<E>::real(e0); d[1]=CNT<E>::real(e2);
00306 l[0]=e1; }
00307
00308 SymMat(const E& e0,
00309 const E& e1,const E& e2,
00310 const E& e3,const E& e4, const E& e5)
00311 { assert(M==3); TDiag& d=updDiag(); TLower& l=updLower();
00312 d[0]=CNT<E>::real(e0);d[1]=CNT<E>::real(e2);d[2]=CNT<E>::real(e5);
00313 l[0]=e1;l[1]=e3;
00314 l[2]=e4; }
00315
00316 SymMat(const E& e0,
00317 const E& e1,const E& e2,
00318 const E& e3,const E& e4, const E& e5,
00319 const E& e6,const E& e7, const E& e8, const E& e9)
00320 { assert(M==4); TDiag& d=updDiag(); TLower& l=updLower();
00321 d[0]=CNT<E>::real(e0);d[1]=CNT<E>::real(e2);d[2]=CNT<E>::real(e5);d[3]=CNT<E>::real(e9);
00322 l[0]=e1;l[1]=e3;l[2]=e6;
00323 l[3]=e4;l[4]=e7;
00324 l[5]=e8; }
00325
00326 SymMat(const E& e0,
00327 const E& e1, const E& e2,
00328 const E& e3, const E& e4, const E& e5,
00329 const E& e6, const E& e7, const E& e8, const E& e9,
00330 const E& e10,const E& e11, const E& e12, const E& e13, const E& e14)
00331 { assert(M==5); TDiag& d=updDiag(); TLower& l=updLower();
00332 d[0]=CNT<E>::real(e0);d[1]=CNT<E>::real(e2);d[2]=CNT<E>::real(e5);d[3]=CNT<E>::real(e9);d[4]=CNT<E>::real(e14);
00333 l[0]=e1;l[1]=e3;l[2]=e6;l[3]=e10;
00334 l[4]=e4;l[5]=e7;l[6]=e11;
00335 l[7]=e8;l[8]=e12;
00336 l[9]=e13; }
00337
00338 SymMat(const E& e0,
00339 const E& e1, const E& e2,
00340 const E& e3, const E& e4, const E& e5,
00341 const E& e6, const E& e7, const E& e8, const E& e9,
00342 const E& e10,const E& e11, const E& e12, const E& e13, const E& e14,
00343 const E& e15,const E& e16, const E& e17, const E& e18, const E& e19, const E& e20)
00344 { assert(M==6); TDiag& d=updDiag(); TLower& l=updLower();
00345 d[0]=CNT<E>::real(e0);d[1]=CNT<E>::real(e2);d[2]=CNT<E>::real(e5);
00346 d[3]=CNT<E>::real(e9);d[4]=CNT<E>::real(e14);d[5]=CNT<E>::real(e20);
00347 l[0] =e1; l[1] =e3; l[2] =e6; l[3]=e10; l[4]=e15;
00348 l[5] =e4; l[6] =e7; l[7] =e11; l[8]=e16;
00349 l[9] =e8; l[10]=e12;l[11]=e17;
00350 l[12]=e13;l[13]=e18;
00351 l[14]=e19; }
00352
00353
00354
00355
00356
00357
00358
00359
00360
00361
00362
00363
00364
00365
00366
00367 template <class EE> explicit SymMat(const EE* p) {
00368 assert(p);
00369 for (int i=0; i<M; ++i) {
00370 const int rowStart = (i*(i+1))/2;
00371 updDiag()[i] = CNT<EE>::real(p[rowStart + i]);
00372 for (int j=0; j<i; ++j)
00373 updEltLower(i,j) = p[rowStart + j];
00374 }
00375 }
00376
00377
00378
00379 template <class EE> SymMat& operator=(const EE* p) {
00380 assert(p);
00381 for (int i=0; i<M; ++i) {
00382 const int rowStart = (i*(i+1))/2;
00383 updDiag()[i] = CNT<EE>::real(p[rowStart + i]);
00384 for (int j=0; j<i; ++j)
00385 updEltLower(i,j) = p[rowStart + j];
00386 }
00387 return *this;
00388 }
00389
00390
00391
00392
00393 template <class EE, int RSS> SymMat& operator=(const SymMat<M,EE,RSS>& mm) {
00394 updAsVec() = mm.getAsVec();
00395 return *this;
00396 }
00397
00398
00399
00400 template <class EE, int RSS> SymMat&
00401 operator+=(const SymMat<M,EE,RSS>& mm) {
00402 updAsVec() += mm.getAsVec();
00403 return *this;
00404 }
00405 template <class EE, int RSS> SymMat&
00406 operator+=(const SymMat<M,negator<EE>,RSS>& mm) {
00407 updAsVec() -= -mm.getAsVec();
00408 return *this;
00409 }
00410
00411 template <class EE, int RSS> SymMat&
00412 operator-=(const SymMat<M,EE,RSS>& mm) {
00413 updAsVec() -= mm.getAsVec();
00414 return *this;
00415 }
00416 template <class EE, int RSS> SymMat&
00417 operator-=(const SymMat<M,negator<EE>,RSS>& mm) {
00418 updAsVec() += -mm.getAsVec();
00419 return *this;
00420 }
00421
00422
00423
00424 template <class EE, int RSS> SymMat&
00425 operator*=(const SymMat<M,EE,RSS>& mm) {
00426 assert(false);
00427 return *this;
00428 }
00429
00430
00431
00432
00433
00434 template <class E2, int RS2>
00435 typename Result<SymMat<M,E2,RS2> >::Add
00436 conformingAdd(const SymMat<M,E2,RS2>& r) const {
00437 return typename Result<SymMat<M,E2,RS2> >::Add
00438 (getAsVec().conformingAdd(r.getAsVec()));
00439 }
00440
00441 template <class E2, int RS2>
00442 typename Result<SymMat<M,E2,RS2> >::Sub
00443 conformingSubtract(const SymMat<M,E2,RS2>& r) const {
00444 return typename Result<SymMat<M,E2,RS2> >::Sub
00445 (getAsVec().conformingSubtract(r.getAsVec()));
00446 }
00447
00448
00449
00450
00451 const E& operator()(int i,int j) const
00452 { return i==j ? getDiag()[i] : getEltLower(i,j); }
00453 E& operator()(int i,int j)
00454 { return i==j ? updDiag()[i] : updEltLower(i,j); }
00455
00456
00457 ScalarNormSq normSqr() const { return scalarNormSqr(); }
00458 typename CNT<ScalarNormSq>::TSqrt
00459 norm() const { return CNT<ScalarNormSq>::sqrt(scalarNormSqr()); }
00460
00461
00462
00463
00464
00465
00466
00467
00468
00469
00470
00471
00472 TNormalize normalize() const {
00473 if (CNT<E>::IsScalar) {
00474 return castAwayNegatorIfAny() / (SignInterpretation*norm());
00475 } else {
00476 TNormalize elementwiseNormalized;
00477
00478 elementwiseNormalized.updAsVec() = getAsVec().normalize();
00479 return elementwiseNormalized;
00480 }
00481 }
00482
00483 TInvert invert() const {assert(false); return TInvert();}
00484
00485 const SymMat& operator+() const { return *this; }
00486 const TNeg& operator-() const { return negate(); }
00487 TNeg& operator-() { return updNegate(); }
00488 const THerm& operator~() const { return transpose(); }
00489 THerm& operator~() { return updTranspose(); }
00490
00491 const TNeg& negate() const { return *reinterpret_cast<const TNeg*>(this); }
00492 TNeg& updNegate() { return *reinterpret_cast<TNeg*>(this); }
00493
00494 const THerm& transpose() const { return *reinterpret_cast<const THerm*>(this); }
00495 THerm& updTranspose() { return *reinterpret_cast<THerm*>(this); }
00496
00497 const TPosTrans& positionalTranspose() const
00498 { return *reinterpret_cast<const TPosTrans*>(this); }
00499 TPosTrans& updPositionalTranspose()
00500 { return *reinterpret_cast<TPosTrans*>(this); }
00501
00502 const TReal& real() const { return *reinterpret_cast<const TReal*>(this); }
00503 TReal& real() { return *reinterpret_cast< TReal*>(this); }
00504
00505
00506 const TImag& imag() const {
00507 const int offs = ImagOffset;
00508 const EImag* p = reinterpret_cast<const EImag*>(this);
00509 return *reinterpret_cast<const TImag*>(p+offs);
00510 }
00511 TImag& imag() {
00512 const int offs = ImagOffset;
00513 EImag* p = reinterpret_cast<EImag*>(this);
00514 return *reinterpret_cast<TImag*>(p+offs);
00515 }
00516
00517 const TWithoutNegator& castAwayNegatorIfAny() const {return *reinterpret_cast<const TWithoutNegator*>(this);}
00518 TWithoutNegator& updCastAwayNegatorIfAny() {return *reinterpret_cast<TWithoutNegator*>(this);}
00519
00520
00521
00522
00523
00524
00525
00526
00527 template <class EE> SymMat<M, typename CNT<E>::template Result<EE>::Mul>
00528 scalarMultiply(const EE& e) const {
00529 SymMat<M, typename CNT<E>::template Result<EE>::Mul> result;
00530 result.updAsVec() = getAsVec().scalarMultiply(e);
00531 return result;
00532 }
00533 template <class EE> SymMat<M, typename CNT<EE>::template Result<E>::Mul>
00534 scalarMultiplyFromLeft(const EE& e) const {
00535 SymMat<M, typename CNT<EE>::template Result<E>::Mul> result;
00536 result.updAsVec() = getAsVec().scalarMultiplyFromLeft(e);
00537 return result;
00538 }
00539
00540
00541
00542 template <class EE> SymMat<M, typename CNT<E>::template Result<EE>::Dvd>
00543 scalarDivide(const EE& e) const {
00544 SymMat<M, typename CNT<E>::template Result<EE>::Dvd> result;
00545 result.updAsVec() = getAsVec().scalarDivide(e);
00546 return result;
00547 }
00548 template <class EE> SymMat<M, typename CNT<EE>::template Result<E>::Dvd>
00549 scalarDivideFromLeft(const EE& e) const {
00550 SymMat<M, typename CNT<EE>::template Result<E>::Dvd> result;
00551 result.updAsVec() = getAsVec().scalarDivideFromLeft(e);
00552 return result;
00553 }
00554
00555
00556 template <class EE> SymMat<M, typename CNT<E>::template Result<EE>::Add>
00557 scalarAdd(const EE& e) const {
00558 SymMat<M, typename CNT<E>::template Result<EE>::Add> result(*this);
00559 result.updDiag() += e;
00560 return result;
00561 }
00562
00563
00564 template <class EE> SymMat<M, typename CNT<E>::template Result<EE>::Sub>
00565 scalarSubtract(const EE& e) const {
00566 SymMat<M, typename CNT<E>::template Result<EE>::Sub> result(*this);
00567 result.diag() -= e;
00568 return result;
00569 }
00570
00571
00572 template <class EE> SymMat<M, typename CNT<EE>::template Result<E>::Sub>
00573 scalarSubtractFromLeft(const EE& e) const {
00574 SymMat<M, typename CNT<EE>::template Result<E>::Sub> result(-(*this));
00575 result.diag() += e;
00576 return result;
00577 }
00578
00579
00580
00581
00582 template <class EE> SymMat& operator =(const EE& e) {return scalarEq(e);}
00583 template <class EE> SymMat& operator+=(const EE& e) {return scalarPlusEq(e);}
00584 template <class EE> SymMat& operator-=(const EE& e) {return scalarMinusEq(e);}
00585 template <class EE> SymMat& operator*=(const EE& e) {return scalarTimesEq(e);}
00586 template <class EE> SymMat& operator/=(const EE& e) {return scalarDivideEq(e);}
00587
00588
00589
00590 template <class EE> SymMat& scalarEq(const EE& ee)
00591 { updDiag() = ee; updLower() = E(0); return *this; }
00592 template <class EE> SymMat& scalarPlusEq(const EE& ee)
00593 { updDiag() += ee; return *this; }
00594 template <class EE> SymMat& scalarMinusEq(const EE& ee)
00595 { updDiag() -= ee; return *this; }
00596
00597
00598 template <class EE> SymMat& scalarMinusEqFromLeft(const EE& ee)
00599 { updLower() *= E(-1); updDiag().scalarMinusEqFromLeft(ee); return *this; }
00600
00601 template <class EE> SymMat& scalarTimesEq(const EE& ee)
00602 { updAsVec() *= ee; return *this; }
00603 template <class EE> SymMat& scalarTimesEqFromLeft(const EE& ee)
00604 { updAsVec().scalarTimesEqFromLeft(ee); return *this; }
00605 template <class EE> SymMat& scalarDivideEq(const EE& ee)
00606 { updAsVec() /= ee; return *this; }
00607 template <class EE> SymMat& scalarDivideEqFromLeft(const EE& ee)
00608 { updAsVec().scalarDivideEqFromLeft(ee); return *this; }
00609
00610 void setToNaN() { updAsVec().setToNaN(); }
00611
00612
00613 static const SymMat& getAs(const ELT* p) {return *reinterpret_cast<const SymMat*>(p);}
00614 static SymMat& updAs(ELT* p) {return *reinterpret_cast<SymMat*>(p);}
00615
00616
00617 static TPacked getNaN() {
00618 return TPacked(CNT<typename TPacked::TDiag>::getNaN(),
00619 CNT<typename TPacked::TLower>::getNaN());
00620 }
00621
00622 const TDiag& getDiag() const {return TDiag::getAs(d);}
00623 TDiag& updDiag() {return TDiag::updAs(d);}
00624
00625
00626 const TDiag& diag() const {return getDiag();}
00627 TDiag& diag() {return updDiag();}
00628
00629 const TLower& getLower() const {return TLower::getAs(&d[M*RS]);}
00630 TLower& updLower() {return TLower::updAs(&d[M*RS]);}
00631
00632 const TUpper& getUpper() const {return reinterpret_cast<const TUpper&>(getLower());}
00633 TUpper& updUpper() {return reinterpret_cast< TUpper&>(updLower());}
00634
00635 const TAsVec& getAsVec() const {return TAsVec::getAs(d);}
00636 TAsVec& updAsVec() {return TAsVec::updAs(d);}
00637
00638
00639 const E& getEltLower(int i, int j) const {return getLower()[lowerIx(i,j)];}
00640 E& updEltLower(int i, int j) {return updLower()[lowerIx(i,j)];}
00641
00642
00643 const EHerm& getEltUpper(int i, int j) const {return getUpper()[lowerIx(j,i)];}
00644 EHerm& updEltUpper(int i, int j) {return updUpper()[lowerIx(j,i)];}
00645
00646 TRow sum() const {
00647 TRow temp(~getDiag());
00648 for (int i = 1; i < M; ++i)
00649 for (int j = 0; j < i; ++j) {
00650 E value = getEltLower(i, j);;
00651 temp[i] += value;
00652 temp[j] += value;
00653 }
00654 return temp;
00655 }
00656
00657 private:
00658 E d[NActualElements];
00659
00660 SymMat(const TDiag& di, const TLower& low) {
00661 updDiag() = di; updLower() = low;
00662 }
00663
00664 explicit SymMat(const TAsVec& v) {
00665 TAsVec::updAs(d) = v;
00666 }
00667
00668
00669
00670 static int lowerIx(int i, int j) {
00671 assert(0 <= j && j < i && i < M);
00672 return (i-j-1) + j*(M-1) - (j*(j-1))/2;
00673 }
00674
00675 template <int MM, class EE, int RSS> friend class SymMat;
00676 };
00677
00679
00680
00682
00683
00684 template <int M, class E1, int S1, class E2, int S2> inline
00685 typename SymMat<M,E1,S1>::template Result< SymMat<M,E2,S2> >::Add
00686 operator+(const SymMat<M,E1,S1>& l, const SymMat<M,E2,S2>& r) {
00687 return SymMat<M,E1,S1>::template Result< SymMat<M,E2,S2> >
00688 ::AddOp::perform(l,r);
00689 }
00690
00691
00692 template <int M, class E1, int S1, class E2, int S2> inline
00693 typename SymMat<M,E1,S1>::template Result< SymMat<M,E2,S2> >::Sub
00694 operator-(const SymMat<M,E1,S1>& l, const SymMat<M,E2,S2>& r) {
00695 return SymMat<M,E1,S1>::template Result< SymMat<M,E2,S2> >
00696 ::SubOp::perform(l,r);
00697 }
00698
00699
00700
00701 template <int M, class E1, int S1, class E2, int S2> inline
00702 typename SymMat<M,E1,S1>::template Result< SymMat<M,E2,S2> >::Mul
00703 operator-(const SymMat<M,E1,S1>& l, const SymMat<M,E2,S2>& r) {
00704 return SymMat<M,E1,S1>::template Result< SymMat<M,E2,S2> >
00705 ::MulOp::perform(l,r);
00706 }
00707
00708
00709 template <int M, class E1, int S1, class E2, int S2> inline bool
00710 operator==(const SymMat<M,E1,S1>& l, const SymMat<M,E2,S2>& r) {
00711 return l.getAsVec() == r.getAsVec();
00712 }
00713
00714
00715 template <int M, class E1, int S1, class E2, int S2> inline bool
00716 operator!=(const SymMat<M,E1,S1>& l, const SymMat<M,E2,S2>& r) {return !(l==r);}
00717
00719
00721
00722
00723
00724
00725 template <int M, class E, int S> inline
00726 typename SymMat<M,E,S>::template Result<float>::Mul
00727 operator*(const SymMat<M,E,S>& l, const float& r)
00728 { return SymMat<M,E,S>::template Result<float>::MulOp::perform(l,r); }
00729 template <int M, class E, int S> inline
00730 typename SymMat<M,E,S>::template Result<float>::Mul
00731 operator*(const float& l, const SymMat<M,E,S>& r) {return r*l;}
00732
00733 template <int M, class E, int S> inline
00734 typename SymMat<M,E,S>::template Result<double>::Mul
00735 operator*(const SymMat<M,E,S>& l, const double& r)
00736 { return SymMat<M,E,S>::template Result<double>::MulOp::perform(l,r); }
00737 template <int M, class E, int S> inline
00738 typename SymMat<M,E,S>::template Result<double>::Mul
00739 operator*(const double& l, const SymMat<M,E,S>& r) {return r*l;}
00740
00741 template <int M, class E, int S> inline
00742 typename SymMat<M,E,S>::template Result<long double>::Mul
00743 operator*(const SymMat<M,E,S>& l, const long double& r)
00744 { return SymMat<M,E,S>::template Result<long double>::MulOp::perform(l,r); }
00745 template <int M, class E, int S> inline
00746 typename SymMat<M,E,S>::template Result<long double>::Mul
00747 operator*(const long double& l, const SymMat<M,E,S>& r) {return r*l;}
00748
00749
00750 template <int M, class E, int S> inline
00751 typename SymMat<M,E,S>::template Result<typename CNT<E>::Precision>::Mul
00752 operator*(const SymMat<M,E,S>& l, int r) {return l * (typename CNT<E>::Precision)r;}
00753 template <int M, class E, int S> inline
00754 typename SymMat<M,E,S>::template Result<typename CNT<E>::Precision>::Mul
00755 operator*(int l, const SymMat<M,E,S>& r) {return r * (typename CNT<E>::Precision)l;}
00756
00757
00758
00759
00760 template <int M, class E, int S, class R> inline
00761 typename SymMat<M,E,S>::template Result<std::complex<R> >::Mul
00762 operator*(const SymMat<M,E,S>& l, const std::complex<R>& r)
00763 { return SymMat<M,E,S>::template Result<std::complex<R> >::MulOp::perform(l,r); }
00764 template <int M, class E, int S, class R> inline
00765 typename SymMat<M,E,S>::template Result<std::complex<R> >::Mul
00766 operator*(const std::complex<R>& l, const SymMat<M,E,S>& r) {return r*l;}
00767
00768
00769 template <int M, class E, int S, class R> inline
00770 typename SymMat<M,E,S>::template Result<std::complex<R> >::Mul
00771 operator*(const SymMat<M,E,S>& l, const conjugate<R>& r) {return l*(std::complex<R>)r;}
00772 template <int M, class E, int S, class R> inline
00773 typename SymMat<M,E,S>::template Result<std::complex<R> >::Mul
00774 operator*(const conjugate<R>& l, const SymMat<M,E,S>& r) {return r*(std::complex<R>)l;}
00775
00776
00777 template <int M, class E, int S, class R> inline
00778 typename SymMat<M,E,S>::template Result<typename negator<R>::StdNumber>::Mul
00779 operator*(const SymMat<M,E,S>& l, const negator<R>& r) {return l * (typename negator<R>::StdNumber)(R)r;}
00780 template <int M, class E, int S, class R> inline
00781 typename SymMat<M,E,S>::template Result<typename negator<R>::StdNumber>::Mul
00782 operator*(const negator<R>& l, const SymMat<M,E,S>& r) {return r * (typename negator<R>::StdNumber)(R)l;}
00783
00784
00785
00786
00787
00788
00789
00790 template <int M, class E, int S> inline
00791 typename SymMat<M,E,S>::template Result<float>::Dvd
00792 operator/(const SymMat<M,E,S>& l, const float& r)
00793 { return SymMat<M,E,S>::template Result<float>::DvdOp::perform(l,r); }
00794 template <int M, class E, int S> inline
00795 typename CNT<float>::template Result<SymMat<M,E,S> >::Dvd
00796 operator/(const float& l, const SymMat<M,E,S>& r)
00797 { return CNT<float>::template Result<SymMat<M,E,S> >::DvdOp::perform(l,r); }
00798
00799 template <int M, class E, int S> inline
00800 typename SymMat<M,E,S>::template Result<double>::Dvd
00801 operator/(const SymMat<M,E,S>& l, const double& r)
00802 { return SymMat<M,E,S>::template Result<double>::DvdOp::perform(l,r); }
00803 template <int M, class E, int S> inline
00804 typename CNT<double>::template Result<SymMat<M,E,S> >::Dvd
00805 operator/(const double& l, const SymMat<M,E,S>& r)
00806 { return CNT<double>::template Result<SymMat<M,E,S> >::DvdOp::perform(l,r); }
00807
00808 template <int M, class E, int S> inline
00809 typename SymMat<M,E,S>::template Result<long double>::Dvd
00810 operator/(const SymMat<M,E,S>& l, const long double& r)
00811 { return SymMat<M,E,S>::template Result<long double>::DvdOp::perform(l,r); }
00812 template <int M, class E, int S> inline
00813 typename CNT<long double>::template Result<SymMat<M,E,S> >::Dvd
00814 operator/(const long double& l, const SymMat<M,E,S>& r)
00815 { return CNT<long double>::template Result<SymMat<M,E,S> >::DvdOp::perform(l,r); }
00816
00817
00818 template <int M, class E, int S> inline
00819 typename SymMat<M,E,S>::template Result<typename CNT<E>::Precision>::Dvd
00820 operator/(const SymMat<M,E,S>& l, int r) {return l / (typename CNT<E>::Precision)r;}
00821 template <int M, class E, int S> inline
00822 typename CNT<typename CNT<E>::Precision>::template Result<SymMat<M,E,S> >::Dvd
00823 operator/(int l, const SymMat<M,E,S>& r) {return (typename CNT<E>::Precision)l / r;}
00824
00825
00826
00827
00828
00829 template <int M, class E, int S, class R> inline
00830 typename SymMat<M,E,S>::template Result<std::complex<R> >::Dvd
00831 operator/(const SymMat<M,E,S>& l, const std::complex<R>& r)
00832 { return SymMat<M,E,S>::template Result<std::complex<R> >::DvdOp::perform(l,r); }
00833 template <int M, class E, int S, class R> inline
00834 typename CNT<std::complex<R> >::template Result<SymMat<M,E,S> >::Dvd
00835 operator/(const std::complex<R>& l, const SymMat<M,E,S>& r)
00836 { return CNT<std::complex<R> >::template Result<SymMat<M,E,S> >::DvdOp::perform(l,r); }
00837
00838
00839 template <int M, class E, int S, class R> inline
00840 typename SymMat<M,E,S>::template Result<std::complex<R> >::Dvd
00841 operator/(const SymMat<M,E,S>& l, const conjugate<R>& r) {return l/(std::complex<R>)r;}
00842 template <int M, class E, int S, class R> inline
00843 typename CNT<std::complex<R> >::template Result<SymMat<M,E,S> >::Dvd
00844 operator/(const conjugate<R>& l, const SymMat<M,E,S>& r) {return (std::complex<R>)l/r;}
00845
00846
00847 template <int M, class E, int S, class R> inline
00848 typename SymMat<M,E,S>::template Result<typename negator<R>::StdNumber>::Dvd
00849 operator/(const SymMat<M,E,S>& l, const negator<R>& r) {return l/(typename negator<R>::StdNumber)(R)r;}
00850 template <int M, class E, int S, class R> inline
00851 typename CNT<R>::template Result<SymMat<M,E,S> >::Dvd
00852 operator/(const negator<R>& l, const SymMat<M,E,S>& r) {return (typename negator<R>::StdNumber)(R)l/r;}
00853
00854
00855
00856
00857
00858
00859
00860
00861
00862 template <int M, class E, int S> inline
00863 typename SymMat<M,E,S>::template Result<float>::Add
00864 operator+(const SymMat<M,E,S>& l, const float& r)
00865 { return SymMat<M,E,S>::template Result<float>::AddOp::perform(l,r); }
00866 template <int M, class E, int S> inline
00867 typename SymMat<M,E,S>::template Result<float>::Add
00868 operator+(const float& l, const SymMat<M,E,S>& r) {return r+l;}
00869
00870 template <int M, class E, int S> inline
00871 typename SymMat<M,E,S>::template Result<double>::Add
00872 operator+(const SymMat<M,E,S>& l, const double& r)
00873 { return SymMat<M,E,S>::template Result<double>::AddOp::perform(l,r); }
00874 template <int M, class E, int S> inline
00875 typename SymMat<M,E,S>::template Result<double>::Add
00876 operator+(const double& l, const SymMat<M,E,S>& r) {return r+l;}
00877
00878 template <int M, class E, int S> inline
00879 typename SymMat<M,E,S>::template Result<long double>::Add
00880 operator+(const SymMat<M,E,S>& l, const long double& r)
00881 { return SymMat<M,E,S>::template Result<long double>::AddOp::perform(l,r); }
00882 template <int M, class E, int S> inline
00883 typename SymMat<M,E,S>::template Result<long double>::Add
00884 operator+(const long double& l, const SymMat<M,E,S>& r) {return r+l;}
00885
00886
00887 template <int M, class E, int S> inline
00888 typename SymMat<M,E,S>::template Result<typename CNT<E>::Precision>::Add
00889 operator+(const SymMat<M,E,S>& l, int r) {return l + (typename CNT<E>::Precision)r;}
00890 template <int M, class E, int S> inline
00891 typename SymMat<M,E,S>::template Result<typename CNT<E>::Precision>::Add
00892 operator+(int l, const SymMat<M,E,S>& r) {return r + (typename CNT<E>::Precision)l;}
00893
00894
00895
00896
00897 template <int M, class E, int S, class R> inline
00898 typename SymMat<M,E,S>::template Result<std::complex<R> >::Add
00899 operator+(const SymMat<M,E,S>& l, const std::complex<R>& r)
00900 { return SymMat<M,E,S>::template Result<std::complex<R> >::AddOp::perform(l,r); }
00901 template <int M, class E, int S, class R> inline
00902 typename SymMat<M,E,S>::template Result<std::complex<R> >::Add
00903 operator+(const std::complex<R>& l, const SymMat<M,E,S>& r) {return r+l;}
00904
00905
00906 template <int M, class E, int S, class R> inline
00907 typename SymMat<M,E,S>::template Result<std::complex<R> >::Add
00908 operator+(const SymMat<M,E,S>& l, const conjugate<R>& r) {return l+(std::complex<R>)r;}
00909 template <int M, class E, int S, class R> inline
00910 typename SymMat<M,E,S>::template Result<std::complex<R> >::Add
00911 operator+(const conjugate<R>& l, const SymMat<M,E,S>& r) {return r+(std::complex<R>)l;}
00912
00913
00914 template <int M, class E, int S, class R> inline
00915 typename SymMat<M,E,S>::template Result<typename negator<R>::StdNumber>::Add
00916 operator+(const SymMat<M,E,S>& l, const negator<R>& r) {return l + (typename negator<R>::StdNumber)(R)r;}
00917 template <int M, class E, int S, class R> inline
00918 typename SymMat<M,E,S>::template Result<typename negator<R>::StdNumber>::Add
00919 operator+(const negator<R>& l, const SymMat<M,E,S>& r) {return r + (typename negator<R>::StdNumber)(R)l;}
00920
00921
00922
00923
00924 template <int M, class E, int S> inline
00925 typename SymMat<M,E,S>::template Result<float>::Sub
00926 operator-(const SymMat<M,E,S>& l, const float& r)
00927 { return SymMat<M,E,S>::template Result<float>::SubOp::perform(l,r); }
00928 template <int M, class E, int S> inline
00929 typename CNT<float>::template Result<SymMat<M,E,S> >::Sub
00930 operator-(const float& l, const SymMat<M,E,S>& r)
00931 { return CNT<float>::template Result<SymMat<M,E,S> >::SubOp::perform(l,r); }
00932
00933 template <int M, class E, int S> inline
00934 typename SymMat<M,E,S>::template Result<double>::Sub
00935 operator-(const SymMat<M,E,S>& l, const double& r)
00936 { return SymMat<M,E,S>::template Result<double>::SubOp::perform(l,r); }
00937 template <int M, class E, int S> inline
00938 typename CNT<double>::template Result<SymMat<M,E,S> >::Sub
00939 operator-(const double& l, const SymMat<M,E,S>& r)
00940 { return CNT<double>::template Result<SymMat<M,E,S> >::SubOp::perform(l,r); }
00941
00942 template <int M, class E, int S> inline
00943 typename SymMat<M,E,S>::template Result<long double>::Sub
00944 operator-(const SymMat<M,E,S>& l, const long double& r)
00945 { return SymMat<M,E,S>::template Result<long double>::SubOp::perform(l,r); }
00946 template <int M, class E, int S> inline
00947 typename CNT<long double>::template Result<SymMat<M,E,S> >::Sub
00948 operator-(const long double& l, const SymMat<M,E,S>& r)
00949 { return CNT<long double>::template Result<SymMat<M,E,S> >::SubOp::perform(l,r); }
00950
00951
00952 template <int M, class E, int S> inline
00953 typename SymMat<M,E,S>::template Result<typename CNT<E>::Precision>::Sub
00954 operator-(const SymMat<M,E,S>& l, int r) {return l - (typename CNT<E>::Precision)r;}
00955 template <int M, class E, int S> inline
00956 typename CNT<typename CNT<E>::Precision>::template Result<SymMat<M,E,S> >::Sub
00957 operator-(int l, const SymMat<M,E,S>& r) {return (typename CNT<E>::Precision)l - r;}
00958
00959
00960
00961
00962
00963 template <int M, class E, int S, class R> inline
00964 typename SymMat<M,E,S>::template Result<std::complex<R> >::Sub
00965 operator-(const SymMat<M,E,S>& l, const std::complex<R>& r)
00966 { return SymMat<M,E,S>::template Result<std::complex<R> >::SubOp::perform(l,r); }
00967 template <int M, class E, int S, class R> inline
00968 typename CNT<std::complex<R> >::template Result<SymMat<M,E,S> >::Sub
00969 operator-(const std::complex<R>& l, const SymMat<M,E,S>& r)
00970 { return CNT<std::complex<R> >::template Result<SymMat<M,E,S> >::SubOp::perform(l,r); }
00971
00972
00973 template <int M, class E, int S, class R> inline
00974 typename SymMat<M,E,S>::template Result<std::complex<R> >::Sub
00975 operator-(const SymMat<M,E,S>& l, const conjugate<R>& r) {return l-(std::complex<R>)r;}
00976 template <int M, class E, int S, class R> inline
00977 typename CNT<std::complex<R> >::template Result<SymMat<M,E,S> >::Sub
00978 operator-(const conjugate<R>& l, const SymMat<M,E,S>& r) {return (std::complex<R>)l-r;}
00979
00980
00981 template <int M, class E, int S, class R> inline
00982 typename SymMat<M,E,S>::template Result<typename negator<R>::StdNumber>::Sub
00983 operator-(const SymMat<M,E,S>& l, const negator<R>& r) {return l-(typename negator<R>::StdNumber)(R)r;}
00984 template <int M, class E, int S, class R> inline
00985 typename CNT<R>::template Result<SymMat<M,E,S> >::Sub
00986 operator-(const negator<R>& l, const SymMat<M,E,S>& r) {return (typename negator<R>::StdNumber)(R)l-r;}
00987
00988
00989
00990 template <int M, class E, int RS, class CHAR, class TRAITS> inline
00991 std::basic_ostream<CHAR,TRAITS>&
00992 operator<<(std::basic_ostream<CHAR,TRAITS>& o, const SymMat<M,E,RS>& m) {
00993 for (int i=0;i<M;++i) {
00994 o << std::endl << "[";
00995 for (int j=0; j<=i; ++j)
00996 o << (j>0?" ":"") << m(i,j);
00997 for (int j=i+1; j<M; ++j)
00998 o << " *";
00999 o << "]";
01000 }
01001 if (M) o << std::endl;
01002 return o;
01003 }
01004
01005 template <int M, class E, int RS, class CHAR, class TRAITS> inline
01006 std::basic_istream<CHAR,TRAITS>&
01007 operator>>(std::basic_istream<CHAR,TRAITS>& is, SymMat<M,E,RS>& m) {
01008
01009 assert(false);
01010 return is;
01011 }
01012
01013 }
01014
01015
01016 #endif //SimTK_SIMMATRIX_SMALLMATRIX_SYMMAT_H_