a0f614b5cbd9d3589ff5e8416d625f5c9dfe5acc
[blender.git] / intern / smoke / intern / tnt / tnt_vec.h
1 /*
2 *
3 * Template Numerical Toolkit (TNT)
4 *
5 * Mathematical and Computational Sciences Division
6 * National Institute of Technology,
7 * Gaithersburg, MD USA
8 *
9 *
10 * This software was developed at the National Institute of Standards and
11 * Technology (NIST) by employees of the Federal Government in the course
12 * of their official duties. Pursuant to title 17 Section 105 of the
13 * United States Code, this software is not subject to copyright protection
14 * and is in the public domain. NIST assumes no responsibility whatsoever for
15 * its use by other parties, and makes no guarantees, expressed or implied,
16 * about its quality, reliability, or any other characteristic.
17 *
18 */
19
20
21
22 #ifndef TNT_VEC_H
23 #define TNT_VEC_H
24
25 #include "tnt_subscript.h"
26 #include <cstdlib>
27 #include <cassert>
28 #include <iostream>
29 #include <sstream>
30
31 namespace TNT
32 {
33
34 /**
35  <b>[Deprecatred]</b>  Value-based vector class from pre-1.0
36         TNT version.  Kept here for backward compatiblity, but should
37         use the newer TNT::Array1D classes instead.
38
39 */
40
41 template <class T>
42 class Vector 
43 {
44
45
46   public:
47
48     typedef Subscript   size_type;
49     typedef         T   value_type;
50     typedef         T   element_type;
51     typedef         T*  pointer;
52     typedef         T*  iterator;
53     typedef         T&  reference;
54     typedef const   T*  const_iterator;
55     typedef const   T&  const_reference;
56
57     Subscript lbound() const { return 1;}
58  
59   protected:
60     T* v_;                  
61     T* vm1_;        // pointer adjustment for optimzied 1-offset indexing
62     Subscript n_;
63
64     // internal helper function to create the array
65     // of row pointers
66
67     void initialize(Subscript N)
68     {
69         // adjust pointers so that they are 1-offset:
70         // v_[] is the internal contiguous array, it is still 0-offset
71         //
72         assert(v_ == NULL);
73         v_ = new T[N];
74         assert(v_  != NULL);
75         vm1_ = v_-1;
76         n_ = N;
77     }
78    
79     void copy(const T*  v)
80     {
81         Subscript N = n_;
82         Subscript i;
83
84 #ifdef TNT_UNROLL_LOOPS
85         Subscript Nmod4 = N & 3;
86         Subscript N4 = N - Nmod4;
87
88         for (i=0; i<N4; i+=4)
89         {
90             v_[i] = v[i];
91             v_[i+1] = v[i+1];
92             v_[i+2] = v[i+2];
93             v_[i+3] = v[i+3];
94         }
95
96         for (i=N4; i< N; i++)
97             v_[i] = v[i];
98 #else
99
100         for (i=0; i< N; i++)
101             v_[i] = v[i];
102 #endif      
103     }
104
105     void set(const T& val)
106     {
107         Subscript N = n_;
108         Subscript i;
109
110 #ifdef TNT_UNROLL_LOOPS
111         Subscript Nmod4 = N & 3;
112         Subscript N4 = N - Nmod4;
113
114         for (i=0; i<N4; i+=4)
115         {
116             v_[i] = val;
117             v_[i+1] = val;
118             v_[i+2] = val;
119             v_[i+3] = val; 
120         }
121
122         for (i=N4; i< N; i++)
123             v_[i] = val;
124 #else
125
126         for (i=0; i< N; i++)
127             v_[i] = val;
128         
129 #endif      
130     }
131     
132
133
134     void destroy()
135     {     
136         /* do nothing, if no memory has been previously allocated */
137         if (v_ == NULL) return ;
138
139         /* if we are here, then matrix was previously allocated */
140         delete [] (v_);     
141
142         v_ = NULL;
143         vm1_ = NULL;
144     }
145
146
147   public:
148
149     // access
150
151     iterator begin() { return v_;}
152     iterator end()   { return v_ + n_; }
153     iterator begin() const { return v_;}
154     iterator end() const  { return v_ + n_; }
155
156     // destructor
157
158     ~Vector() 
159     {
160         destroy();
161     }
162
163     // constructors
164
165     Vector() : v_(0), vm1_(0), n_(0)  {};
166
167     Vector(const Vector<T> &A) : v_(0), vm1_(0), n_(0)
168     {
169         initialize(A.n_);
170         copy(A.v_);
171     }
172
173     Vector(Subscript N, const T& value = T()) :  v_(0), vm1_(0), n_(0)
174     {
175         initialize(N);
176         set(value);
177     }
178
179     Vector(Subscript N, const T* v) :  v_(0), vm1_(0), n_(0)
180     {
181         initialize(N);
182         copy(v);
183     }
184
185     Vector(Subscript N, char *s) :  v_(0), vm1_(0), n_(0)
186     {
187         initialize(N);
188         std::istringstream ins(s);
189
190         Subscript i;
191
192         for (i=0; i<N; i++)
193                 ins >> v_[i];
194     }
195
196
197     // methods
198     // 
199     Vector<T>& newsize(Subscript N)
200     {
201         if (n_ == N) return *this;
202
203         destroy();
204         initialize(N);
205
206         return *this;
207     }
208
209
210     // assignments
211     //
212     Vector<T>& operator=(const Vector<T> &A)
213     {
214         if (v_ == A.v_)
215             return *this;
216
217         if (n_ == A.n_)         // no need to re-alloc
218             copy(A.v_);
219
220         else
221         {
222             destroy();
223             initialize(A.n_);
224             copy(A.v_);
225         }
226
227         return *this;
228     }
229         
230     Vector<T>& operator=(const T& scalar)
231     { 
232         set(scalar);  
233         return *this;
234     }
235
236     inline Subscript dim() const 
237     {
238         return  n_; 
239     }
240
241     inline Subscript size() const 
242     {
243         return  n_; 
244     }
245
246
247     inline reference operator()(Subscript i)
248     { 
249 #ifdef TNT_BOUNDS_CHECK
250         assert(1<=i);
251         assert(i <= n_) ;
252 #endif
253         return vm1_[i]; 
254     }
255
256     inline const_reference operator() (Subscript i) const
257     {
258 #ifdef TNT_BOUNDS_CHECK
259         assert(1<=i);
260         assert(i <= n_) ;
261 #endif
262         return vm1_[i]; 
263     }
264
265     inline reference operator[](Subscript i)
266     { 
267 #ifdef TNT_BOUNDS_CHECK
268         assert(0<=i);
269         assert(i < n_) ;
270 #endif
271         return v_[i]; 
272     }
273
274     inline const_reference operator[](Subscript i) const
275     {
276 #ifdef TNT_BOUNDS_CHECK
277         assert(0<=i);
278
279
280
281
282
283
284         assert(i < n_) ;
285 #endif
286         return v_[i]; 
287     }
288
289
290
291 };
292
293
294 /* ***************************  I/O  ********************************/
295
296 template <class T>
297 std::ostream& operator<<(std::ostream &s, const Vector<T> &A)
298 {
299     Subscript N=A.dim();
300
301     s <<  N << "\n";
302
303     for (Subscript i=0; i<N; i++)
304         s   << A[i] << " " << "\n";
305     s << "\n";
306
307     return s;
308 }
309
310 template <class T>
311 std::istream & operator>>(std::istream &s, Vector<T> &A)
312 {
313
314     Subscript N;
315
316     s >> N;
317
318     if ( !(N == A.size() ))
319     {
320         A.newsize(N);
321     }
322
323
324     for (Subscript i=0; i<N; i++)
325             s >>  A[i];
326
327
328     return s;
329 }
330
331 // *******************[ basic matrix algorithms ]***************************
332
333
334 template <class T>
335 Vector<T> operator+(const Vector<T> &A, 
336     const Vector<T> &B)
337 {
338     Subscript N = A.dim();
339
340     assert(N==B.dim());
341
342     Vector<T> tmp(N);
343     Subscript i;
344
345     for (i=0; i<N; i++)
346             tmp[i] = A[i] + B[i];
347
348     return tmp;
349 }
350
351 template <class T>
352 Vector<T> operator-(const Vector<T> &A, 
353     const Vector<T> &B)
354 {
355     Subscript N = A.dim();
356
357     assert(N==B.dim());
358
359     Vector<T> tmp(N);
360     Subscript i;
361
362     for (i=0; i<N; i++)
363             tmp[i] = A[i] - B[i];
364
365     return tmp;
366 }
367
368 template <class T>
369 Vector<T> operator*(const Vector<T> &A, 
370     const Vector<T> &B)
371 {
372     Subscript N = A.dim();
373
374     assert(N==B.dim());
375
376     Vector<T> tmp(N);
377     Subscript i;
378
379     for (i=0; i<N; i++)
380             tmp[i] = A[i] * B[i];
381
382     return tmp;
383 }
384
385
386 template <class T>
387 T dot_prod(const Vector<T> &A, const Vector<T> &B)
388 {
389     Subscript N = A.dim();
390     assert(N == B.dim());
391
392     Subscript i;
393     T sum = 0;
394
395     for (i=0; i<N; i++)
396         sum += A[i] * B[i];
397
398     return sum;
399 }
400
401 }   /* namespace TNT */
402
403 #endif
404 // TNT_VEC_H