SkylineStorage.h
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 // Eigen is free software; you can redistribute it and/or
00007 // modify it under the terms of the GNU Lesser General Public
00008 // License as published by the Free Software Foundation; either
00009 // version 3 of the License, or (at your option) any later version.
00010 //
00011 // Alternatively, you can redistribute it and/or
00012 // modify it under the terms of the GNU General Public License as
00013 // published by the Free Software Foundation; either version 2 of
00014 // the License, or (at your option) any later version.
00015 //
00016 // Eigen is distributed in the hope that it will be useful, but WITHOUT ANY
00017 // WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS
00018 // FOR A PARTICULAR PURPOSE. See the GNU Lesser General Public License or the
00019 // GNU General Public License for more details.
00020 //
00021 // You should have received a copy of the GNU Lesser General Public
00022 // License and a copy of the GNU General Public License along with
00023 // Eigen. If not, see <http://www.gnu.org/licenses/>.
00024 
00025 #ifndef EIGEN_SKYLINE_STORAGE_H
00026 #define EIGEN_SKYLINE_STORAGE_H
00027 
00028 namespace Eigen { 
00029 
00036 template<typename Scalar>
00037 class SkylineStorage {
00038     typedef typename NumTraits<Scalar>::Real RealScalar;
00039     typedef SparseIndex Index;
00040 public:
00041 
00042     SkylineStorage()
00043     : m_diag(0),
00044     m_lower(0),
00045     m_upper(0),
00046     m_lowerProfile(0),
00047     m_upperProfile(0),
00048     m_diagSize(0),
00049     m_upperSize(0),
00050     m_lowerSize(0),
00051     m_upperProfileSize(0),
00052     m_lowerProfileSize(0),
00053     m_allocatedSize(0) {
00054     }
00055 
00056     SkylineStorage(const SkylineStorage& other)
00057     : m_diag(0),
00058     m_lower(0),
00059     m_upper(0),
00060     m_lowerProfile(0),
00061     m_upperProfile(0),
00062     m_diagSize(0),
00063     m_upperSize(0),
00064     m_lowerSize(0),
00065     m_upperProfileSize(0),
00066     m_lowerProfileSize(0),
00067     m_allocatedSize(0) {
00068         *this = other;
00069     }
00070 
00071     SkylineStorage & operator=(const SkylineStorage& other) {
00072         resize(other.diagSize(), other.m_upperProfileSize, other.m_lowerProfileSize, other.upperSize(), other.lowerSize());
00073         memcpy(m_diag, other.m_diag, m_diagSize * sizeof (Scalar));
00074         memcpy(m_upper, other.m_upper, other.upperSize() * sizeof (Scalar));
00075         memcpy(m_lower, other.m_lower, other.lowerSize() * sizeof (Scalar));
00076         memcpy(m_upperProfile, other.m_upperProfile, m_upperProfileSize * sizeof (Index));
00077         memcpy(m_lowerProfile, other.m_lowerProfile, m_lowerProfileSize * sizeof (Index));
00078         return *this;
00079     }
00080 
00081     void swap(SkylineStorage& other) {
00082         std::swap(m_diag, other.m_diag);
00083         std::swap(m_upper, other.m_upper);
00084         std::swap(m_lower, other.m_lower);
00085         std::swap(m_upperProfile, other.m_upperProfile);
00086         std::swap(m_lowerProfile, other.m_lowerProfile);
00087         std::swap(m_diagSize, other.m_diagSize);
00088         std::swap(m_upperSize, other.m_upperSize);
00089         std::swap(m_lowerSize, other.m_lowerSize);
00090         std::swap(m_allocatedSize, other.m_allocatedSize);
00091     }
00092 
00093     ~SkylineStorage() {
00094         delete[] m_diag;
00095         delete[] m_upper;
00096         if (m_upper != m_lower)
00097             delete[] m_lower;
00098         delete[] m_upperProfile;
00099         delete[] m_lowerProfile;
00100     }
00101 
00102     void reserve(Index size, Index upperProfileSize, Index lowerProfileSize, Index upperSize, Index lowerSize) {
00103         Index newAllocatedSize = size + upperSize + lowerSize;
00104         if (newAllocatedSize > m_allocatedSize)
00105             reallocate(size, upperProfileSize, lowerProfileSize, upperSize, lowerSize);
00106     }
00107 
00108     void squeeze() {
00109         if (m_allocatedSize > m_diagSize + m_upperSize + m_lowerSize)
00110             reallocate(m_diagSize, m_upperProfileSize, m_lowerProfileSize, m_upperSize, m_lowerSize);
00111     }
00112 
00113     void resize(Index diagSize, Index upperProfileSize, Index lowerProfileSize, Index upperSize, Index lowerSize, float reserveSizeFactor = 0) {
00114         if (m_allocatedSize < diagSize + upperSize + lowerSize)
00115             reallocate(diagSize, upperProfileSize, lowerProfileSize, upperSize + Index(reserveSizeFactor * upperSize), lowerSize + Index(reserveSizeFactor * lowerSize));
00116         m_diagSize = diagSize;
00117         m_upperSize = upperSize;
00118         m_lowerSize = lowerSize;
00119         m_upperProfileSize = upperProfileSize;
00120         m_lowerProfileSize = lowerProfileSize;
00121     }
00122 
00123     inline Index diagSize() const {
00124         return m_diagSize;
00125     }
00126 
00127     inline Index upperSize() const {
00128         return m_upperSize;
00129     }
00130 
00131     inline Index lowerSize() const {
00132         return m_lowerSize;
00133     }
00134 
00135     inline Index upperProfileSize() const {
00136         return m_upperProfileSize;
00137     }
00138 
00139     inline Index lowerProfileSize() const {
00140         return m_lowerProfileSize;
00141     }
00142 
00143     inline Index allocatedSize() const {
00144         return m_allocatedSize;
00145     }
00146 
00147     inline void clear() {
00148         m_diagSize = 0;
00149     }
00150 
00151     inline Scalar& diag(Index i) {
00152         return m_diag[i];
00153     }
00154 
00155     inline const Scalar& diag(Index i) const {
00156         return m_diag[i];
00157     }
00158 
00159     inline Scalar& upper(Index i) {
00160         return m_upper[i];
00161     }
00162 
00163     inline const Scalar& upper(Index i) const {
00164         return m_upper[i];
00165     }
00166 
00167     inline Scalar& lower(Index i) {
00168         return m_lower[i];
00169     }
00170 
00171     inline const Scalar& lower(Index i) const {
00172         return m_lower[i];
00173     }
00174 
00175     inline Index& upperProfile(Index i) {
00176         return m_upperProfile[i];
00177     }
00178 
00179     inline const Index& upperProfile(Index i) const {
00180         return m_upperProfile[i];
00181     }
00182 
00183     inline Index& lowerProfile(Index i) {
00184         return m_lowerProfile[i];
00185     }
00186 
00187     inline const Index& lowerProfile(Index i) const {
00188         return m_lowerProfile[i];
00189     }
00190 
00191     static SkylineStorage Map(Index* upperProfile, Index* lowerProfile, Scalar* diag, Scalar* upper, Scalar* lower, Index size, Index upperSize, Index lowerSize) {
00192         SkylineStorage res;
00193         res.m_upperProfile = upperProfile;
00194         res.m_lowerProfile = lowerProfile;
00195         res.m_diag = diag;
00196         res.m_upper = upper;
00197         res.m_lower = lower;
00198         res.m_allocatedSize = res.m_diagSize = size;
00199         res.m_upperSize = upperSize;
00200         res.m_lowerSize = lowerSize;
00201         return res;
00202     }
00203 
00204     inline void reset() {
00205         memset(m_diag, 0, m_diagSize * sizeof (Scalar));
00206         memset(m_upper, 0, m_upperSize * sizeof (Scalar));
00207         memset(m_lower, 0, m_lowerSize * sizeof (Scalar));
00208         memset(m_upperProfile, 0, m_diagSize * sizeof (Index));
00209         memset(m_lowerProfile, 0, m_diagSize * sizeof (Index));
00210     }
00211 
00212     void prune(Scalar reference, RealScalar epsilon = dummy_precision<RealScalar>()) {
00213         //TODO
00214     }
00215 
00216 protected:
00217 
00218     inline void reallocate(Index diagSize, Index upperProfileSize, Index lowerProfileSize, Index upperSize, Index lowerSize) {
00219 
00220         Scalar* diag = new Scalar[diagSize];
00221         Scalar* upper = new Scalar[upperSize];
00222         Scalar* lower = new Scalar[lowerSize];
00223         Index* upperProfile = new Index[upperProfileSize];
00224         Index* lowerProfile = new Index[lowerProfileSize];
00225 
00226         Index copyDiagSize = (std::min)(diagSize, m_diagSize);
00227         Index copyUpperSize = (std::min)(upperSize, m_upperSize);
00228         Index copyLowerSize = (std::min)(lowerSize, m_lowerSize);
00229         Index copyUpperProfileSize = (std::min)(upperProfileSize, m_upperProfileSize);
00230         Index copyLowerProfileSize = (std::min)(lowerProfileSize, m_lowerProfileSize);
00231 
00232         // copy
00233         memcpy(diag, m_diag, copyDiagSize * sizeof (Scalar));
00234         memcpy(upper, m_upper, copyUpperSize * sizeof (Scalar));
00235         memcpy(lower, m_lower, copyLowerSize * sizeof (Scalar));
00236         memcpy(upperProfile, m_upperProfile, copyUpperProfileSize * sizeof (Index));
00237         memcpy(lowerProfile, m_lowerProfile, copyLowerProfileSize * sizeof (Index));
00238 
00239 
00240 
00241         // delete old stuff
00242         delete[] m_diag;
00243         delete[] m_upper;
00244         delete[] m_lower;
00245         delete[] m_upperProfile;
00246         delete[] m_lowerProfile;
00247         m_diag = diag;
00248         m_upper = upper;
00249         m_lower = lower;
00250         m_upperProfile = upperProfile;
00251         m_lowerProfile = lowerProfile;
00252         m_allocatedSize = diagSize + upperSize + lowerSize;
00253         m_upperSize = upperSize;
00254         m_lowerSize = lowerSize;
00255     }
00256 
00257 public:
00258     Scalar* m_diag;
00259     Scalar* m_upper;
00260     Scalar* m_lower;
00261     Index* m_upperProfile;
00262     Index* m_lowerProfile;
00263     Index m_diagSize;
00264     Index m_upperSize;
00265     Index m_lowerSize;
00266     Index m_upperProfileSize;
00267     Index m_lowerProfileSize;
00268     Index m_allocatedSize;
00269 
00270 };
00271 
00272 } // end namespace Eigen
00273 
00274 #endif // EIGEN_COMPRESSED_STORAGE_H