11#include "CLHEP/Matrix/defs.h"
12#include "CLHEP/Random/Random.h"
13#include "CLHEP/Matrix/DiagMatrix.h"
14#include "CLHEP/Matrix/Matrix.h"
15#include "CLHEP/Matrix/SymMatrix.h"
16#include "CLHEP/Matrix/Vector.h"
18#ifdef HEP_DEBUG_INLINE
19#include "CLHEP/Matrix/DiagMatrix.icc"
26#define SIMPLE_UOP(OPER) \
27 HepMatrix::mIter a=m.begin(); \
28 HepMatrix::mIter e=m.begin()+num_size(); \
29 for(;a<e; a++) (*a) OPER t;
31#define SIMPLE_BOP(OPER) \
32 HepMatrix::mIter a=m.begin(); \
33 HepMatrix::mcIter b=hm2.m.begin(); \
34 HepMatrix::mIter e=m.begin()+num_size(); \
35 for(;a<e; a++, b++) (*a) OPER (*b);
37#define SIMPLE_TOP(OPER) \
38 HepMatrix::mcIter a=hm1.m.begin(); \
39 HepMatrix::mcIter b=hm2.m.begin(); \
40 HepMatrix::mIter t=mret.m.begin(); \
41 HepMatrix::mcIter e=hm1.m.begin()+hm1.nrow; \
42 for( ;a<e; a++, b++, t++) (*t) = (*a) OPER (*b);
44#define CHK_DIM_2(r1,r2,c1,c2,fun) \
45 if (r1!=r2 || c1!=c2) { \
46 HepGenMatrix::error("Range error in DiagMatrix function " #fun "(1)."); \
49#define CHK_DIM_1(c1,r2,fun) \
51 HepGenMatrix::error("Range error in DiagMatrix function " #fun "(2)."); \
56#if defined(__sun) || !defined(__GNUG__)
60double HepDiagMatrix::zero = 0;
62const double HepDiagMatrix::zero = 0;
85 for( ; a<b; a++) *a = 1.0;
89 error(
"DiagMatrix: initialization must be either 0 or 1.");
98 for(;a<b;a++) *a = r();
119#ifdef HEP_GNU_OPTIMIZED_RETURN
120return mret(max_row-min_row+1);
126 if(max_row > num_row())
127 error(
"HepDiagMatrix::sub: Index out of range");
131 for(;a<e;) *(a++) = *(b++);
139 error(
"HepDiagMatrix::sub: Index out of range");
143 for(;a<e;) *(a++) = *(b++);
150 error(
"HepDiagMatrix::sub: Index out of range");
154 for(;a<e;) *(b++) = *(a++);
163#ifdef HEP_GNU_OPTIMIZED_RETURN
177#ifdef HEP_GNU_OPTIMIZED_RETURN
187 for(;a<e; a++, b++) (*b) = -(*a);
194#ifdef HEP_GNU_OPTIMIZED_RETURN
208#ifdef HEP_GNU_OPTIMIZED_RETURN
222#ifdef HEP_GNU_OPTIMIZED_RETURN
223 return mret(hm1.nrow);
235#ifdef HEP_GNU_OPTIMIZED_RETURN
248#ifdef HEP_GNU_OPTIMIZED_RETURN
265#ifdef HEP_GNU_OPTIMIZED_RETURN
278#ifdef HEP_GNU_OPTIMIZED_RETURN
292#ifdef HEP_GNU_OPTIMIZED_RETURN
293 return mret(hm1.nrow);
304#ifdef HEP_GNU_OPTIMIZED_RETURN
317#ifdef HEP_GNU_OPTIMIZED_RETURN
337#ifdef HEP_GNU_OPTIMIZED_RETURN
349#ifdef HEP_GNU_OPTIMIZED_RETURN
361#ifdef HEP_GNU_OPTIMIZED_RETURN
373#ifdef HEP_GNU_OPTIMIZED_RETURN
383 for(
int irow=1;irow<=hm1.
num_row();irow++) {
385 for(
int icol=1;icol<=hm1.
num_col();icol++) {
386 *(mir++) = *(mit1++) * (*(mcc++));
393#ifdef HEP_GNU_OPTIMIZED_RETURN
404 for(
int irow=1;irow<=hm2.
num_row();irow++) {
405 for(
int icol=1;icol<=hm2.
num_col();icol++) {
406 *(mir++) = *(mit1++) * (*mrr);
414#ifdef HEP_GNU_OPTIMIZED_RETURN
426 for(;a<e;) *(a++) = *(b++) * (*(c++));
431#ifdef HEP_GNU_OPTIMIZED_RETURN
441 for(
int icol=1;icol<=hm1.
num_col();icol++) {
442 *(mir++) = *(mi1++) * *(mi2++);
455 mIter mrr = m.begin();
457 for(
int r=1;r<=n;r++) {
459 if(r<n) mrr += (n+1);
469 for(
int i=1;i<=
num_row();i++) {
487 mIter mrr = m.begin();
489 for(
int r=1;r<=n;r++) {
491 if(r<n) mrr += (n+1);
501 for(
int i=1;i<=
num_row();i++) {
529 if(hm1.nrow*hm1.nrow != size_)
531 size_ = hm1.nrow * hm1.nrow;
538 mIter mrr = m.begin();
540 for(
int r=1;r<=n;r++) {
542 if(r<n) mrr += (n+1);
565 if(os.flags() & std::ios::fixed)
566 width = os.precision()+3;
568 width = os.precision()+7;
569 for(
int irow = 1; irow<= q.
num_row(); irow++)
571 for(
int icol = 1; icol <= q.
num_col(); icol++)
574 os << q(irow,icol) <<
" ";
582apply(
double (*
f)(
double,
int,
int))
const
583#ifdef HEP_GNU_OPTIMIZED_RETURN
584return mret(num_row());
592 for(
int ir=1;ir<=num_row();ir++) {
593 *(b++) = (*
f)(*(a++), ir, ir);
607 for(
int r=1;r<=nrow;r++) {
609 if(r<nrow) a += (nrow+1);
622 for(
int r=1;r<=nrow;r++) {
624 if(r<nrow) a += (r+1);
629#ifdef HEP_GNU_OPTIMIZED_RETURN
641 for(
int r=1;r<=mret.num_row();r++) {
644 for(
int c=1;c<=r;c++) {
648 for(
int i=0;i<hm1.
num_col();i++)
649 tmp+=*(mr++) * *(mc++) * *(mi++);
662 mret = *(mv)* *(mv)* *(mi++);
664 for(
int i=2;i<=hm1.
num_row();i++) {
665 mret+=*(mv)* *(mv)* *(mi++);
672#ifdef HEP_GNU_OPTIMIZED_RETURN
683 for(
int r=1;r<=mret.num_row();r++)
684 for(
int c=1;c<=r;c++)
687 double tmp = hm1(1,r)*hm1(1,c)* *(mi++);
688 for(
int i=2;i<=hm1.
num_row();i++)
689 tmp+=hm1(i,r)*hm1(i,c)* *(mi++);
690 mret.
fast(r,c) = tmp;
701 if(*(hmm++)==0)
return;
#define CHK_DIM_2(r1, r2, c1, c2, fun)
#define CHK_DIM_1(c1, r2, fun)
double & fast(int row, int col)
double determinant() const
HepDiagMatrix & operator*=(double t)
HepDiagMatrix sub(int min_row, int max_row) const
HepDiagMatrix & operator/=(double t)
HepDiagMatrix & operator+=(const HepDiagMatrix &hm2)
void assign(const HepMatrix &hm2)
HepDiagMatrix operator-() const
HepSymMatrix similarityT(const HepMatrix &hm1) const
HepSymMatrix similarity(const HepMatrix &hm1) const
HepDiagMatrix & operator-=(const HepDiagMatrix &hm2)
HepDiagMatrix & operator=(const HepDiagMatrix &hm2)
HepDiagMatrix apply(double(*f)(double, int, int)) const
std::vector< double, Alloc< double, 25 > >::const_iterator mcIter
std::vector< double, Alloc< double, 25 > >::iterator mIter
static void error(const char *s)
HepMatrix & operator=(const HepMatrix &)
virtual int num_col() const
virtual int num_row() const
HepMatrix & operator+=(const HepMatrix &)
HepMatrix & operator-=(const HepMatrix &)
HepSymMatrix & operator+=(const HepSymMatrix &hm2)
HepSymMatrix & operator-=(const HepSymMatrix &hm2)
virtual int num_row() const
std::ostream & operator<<(std::ostream &s, const HepDiagMatrix &q)
HepMatrix operator+(const HepMatrix &hm1, const HepDiagMatrix &d2)
HepMatrix operator-(const HepMatrix &hm1, const HepDiagMatrix &d2)
HepDiagMatrix operator/(const HepDiagMatrix &hm1, double t)
HepMatrix operator*(const HepMatrix &hm1, const HepDiagMatrix &hm2)
HepDiagMatrix dsum(const HepDiagMatrix &s1, const HepDiagMatrix &s2)