1 // This file is part of Eigen, a lightweight C++ template library
2 // for linear algebra.
3 //
4 // Copyright (C) 2008-2009 Guillaume Saupin <guillaume.saupin@cea.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_SKYLINE_STORAGE_H
11 #define EIGEN_SKYLINE_STORAGE_H
12 
13 namespace Eigen {
14 
15 /** Stores a skyline set of values in three structures :
16  * The diagonal elements
17  * The upper elements
18  * The lower elements
19  *
20  */
21 template<typename Scalar>
22 class SkylineStorage {
23     typedef typename NumTraits<Scalar>::Real RealScalar;
24     typedef SparseIndex Index;
25 public:
26 
SkylineStorage()27     SkylineStorage()
28     : m_diag(0),
29     m_lower(0),
30     m_upper(0),
31     m_lowerProfile(0),
32     m_upperProfile(0),
33     m_diagSize(0),
34     m_upperSize(0),
35     m_lowerSize(0),
36     m_upperProfileSize(0),
37     m_lowerProfileSize(0),
38     m_allocatedSize(0) {
39     }
40 
SkylineStorage(const SkylineStorage & other)41     SkylineStorage(const SkylineStorage& other)
42     : m_diag(0),
43     m_lower(0),
44     m_upper(0),
45     m_lowerProfile(0),
46     m_upperProfile(0),
47     m_diagSize(0),
48     m_upperSize(0),
49     m_lowerSize(0),
50     m_upperProfileSize(0),
51     m_lowerProfileSize(0),
52     m_allocatedSize(0) {
53         *this = other;
54     }
55 
56     SkylineStorage & operator=(const SkylineStorage& other) {
57         resize(other.diagSize(), other.m_upperProfileSize, other.m_lowerProfileSize, other.upperSize(), other.lowerSize());
58         memcpy(m_diag, other.m_diag, m_diagSize * sizeof (Scalar));
59         memcpy(m_upper, other.m_upper, other.upperSize() * sizeof (Scalar));
60         memcpy(m_lower, other.m_lower, other.lowerSize() * sizeof (Scalar));
61         memcpy(m_upperProfile, other.m_upperProfile, m_upperProfileSize * sizeof (Index));
62         memcpy(m_lowerProfile, other.m_lowerProfile, m_lowerProfileSize * sizeof (Index));
63         return *this;
64     }
65 
swap(SkylineStorage & other)66     void swap(SkylineStorage& other) {
67         std::swap(m_diag, other.m_diag);
68         std::swap(m_upper, other.m_upper);
69         std::swap(m_lower, other.m_lower);
70         std::swap(m_upperProfile, other.m_upperProfile);
71         std::swap(m_lowerProfile, other.m_lowerProfile);
72         std::swap(m_diagSize, other.m_diagSize);
73         std::swap(m_upperSize, other.m_upperSize);
74         std::swap(m_lowerSize, other.m_lowerSize);
75         std::swap(m_allocatedSize, other.m_allocatedSize);
76     }
77 
~SkylineStorage()78     ~SkylineStorage() {
79         delete[] m_diag;
80         delete[] m_upper;
81         if (m_upper != m_lower)
82             delete[] m_lower;
83         delete[] m_upperProfile;
84         delete[] m_lowerProfile;
85     }
86 
reserve(Index size,Index upperProfileSize,Index lowerProfileSize,Index upperSize,Index lowerSize)87     void reserve(Index size, Index upperProfileSize, Index lowerProfileSize, Index upperSize, Index lowerSize) {
88         Index newAllocatedSize = size + upperSize + lowerSize;
89         if (newAllocatedSize > m_allocatedSize)
90             reallocate(size, upperProfileSize, lowerProfileSize, upperSize, lowerSize);
91     }
92 
squeeze()93     void squeeze() {
94         if (m_allocatedSize > m_diagSize + m_upperSize + m_lowerSize)
95             reallocate(m_diagSize, m_upperProfileSize, m_lowerProfileSize, m_upperSize, m_lowerSize);
96     }
97 
98     void resize(Index diagSize, Index upperProfileSize, Index lowerProfileSize, Index upperSize, Index lowerSize, float reserveSizeFactor = 0) {
99         if (m_allocatedSize < diagSize + upperSize + lowerSize)
100             reallocate(diagSize, upperProfileSize, lowerProfileSize, upperSize + Index(reserveSizeFactor * upperSize), lowerSize + Index(reserveSizeFactor * lowerSize));
101         m_diagSize = diagSize;
102         m_upperSize = upperSize;
103         m_lowerSize = lowerSize;
104         m_upperProfileSize = upperProfileSize;
105         m_lowerProfileSize = lowerProfileSize;
106     }
107 
diagSize()108     inline Index diagSize() const {
109         return m_diagSize;
110     }
111 
upperSize()112     inline Index upperSize() const {
113         return m_upperSize;
114     }
115 
lowerSize()116     inline Index lowerSize() const {
117         return m_lowerSize;
118     }
119 
upperProfileSize()120     inline Index upperProfileSize() const {
121         return m_upperProfileSize;
122     }
123 
lowerProfileSize()124     inline Index lowerProfileSize() const {
125         return m_lowerProfileSize;
126     }
127 
allocatedSize()128     inline Index allocatedSize() const {
129         return m_allocatedSize;
130     }
131 
clear()132     inline void clear() {
133         m_diagSize = 0;
134     }
135 
diag(Index i)136     inline Scalar& diag(Index i) {
137         return m_diag[i];
138     }
139 
diag(Index i)140     inline const Scalar& diag(Index i) const {
141         return m_diag[i];
142     }
143 
upper(Index i)144     inline Scalar& upper(Index i) {
145         return m_upper[i];
146     }
147 
upper(Index i)148     inline const Scalar& upper(Index i) const {
149         return m_upper[i];
150     }
151 
lower(Index i)152     inline Scalar& lower(Index i) {
153         return m_lower[i];
154     }
155 
lower(Index i)156     inline const Scalar& lower(Index i) const {
157         return m_lower[i];
158     }
159 
upperProfile(Index i)160     inline Index& upperProfile(Index i) {
161         return m_upperProfile[i];
162     }
163 
upperProfile(Index i)164     inline const Index& upperProfile(Index i) const {
165         return m_upperProfile[i];
166     }
167 
lowerProfile(Index i)168     inline Index& lowerProfile(Index i) {
169         return m_lowerProfile[i];
170     }
171 
lowerProfile(Index i)172     inline const Index& lowerProfile(Index i) const {
173         return m_lowerProfile[i];
174     }
175 
Map(Index * upperProfile,Index * lowerProfile,Scalar * diag,Scalar * upper,Scalar * lower,Index size,Index upperSize,Index lowerSize)176     static SkylineStorage Map(Index* upperProfile, Index* lowerProfile, Scalar* diag, Scalar* upper, Scalar* lower, Index size, Index upperSize, Index lowerSize) {
177         SkylineStorage res;
178         res.m_upperProfile = upperProfile;
179         res.m_lowerProfile = lowerProfile;
180         res.m_diag = diag;
181         res.m_upper = upper;
182         res.m_lower = lower;
183         res.m_allocatedSize = res.m_diagSize = size;
184         res.m_upperSize = upperSize;
185         res.m_lowerSize = lowerSize;
186         return res;
187     }
188 
reset()189     inline void reset() {
190         memset(m_diag, 0, m_diagSize * sizeof (Scalar));
191         memset(m_upper, 0, m_upperSize * sizeof (Scalar));
192         memset(m_lower, 0, m_lowerSize * sizeof (Scalar));
193         memset(m_upperProfile, 0, m_diagSize * sizeof (Index));
194         memset(m_lowerProfile, 0, m_diagSize * sizeof (Index));
195     }
196 
197     void prune(Scalar reference, RealScalar epsilon = dummy_precision<RealScalar>()) {
198         //TODO
199     }
200 
201 protected:
202 
reallocate(Index diagSize,Index upperProfileSize,Index lowerProfileSize,Index upperSize,Index lowerSize)203     inline void reallocate(Index diagSize, Index upperProfileSize, Index lowerProfileSize, Index upperSize, Index lowerSize) {
204 
205         Scalar* diag = new Scalar[diagSize];
206         Scalar* upper = new Scalar[upperSize];
207         Scalar* lower = new Scalar[lowerSize];
208         Index* upperProfile = new Index[upperProfileSize];
209         Index* lowerProfile = new Index[lowerProfileSize];
210 
211         Index copyDiagSize = (std::min)(diagSize, m_diagSize);
212         Index copyUpperSize = (std::min)(upperSize, m_upperSize);
213         Index copyLowerSize = (std::min)(lowerSize, m_lowerSize);
214         Index copyUpperProfileSize = (std::min)(upperProfileSize, m_upperProfileSize);
215         Index copyLowerProfileSize = (std::min)(lowerProfileSize, m_lowerProfileSize);
216 
217         // copy
218         memcpy(diag, m_diag, copyDiagSize * sizeof (Scalar));
219         memcpy(upper, m_upper, copyUpperSize * sizeof (Scalar));
220         memcpy(lower, m_lower, copyLowerSize * sizeof (Scalar));
221         memcpy(upperProfile, m_upperProfile, copyUpperProfileSize * sizeof (Index));
222         memcpy(lowerProfile, m_lowerProfile, copyLowerProfileSize * sizeof (Index));
223 
224 
225 
226         // delete old stuff
227         delete[] m_diag;
228         delete[] m_upper;
229         delete[] m_lower;
230         delete[] m_upperProfile;
231         delete[] m_lowerProfile;
232         m_diag = diag;
233         m_upper = upper;
234         m_lower = lower;
235         m_upperProfile = upperProfile;
236         m_lowerProfile = lowerProfile;
237         m_allocatedSize = diagSize + upperSize + lowerSize;
238         m_upperSize = upperSize;
239         m_lowerSize = lowerSize;
240     }
241 
242 public:
243     Scalar* m_diag;
244     Scalar* m_upper;
245     Scalar* m_lower;
246     Index* m_upperProfile;
247     Index* m_lowerProfile;
248     Index m_diagSize;
249     Index m_upperSize;
250     Index m_lowerSize;
251     Index m_upperProfileSize;
252     Index m_lowerProfileSize;
253     Index m_allocatedSize;
254 
255 };
256 
257 } // end namespace Eigen
258 
259 #endif // EIGEN_COMPRESSED_STORAGE_H
260