~ubuntu-branches/ubuntu/raring/blitz++/raring

« back to all changes in this revision

Viewing changes to blitz/array/methods.cc

  • Committer: Bazaar Package Importer
  • Author(s): Konstantinos Margaritis
  • Date: 2005-02-28 20:25:01 UTC
  • mfrom: (2.1.2 hoary)
  • Revision ID: james.westby@ubuntu.com-20050228202501-3i4f2sknnprsqfhz
Tags: 1:0.8-4
Added missing build-depends (Closes: #297323)

Show diffs side-by-side

added added

removed removed

Lines of Context:
1
 
/*
2
 
 * $Id: methods.cc,v 1.1.1.1 2000/06/19 12:26:13 tveldhui Exp $
3
 
 *
4
 
 * $Log: methods.cc,v $
5
 
 * Revision 1.1.1.1  2000/06/19 12:26:13  tveldhui
6
 
 * Imported sources
7
 
 *
8
 
 * Revision 1.4  1998/03/14 00:04:47  tveldhui
9
 
 * 0.2-alpha-05
10
 
 *
11
 
 * Revision 1.3  1997/08/18 19:13:08  tveldhui
12
 
 * Just prior to implementing fastRead() optimization for array
13
 
 * expression evaluation.
14
 
 *
15
 
 * Revision 1.2  1997/08/15 21:14:10  tveldhui
16
 
 * Just prior to loop-collapse change
17
 
 *
18
 
 */
19
 
 
20
1
#ifndef BZ_ARRAYMETHODS_CC
21
2
#define BZ_ARRAYMETHODS_CC
22
3
 
26
7
 
27
8
BZ_NAMESPACE(blitz)
28
9
 
29
 
template<class P_numtype, int N_rank> template<class T_expr>
 
10
template<typename P_numtype, int N_rank> template<typename T_expr>
30
11
Array<P_numtype,N_rank>::Array(_bz_ArrayExpr<T_expr> expr)
31
12
{
32
13
    // Determine extent of the array expression
33
14
 
34
15
    TinyVector<int,N_rank> lbound, extent, ordering;
35
16
    TinyVector<bool,N_rank> ascendingFlag;
 
17
    TinyVector<bool,N_rank> in_ordering;
 
18
    in_ordering = false;
36
19
 
 
20
    int j = 0;
37
21
    for (int i=0; i < N_rank; ++i)
38
22
    {
39
23
        lbound(i) = expr.lbound(i);
40
24
        int ubound = expr.ubound(i);
41
25
        extent(i) = ubound - lbound(i) + 1;
42
 
        ordering(i) = expr.ordering(i);
 
26
        int orderingj = expr.ordering(i);
 
27
        if (orderingj != INT_MIN && orderingj < N_rank &&
 
28
            !in_ordering( orderingj )) { // unique value in ordering array
 
29
            in_ordering( orderingj ) = true;
 
30
            ordering(j++) = orderingj;
 
31
        }
43
32
        int ascending = expr.ascending(i);
44
33
        ascendingFlag(i) = (ascending == 1);
45
34
 
46
35
#ifdef BZ_DEBUG
47
36
        if ((lbound(i) == INT_MIN) || (ubound == INT_MAX) 
48
 
          || (ordering(i) == INT_MAX) || (ascending == INT_MAX))
 
37
          || (ordering(i) == INT_MIN) || (ascending == INT_MIN))
49
38
        {
50
39
          BZPRECHECK(0,
51
40
           "Attempted to construct an array from an expression " << endl
57
46
#endif
58
47
    }
59
48
 
60
 
    Array<P_numtype,N_rank> A(lbound,extent,
 
49
    // It is possible that ordering is not a permutation of 0,...,N_rank-1.
 
50
    // In that case j will be less than N_rank. We fill in ordering with the
 
51
    // usused values in decreasing order.
 
52
    for (int i = N_rank-1; j < N_rank; ++j) {
 
53
        while (in_ordering(i))
 
54
          --i;
 
55
        ordering(j) = i--;
 
56
    }
 
57
 
 
58
    Array<T_numtype,N_rank> A(lbound,extent,
61
59
        GeneralArrayStorage<N_rank>(ordering,ascendingFlag));
62
60
    A = expr;
63
61
    reference(A);
64
62
}
65
63
 
66
 
template<class T_numtype, int N_rank>
67
 
Array<T_numtype,N_rank>::Array(const TinyVector<int, N_rank>& lbounds,
 
64
template<typename P_numtype, int N_rank>
 
65
Array<P_numtype,N_rank>::Array(const TinyVector<int, N_rank>& lbounds,
68
66
    const TinyVector<int, N_rank>& extent,
69
67
    const GeneralArrayStorage<N_rank>& storage)
70
68
    : storage_(storage)
81
79
 * of the array (length_[]) and computes the stride vector
82
80
 * (stride_[]) and the zero offset (see explanation in array.h).
83
81
 */
84
 
template<class P_numtype, int N_rank>
 
82
template<typename P_numtype, int N_rank>
85
83
_bz_inline2 void Array<P_numtype, N_rank>::computeStrides()
86
84
{
87
85
    if (N_rank > 1)
90
88
 
91
89
      // This flag simplifies the code in the loop, encouraging
92
90
      // compile-time computation of strides through constant folding.
93
 
      _bz_bool allAscending = storage_.allRanksStoredAscending();
 
91
      bool allAscending = storage_.allRanksStoredAscending();
94
92
 
95
93
      // BZ_OLD_FOR_SCOPING
96
94
      int n;
127
125
    calculateZeroOffset();
128
126
}
129
127
 
130
 
template<class T_numtype, int N_rank>
131
 
void Array<T_numtype, N_rank>::calculateZeroOffset()
 
128
template<typename P_numtype, int N_rank>
 
129
void Array<P_numtype, N_rank>::calculateZeroOffset()
132
130
{
133
131
    // Calculate the offset of (0,0,...,0)
134
132
    zeroOffset_ = 0;
144
142
    }
145
143
}
146
144
 
147
 
template<class T_numtype, int N_rank>
148
 
_bz_bool Array<T_numtype, N_rank>::isStorageContiguous() const
 
145
template<typename P_numtype, int N_rank>
 
146
bool Array<P_numtype, N_rank>::isStorageContiguous() const
149
147
{
150
148
    // The storage is contiguous if for the set
151
149
    // { | stride[i] * extent[i] | }, i = 0..N_rank-1,
156
154
    // to imagine this being a serious problem.
157
155
 
158
156
    int numStridesMissing = 0;
159
 
    bool haveUnitStride = _bz_false;
 
157
    bool haveUnitStride = false;
160
158
 
161
159
    for (int i=0; i < N_rank; ++i)
162
160
    {
163
161
        int stride = BZ_MATHFN_SCOPE(abs)(stride_[i]);
164
162
        if (stride == 1)
165
 
            haveUnitStride = _bz_true;
 
163
            haveUnitStride = true;
166
164
 
167
165
        int vi = stride * length_[i];
168
166
 
175
173
        {
176
174
            ++numStridesMissing;
177
175
            if (numStridesMissing == 2)
178
 
                return _bz_false;
 
176
                return false;
179
177
        }
180
178
    }
181
179
 
182
180
    return haveUnitStride;
183
181
}
184
182
 
185
 
template<class P_numtype, int N_rank>
 
183
template<typename P_numtype, int N_rank>
186
184
void Array<P_numtype, N_rank>::dumpStructureInformation(ostream& os) const
187
185
{
188
186
    os << "Dump of Array<" << BZ_DEBUG_TEMPLATE_AS_STRING_LITERAL(P_numtype) 
200
198
/*
201
199
 * Make this array a view of another array's data.
202
200
 */
203
 
template<class P_numtype, int N_rank>
 
201
template<typename P_numtype, int N_rank>
204
202
void Array<P_numtype, N_rank>::reference(const Array<P_numtype, N_rank>& array)
205
203
{
206
204
    storage_ = array.storage_;
208
206
    stride_ = array.stride_;
209
207
    zeroOffset_ = array.zeroOffset_;
210
208
 
211
 
    MemoryBlockReference<P_numtype>::changeBlock(array.noConst(),
212
 
        array.zeroOffset_);
 
209
    MemoryBlockReference<P_numtype>::changeBlock(array.noConst());
 
210
}
213
211
 
214
 
    data_ = const_cast<P_numtype*>(array.data_);
 
212
/*
 
213
 * Modify the Array storage.  Array must be unallocated.
 
214
 */
 
215
template<typename P_numtype, int N_rank>
 
216
void Array<P_numtype, N_rank>::setStorage(GeneralArrayStorage<N_rank> x)
 
217
{
 
218
#ifdef BZ_DEBUG
 
219
    if (size() != 0) {
 
220
        BZPRECHECK(0,"Cannot modify storage format of an Array that has already been allocated!" << endl);
 
221
        return;
 
222
    }
 
223
#endif
 
224
    storage_ = x;
 
225
    return;
215
226
}
216
227
 
217
228
/*
218
229
 * This method is called to allocate memory for a new array.  
219
230
 */
220
 
template<class P_numtype, int N_rank>
 
231
template<typename P_numtype, int N_rank>
221
232
_bz_inline2 void Array<P_numtype, N_rank>::setupStorage(int lastRankInitialized)
222
233
{
223
234
    TAU_TYPE_STRING(p1, "Array<T,N>::setupStorage() [T="
240
251
    computeStrides();
241
252
 
242
253
    // Allocate a block of memory
243
 
    MemoryBlockReference<P_numtype>::newBlock(numElements());
 
254
    int numElem = numElements();
 
255
    if (numElem==0)
 
256
        MemoryBlockReference<P_numtype>::changeToNullBlock();
 
257
    else
 
258
        MemoryBlockReference<P_numtype>::newBlock(numElem);
244
259
 
245
260
    // Adjust the base of the array to account for non-zero base
246
261
    // indices and reversals
247
262
    data_ += zeroOffset_;
248
263
}
249
264
 
250
 
template<class T_numtype, int N_rank>
251
 
Array<T_numtype, N_rank> Array<T_numtype, N_rank>::copy() const
 
265
template<typename P_numtype, int N_rank>
 
266
Array<P_numtype, N_rank> Array<P_numtype, N_rank>::copy() const
252
267
{
253
268
    if (numElements())
254
269
    {
262
277
    }
263
278
}
264
279
 
265
 
template<class T_numtype, int N_rank>
266
 
void Array<T_numtype, N_rank>::makeUnique()
 
280
template<typename P_numtype, int N_rank>
 
281
void Array<P_numtype, N_rank>::makeUnique()
267
282
{
268
283
    if (numReferences() > 1)
269
284
    {
272
287
    }
273
288
}
274
289
 
275
 
template<class T_numtype, int N_rank>
276
 
Array<T_numtype, N_rank> Array<T_numtype, N_rank>::transpose(int r0, int r1, 
 
290
template<typename P_numtype, int N_rank>
 
291
Array<P_numtype, N_rank> Array<P_numtype, N_rank>::transpose(int r0, int r1, 
277
292
    int r2, int r3, int r4, int r5, int r6, int r7, int r8, int r9, int r10)
278
293
{
279
294
    T_array B(*this);
281
296
    return B;
282
297
}
283
298
 
284
 
template<class T_numtype, int N_rank>
285
 
void Array<T_numtype, N_rank>::transposeSelf(int r0, int r1, int r2, int r3,
 
299
template<typename P_numtype, int N_rank>
 
300
void Array<P_numtype, N_rank>::transposeSelf(int r0, int r1, int r2, int r3,
286
301
    int r4, int r5, int r6, int r7, int r8, int r9, int r10)
287
302
{
288
303
    BZPRECHECK(r0+r1+r2+r3+r4+r5+r6+r7+r8+r9+r10 == N_rank * (N_rank-1) / 2,
307
322
    doTranspose(10, r10, x);
308
323
}
309
324
 
310
 
template<class T_numtype, int N_rank>
311
 
void Array<T_numtype, N_rank>::doTranspose(int destRank, int sourceRank,
 
325
template<typename P_numtype, int N_rank>
 
326
void Array<P_numtype, N_rank>::doTranspose(int destRank, int sourceRank,
312
327
    Array<T_numtype, N_rank>& array)
313
328
{
314
329
    // BZ_NEEDS_WORK: precondition check
335
350
    storage_.setOrdering(i, destRank);
336
351
}
337
352
 
338
 
template<class T_numtype, int N_rank>
339
 
void Array<T_numtype, N_rank>::reverseSelf(int rank)
 
353
template<typename P_numtype, int N_rank>
 
354
void Array<P_numtype, N_rank>::reverseSelf(int rank)
340
355
{
341
356
    BZPRECONDITION(rank < N_rank);
342
357
 
348
363
    stride_[rank] *= -1;
349
364
}
350
365
 
351
 
template<class T_numtype, int N_rank>
352
 
Array<T_numtype, N_rank> Array<T_numtype,N_rank>::reverse(int rank)
 
366
template<typename P_numtype, int N_rank>
 
367
Array<P_numtype, N_rank> Array<P_numtype,N_rank>::reverse(int rank)
353
368
{
354
369
    T_array B(*this);
355
370
    B.reverseSelf(rank);
356
371
    return B;
357
372
}
358
373
 
359
 
template<class T_numtype, int N_rank> template<class T_numtype2>
360
 
Array<T_numtype2,N_rank> Array<T_numtype,N_rank>::extractComponent(T_numtype2, 
 
374
template<typename P_numtype, int N_rank> template<typename P_numtype2>
 
375
Array<P_numtype2,N_rank> Array<P_numtype,N_rank>::extractComponent(P_numtype2, 
361
376
    int componentNumber, int numComponents) const
362
377
{
363
378
    BZPRECONDITION((componentNumber >= 0) 
364
379
        && (componentNumber < numComponents));
365
380
 
366
381
    TinyVector<int,N_rank> stride2;
367
 
    stride2 = stride_ * numComponents;
368
 
    const T_numtype2* dataFirst2 = 
369
 
        ((const T_numtype2*)dataFirst()) + componentNumber;
370
 
    return Array<T_numtype2,N_rank>(const_cast<T_numtype2*>(dataFirst2), 
 
382
    for (int i=0; i < N_rank; ++i)
 
383
      stride2(i) = stride_(i) * numComponents;
 
384
    const P_numtype2* dataFirst2 = 
 
385
        ((const P_numtype2*)dataFirst()) + componentNumber;
 
386
    return Array<P_numtype2,N_rank>(const_cast<P_numtype2*>(dataFirst2), 
371
387
        length_, stride2, storage_);
372
388
}
373
389
 
377
393
 * of the current array, leaving the current array unmodified.
378
394
 * (Contributed by Derrick Bass)
379
395
 */
380
 
template<class P_numtype, int N_rank>
 
396
template<typename P_numtype, int N_rank>
381
397
_bz_inline2 void Array<P_numtype, N_rank>::reindexSelf(const 
382
398
    TinyVector<int, N_rank>& newBase) 
383
399
{
384
 
    data_ += dot(base() - newBase, stride_);
 
400
    int delta = 0;
 
401
    for (int i=0; i < N_rank; ++i)
 
402
      delta += (base(i) - newBase(i)) * stride_(i);
 
403
 
 
404
    data_ += delta;
 
405
 
 
406
    // WAS: dot(base() - newBase, stride_);
 
407
 
385
408
    storage_.setBase(newBase);
386
409
    calculateZeroOffset();
387
410
}
388
411
 
389
 
template<class P_numtype, int N_rank>
 
412
template<typename P_numtype, int N_rank>
390
413
_bz_inline2 Array<P_numtype, N_rank> 
391
414
Array<P_numtype, N_rank>::reindex(const TinyVector<int, N_rank>& newBase) 
392
415
{