$treeview $search $mathjax
Eigen-unsupported
3.2.5
$projectbrief
|
$projectbrief
|
$searchbox |
00001 // This file is part of Eigen, a lightweight C++ template library 00002 // for linear algebra. 00003 // 00004 // Copyright (C) 2008-2009 Guillaume Saupin <guillaume.saupin@cea.fr> 00005 // 00006 // This Source Code Form is subject to the terms of the Mozilla 00007 // Public License v. 2.0. If a copy of the MPL was not distributed 00008 // with this file, You can obtain one at http://mozilla.org/MPL/2.0/. 00009 00010 #ifndef EIGEN_SKYLINEMATRIX_H 00011 #define EIGEN_SKYLINEMATRIX_H 00012 00013 #include "SkylineStorage.h" 00014 #include "SkylineMatrixBase.h" 00015 00016 namespace Eigen { 00017 00033 namespace internal { 00034 template<typename _Scalar, int _Options> 00035 struct traits<SkylineMatrix<_Scalar, _Options> > { 00036 typedef _Scalar Scalar; 00037 typedef Sparse StorageKind; 00038 00039 enum { 00040 RowsAtCompileTime = Dynamic, 00041 ColsAtCompileTime = Dynamic, 00042 MaxRowsAtCompileTime = Dynamic, 00043 MaxColsAtCompileTime = Dynamic, 00044 Flags = SkylineBit | _Options, 00045 CoeffReadCost = NumTraits<Scalar>::ReadCost, 00046 }; 00047 }; 00048 } 00049 00050 template<typename _Scalar, int _Options> 00051 class SkylineMatrix 00052 : public SkylineMatrixBase<SkylineMatrix<_Scalar, _Options> > { 00053 public: 00054 EIGEN_SKYLINE_GENERIC_PUBLIC_INTERFACE(SkylineMatrix) 00055 EIGEN_SKYLINE_INHERIT_ASSIGNMENT_OPERATOR(SkylineMatrix, +=) 00056 EIGEN_SKYLINE_INHERIT_ASSIGNMENT_OPERATOR(SkylineMatrix, -=) 00057 00058 using Base::IsRowMajor; 00059 00060 protected: 00061 00062 typedef SkylineMatrix<Scalar, (Flags&~RowMajorBit) | (IsRowMajor ? RowMajorBit : 0) > TransposedSkylineMatrix; 00063 00064 Index m_outerSize; 00065 Index m_innerSize; 00066 00067 public: 00068 Index* m_colStartIndex; 00069 Index* m_rowStartIndex; 00070 SkylineStorage<Scalar> m_data; 00071 00072 public: 00073 00074 inline Index rows() const { 00075 return IsRowMajor ? m_outerSize : m_innerSize; 00076 } 00077 00078 inline Index cols() const { 00079 return IsRowMajor ? m_innerSize : m_outerSize; 00080 } 00081 00082 inline Index innerSize() const { 00083 return m_innerSize; 00084 } 00085 00086 inline Index outerSize() const { 00087 return m_outerSize; 00088 } 00089 00090 inline Index upperNonZeros() const { 00091 return m_data.upperSize(); 00092 } 00093 00094 inline Index lowerNonZeros() const { 00095 return m_data.lowerSize(); 00096 } 00097 00098 inline Index upperNonZeros(Index j) const { 00099 return m_colStartIndex[j + 1] - m_colStartIndex[j]; 00100 } 00101 00102 inline Index lowerNonZeros(Index j) const { 00103 return m_rowStartIndex[j + 1] - m_rowStartIndex[j]; 00104 } 00105 00106 inline const Scalar* _diagPtr() const { 00107 return &m_data.diag(0); 00108 } 00109 00110 inline Scalar* _diagPtr() { 00111 return &m_data.diag(0); 00112 } 00113 00114 inline const Scalar* _upperPtr() const { 00115 return &m_data.upper(0); 00116 } 00117 00118 inline Scalar* _upperPtr() { 00119 return &m_data.upper(0); 00120 } 00121 00122 inline const Scalar* _lowerPtr() const { 00123 return &m_data.lower(0); 00124 } 00125 00126 inline Scalar* _lowerPtr() { 00127 return &m_data.lower(0); 00128 } 00129 00130 inline const Index* _upperProfilePtr() const { 00131 return &m_data.upperProfile(0); 00132 } 00133 00134 inline Index* _upperProfilePtr() { 00135 return &m_data.upperProfile(0); 00136 } 00137 00138 inline const Index* _lowerProfilePtr() const { 00139 return &m_data.lowerProfile(0); 00140 } 00141 00142 inline Index* _lowerProfilePtr() { 00143 return &m_data.lowerProfile(0); 00144 } 00145 00146 inline Scalar coeff(Index row, Index col) const { 00147 const Index outer = IsRowMajor ? row : col; 00148 const Index inner = IsRowMajor ? col : row; 00149 00150 eigen_assert(outer < outerSize()); 00151 eigen_assert(inner < innerSize()); 00152 00153 if (outer == inner) 00154 return this->m_data.diag(outer); 00155 00156 if (IsRowMajor) { 00157 if (inner > outer) //upper matrix 00158 { 00159 const Index minOuterIndex = inner - m_data.upperProfile(inner); 00160 if (outer >= minOuterIndex) 00161 return this->m_data.upper(m_colStartIndex[inner] + outer - (inner - m_data.upperProfile(inner))); 00162 else 00163 return Scalar(0); 00164 } 00165 if (inner < outer) //lower matrix 00166 { 00167 const Index minInnerIndex = outer - m_data.lowerProfile(outer); 00168 if (inner >= minInnerIndex) 00169 return this->m_data.lower(m_rowStartIndex[outer] + inner - (outer - m_data.lowerProfile(outer))); 00170 else 00171 return Scalar(0); 00172 } 00173 return m_data.upper(m_colStartIndex[inner] + outer - inner); 00174 } else { 00175 if (outer > inner) //upper matrix 00176 { 00177 const Index maxOuterIndex = inner + m_data.upperProfile(inner); 00178 if (outer <= maxOuterIndex) 00179 return this->m_data.upper(m_colStartIndex[inner] + (outer - inner)); 00180 else 00181 return Scalar(0); 00182 } 00183 if (outer < inner) //lower matrix 00184 { 00185 const Index maxInnerIndex = outer + m_data.lowerProfile(outer); 00186 00187 if (inner <= maxInnerIndex) 00188 return this->m_data.lower(m_rowStartIndex[outer] + (inner - outer)); 00189 else 00190 return Scalar(0); 00191 } 00192 } 00193 } 00194 00195 inline Scalar& coeffRef(Index row, Index col) { 00196 const Index outer = IsRowMajor ? row : col; 00197 const Index inner = IsRowMajor ? col : row; 00198 00199 eigen_assert(outer < outerSize()); 00200 eigen_assert(inner < innerSize()); 00201 00202 if (outer == inner) 00203 return this->m_data.diag(outer); 00204 00205 if (IsRowMajor) { 00206 if (col > row) //upper matrix 00207 { 00208 const Index minOuterIndex = inner - m_data.upperProfile(inner); 00209 eigen_assert(outer >= minOuterIndex && "you try to acces a coeff that do not exist in the storage"); 00210 return this->m_data.upper(m_colStartIndex[inner] + outer - (inner - m_data.upperProfile(inner))); 00211 } 00212 if (col < row) //lower matrix 00213 { 00214 const Index minInnerIndex = outer - m_data.lowerProfile(outer); 00215 eigen_assert(inner >= minInnerIndex && "you try to acces a coeff that do not exist in the storage"); 00216 return this->m_data.lower(m_rowStartIndex[outer] + inner - (outer - m_data.lowerProfile(outer))); 00217 } 00218 } else { 00219 if (outer > inner) //upper matrix 00220 { 00221 const Index maxOuterIndex = inner + m_data.upperProfile(inner); 00222 eigen_assert(outer <= maxOuterIndex && "you try to acces a coeff that do not exist in the storage"); 00223 return this->m_data.upper(m_colStartIndex[inner] + (outer - inner)); 00224 } 00225 if (outer < inner) //lower matrix 00226 { 00227 const Index maxInnerIndex = outer + m_data.lowerProfile(outer); 00228 eigen_assert(inner <= maxInnerIndex && "you try to acces a coeff that do not exist in the storage"); 00229 return this->m_data.lower(m_rowStartIndex[outer] + (inner - outer)); 00230 } 00231 } 00232 } 00233 00234 inline Scalar coeffDiag(Index idx) const { 00235 eigen_assert(idx < outerSize()); 00236 eigen_assert(idx < innerSize()); 00237 return this->m_data.diag(idx); 00238 } 00239 00240 inline Scalar coeffLower(Index row, Index col) const { 00241 const Index outer = IsRowMajor ? row : col; 00242 const Index inner = IsRowMajor ? col : row; 00243 00244 eigen_assert(outer < outerSize()); 00245 eigen_assert(inner < innerSize()); 00246 eigen_assert(inner != outer); 00247 00248 if (IsRowMajor) { 00249 const Index minInnerIndex = outer - m_data.lowerProfile(outer); 00250 if (inner >= minInnerIndex) 00251 return this->m_data.lower(m_rowStartIndex[outer] + inner - (outer - m_data.lowerProfile(outer))); 00252 else 00253 return Scalar(0); 00254 00255 } else { 00256 const Index maxInnerIndex = outer + m_data.lowerProfile(outer); 00257 if (inner <= maxInnerIndex) 00258 return this->m_data.lower(m_rowStartIndex[outer] + (inner - outer)); 00259 else 00260 return Scalar(0); 00261 } 00262 } 00263 00264 inline Scalar coeffUpper(Index row, Index col) const { 00265 const Index outer = IsRowMajor ? row : col; 00266 const Index inner = IsRowMajor ? col : row; 00267 00268 eigen_assert(outer < outerSize()); 00269 eigen_assert(inner < innerSize()); 00270 eigen_assert(inner != outer); 00271 00272 if (IsRowMajor) { 00273 const Index minOuterIndex = inner - m_data.upperProfile(inner); 00274 if (outer >= minOuterIndex) 00275 return this->m_data.upper(m_colStartIndex[inner] + outer - (inner - m_data.upperProfile(inner))); 00276 else 00277 return Scalar(0); 00278 } else { 00279 const Index maxOuterIndex = inner + m_data.upperProfile(inner); 00280 if (outer <= maxOuterIndex) 00281 return this->m_data.upper(m_colStartIndex[inner] + (outer - inner)); 00282 else 00283 return Scalar(0); 00284 } 00285 } 00286 00287 inline Scalar& coeffRefDiag(Index idx) { 00288 eigen_assert(idx < outerSize()); 00289 eigen_assert(idx < innerSize()); 00290 return this->m_data.diag(idx); 00291 } 00292 00293 inline Scalar& coeffRefLower(Index row, Index col) { 00294 const Index outer = IsRowMajor ? row : col; 00295 const Index inner = IsRowMajor ? col : row; 00296 00297 eigen_assert(outer < outerSize()); 00298 eigen_assert(inner < innerSize()); 00299 eigen_assert(inner != outer); 00300 00301 if (IsRowMajor) { 00302 const Index minInnerIndex = outer - m_data.lowerProfile(outer); 00303 eigen_assert(inner >= minInnerIndex && "you try to acces a coeff that do not exist in the storage"); 00304 return this->m_data.lower(m_rowStartIndex[outer] + inner - (outer - m_data.lowerProfile(outer))); 00305 } else { 00306 const Index maxInnerIndex = outer + m_data.lowerProfile(outer); 00307 eigen_assert(inner <= maxInnerIndex && "you try to acces a coeff that do not exist in the storage"); 00308 return this->m_data.lower(m_rowStartIndex[outer] + (inner - outer)); 00309 } 00310 } 00311 00312 inline bool coeffExistLower(Index row, Index col) { 00313 const Index outer = IsRowMajor ? row : col; 00314 const Index inner = IsRowMajor ? col : row; 00315 00316 eigen_assert(outer < outerSize()); 00317 eigen_assert(inner < innerSize()); 00318 eigen_assert(inner != outer); 00319 00320 if (IsRowMajor) { 00321 const Index minInnerIndex = outer - m_data.lowerProfile(outer); 00322 return inner >= minInnerIndex; 00323 } else { 00324 const Index maxInnerIndex = outer + m_data.lowerProfile(outer); 00325 return inner <= maxInnerIndex; 00326 } 00327 } 00328 00329 inline Scalar& coeffRefUpper(Index row, Index col) { 00330 const Index outer = IsRowMajor ? row : col; 00331 const Index inner = IsRowMajor ? col : row; 00332 00333 eigen_assert(outer < outerSize()); 00334 eigen_assert(inner < innerSize()); 00335 eigen_assert(inner != outer); 00336 00337 if (IsRowMajor) { 00338 const Index minOuterIndex = inner - m_data.upperProfile(inner); 00339 eigen_assert(outer >= minOuterIndex && "you try to acces a coeff that do not exist in the storage"); 00340 return this->m_data.upper(m_colStartIndex[inner] + outer - (inner - m_data.upperProfile(inner))); 00341 } else { 00342 const Index maxOuterIndex = inner + m_data.upperProfile(inner); 00343 eigen_assert(outer <= maxOuterIndex && "you try to acces a coeff that do not exist in the storage"); 00344 return this->m_data.upper(m_colStartIndex[inner] + (outer - inner)); 00345 } 00346 } 00347 00348 inline bool coeffExistUpper(Index row, Index col) { 00349 const Index outer = IsRowMajor ? row : col; 00350 const Index inner = IsRowMajor ? col : row; 00351 00352 eigen_assert(outer < outerSize()); 00353 eigen_assert(inner < innerSize()); 00354 eigen_assert(inner != outer); 00355 00356 if (IsRowMajor) { 00357 const Index minOuterIndex = inner - m_data.upperProfile(inner); 00358 return outer >= minOuterIndex; 00359 } else { 00360 const Index maxOuterIndex = inner + m_data.upperProfile(inner); 00361 return outer <= maxOuterIndex; 00362 } 00363 } 00364 00365 00366 protected: 00367 00368 public: 00369 class InnerUpperIterator; 00370 class InnerLowerIterator; 00371 00372 class OuterUpperIterator; 00373 class OuterLowerIterator; 00374 00376 inline void setZero() { 00377 m_data.clear(); 00378 memset(m_colStartIndex, 0, (m_outerSize + 1) * sizeof (Index)); 00379 memset(m_rowStartIndex, 0, (m_outerSize + 1) * sizeof (Index)); 00380 } 00381 00383 inline Index nonZeros() const { 00384 return m_data.diagSize() + m_data.upperSize() + m_data.lowerSize(); 00385 } 00386 00388 inline void reserve(Index reserveSize, Index reserveUpperSize, Index reserveLowerSize) { 00389 m_data.reserve(reserveSize, reserveUpperSize, reserveLowerSize); 00390 } 00391 00400 EIGEN_DONT_INLINE Scalar & insert(Index row, Index col) { 00401 const Index outer = IsRowMajor ? row : col; 00402 const Index inner = IsRowMajor ? col : row; 00403 00404 eigen_assert(outer < outerSize()); 00405 eigen_assert(inner < innerSize()); 00406 00407 if (outer == inner) 00408 return m_data.diag(col); 00409 00410 if (IsRowMajor) { 00411 if (outer < inner) //upper matrix 00412 { 00413 Index minOuterIndex = 0; 00414 minOuterIndex = inner - m_data.upperProfile(inner); 00415 00416 if (outer < minOuterIndex) //The value does not yet exist 00417 { 00418 const Index previousProfile = m_data.upperProfile(inner); 00419 00420 m_data.upperProfile(inner) = inner - outer; 00421 00422 00423 const Index bandIncrement = m_data.upperProfile(inner) - previousProfile; 00424 //shift data stored after this new one 00425 const Index stop = m_colStartIndex[cols()]; 00426 const Index start = m_colStartIndex[inner]; 00427 00428 00429 for (Index innerIdx = stop; innerIdx >= start; innerIdx--) { 00430 m_data.upper(innerIdx + bandIncrement) = m_data.upper(innerIdx); 00431 } 00432 00433 for (Index innerIdx = cols(); innerIdx > inner; innerIdx--) { 00434 m_colStartIndex[innerIdx] += bandIncrement; 00435 } 00436 00437 //zeros new data 00438 memset(this->_upperPtr() + start, 0, (bandIncrement - 1) * sizeof (Scalar)); 00439 00440 return m_data.upper(m_colStartIndex[inner]); 00441 } else { 00442 return m_data.upper(m_colStartIndex[inner] + outer - (inner - m_data.upperProfile(inner))); 00443 } 00444 } 00445 00446 if (outer > inner) //lower matrix 00447 { 00448 const Index minInnerIndex = outer - m_data.lowerProfile(outer); 00449 if (inner < minInnerIndex) //The value does not yet exist 00450 { 00451 const Index previousProfile = m_data.lowerProfile(outer); 00452 m_data.lowerProfile(outer) = outer - inner; 00453 00454 const Index bandIncrement = m_data.lowerProfile(outer) - previousProfile; 00455 //shift data stored after this new one 00456 const Index stop = m_rowStartIndex[rows()]; 00457 const Index start = m_rowStartIndex[outer]; 00458 00459 00460 for (Index innerIdx = stop; innerIdx >= start; innerIdx--) { 00461 m_data.lower(innerIdx + bandIncrement) = m_data.lower(innerIdx); 00462 } 00463 00464 for (Index innerIdx = rows(); innerIdx > outer; innerIdx--) { 00465 m_rowStartIndex[innerIdx] += bandIncrement; 00466 } 00467 00468 //zeros new data 00469 memset(this->_lowerPtr() + start, 0, (bandIncrement - 1) * sizeof (Scalar)); 00470 return m_data.lower(m_rowStartIndex[outer]); 00471 } else { 00472 return m_data.lower(m_rowStartIndex[outer] + inner - (outer - m_data.lowerProfile(outer))); 00473 } 00474 } 00475 } else { 00476 if (outer > inner) //upper matrix 00477 { 00478 const Index maxOuterIndex = inner + m_data.upperProfile(inner); 00479 if (outer > maxOuterIndex) //The value does not yet exist 00480 { 00481 const Index previousProfile = m_data.upperProfile(inner); 00482 m_data.upperProfile(inner) = outer - inner; 00483 00484 const Index bandIncrement = m_data.upperProfile(inner) - previousProfile; 00485 //shift data stored after this new one 00486 const Index stop = m_rowStartIndex[rows()]; 00487 const Index start = m_rowStartIndex[inner + 1]; 00488 00489 for (Index innerIdx = stop; innerIdx >= start; innerIdx--) { 00490 m_data.upper(innerIdx + bandIncrement) = m_data.upper(innerIdx); 00491 } 00492 00493 for (Index innerIdx = inner + 1; innerIdx < outerSize() + 1; innerIdx++) { 00494 m_rowStartIndex[innerIdx] += bandIncrement; 00495 } 00496 memset(this->_upperPtr() + m_rowStartIndex[inner] + previousProfile + 1, 0, (bandIncrement - 1) * sizeof (Scalar)); 00497 return m_data.upper(m_rowStartIndex[inner] + m_data.upperProfile(inner)); 00498 } else { 00499 return m_data.upper(m_rowStartIndex[inner] + (outer - inner)); 00500 } 00501 } 00502 00503 if (outer < inner) //lower matrix 00504 { 00505 const Index maxInnerIndex = outer + m_data.lowerProfile(outer); 00506 if (inner > maxInnerIndex) //The value does not yet exist 00507 { 00508 const Index previousProfile = m_data.lowerProfile(outer); 00509 m_data.lowerProfile(outer) = inner - outer; 00510 00511 const Index bandIncrement = m_data.lowerProfile(outer) - previousProfile; 00512 //shift data stored after this new one 00513 const Index stop = m_colStartIndex[cols()]; 00514 const Index start = m_colStartIndex[outer + 1]; 00515 00516 for (Index innerIdx = stop; innerIdx >= start; innerIdx--) { 00517 m_data.lower(innerIdx + bandIncrement) = m_data.lower(innerIdx); 00518 } 00519 00520 for (Index innerIdx = outer + 1; innerIdx < outerSize() + 1; innerIdx++) { 00521 m_colStartIndex[innerIdx] += bandIncrement; 00522 } 00523 memset(this->_lowerPtr() + m_colStartIndex[outer] + previousProfile + 1, 0, (bandIncrement - 1) * sizeof (Scalar)); 00524 return m_data.lower(m_colStartIndex[outer] + m_data.lowerProfile(outer)); 00525 } else { 00526 return m_data.lower(m_colStartIndex[outer] + (inner - outer)); 00527 } 00528 } 00529 } 00530 } 00531 00534 inline void finalize() { 00535 if (IsRowMajor) { 00536 if (rows() > cols()) 00537 m_data.resize(cols(), cols(), rows(), m_colStartIndex[cols()] + 1, m_rowStartIndex[rows()] + 1); 00538 else 00539 m_data.resize(rows(), cols(), rows(), m_colStartIndex[cols()] + 1, m_rowStartIndex[rows()] + 1); 00540 00541 // eigen_assert(rows() == cols() && "memory reorganisatrion only works with suare matrix"); 00542 // 00543 // Scalar* newArray = new Scalar[m_colStartIndex[cols()] + 1 + m_rowStartIndex[rows()] + 1]; 00544 // Index dataIdx = 0; 00545 // for (Index row = 0; row < rows(); row++) { 00546 // 00547 // const Index nbLowerElts = m_rowStartIndex[row + 1] - m_rowStartIndex[row]; 00548 // // std::cout << "nbLowerElts" << nbLowerElts << std::endl; 00549 // memcpy(newArray + dataIdx, m_data.m_lower + m_rowStartIndex[row], nbLowerElts * sizeof (Scalar)); 00550 // m_rowStartIndex[row] = dataIdx; 00551 // dataIdx += nbLowerElts; 00552 // 00553 // const Index nbUpperElts = m_colStartIndex[row + 1] - m_colStartIndex[row]; 00554 // memcpy(newArray + dataIdx, m_data.m_upper + m_colStartIndex[row], nbUpperElts * sizeof (Scalar)); 00555 // m_colStartIndex[row] = dataIdx; 00556 // dataIdx += nbUpperElts; 00557 // 00558 // 00559 // } 00560 // //todo : don't access m_data profile directly : add an accessor from SkylineMatrix 00561 // m_rowStartIndex[rows()] = m_rowStartIndex[rows()-1] + m_data.lowerProfile(rows()-1); 00562 // m_colStartIndex[cols()] = m_colStartIndex[cols()-1] + m_data.upperProfile(cols()-1); 00563 // 00564 // delete[] m_data.m_lower; 00565 // delete[] m_data.m_upper; 00566 // 00567 // m_data.m_lower = newArray; 00568 // m_data.m_upper = newArray; 00569 } else { 00570 if (rows() > cols()) 00571 m_data.resize(cols(), rows(), cols(), m_rowStartIndex[cols()] + 1, m_colStartIndex[cols()] + 1); 00572 else 00573 m_data.resize(rows(), rows(), cols(), m_rowStartIndex[rows()] + 1, m_colStartIndex[rows()] + 1); 00574 } 00575 } 00576 00577 inline void squeeze() { 00578 finalize(); 00579 m_data.squeeze(); 00580 } 00581 00582 void prune(Scalar reference, RealScalar epsilon = dummy_precision<RealScalar > ()) { 00583 //TODO 00584 } 00585 00589 void resize(size_t rows, size_t cols) { 00590 const Index diagSize = rows > cols ? cols : rows; 00591 m_innerSize = IsRowMajor ? cols : rows; 00592 00593 eigen_assert(rows == cols && "Skyline matrix must be square matrix"); 00594 00595 if (diagSize % 2) { // diagSize is odd 00596 const Index k = (diagSize - 1) / 2; 00597 00598 m_data.resize(diagSize, IsRowMajor ? cols : rows, IsRowMajor ? rows : cols, 00599 2 * k * k + k + 1, 00600 2 * k * k + k + 1); 00601 00602 } else // diagSize is even 00603 { 00604 const Index k = diagSize / 2; 00605 m_data.resize(diagSize, IsRowMajor ? cols : rows, IsRowMajor ? rows : cols, 00606 2 * k * k - k + 1, 00607 2 * k * k - k + 1); 00608 } 00609 00610 if (m_colStartIndex && m_rowStartIndex) { 00611 delete[] m_colStartIndex; 00612 delete[] m_rowStartIndex; 00613 } 00614 m_colStartIndex = new Index [cols + 1]; 00615 m_rowStartIndex = new Index [rows + 1]; 00616 m_outerSize = diagSize; 00617 00618 m_data.reset(); 00619 m_data.clear(); 00620 00621 m_outerSize = diagSize; 00622 memset(m_colStartIndex, 0, (cols + 1) * sizeof (Index)); 00623 memset(m_rowStartIndex, 0, (rows + 1) * sizeof (Index)); 00624 } 00625 00626 void resizeNonZeros(Index size) { 00627 m_data.resize(size); 00628 } 00629 00630 inline SkylineMatrix() 00631 : m_outerSize(-1), m_innerSize(0), m_colStartIndex(0), m_rowStartIndex(0) { 00632 resize(0, 0); 00633 } 00634 00635 inline SkylineMatrix(size_t rows, size_t cols) 00636 : m_outerSize(0), m_innerSize(0), m_colStartIndex(0), m_rowStartIndex(0) { 00637 resize(rows, cols); 00638 } 00639 00640 template<typename OtherDerived> 00641 inline SkylineMatrix(const SkylineMatrixBase<OtherDerived>& other) 00642 : m_outerSize(0), m_innerSize(0), m_colStartIndex(0), m_rowStartIndex(0) { 00643 *this = other.derived(); 00644 } 00645 00646 inline SkylineMatrix(const SkylineMatrix & other) 00647 : Base(), m_outerSize(0), m_innerSize(0), m_colStartIndex(0), m_rowStartIndex(0) { 00648 *this = other.derived(); 00649 } 00650 00651 inline void swap(SkylineMatrix & other) { 00652 //EIGEN_DBG_SKYLINE(std::cout << "SkylineMatrix:: swap\n"); 00653 std::swap(m_colStartIndex, other.m_colStartIndex); 00654 std::swap(m_rowStartIndex, other.m_rowStartIndex); 00655 std::swap(m_innerSize, other.m_innerSize); 00656 std::swap(m_outerSize, other.m_outerSize); 00657 m_data.swap(other.m_data); 00658 } 00659 00660 inline SkylineMatrix & operator=(const SkylineMatrix & other) { 00661 std::cout << "SkylineMatrix& operator=(const SkylineMatrix& other)\n"; 00662 if (other.isRValue()) { 00663 swap(other.const_cast_derived()); 00664 } else { 00665 resize(other.rows(), other.cols()); 00666 memcpy(m_colStartIndex, other.m_colStartIndex, (m_outerSize + 1) * sizeof (Index)); 00667 memcpy(m_rowStartIndex, other.m_rowStartIndex, (m_outerSize + 1) * sizeof (Index)); 00668 m_data = other.m_data; 00669 } 00670 return *this; 00671 } 00672 00673 template<typename OtherDerived> 00674 inline SkylineMatrix & operator=(const SkylineMatrixBase<OtherDerived>& other) { 00675 const bool needToTranspose = (Flags & RowMajorBit) != (OtherDerived::Flags & RowMajorBit); 00676 if (needToTranspose) { 00677 // TODO 00678 // return *this; 00679 } else { 00680 // there is no special optimization 00681 return SkylineMatrixBase<SkylineMatrix>::operator=(other.derived()); 00682 } 00683 } 00684 00685 friend std::ostream & operator <<(std::ostream & s, const SkylineMatrix & m) { 00686 00687 EIGEN_DBG_SKYLINE( 00688 std::cout << "upper elements : " << std::endl; 00689 for (Index i = 0; i < m.m_data.upperSize(); i++) 00690 std::cout << m.m_data.upper(i) << "\t"; 00691 std::cout << std::endl; 00692 std::cout << "upper profile : " << std::endl; 00693 for (Index i = 0; i < m.m_data.upperProfileSize(); i++) 00694 std::cout << m.m_data.upperProfile(i) << "\t"; 00695 std::cout << std::endl; 00696 std::cout << "lower startIdx : " << std::endl; 00697 for (Index i = 0; i < m.m_data.upperProfileSize(); i++) 00698 std::cout << (IsRowMajor ? m.m_colStartIndex[i] : m.m_rowStartIndex[i]) << "\t"; 00699 std::cout << std::endl; 00700 00701 00702 std::cout << "lower elements : " << std::endl; 00703 for (Index i = 0; i < m.m_data.lowerSize(); i++) 00704 std::cout << m.m_data.lower(i) << "\t"; 00705 std::cout << std::endl; 00706 std::cout << "lower profile : " << std::endl; 00707 for (Index i = 0; i < m.m_data.lowerProfileSize(); i++) 00708 std::cout << m.m_data.lowerProfile(i) << "\t"; 00709 std::cout << std::endl; 00710 std::cout << "lower startIdx : " << std::endl; 00711 for (Index i = 0; i < m.m_data.lowerProfileSize(); i++) 00712 std::cout << (IsRowMajor ? m.m_rowStartIndex[i] : m.m_colStartIndex[i]) << "\t"; 00713 std::cout << std::endl; 00714 ); 00715 for (Index rowIdx = 0; rowIdx < m.rows(); rowIdx++) { 00716 for (Index colIdx = 0; colIdx < m.cols(); colIdx++) { 00717 s << m.coeff(rowIdx, colIdx) << "\t"; 00718 } 00719 s << std::endl; 00720 } 00721 return s; 00722 } 00723 00725 inline ~SkylineMatrix() { 00726 delete[] m_colStartIndex; 00727 delete[] m_rowStartIndex; 00728 } 00729 00731 Scalar sum() const; 00732 }; 00733 00734 template<typename Scalar, int _Options> 00735 class SkylineMatrix<Scalar, _Options>::InnerUpperIterator { 00736 public: 00737 00738 InnerUpperIterator(const SkylineMatrix& mat, Index outer) 00739 : m_matrix(mat), m_outer(outer), 00740 m_id(_Options == RowMajor ? mat.m_colStartIndex[outer] : mat.m_rowStartIndex[outer] + 1), 00741 m_start(m_id), 00742 m_end(_Options == RowMajor ? mat.m_colStartIndex[outer + 1] : mat.m_rowStartIndex[outer + 1] + 1) { 00743 } 00744 00745 inline InnerUpperIterator & operator++() { 00746 m_id++; 00747 return *this; 00748 } 00749 00750 inline InnerUpperIterator & operator+=(Index shift) { 00751 m_id += shift; 00752 return *this; 00753 } 00754 00755 inline Scalar value() const { 00756 return m_matrix.m_data.upper(m_id); 00757 } 00758 00759 inline Scalar* valuePtr() { 00760 return const_cast<Scalar*> (&(m_matrix.m_data.upper(m_id))); 00761 } 00762 00763 inline Scalar& valueRef() { 00764 return const_cast<Scalar&> (m_matrix.m_data.upper(m_id)); 00765 } 00766 00767 inline Index index() const { 00768 return IsRowMajor ? m_outer - m_matrix.m_data.upperProfile(m_outer) + (m_id - m_start) : 00769 m_outer + (m_id - m_start) + 1; 00770 } 00771 00772 inline Index row() const { 00773 return IsRowMajor ? index() : m_outer; 00774 } 00775 00776 inline Index col() const { 00777 return IsRowMajor ? m_outer : index(); 00778 } 00779 00780 inline size_t size() const { 00781 return m_matrix.m_data.upperProfile(m_outer); 00782 } 00783 00784 inline operator bool() const { 00785 return (m_id < m_end) && (m_id >= m_start); 00786 } 00787 00788 protected: 00789 const SkylineMatrix& m_matrix; 00790 const Index m_outer; 00791 Index m_id; 00792 const Index m_start; 00793 const Index m_end; 00794 }; 00795 00796 template<typename Scalar, int _Options> 00797 class SkylineMatrix<Scalar, _Options>::InnerLowerIterator { 00798 public: 00799 00800 InnerLowerIterator(const SkylineMatrix& mat, Index outer) 00801 : m_matrix(mat), 00802 m_outer(outer), 00803 m_id(_Options == RowMajor ? mat.m_rowStartIndex[outer] : mat.m_colStartIndex[outer] + 1), 00804 m_start(m_id), 00805 m_end(_Options == RowMajor ? mat.m_rowStartIndex[outer + 1] : mat.m_colStartIndex[outer + 1] + 1) { 00806 } 00807 00808 inline InnerLowerIterator & operator++() { 00809 m_id++; 00810 return *this; 00811 } 00812 00813 inline InnerLowerIterator & operator+=(Index shift) { 00814 m_id += shift; 00815 return *this; 00816 } 00817 00818 inline Scalar value() const { 00819 return m_matrix.m_data.lower(m_id); 00820 } 00821 00822 inline Scalar* valuePtr() { 00823 return const_cast<Scalar*> (&(m_matrix.m_data.lower(m_id))); 00824 } 00825 00826 inline Scalar& valueRef() { 00827 return const_cast<Scalar&> (m_matrix.m_data.lower(m_id)); 00828 } 00829 00830 inline Index index() const { 00831 return IsRowMajor ? m_outer - m_matrix.m_data.lowerProfile(m_outer) + (m_id - m_start) : 00832 m_outer + (m_id - m_start) + 1; 00833 ; 00834 } 00835 00836 inline Index row() const { 00837 return IsRowMajor ? m_outer : index(); 00838 } 00839 00840 inline Index col() const { 00841 return IsRowMajor ? index() : m_outer; 00842 } 00843 00844 inline size_t size() const { 00845 return m_matrix.m_data.lowerProfile(m_outer); 00846 } 00847 00848 inline operator bool() const { 00849 return (m_id < m_end) && (m_id >= m_start); 00850 } 00851 00852 protected: 00853 const SkylineMatrix& m_matrix; 00854 const Index m_outer; 00855 Index m_id; 00856 const Index m_start; 00857 const Index m_end; 00858 }; 00859 00860 } // end namespace Eigen 00861 00862 #endif // EIGEN_SkylineMatrix_H