CompressedStorage.h
Go to the documentation of this file.
00001 // This file is part of Eigen, a lightweight C++ template library
00002 // for linear algebra.
00003 //
00004 // Copyright (C) 2008 Gael Guennebaud <gael.guennebaud@inria.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_COMPRESSED_STORAGE_H
00026 #define EIGEN_COMPRESSED_STORAGE_H
00027 
00028 namespace Eigen { 
00029 
00030 namespace internal {
00031 
00036 template<typename _Scalar,typename _Index>
00037 class CompressedStorage
00038 {
00039   public:
00040 
00041     typedef _Scalar Scalar;
00042     typedef _Index Index;
00043 
00044   protected:
00045 
00046     typedef typename NumTraits<Scalar>::Real RealScalar;
00047 
00048   public:
00049 
00050     CompressedStorage()
00051       : m_values(0), m_indices(0), m_size(0), m_allocatedSize(0)
00052     {}
00053 
00054     CompressedStorage(size_t size)
00055       : m_values(0), m_indices(0), m_size(0), m_allocatedSize(0)
00056     {
00057       resize(size);
00058     }
00059 
00060     CompressedStorage(const CompressedStorage& other)
00061       : m_values(0), m_indices(0), m_size(0), m_allocatedSize(0)
00062     {
00063       *this = other;
00064     }
00065 
00066     CompressedStorage& operator=(const CompressedStorage& other)
00067     {
00068       resize(other.size());
00069       memcpy(m_values, other.m_values, m_size * sizeof(Scalar));
00070       memcpy(m_indices, other.m_indices, m_size * sizeof(Index));
00071       return *this;
00072     }
00073 
00074     void swap(CompressedStorage& other)
00075     {
00076       std::swap(m_values, other.m_values);
00077       std::swap(m_indices, other.m_indices);
00078       std::swap(m_size, other.m_size);
00079       std::swap(m_allocatedSize, other.m_allocatedSize);
00080     }
00081 
00082     ~CompressedStorage()
00083     {
00084       delete[] m_values;
00085       delete[] m_indices;
00086     }
00087 
00088     void reserve(size_t size)
00089     {
00090       size_t newAllocatedSize = m_size + size;
00091       if (newAllocatedSize > m_allocatedSize)
00092         reallocate(newAllocatedSize);
00093     }
00094 
00095     void squeeze()
00096     {
00097       if (m_allocatedSize>m_size)
00098         reallocate(m_size);
00099     }
00100 
00101     void resize(size_t size, float reserveSizeFactor = 0)
00102     {
00103       if (m_allocatedSize<size)
00104         reallocate(size + size_t(reserveSizeFactor*size));
00105       m_size = size;
00106     }
00107 
00108     void append(const Scalar& v, Index i)
00109     {
00110       Index id = static_cast<Index>(m_size);
00111       resize(m_size+1, 1);
00112       m_values[id] = v;
00113       m_indices[id] = i;
00114     }
00115 
00116     inline size_t size() const { return m_size; }
00117     inline size_t allocatedSize() const { return m_allocatedSize; }
00118     inline void clear() { m_size = 0; }
00119 
00120     inline Scalar& value(size_t i) { return m_values[i]; }
00121     inline const Scalar& value(size_t i) const { return m_values[i]; }
00122 
00123     inline Index& index(size_t i) { return m_indices[i]; }
00124     inline const Index& index(size_t i) const { return m_indices[i]; }
00125 
00126     static CompressedStorage Map(Index* indices, Scalar* values, size_t size)
00127     {
00128       CompressedStorage res;
00129       res.m_indices = indices;
00130       res.m_values = values;
00131       res.m_allocatedSize = res.m_size = size;
00132       return res;
00133     }
00134 
00136     inline Index searchLowerIndex(Index key) const
00137     {
00138       return searchLowerIndex(0, m_size, key);
00139     }
00140 
00142     inline Index searchLowerIndex(size_t start, size_t end, Index key) const
00143     {
00144       while(end>start)
00145       {
00146         size_t mid = (end+start)>>1;
00147         if (m_indices[mid]<key)
00148           start = mid+1;
00149         else
00150           end = mid;
00151       }
00152       return static_cast<Index>(start);
00153     }
00154 
00157     inline Scalar at(Index key, Scalar defaultValue = Scalar(0)) const
00158     {
00159       if (m_size==0)
00160         return defaultValue;
00161       else if (key==m_indices[m_size-1])
00162         return m_values[m_size-1];
00163       // ^^  optimization: let's first check if it is the last coefficient
00164       // (very common in high level algorithms)
00165       const size_t id = searchLowerIndex(0,m_size-1,key);
00166       return ((id<m_size) && (m_indices[id]==key)) ? m_values[id] : defaultValue;
00167     }
00168 
00170     inline Scalar atInRange(size_t start, size_t end, Index key, Scalar defaultValue = Scalar(0)) const
00171     {
00172       if (start>=end)
00173         return Scalar(0);
00174       else if (end>start && key==m_indices[end-1])
00175         return m_values[end-1];
00176       // ^^  optimization: let's first check if it is the last coefficient
00177       // (very common in high level algorithms)
00178       const size_t id = searchLowerIndex(start,end-1,key);
00179       return ((id<end) && (m_indices[id]==key)) ? m_values[id] : defaultValue;
00180     }
00181 
00185     inline Scalar& atWithInsertion(Index key, Scalar defaultValue = Scalar(0))
00186     {
00187       size_t id = searchLowerIndex(0,m_size,key);
00188       if (id>=m_size || m_indices[id]!=key)
00189       {
00190         resize(m_size+1,1);
00191         for (size_t j=m_size-1; j>id; --j)
00192         {
00193           m_indices[j] = m_indices[j-1];
00194           m_values[j] = m_values[j-1];
00195         }
00196         m_indices[id] = key;
00197         m_values[id] = defaultValue;
00198       }
00199       return m_values[id];
00200     }
00201 
00202     void prune(Scalar reference, RealScalar epsilon = NumTraits<RealScalar>::dummy_precision())
00203     {
00204       size_t k = 0;
00205       size_t n = size();
00206       for (size_t i=0; i<n; ++i)
00207       {
00208         if (!internal::isMuchSmallerThan(value(i), reference, epsilon))
00209         {
00210           value(k) = value(i);
00211           index(k) = index(i);
00212           ++k;
00213         }
00214       }
00215       resize(k,0);
00216     }
00217 
00218   protected:
00219 
00220     inline void reallocate(size_t size)
00221     {
00222       Scalar* newValues  = new Scalar[size];
00223       Index* newIndices = new Index[size];
00224       size_t copySize = (std::min)(size, m_size);
00225       // copy
00226       internal::smart_copy(m_values, m_values+copySize, newValues);
00227       internal::smart_copy(m_indices, m_indices+copySize, newIndices);
00228       // delete old stuff
00229       delete[] m_values;
00230       delete[] m_indices;
00231       m_values = newValues;
00232       m_indices = newIndices;
00233       m_allocatedSize = size;
00234     }
00235 
00236   protected:
00237     Scalar* m_values;
00238     Index* m_indices;
00239     size_t m_size;
00240     size_t m_allocatedSize;
00241 
00242 };
00243 
00244 } // end namespace internal
00245 
00246 } // end namespace Eigen
00247 
00248 #endif // EIGEN_COMPRESSED_STORAGE_H