1 | // This file is part of Eigen, a lightweight C++ template library
|
---|
2 | // for linear algebra.
|
---|
3 | //
|
---|
4 | // Copyright (C) 2008 Gael Guennebaud <gael.guennebaud@inria.fr>
|
---|
5 | //
|
---|
6 | // This Source Code Form is subject to the terms of the Mozilla
|
---|
7 | // Public License v. 2.0. If a copy of the MPL was not distributed
|
---|
8 | // with this file, You can obtain one at http://mozilla.org/MPL/2.0/.
|
---|
9 |
|
---|
10 | #ifndef EIGEN_COMPRESSED_STORAGE_H
|
---|
11 | #define EIGEN_COMPRESSED_STORAGE_H
|
---|
12 |
|
---|
13 | namespace Eigen {
|
---|
14 |
|
---|
15 | namespace internal {
|
---|
16 |
|
---|
17 | /** \internal
|
---|
18 | * Stores a sparse set of values as a list of values and a list of indices.
|
---|
19 | *
|
---|
20 | */
|
---|
21 | template<typename _Scalar,typename _Index>
|
---|
22 | class CompressedStorage
|
---|
23 | {
|
---|
24 | public:
|
---|
25 |
|
---|
26 | typedef _Scalar Scalar;
|
---|
27 | typedef _Index Index;
|
---|
28 |
|
---|
29 | protected:
|
---|
30 |
|
---|
31 | typedef typename NumTraits<Scalar>::Real RealScalar;
|
---|
32 |
|
---|
33 | public:
|
---|
34 |
|
---|
35 | CompressedStorage()
|
---|
36 | : m_values(0), m_indices(0), m_size(0), m_allocatedSize(0)
|
---|
37 | {}
|
---|
38 |
|
---|
39 | CompressedStorage(size_t size)
|
---|
40 | : m_values(0), m_indices(0), m_size(0), m_allocatedSize(0)
|
---|
41 | {
|
---|
42 | resize(size);
|
---|
43 | }
|
---|
44 |
|
---|
45 | CompressedStorage(const CompressedStorage& other)
|
---|
46 | : m_values(0), m_indices(0), m_size(0), m_allocatedSize(0)
|
---|
47 | {
|
---|
48 | *this = other;
|
---|
49 | }
|
---|
50 |
|
---|
51 | CompressedStorage& operator=(const CompressedStorage& other)
|
---|
52 | {
|
---|
53 | resize(other.size());
|
---|
54 | internal::smart_copy(other.m_values, other.m_values + m_size, m_values);
|
---|
55 | internal::smart_copy(other.m_indices, other.m_indices + m_size, m_indices);
|
---|
56 | return *this;
|
---|
57 | }
|
---|
58 |
|
---|
59 | void swap(CompressedStorage& other)
|
---|
60 | {
|
---|
61 | std::swap(m_values, other.m_values);
|
---|
62 | std::swap(m_indices, other.m_indices);
|
---|
63 | std::swap(m_size, other.m_size);
|
---|
64 | std::swap(m_allocatedSize, other.m_allocatedSize);
|
---|
65 | }
|
---|
66 |
|
---|
67 | ~CompressedStorage()
|
---|
68 | {
|
---|
69 | delete[] m_values;
|
---|
70 | delete[] m_indices;
|
---|
71 | }
|
---|
72 |
|
---|
73 | void reserve(size_t size)
|
---|
74 | {
|
---|
75 | size_t newAllocatedSize = m_size + size;
|
---|
76 | if (newAllocatedSize > m_allocatedSize)
|
---|
77 | reallocate(newAllocatedSize);
|
---|
78 | }
|
---|
79 |
|
---|
80 | void squeeze()
|
---|
81 | {
|
---|
82 | if (m_allocatedSize>m_size)
|
---|
83 | reallocate(m_size);
|
---|
84 | }
|
---|
85 |
|
---|
86 | void resize(size_t size, double reserveSizeFactor = 0)
|
---|
87 | {
|
---|
88 | if (m_allocatedSize<size)
|
---|
89 | reallocate(size + size_t(reserveSizeFactor*double(size)));
|
---|
90 | m_size = size;
|
---|
91 | }
|
---|
92 |
|
---|
93 | void append(const Scalar& v, Index i)
|
---|
94 | {
|
---|
95 | Index id = static_cast<Index>(m_size);
|
---|
96 | resize(m_size+1, 1);
|
---|
97 | m_values[id] = v;
|
---|
98 | m_indices[id] = i;
|
---|
99 | }
|
---|
100 |
|
---|
101 | inline size_t size() const { return m_size; }
|
---|
102 | inline size_t allocatedSize() const { return m_allocatedSize; }
|
---|
103 | inline void clear() { m_size = 0; }
|
---|
104 |
|
---|
105 | const Scalar* valuePtr() const { return m_values; }
|
---|
106 | Scalar* valuePtr() { return m_values; }
|
---|
107 | const Index* indexPtr() const { return m_indices; }
|
---|
108 | Index* indexPtr() { return m_indices; }
|
---|
109 |
|
---|
110 | inline Scalar& value(size_t i) { return m_values[i]; }
|
---|
111 | inline const Scalar& value(size_t i) const { return m_values[i]; }
|
---|
112 |
|
---|
113 | inline Index& index(size_t i) { return m_indices[i]; }
|
---|
114 | inline const Index& index(size_t i) const { return m_indices[i]; }
|
---|
115 |
|
---|
116 | static CompressedStorage Map(Index* indices, Scalar* values, size_t size)
|
---|
117 | {
|
---|
118 | CompressedStorage res;
|
---|
119 | res.m_indices = indices;
|
---|
120 | res.m_values = values;
|
---|
121 | res.m_allocatedSize = res.m_size = size;
|
---|
122 | return res;
|
---|
123 | }
|
---|
124 |
|
---|
125 | /** \returns the largest \c k such that for all \c j in [0,k) index[\c j]\<\a key */
|
---|
126 | inline Index searchLowerIndex(Index key) const
|
---|
127 | {
|
---|
128 | return searchLowerIndex(0, m_size, key);
|
---|
129 | }
|
---|
130 |
|
---|
131 | /** \returns the largest \c k in [start,end) such that for all \c j in [start,k) index[\c j]\<\a key */
|
---|
132 | inline Index searchLowerIndex(size_t start, size_t end, Index key) const
|
---|
133 | {
|
---|
134 | while(end>start)
|
---|
135 | {
|
---|
136 | size_t mid = (end+start)>>1;
|
---|
137 | if (m_indices[mid]<key)
|
---|
138 | start = mid+1;
|
---|
139 | else
|
---|
140 | end = mid;
|
---|
141 | }
|
---|
142 | return static_cast<Index>(start);
|
---|
143 | }
|
---|
144 |
|
---|
145 | /** \returns the stored value at index \a key
|
---|
146 | * If the value does not exist, then the value \a defaultValue is returned without any insertion. */
|
---|
147 | inline Scalar at(Index key, const Scalar& defaultValue = Scalar(0)) const
|
---|
148 | {
|
---|
149 | if (m_size==0)
|
---|
150 | return defaultValue;
|
---|
151 | else if (key==m_indices[m_size-1])
|
---|
152 | return m_values[m_size-1];
|
---|
153 | // ^^ optimization: let's first check if it is the last coefficient
|
---|
154 | // (very common in high level algorithms)
|
---|
155 | const size_t id = searchLowerIndex(0,m_size-1,key);
|
---|
156 | return ((id<m_size) && (m_indices[id]==key)) ? m_values[id] : defaultValue;
|
---|
157 | }
|
---|
158 |
|
---|
159 | /** Like at(), but the search is performed in the range [start,end) */
|
---|
160 | inline Scalar atInRange(size_t start, size_t end, Index key, const Scalar& defaultValue = Scalar(0)) const
|
---|
161 | {
|
---|
162 | if (start>=end)
|
---|
163 | return Scalar(0);
|
---|
164 | else if (end>start && key==m_indices[end-1])
|
---|
165 | return m_values[end-1];
|
---|
166 | // ^^ optimization: let's first check if it is the last coefficient
|
---|
167 | // (very common in high level algorithms)
|
---|
168 | const size_t id = searchLowerIndex(start,end-1,key);
|
---|
169 | return ((id<end) && (m_indices[id]==key)) ? m_values[id] : defaultValue;
|
---|
170 | }
|
---|
171 |
|
---|
172 | /** \returns a reference to the value at index \a key
|
---|
173 | * If the value does not exist, then the value \a defaultValue is inserted
|
---|
174 | * such that the keys are sorted. */
|
---|
175 | inline Scalar& atWithInsertion(Index key, const Scalar& defaultValue = Scalar(0))
|
---|
176 | {
|
---|
177 | size_t id = searchLowerIndex(0,m_size,key);
|
---|
178 | if (id>=m_size || m_indices[id]!=key)
|
---|
179 | {
|
---|
180 | resize(m_size+1,1);
|
---|
181 | for (size_t j=m_size-1; j>id; --j)
|
---|
182 | {
|
---|
183 | m_indices[j] = m_indices[j-1];
|
---|
184 | m_values[j] = m_values[j-1];
|
---|
185 | }
|
---|
186 | m_indices[id] = key;
|
---|
187 | m_values[id] = defaultValue;
|
---|
188 | }
|
---|
189 | return m_values[id];
|
---|
190 | }
|
---|
191 |
|
---|
192 | void prune(const Scalar& reference, const RealScalar& epsilon = NumTraits<RealScalar>::dummy_precision())
|
---|
193 | {
|
---|
194 | size_t k = 0;
|
---|
195 | size_t n = size();
|
---|
196 | for (size_t i=0; i<n; ++i)
|
---|
197 | {
|
---|
198 | if (!internal::isMuchSmallerThan(value(i), reference, epsilon))
|
---|
199 | {
|
---|
200 | value(k) = value(i);
|
---|
201 | index(k) = index(i);
|
---|
202 | ++k;
|
---|
203 | }
|
---|
204 | }
|
---|
205 | resize(k,0);
|
---|
206 | }
|
---|
207 |
|
---|
208 | protected:
|
---|
209 |
|
---|
210 | inline void reallocate(size_t size)
|
---|
211 | {
|
---|
212 | Scalar* newValues = new Scalar[size];
|
---|
213 | Index* newIndices = new Index[size];
|
---|
214 | size_t copySize = (std::min)(size, m_size);
|
---|
215 | // copy
|
---|
216 | if (copySize>0) {
|
---|
217 | internal::smart_copy(m_values, m_values+copySize, newValues);
|
---|
218 | internal::smart_copy(m_indices, m_indices+copySize, newIndices);
|
---|
219 | }
|
---|
220 | // delete old stuff
|
---|
221 | delete[] m_values;
|
---|
222 | delete[] m_indices;
|
---|
223 | m_values = newValues;
|
---|
224 | m_indices = newIndices;
|
---|
225 | m_allocatedSize = size;
|
---|
226 | }
|
---|
227 |
|
---|
228 | protected:
|
---|
229 | Scalar* m_values;
|
---|
230 | Index* m_indices;
|
---|
231 | size_t m_size;
|
---|
232 | size_t m_allocatedSize;
|
---|
233 |
|
---|
234 | };
|
---|
235 |
|
---|
236 | } // end namespace internal
|
---|
237 |
|
---|
238 | } // end namespace Eigen
|
---|
239 |
|
---|
240 | #endif // EIGEN_COMPRESSED_STORAGE_H
|
---|