← Back to branch summary

~chaffra/+junk/trilinos

« back to all changes in this revision

Viewing changes to packages/phalanx/src/multidimensional_array/version_2/Array.hpp

  • Committer: Bazaar Package Importer
  • Author(s): Christophe Prud'homme, Christophe Prud'homme, Johannes Ring
  • Date: 2009-12-13 12:53:22 UTC
  • mfrom: (5.1.2 sid)
  • Revision ID: james.westby@ubuntu.com-20091213125322-in0nrdjc55deqsw9
Tags: 10.0.3.dfsg-1
[Christophe Prud'homme]
* New upstream release

[Johannes Ring]
* debian/patches/libname.patch: Add prefix 'libtrilinos_' to all
  libraries. 
* debian/patches/soname.patch: Add soversion to libraries.
* debian/watch: Update download URL.
* debian/control:
  - Remove python-numeric from Build-Depends (virtual package).
  - Remove automake and autotools from Build-Depends and add cmake to
    reflect switch to CMake.
  - Add python-support to Build-Depends.
* debian/rules: 
  - Cleanup and updates for switch to CMake.

Show diffs side-by-side

added added

removed removed

Lines of Context:
packages/phalanx/src/multidimensional_array/version_2/Array.hpp
 
1
/*------------------------------------------------------------------------*/
 
2
/*      phdMesh : Parallel Heterogneous Dynamic unstructured Mesh         */
 
3
/*                Copyright (2008) Sandia Corporation                     */
 
4
/*                                                                        */
 
5
/*  Under terms of Contract DE-AC04-94AL85000, there is a non-exclusive   */
 
6
/*  license for use of this work by or on behalf of the U.S. Government.  */
 
7
/*                                                                        */
 
8
/*  This library is free software; you can redistribute it and/or modify  */
 
9
/*  it under the terms of the GNU Lesser General Public License as        */
 
10
/*  published by the Free Software Foundation; either version 2.1 of the  */
 
11
/*  License, or (at your option) any later version.                       */
 
12
/*                                                                        */
 
13
/*  This library is distributed in the hope that it will be useful,       */
 
14
/*  but WITHOUT ANY WARRANTY; without even the implied warranty of        */
 
15
/*  MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU     */
 
16
/*  Lesser General Public License for more details.                       */
 
17
/*                                                                        */
 
18
/*  You should have received a copy of the GNU Lesser General Public      */
 
19
/*  License along with this library; if not, write to the Free Software   */
 
20
/*  Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307   */
 
21
/*  USA                                                                   */
 
22
/*------------------------------------------------------------------------*/
 
23
/**
 
24
 * @author H. Carter Edwards  <hcedwar@sandia.gov>
 
25
 * @date   June 2008
 
26
 */
 
27
 
 
28
#ifndef util_Array_hpp
 
29
#define util_Array_hpp
 
30
 
 
31
//----------------------------------------------------------------------
 
32
 
 
33
namespace phdmesh {
 
34
 
 
35
/** \class ArrayNatural
 
36
 *  \brief A multidimensional array index mapping into a contiguous storage.
 
37
 *
 
38
 *  The Scalar template parameter defines the member type of the array.
 
39
 *  The Tag? template parameters define the ordinates.
 
40
 *
 
41
 *  Template parameter for a multi-index mapping class.
 
42
 *  This class is required to provide:
 
43
 *  -  enum { Rank = ... };
 
44
 *     which is the rank of the dimension and
 
45
 *  -  unsigned operator()( unsigned i1 , ... ) const ;
 
46
 *     which is the multi-index to offset map.
 
47
 *
 
48
 */
 
49
template< typename Scalar ,
 
50
          class Tag1 = void , class Tag2 = void ,
 
51
          class Tag3 = void , class Tag4 = void ,
 
52
          class Tag5 = void , class Tag6 = void ,
 
53
          class Tag7 = void , class Tag8 = void >
 
54
class ArrayNatural ;
 
55
 
 
56
template< typename Scalar ,
 
57
          class Tag1 = void , class Tag2 = void ,
 
58
          class Tag3 = void , class Tag4 = void ,
 
59
          class Tag5 = void , class Tag6 = void ,
 
60
          class Tag7 = void , class Tag8 = void >
 
61
class ArrayFortran ;
 
62
 
 
63
template< class Tag1 = void , class Tag2 = void ,
 
64
          class Tag3 = void , class Tag4 = void ,
 
65
          class Tag5 = void , class Tag6 = void ,
 
66
          class Tag7 = void , class Tag8 = void >
 
67
class ArrayDimension ;
 
68
 
 
69
class ArrayDimTag ;
 
70
 
 
71
}
 
72
 
 
73
//----------------------------------------------------------------------
 
74
 
 
75
#include <ArrayHelper.hpp>
 
76
 
 
77
//----------------------------------------------------------------------
 
78
 
 
79
namespace phdmesh {
 
80
 
 
81
struct ArrayDimTag {
 
82
  virtual const char * name() const = 0 ;
 
83
 
 
84
  virtual std::string to_string( size_t , int ) const ;
 
85
  virtual int          to_index( size_t , const std::string & ) const ;
 
86
 
 
87
  // A derived type must provide the following static method:
 
88
  //
 
89
  // static const Derived_ArrayDimTag_Type & descriptor();
 
90
 
 
91
protected:
 
92
  virtual ~ArrayDimTag();
 
93
  ArrayDimTag() {}
 
94
private:
 
95
  ArrayDimTag( const ArrayDimTag & );
 
96
  ArrayDimTag & operator = ( const ArrayDimTag & );
 
97
};
 
98
 
 
99
//----------------------------------------------------------------------
 
100
 
 
101
/** \cond */
 
102
 
 
103
template<>
 
104
class ArrayDimension<void,void,void,void,void,void,void,void> {
 
105
public:
 
106
  size_t              rank ;
 
107
  size_t              dimension[8];
 
108
  const ArrayDimTag * tags[8];
 
109
};
 
110
 
 
111
template< class Tag1 >
 
112
class ArrayDimension<Tag1,void,void,void,void,void,void,void> {
 
113
public:
 
114
  enum { Rank = 1 };
 
115
  size_t dimension[ Rank ];
 
116
};
 
117
 
 
118
template< class Tag1 , class Tag2 >
 
119
class ArrayDimension<Tag1,Tag2,void,void,void,void,void,void> {
 
120
public:
 
121
  enum { Rank = 2 };
 
122
  size_t dimension[ Rank ];
 
123
};
 
124
 
 
125
template< class Tag1 , class Tag2 , class Tag3 >
 
126
class ArrayDimension<Tag1,Tag2,Tag3,void,void,void,void,void> {
 
127
public:
 
128
  enum { Rank = 3 };
 
129
  size_t dimension[ Rank ];
 
130
};
 
131
 
 
132
template< class Tag1 , class Tag2 , class Tag3 , class Tag4 >
 
133
class ArrayDimension<Tag1,Tag2,Tag3,Tag4,void,void,void,void> {
 
134
public:
 
135
  enum { Rank = 4 };
 
136
  size_t dimension[ Rank ];
 
137
};
 
138
 
 
139
template< class Tag1 , class Tag2 , class Tag3 , class Tag4 ,
 
140
          class Tag5 >
 
141
class ArrayDimension<Tag1,Tag2,Tag3,Tag4,Tag5,void,void,void> {
 
142
public:
 
143
  enum { Rank = 5 };
 
144
  size_t dimension[ Rank ];
 
145
};
 
146
 
 
147
template< class Tag1 , class Tag2 , class Tag3 , class Tag4 ,
 
148
          class Tag5 , class Tag6 >
 
149
class ArrayDimension<Tag1,Tag2,Tag3,Tag4,Tag5,Tag6,void,void> {
 
150
public:
 
151
  enum { Rank = 6 };
 
152
  size_t dimension[ Rank ];
 
153
};
 
154
 
 
155
template< class Tag1 , class Tag2 , class Tag3 , class Tag4 ,
 
156
          class Tag5 , class Tag6 , class Tag7 >
 
157
class ArrayDimension<Tag1,Tag2,Tag3,Tag4,Tag5,Tag6,Tag7,void> {
 
158
public:
 
159
  enum { Rank = 7 };
 
160
  size_t dimension[ Rank ];
 
161
};
 
162
 
 
163
template< class Tag1 , class Tag2 , class Tag3 , class Tag4 ,
 
164
          class Tag5 , class Tag6 , class Tag7 , class Tag8 >
 
165
class ArrayDimension {
 
166
public:
 
167
  enum { Rank = 8 };
 
168
  size_t dimension[ Rank ];
 
169
};
 
170
 
 
171
//----------------------------------------------------------------------
 
172
 
 
173
template< typename Scalar ,
 
174
          class Tag1 , class Tag2 , class Tag3 , class Tag4 ,
 
175
          class Tag5 , class Tag6 , class Tag7 , class Tag8 >
 
176
class ArrayNatural
 
177
{
 
178
public:
 
179
  // Required by all arrays:
 
180
  typedef Scalar  value_type ;
 
181
  typedef int     index_type ;
 
182
  typedef size_t  size_type ;
 
183
  typedef const ArrayDimTag * tag_type ;
 
184
 
 
185
  // Required by compile-time knowledge arrays:
 
186
 
 
187
  typedef ArrayDimTagList<Tag1,Tag2,Tag3,Tag4,Tag5,Tag6,Tag7,Tag8> TagList ;
 
188
 
 
189
  enum { Rank       = TagList::Rank };
 
190
  enum { Natural    = true };
 
191
  enum { Reverse    = false };
 
192
  enum { Contiguous = true };
 
193
 
 
194
  template < unsigned Ord >
 
195
  struct Tag {
 
196
    typedef typename ArrayDimTagListAt<TagList,Ord>::type type ;
 
197
  };
 
198
 
 
199
private:
 
200
 
 
201
  typedef ArrayHelper< TagList::Rank > Helper ;
 
202
 
 
203
  Scalar * m_ptr ;
 
204
  size_t   m_stride[ TagList::Rank ];
 
205
 
 
206
  template< typename , class , class , class , class ,
 
207
                       class , class , class , class >
 
208
  friend class ArrayNatural ;
 
209
 
 
210
  template< typename , class , class , class , class ,
 
211
                       class , class , class , class >
 
212
  friend class ArrayFortran ;
 
213
 
 
214
public:
 
215
 
 
216
  //----------------------------------
 
217
 
 
218
  unsigned rank() const { return Rank ; }
 
219
  bool natural() const { return Natural ; }
 
220
  bool reverse() const { return Reverse ; }
 
221
  bool contiguous() const { return Contiguous ; }
 
222
 
 
223
  tag_type tag( unsigned ord ) const
 
224
    {
 
225
      array_bounds_checking( Rank , ord );
 
226
      return array_tags<TagList>()[ord] ;
 
227
    }
 
228
 
 
229
  //----------------------------------
 
230
 
 
231
  size_type size() const { return m_stride[ Rank - 1 ]; }
 
232
 
 
233
  template < unsigned Ord > size_type dimension() const
 
234
    {
 
235
      array_bounds_check_ordinal_is_less<Rank,Ord>();
 
236
      enum { I = ( Rank - 1 ) - Ord };
 
237
      return I ? m_stride[I] / m_stride[I-1] : m_stride[I] ;
 
238
    }
 
239
 
 
240
  size_type dimension( index_type ord ) const
 
241
    {
 
242
      array_bounds_checking( Rank , ord );
 
243
      const int i = ( Rank - 1 ) - ord ;
 
244
      return i ? m_stride[i] / m_stride[i-1] : m_stride[i] ;
 
245
    }
 
246
 
 
247
  void dimensions( std::vector<size_type> & n )
 
248
    {
 
249
      n.resize( Rank );
 
250
      n[Rank-1] = m_stride[0] ;
 
251
      for ( int i = 1 ; i < Rank ; ++i ) {
 
252
        n[ ( Rank - 1 ) - i ] = m_stride[i] / m_stride[i-1] ;
 
253
      }
 
254
    }
 
255
  
 
256
  //----------------------------------
 
257
  /** \brief Access member data */
 
258
  value_type * contiguous_data() const { return m_ptr ; }
 
259
 
 
260
  /** \brief Access member via full ordering of members. */
 
261
  value_type & operator[]( size_type i ) const
 
262
    {
 
263
      ARRAY_BOUNDS_CHECKING_ORDINAL( size() , i );
 
264
      return m_ptr[ i ];
 
265
    }
 
266
 
 
267
  /** \brief Access member via Rank 8 multi-index */
 
268
  value_type & operator()( index_type i1 , index_type i2 ,
 
269
                           index_type i3 , index_type i4 ,
 
270
                           index_type i5 , index_type i6 ,
 
271
                           index_type i7 , index_type i8 ) const
 
272
    {
 
273
      array_bounds_check_rank_is_equal<Rank,8>();
 
274
      ARRAY_BOUNDS_CHECKING_8(i1,i2,i3,i4,i5,i6,i7,i8);
 
275
      return m_ptr[ ARRAY_INDEX_NATURAL_8(m_stride,i1,i2,i3,i4,i5,i6,i7,i8) ];
 
276
    }
 
277
 
 
278
  /** \brief Access member via Rank 7 multi-index */
 
279
  value_type & operator()( index_type i1 , index_type i2 ,
 
280
                           index_type i3 , index_type i4 ,
 
281
                           index_type i5 , index_type i6 ,
 
282
                           index_type i7 ) const
 
283
    {
 
284
      array_bounds_check_rank_is_equal<Rank,7>();
 
285
      ARRAY_BOUNDS_CHECKING_7(i1,i2,i3,i4,i5,i6,i7);
 
286
      return m_ptr[ ARRAY_INDEX_NATURAL_7(m_stride,i1,i2,i3,i4,i5,i6,i7) ];
 
287
    }
 
288
 
 
289
  /** \brief Access member via Rank 6 multi-index */
 
290
  value_type & operator()( index_type i1 , index_type i2 ,
 
291
                           index_type i3 , index_type i4 ,
 
292
                           index_type i5 , index_type i6 ) const
 
293
    {
 
294
      array_bounds_check_rank_is_equal<Rank,6>();
 
295
      ARRAY_BOUNDS_CHECKING_6(i1,i2,i3,i4,i5,i6);
 
296
      return m_ptr[ ARRAY_INDEX_NATURAL_6(m_stride,i1,i2,i3,i4,i5,i6) ];
 
297
    }
 
298
 
 
299
  /** \brief Access member via Rank 5 multi-index */
 
300
  value_type & operator()( index_type i1 , index_type i2 ,
 
301
                           index_type i3 , index_type i4 ,
 
302
                           index_type i5 ) const
 
303
    {
 
304
      array_bounds_check_rank_is_equal<Rank,5>();
 
305
      ARRAY_BOUNDS_CHECKING_5(i1,i2,i3,i4,i5);
 
306
      return m_ptr[ ARRAY_INDEX_NATURAL_5(m_stride,i1,i2,i3,i4,i5) ];
 
307
    }
 
308
 
 
309
  /** \brief Access member via Rank 4 multi-index */
 
310
  value_type & operator()( index_type i1 , index_type i2 ,
 
311
                           index_type i3 , index_type i4 ) const
 
312
    {
 
313
      array_bounds_check_rank_is_equal<Rank,4>();
 
314
      ARRAY_BOUNDS_CHECKING_4(i1,i2,i3,i4);
 
315
      return m_ptr[ ARRAY_INDEX_NATURAL_4(m_stride,i1,i2,i3,i4) ];
 
316
    }
 
317
 
 
318
  /** \brief Access member via Rank 3 multi-index */
 
319
  value_type & operator()( index_type i1 , index_type i2 ,
 
320
                           index_type i3 ) const
 
321
    {
 
322
      array_bounds_check_rank_is_equal<Rank,3>();
 
323
      ARRAY_BOUNDS_CHECKING_3(i1,i2,i3);
 
324
      return m_ptr[ ARRAY_INDEX_NATURAL_3(m_stride,i1,i2,i3) ];
 
325
    }
 
326
 
 
327
  /** \brief Access member via Rank 2 multi-index */
 
328
  value_type & operator()( index_type i1 , index_type i2 ) const
 
329
    {
 
330
      array_bounds_check_rank_is_equal<Rank,2>();
 
331
      ARRAY_BOUNDS_CHECKING_2(i1,i2);
 
332
      return m_ptr[ ARRAY_INDEX_NATURAL_2(m_stride,i1,i2) ];
 
333
    }
 
334
 
 
335
  /** \brief Access member via Rank 1 multi-index */
 
336
  value_type & operator()( index_type i1 ) const
 
337
    {
 
338
      array_bounds_check_rank_is_equal<Rank,1>();
 
339
      ARRAY_BOUNDS_CHECKING_1(i1);
 
340
      return m_ptr[ i1 ];
 
341
    }
 
342
 
 
343
  //----------------------------------
 
344
  // Required constructors and assignment operators:
 
345
 
 
346
  ArrayNatural() : m_ptr(NULL) { Helper::zero( m_stride ); }
 
347
 
 
348
  ArrayNatural( const ArrayNatural & rhs )
 
349
    : m_ptr( rhs.m_ptr ) { Helper::copy( m_stride , rhs.m_stride ); }
 
350
 
 
351
  ArrayNatural & operator = ( const ArrayNatural & rhs )
 
352
    {
 
353
      m_ptr = rhs.m_ptr ;
 
354
      Helper::copy( m_stride , rhs.m_stride );
 
355
      return *this ;
 
356
    }
 
357
 
 
358
  //----------------------------------
 
359
  // Copy a reverse-map array
 
360
 
 
361
  typedef typename
 
362
    ArrayReverse<Scalar,Tag1,Tag2,Tag3,Tag4,Tag5,Tag6,Tag7,Tag8>::fortran_type
 
363
      ReverseType ;
 
364
 
 
365
  ArrayNatural( const ReverseType & rhs )
 
366
    : m_ptr( rhs.m_ptr ) { Helper::copy( m_stride , rhs.m_stride ); }
 
367
 
 
368
  ArrayNatural & operator = ( const ReverseType & rhs )
 
369
    {
 
370
      m_ptr = rhs.m_ptr ;
 
371
      Helper::copy( m_stride , rhs.m_stride );
 
372
      return *this ;
 
373
    }
 
374
 
 
375
  //----------------------------------
 
376
  // Truncated view
 
377
 
 
378
  typedef typename
 
379
    ArrayTruncate<Scalar,Tag1,Tag2,Tag3,Tag4,Tag5,Tag6,Tag7,Tag8>::natural_type
 
380
      TruncatedViewType ;
 
381
 
 
382
  TruncatedViewType truncated_view( index_type i ) const
 
383
    {
 
384
      return TruncatedViewType( m_ptr + m_stride[ Rank - 2 ] * i ,
 
385
                                ArrayStride(), m_stride );
 
386
    }
 
387
 
 
388
  // "Power user" constructor, promising a properly strided input
 
389
  ArrayNatural( value_type * arg_ptr ,
 
390
                const ArrayStride & , const size_type * arg_stride )
 
391
    : m_ptr( arg_ptr ) { Helper::copy( m_stride , arg_stride ); }
 
392
 
 
393
  // "Power user" constructor, promising a properly strided input
 
394
  ArrayNatural( value_type * arg_ptr ,
 
395
                const ArrayStride & ,
 
396
                const size_type * arg_stride ,
 
397
                const size_type   arg_append )
 
398
    : m_ptr( arg_ptr )
 
399
    { Helper::stride_append( m_stride , arg_stride , arg_append ); }
 
400
 
 
401
  //----------------------------------
 
402
  // Take pointer to other compatible array types.
 
403
 
 
404
  template< class OtherArrayType >
 
405
  explicit ArrayNatural( const OtherArrayType & rhs )
 
406
    : m_ptr( rhs.contiguous_data() )
 
407
    { Helper::stride( m_stride , rhs ); }
 
408
 
 
409
  template< class OtherArrayType >
 
410
  ArrayNatural & operator = ( const OtherArrayType & rhs )
 
411
    {
 
412
      m_ptr = rhs.contiguous_data();
 
413
      Helper::stride( m_stride , rhs );
 
414
      return *this ;
 
415
    }
 
416
 
 
417
  //----------------------------------
 
418
  // Class specific constructors:
 
419
 
 
420
  ArrayNatural( value_type * arg_ptr ,
 
421
         size_type n1 , size_type n2 , size_type n3 , size_type n4 ,
 
422
         size_type n5 , size_type n6 , size_type n7 , size_type n8 )
 
423
    : m_ptr( arg_ptr )
 
424
    {
 
425
      array_bounds_check_rank_is_equal<Rank,8>();
 
426
      ARRAY_STRIDE_NATURAL_8( m_stride, n1, n2, n3, n4, n5, n6, n7, n8 );
 
427
    }
 
428
 
 
429
  ArrayNatural( value_type * arg_ptr ,
 
430
         size_type n1 , size_type n2 , size_type n3 , size_type n4 ,
 
431
         size_type n5 , size_type n6 , size_type n7 )
 
432
    : m_ptr( arg_ptr )
 
433
    {
 
434
      array_bounds_check_rank_is_equal<Rank,7>();
 
435
      ARRAY_STRIDE_NATURAL_7( m_stride, n1, n2, n3, n4, n5, n6, n7 );
 
436
    }
 
437
 
 
438
  ArrayNatural( value_type * arg_ptr ,
 
439
         size_type n1 , size_type n2 , size_type n3 , size_type n4 ,
 
440
         size_type n5 , size_type n6 )
 
441
    : m_ptr( arg_ptr )
 
442
    {
 
443
      array_bounds_check_rank_is_equal<Rank,6>();
 
444
      ARRAY_STRIDE_NATURAL_6( m_stride, n1, n2, n3, n4, n5, n6 );
 
445
    }
 
446
 
 
447
  ArrayNatural( value_type * arg_ptr ,
 
448
         size_type n1 , size_type n2 , size_type n3 , size_type n4 ,
 
449
         size_type n5 )
 
450
    : m_ptr( arg_ptr )
 
451
    {
 
452
      array_bounds_check_rank_is_equal<Rank,5>();
 
453
      ARRAY_STRIDE_NATURAL_5( m_stride, n1, n2, n3, n4, n5 );
 
454
    }
 
455
 
 
456
  ArrayNatural( value_type * arg_ptr ,
 
457
         size_type n1 , size_type n2 , size_type n3 , size_type n4 )
 
458
    : m_ptr( arg_ptr )
 
459
    {
 
460
      array_bounds_check_rank_is_equal<Rank,4>();
 
461
      ARRAY_STRIDE_NATURAL_4( m_stride, n1, n2, n3, n4 );
 
462
    }
 
463
 
 
464
  ArrayNatural( value_type * arg_ptr ,
 
465
         size_type n1 , size_type n2 , size_type n3 )
 
466
    : m_ptr( arg_ptr )
 
467
    {
 
468
      array_bounds_check_rank_is_equal<Rank,3>();
 
469
      ARRAY_STRIDE_NATURAL_3( m_stride, n1, n2, n3 );
 
470
    }
 
471
 
 
472
  ArrayNatural( value_type * arg_ptr ,
 
473
         size_type n1 , size_type n2 )
 
474
    : m_ptr( arg_ptr )
 
475
    {
 
476
      array_bounds_check_rank_is_equal<Rank,2>();
 
477
      ARRAY_STRIDE_NATURAL_2( m_stride, n1, n2 );
 
478
    }
 
479
 
 
480
  ArrayNatural( value_type * arg_ptr , size_type n1 )
 
481
    : m_ptr( arg_ptr )
 
482
    {
 
483
      array_bounds_check_rank_is_equal<Rank,1>();
 
484
      m_stride[0] = n1 ;
 
485
    }
 
486
};
 
487
 
 
488
//----------------------------------------------------------------------
 
489
 
 
490
template< typename Scalar ,
 
491
          class Tag1 , class Tag2 , class Tag3 , class Tag4 ,
 
492
          class Tag5 , class Tag6 , class Tag7 , class Tag8 >
 
493
class ArrayFortran
 
494
{
 
495
public:
 
496
  // Required by all arrays:
 
497
  typedef Scalar           value_type ;
 
498
  typedef int              index_type ;
 
499
  typedef size_t           size_type ;
 
500
  typedef const ArrayDimTag * tag_type ;
 
501
 
 
502
  // Required by compile-time knowledge arrays:
 
503
 
 
504
  typedef ArrayDimTagList<Tag1,Tag2,Tag3,Tag4,Tag5,Tag6,Tag7,Tag8> TagList ;
 
505
 
 
506
  enum { Rank       = TagList::Rank };
 
507
  enum { Natural    = false };
 
508
  enum { Reverse    = true };
 
509
  enum { Contiguous = true };
 
510
 
 
511
  template < unsigned Ord >
 
512
  struct Tag {
 
513
    typedef typename ArrayDimTagListAt<TagList,Ord>::type type ;
 
514
  };
 
515
 
 
516
private:
 
517
 
 
518
  typedef ArrayHelper< TagList::Rank > Helper ;
 
519
 
 
520
  Scalar * m_ptr ;
 
521
  size_t m_stride[ TagList::Rank ];
 
522
 
 
523
  template< typename ,
 
524
            class , class , class , class ,
 
525
            class , class , class , class >
 
526
  friend class ArrayNatural ;
 
527
 
 
528
  template< typename ,
 
529
            class , class , class , class ,
 
530
            class , class , class , class >
 
531
  friend class ArrayFortran ;
 
532
 
 
533
public:
 
534
 
 
535
  //----------------------------------
 
536
 
 
537
  unsigned rank() const { return Rank ; }
 
538
  bool natural() const { return Natural ; }
 
539
  bool reverse() const { return Reverse ; }
 
540
  bool contiguous() const { return Contiguous ; }
 
541
 
 
542
  tag_type tag( unsigned ord ) const
 
543
    {
 
544
      array_bounds_checking( Rank , ord );
 
545
      return array_tags<TagList>()[ord] ;
 
546
    }
 
547
 
 
548
  //----------------------------------
 
549
 
 
550
  size_type size() const { return m_stride[ Rank - 1 ]; }
 
551
 
 
552
  template < unsigned Ord > size_type dimension() const
 
553
    {
 
554
      array_bounds_check_ordinal_is_less<Rank,Ord>();
 
555
      return Ord ? m_stride[Ord] / m_stride[Ord-1] : m_stride[Ord] ;
 
556
    }
 
557
 
 
558
  size_type dimension( index_type ord ) const
 
559
    {
 
560
      array_bounds_checking( Rank , ord );
 
561
      return ord ? m_stride[ord] / m_stride[ord-1] : m_stride[ord] ;
 
562
    }
 
563
 
 
564
  void dimensions( std::vector<size_type> & n )
 
565
    {
 
566
      n.resize( Rank );
 
567
      n[0] = m_stride[0] ;
 
568
      for ( int i = 1 ; i < Rank ; ++i ) {
 
569
        n[i] = m_stride[i] / m_stride[i-1] ;
 
570
      }
 
571
    }
 
572
  
 
573
  //----------------------------------
 
574
  /** \brief Access member data */
 
575
  value_type * contiguous_data() const { return m_ptr ; }
 
576
 
 
577
  /** \brief Access member via full ordering of members. */
 
578
  value_type & operator[]( size_type i ) const
 
579
    {
 
580
      ARRAY_BOUNDS_CHECKING_ORDINAL( size() , i );
 
581
      return m_ptr[ i ];
 
582
    }
 
583
 
 
584
  /** \brief Access member via Rank 8 multi-index */
 
585
  value_type & operator()( index_type i1 , index_type i2 ,
 
586
                           index_type i3 , index_type i4 ,
 
587
                           index_type i5 , index_type i6 ,
 
588
                           index_type i7 , index_type i8 ) const
 
589
    {
 
590
      array_bounds_check_rank_is_equal<Rank,8>();
 
591
      ARRAY_BOUNDS_CHECKING_8(i1,i2,i3,i4,i5,i6,i7,i8);
 
592
      return m_ptr[ ARRAY_INDEX_FORTRAN_8(m_stride,i1,i2,i3,i4,i5,i6,i7,i8) ];
 
593
    }
 
594
 
 
595
  /** \brief Access member via Rank 7 multi-index */
 
596
  value_type & operator()( index_type i1 , index_type i2 ,
 
597
                           index_type i3 , index_type i4 ,
 
598
                           index_type i5 , index_type i6 ,
 
599
                           index_type i7 ) const
 
600
    {
 
601
      array_bounds_check_rank_is_equal<Rank,7>();
 
602
      ARRAY_BOUNDS_CHECKING_7(i1,i2,i3,i4,i5,i6,i7);
 
603
      return m_ptr[ ARRAY_INDEX_FORTRAN_7(m_stride,i1,i2,i3,i4,i5,i6,i7) ];
 
604
    }
 
605
 
 
606
  /** \brief Access member via Rank 6 multi-index */
 
607
  value_type & operator()( index_type i1 , index_type i2 ,
 
608
                           index_type i3 , index_type i4 ,
 
609
                           index_type i5 , index_type i6 ) const
 
610
    {
 
611
      array_bounds_check_rank_is_equal<Rank,6>();
 
612
      ARRAY_BOUNDS_CHECKING_6(i1,i2,i3,i4,i5,i6);
 
613
      return m_ptr[ ARRAY_INDEX_FORTRAN_6(m_stride,i1,i2,i3,i4,i5,i6) ];
 
614
    }
 
615
 
 
616
  /** \brief Access member via Rank 5 multi-index */
 
617
  value_type & operator()( index_type i1 , index_type i2 ,
 
618
                           index_type i3 , index_type i4 ,
 
619
                           index_type i5 ) const
 
620
    {
 
621
      array_bounds_check_rank_is_equal<Rank,5>();
 
622
      ARRAY_BOUNDS_CHECKING_5(i1,i2,i3,i4,i5);
 
623
      return m_ptr[ ARRAY_INDEX_FORTRAN_5(m_stride,i1,i2,i3,i4,i5) ];
 
624
    }
 
625
 
 
626
  /** \brief Access member via Rank 4 multi-index */
 
627
  value_type & operator()( index_type i1 , index_type i2 ,
 
628
                           index_type i3 , index_type i4 ) const
 
629
    {
 
630
      array_bounds_check_rank_is_equal<Rank,4>();
 
631
      ARRAY_BOUNDS_CHECKING_4(i1,i2,i3,i4);
 
632
      return m_ptr[ ARRAY_INDEX_FORTRAN_4(m_stride,i1,i2,i3,i4) ];
 
633
    }
 
634
 
 
635
  /** \brief Access member via Rank 3 multi-index */
 
636
  value_type & operator()( index_type i1 , index_type i2 ,
 
637
                           index_type i3 ) const
 
638
    {
 
639
      array_bounds_check_rank_is_equal<Rank,3>();
 
640
      ARRAY_BOUNDS_CHECKING_3(i1,i2,i3);
 
641
      return m_ptr[ ARRAY_INDEX_FORTRAN_3(m_stride,i1,i2,i3) ];
 
642
    }
 
643
 
 
644
  /** \brief Access member via Rank 2 multi-index */
 
645
  value_type & operator()( index_type i1 , index_type i2 ) const
 
646
    {
 
647
      array_bounds_check_rank_is_equal<Rank,2>();
 
648
      ARRAY_BOUNDS_CHECKING_2(i1,i2);
 
649
      return m_ptr[ ARRAY_INDEX_FORTRAN_2(m_stride,i1,i2) ];
 
650
    }
 
651
 
 
652
  /** \brief Access member via Rank 1 multi-index */
 
653
  value_type & operator()( index_type i1 ) const
 
654
    {
 
655
      array_bounds_check_rank_is_equal<Rank,1>();
 
656
      ARRAY_BOUNDS_CHECKING_1(i1);
 
657
      return m_ptr[ i1 ];
 
658
    }
 
659
 
 
660
  //----------------------------------
 
661
  // Required constructors and assignment operators:
 
662
 
 
663
  ArrayFortran() : m_ptr(NULL) { Helper::zero( m_stride ); }
 
664
 
 
665
  ArrayFortran( const ArrayFortran & rhs )
 
666
    : m_ptr( rhs.m_ptr ) { Helper::copy( m_stride , rhs.m_stride ); }
 
667
 
 
668
  ArrayFortran & operator = ( const ArrayFortran & rhs )
 
669
    {
 
670
      m_ptr = rhs.m_ptr ;
 
671
      Helper::copy( m_stride , rhs.m_stride );
 
672
      return *this ;
 
673
    }
 
674
 
 
675
  //----------------------------------
 
676
  // Copy a reverse-map array
 
677
 
 
678
  typedef typename
 
679
    ArrayReverse<Scalar,Tag1,Tag2,Tag3,Tag4,Tag5,Tag6,Tag7,Tag8>::natural_type
 
680
      ReverseType ;
 
681
 
 
682
  ArrayFortran( const ReverseType & rhs )
 
683
    : m_ptr( rhs.m_ptr ) { Helper::copy( m_stride , rhs.m_stride ); }
 
684
 
 
685
  ArrayFortran & operator = ( const ReverseType & rhs )
 
686
    {
 
687
      m_ptr = rhs.m_ptr ;
 
688
      Helper::copy( m_stride , rhs.m_stride );
 
689
      return *this ;
 
690
    }
 
691
 
 
692
  //----------------------------------
 
693
  // Truncated view
 
694
 
 
695
  typedef typename
 
696
    ArrayTruncate<Scalar,Tag1,Tag2,Tag3,Tag4,Tag5,Tag6,Tag7,Tag8>::fortran_type
 
697
      TruncatedViewType ;
 
698
 
 
699
  TruncatedViewType truncated_view( index_type i ) const
 
700
    {
 
701
      return TruncatedViewType( m_ptr + m_stride[ Rank - 2 ] * i ,
 
702
                                ArrayStride() , m_stride );
 
703
    }
 
704
 
 
705
  // "Power user" constructor, promising a properly strided input
 
706
  ArrayFortran( value_type * arg_ptr ,
 
707
                const ArrayStride & , const size_type * arg_stride )
 
708
    : m_ptr( arg_ptr ) { Helper::copy( m_stride , arg_stride ); }
 
709
 
 
710
  // "Power user" constructor, promising a properly strided input
 
711
  ArrayFortran( value_type * arg_ptr ,
 
712
                const ArrayStride & ,
 
713
                const size_type * arg_stride ,
 
714
                const size_type   arg_append )
 
715
    : m_ptr( arg_ptr )
 
716
    { Helper::stride_append( m_stride , arg_stride , arg_append ); }
 
717
 
 
718
  //----------------------------------
 
719
  // Take pointer to other compatible array types.
 
720
 
 
721
  template< class OtherArrayType >
 
722
  explicit ArrayFortran( const OtherArrayType & rhs )
 
723
    : m_ptr( rhs.contiguous_data() )
 
724
    { Helper::stride( m_stride , rhs ); }
 
725
 
 
726
  template< class OtherArrayType >
 
727
  ArrayFortran & operator = ( const OtherArrayType & rhs )
 
728
    {
 
729
      m_ptr = rhs.contiguous_data();
 
730
      Helper::stride( m_stride , rhs );
 
731
      return *this ;
 
732
    }
 
733
 
 
734
  //----------------------------------
 
735
  // Class specific constructors:
 
736
 
 
737
  ArrayFortran( value_type * arg_ptr ,
 
738
         size_type n1 , size_type n2 , size_type n3 , size_type n4 ,
 
739
         size_type n5 , size_type n6 , size_type n7 , size_type n8 )
 
740
    : m_ptr( arg_ptr )
 
741
    {
 
742
      array_bounds_check_rank_is_equal<Rank,8>();
 
743
      ARRAY_STRIDE_FORTRAN_8( m_stride, n1, n2, n3, n4, n5, n6, n7, n8 );
 
744
    }
 
745
 
 
746
  ArrayFortran( value_type * arg_ptr ,
 
747
         size_type n1 , size_type n2 , size_type n3 , size_type n4 ,
 
748
         size_type n5 , size_type n6 , size_type n7 )
 
749
    : m_ptr( arg_ptr )
 
750
    {
 
751
      array_bounds_check_rank_is_equal<Rank,7>();
 
752
      ARRAY_STRIDE_FORTRAN_7( m_stride, n1, n2, n3, n4, n5, n6, n7 );
 
753
    }
 
754
 
 
755
  ArrayFortran( value_type * arg_ptr ,
 
756
         size_type n1 , size_type n2 , size_type n3 , size_type n4 ,
 
757
         size_type n5 , size_type n6 )
 
758
    : m_ptr( arg_ptr )
 
759
    {
 
760
      array_bounds_check_rank_is_equal<Rank,6>();
 
761
      ARRAY_STRIDE_FORTRAN_6( m_stride, n1, n2, n3, n4, n5, n6 );
 
762
    }
 
763
 
 
764
  ArrayFortran( value_type * arg_ptr ,
 
765
         size_type n1 , size_type n2 , size_type n3 , size_type n4 ,
 
766
         size_type n5 )
 
767
    : m_ptr( arg_ptr )
 
768
    {
 
769
      array_bounds_check_rank_is_equal<Rank,5>();
 
770
      ARRAY_STRIDE_FORTRAN_5( m_stride, n1, n2, n3, n4, n5 );
 
771
    }
 
772
 
 
773
  ArrayFortran( value_type * arg_ptr ,
 
774
         size_type n1 , size_type n2 , size_type n3 , size_type n4 )
 
775
    : m_ptr( arg_ptr )
 
776
    {
 
777
      array_bounds_check_rank_is_equal<Rank,4>();
 
778
      ARRAY_STRIDE_FORTRAN_4( m_stride, n1, n2, n3, n4 );
 
779
    }
 
780
 
 
781
  ArrayFortran( value_type * arg_ptr ,
 
782
         size_type n1 , size_type n2 , size_type n3 )
 
783
    : m_ptr( arg_ptr )
 
784
    {
 
785
      array_bounds_check_rank_is_equal<Rank,3>();
 
786
      ARRAY_STRIDE_FORTRAN_3( m_stride, n1, n2, n3 );
 
787
    }
 
788
 
 
789
  ArrayFortran( value_type * arg_ptr ,
 
790
         size_type n1 , size_type n2 )
 
791
    : m_ptr( arg_ptr )
 
792
    {
 
793
      array_bounds_check_rank_is_equal<Rank,2>();
 
794
      ARRAY_STRIDE_FORTRAN_2( m_stride, n1, n2 );
 
795
    }
 
796
 
 
797
  ArrayFortran( value_type * arg_ptr ,
 
798
         size_type n1 )
 
799
    : m_ptr( arg_ptr )
 
800
    {
 
801
      array_bounds_check_rank_is_equal<Rank,1>();
 
802
      m_stride[0] = n1 ;
 
803
    }
 
804
};
 
805
 
 
806
//----------------------------------------------------------------------
 
807
// Specialization for arrays with rank and dimension tags as runtime info.
 
808
 
 
809
template< typename Scalar >
 
810
class ArrayNatural< Scalar, void,void,void,void,void,void,void,void>
 
811
{
 
812
public:
 
813
 
 
814
  typedef Scalar           value_type ;
 
815
  typedef int              index_type ;
 
816
  typedef size_t           size_type ;
 
817
  typedef const ArrayDimTag * tag_type ;
 
818
 
 
819
  enum { Natural    = true };
 
820
  enum { Reverse    = false };
 
821
  enum { Contiguous = true };
 
822
 
 
823
private:
 
824
 
 
825
  Scalar  * m_ptr ;
 
826
  size_type m_rank ;
 
827
  size_type m_stride[8];
 
828
  tag_type  m_tag[8] ;
 
829
 
 
830
  typedef ArrayHelper<8> Helper ;
 
831
 
 
832
  template< typename ,
 
833
            class , class , class , class ,
 
834
            class , class , class , class >
 
835
  friend class ArrayNatural ;
 
836
 
 
837
  template< typename ,
 
838
            class , class , class , class ,
 
839
            class , class , class , class >
 
840
  friend class ArrayFortran ;
 
841
 
 
842
public:
 
843
 
 
844
  //----------------------------------
 
845
 
 
846
  unsigned rank() const { return m_rank ; }
 
847
  bool natural() const { return Natural ; }
 
848
  bool reverse() const { return Reverse ; }
 
849
  bool contiguous() const { return Contiguous ; }
 
850
 
 
851
  tag_type tag( unsigned ord ) const
 
852
    {
 
853
      array_bounds_checking( m_rank , ord );
 
854
      return m_tag[ ( m_rank - 1 ) - ord ];
 
855
    }
 
856
 
 
857
  //----------------------------------
 
858
 
 
859
  size_type size() const { return m_stride[ m_rank - 1 ]; }
 
860
 
 
861
  size_type dimension( index_type ord ) const
 
862
    {
 
863
      array_bounds_checking( m_rank , ord );
 
864
      const int i = ( m_rank - 1 ) - ord ;
 
865
      return i ? m_stride[i] / m_stride[i-1] : m_stride[i] ;
 
866
    }
 
867
 
 
868
  void dimensions( std::vector<size_type> & n )
 
869
    {
 
870
      n.resize( m_rank );
 
871
      n[m_rank-1] = m_stride[0] ;
 
872
      for ( int i = 1 ; i < m_rank ; ++i ) {
 
873
        n[ ( m_rank - 1 ) - i ] = m_stride[i] / m_stride[i-1] ;
 
874
      }
 
875
    }
 
876
  
 
877
  //----------------------------------
 
878
  /** \brief Access member data */
 
879
  value_type * contiguous_data() const { return m_ptr ; }
 
880
 
 
881
  /** \brief Access member via full ordering of members. */
 
882
  value_type & operator[]( size_type i ) const
 
883
    {
 
884
      ARRAY_BOUNDS_CHECKING_ORDINAL( size() , i );
 
885
      return m_ptr[ i ];
 
886
    }
 
887
 
 
888
  /** \brief Access member via Rank 8 multi-index */
 
889
  value_type & operator()( index_type i1 , index_type i2 ,
 
890
                           index_type i3 , index_type i4 ,
 
891
                           index_type i5 , index_type i6 ,
 
892
                           index_type i7 , index_type i8 ) const
 
893
    {
 
894
      ARRAY_BOUNDS_CHECKING_8(i1,i2,i3,i4,i5,i6,i7,i8);
 
895
      return m_ptr[ ARRAY_INDEX_NATURAL_8(m_stride,i1,i2,i3,i4,i5,i6,i7,i8) ];
 
896
    }
 
897
 
 
898
  /** \brief Access member via Rank 7 multi-index */
 
899
  value_type & operator()( index_type i1 , index_type i2 ,
 
900
                           index_type i3 , index_type i4 ,
 
901
                           index_type i5 , index_type i6 ,
 
902
                           index_type i7 ) const
 
903
    {
 
904
      ARRAY_BOUNDS_CHECKING_7(i1,i2,i3,i4,i5,i6,i7);
 
905
      return m_ptr[ ARRAY_INDEX_NATURAL_7(m_stride,i1,i2,i3,i4,i5,i6,i7) ];
 
906
    }
 
907
 
 
908
  /** \brief Access member via Rank 6 multi-index */
 
909
  value_type & operator()( index_type i1 , index_type i2 ,
 
910
                           index_type i3 , index_type i4 ,
 
911
                           index_type i5 , index_type i6 ) const
 
912
    {
 
913
      ARRAY_BOUNDS_CHECKING_6(i1,i2,i3,i4,i5,i6);
 
914
      return m_ptr[ ARRAY_INDEX_NATURAL_6(m_stride,i1,i2,i3,i4,i5,i6) ];
 
915
    }
 
916
 
 
917
  /** \brief Access member via Rank 5 multi-index */
 
918
  value_type & operator()( index_type i1 , index_type i2 ,
 
919
                           index_type i3 , index_type i4 ,
 
920
                           index_type i5 ) const
 
921
    {
 
922
      ARRAY_BOUNDS_CHECKING_5(i1,i2,i3,i4,i5);
 
923
      return m_ptr[ ARRAY_INDEX_NATURAL_5(m_stride,i1,i2,i3,i4,i5) ];
 
924
    }
 
925
 
 
926
  /** \brief Access member via Rank 4 multi-index */
 
927
  value_type & operator()( index_type i1 , index_type i2 ,
 
928
                           index_type i3 , index_type i4 ) const
 
929
    {
 
930
      ARRAY_BOUNDS_CHECKING_4(i1,i2,i3,i4);
 
931
      return m_ptr[ ARRAY_INDEX_NATURAL_4(m_stride,i1,i2,i3,i4) ];
 
932
    }
 
933
 
 
934
  /** \brief Access member via Rank 3 multi-index */
 
935
  value_type & operator()( index_type i1 , index_type i2 ,
 
936
                           index_type i3 ) const
 
937
    {
 
938
      ARRAY_BOUNDS_CHECKING_3(i1,i2,i3);
 
939
      return m_ptr[ ARRAY_INDEX_NATURAL_3(m_stride,i1,i2,i3) ];
 
940
    }
 
941
 
 
942
  /** \brief Access member via Rank 2 multi-index */
 
943
  value_type & operator()( index_type i1 , index_type i2 ) const
 
944
    {
 
945
      ARRAY_BOUNDS_CHECKING_2(i1,i2);
 
946
      return m_ptr[ ARRAY_INDEX_NATURAL_2(m_stride,i1,i2) ];
 
947
    }
 
948
 
 
949
  /** \brief Access member via Rank 1 multi-index */
 
950
  value_type & operator()( index_type i1 ) const
 
951
    {
 
952
      ARRAY_BOUNDS_CHECKING_1(i1);
 
953
      return m_ptr[ i1 ];
 
954
    }
 
955
 
 
956
  //----------------------------------
 
957
  /** \brief Truncate view of the array */
 
958
 
 
959
  typedef ArrayNatural<Scalar> TruncateViewType ;
 
960
 
 
961
  TruncateViewType truncated_view( index_type i ) const
 
962
    {
 
963
      return TruncateViewType( m_ptr + m_stride[ m_rank - 2 ] * i ,
 
964
                               m_rank - 1 ,
 
965
                               ArrayStride() ,
 
966
                               m_stride ,
 
967
                               m_tag );
 
968
    }
 
969
 
 
970
  ArrayNatural( value_type      * arg_ptr ,
 
971
                size_type         arg_rank ,
 
972
                const ArrayStride & ,
 
973
                const size_type * arg_stride ,
 
974
                tag_type        * arg_tag )
 
975
    : m_ptr( arg_ptr ), m_rank( arg_rank )
 
976
    {
 
977
      Helper::copy( m_stride , arg_stride );
 
978
      Helper::copy( m_tag ,    arg_tag );
 
979
    }
 
980
 
 
981
  //----------------------------------
 
982
  // Required constructors and assignment operators:
 
983
 
 
984
  ArrayNatural()
 
985
    : m_ptr(NULL), m_rank(0)
 
986
    {
 
987
      Helper::zero( m_stride );
 
988
      Helper::zero( m_tag );
 
989
    }
 
990
 
 
991
  ArrayNatural( const ArrayNatural & rhs )
 
992
    : m_ptr( rhs.m_ptr ), m_rank( rhs.m_rank )
 
993
    {
 
994
      Helper::copy( m_stride , rhs.m_stride );
 
995
      Helper::copy( m_tag ,    rhs.m_tag );
 
996
    }
 
997
 
 
998
  ArrayNatural & operator = ( const ArrayNatural & rhs )
 
999
    {
 
1000
      m_ptr  = rhs.m_ptr ;
 
1001
      m_rank = rhs.m_rank ;
 
1002
      Helper::copy( m_stride , rhs.m_stride );
 
1003
      Helper::copy( m_tag ,    rhs.m_tag );
 
1004
      return *this ;
 
1005
    }
 
1006
 
 
1007
  //----------------------------------
 
1008
 
 
1009
  ArrayNatural( const ArrayFortran<Scalar> & rhs )
 
1010
    : m_ptr( rhs.m_ptr ), m_rank( rhs.m_rank )
 
1011
    {
 
1012
      Helper::copy( m_stride , rhs.m_stride );
 
1013
      Helper::copy( m_tag ,    rhs.m_tag );
 
1014
    }
 
1015
 
 
1016
  ArrayNatural & operator = ( const ArrayFortran<Scalar> & rhs )
 
1017
    {
 
1018
      m_ptr  = rhs.m_ptr ;
 
1019
      m_rank = rhs.m_rank ;
 
1020
      Helper::copy( m_stride , rhs.m_stride );
 
1021
      Helper::copy( m_tag ,    rhs.m_tag );
 
1022
      return *this ;
 
1023
    }
 
1024
 
 
1025
  //----------------------------------
 
1026
  // Take pointer to other compatible array types.
 
1027
 
 
1028
  template< class OtherArrayType >
 
1029
  ArrayNatural( const OtherArrayType & rhs )
 
1030
    : m_ptr( rhs.contiguous_data() ), m_rank( rhs.rank() )
 
1031
    {
 
1032
      Helper::stride( m_stride , rhs );
 
1033
      Helper::tags( m_tag , rhs );
 
1034
    }
 
1035
 
 
1036
  template< class OtherArrayType >
 
1037
  ArrayNatural & operator = ( const OtherArrayType & rhs )
 
1038
    {
 
1039
      m_ptr  = rhs.contiguous_data();
 
1040
      m_rank = rhs.rank();
 
1041
      Helper::stride( m_stride , rhs );
 
1042
      Helper::tags(   m_tag ,    rhs );
 
1043
      return *this ;
 
1044
    }
 
1045
 
 
1046
  //----------------------------------
 
1047
  // Class specific constructors:
 
1048
 
 
1049
  ArrayNatural( value_type * arg_ptr ,
 
1050
                const std::vector<size_type> & arg_size ,
 
1051
                const std::vector<tag_type>  & arg_tag )
 
1052
    : m_ptr( arg_ptr ),
 
1053
      m_rank( arg_size.size() == arg_tag.size() ? arg_size.size() : 0 )
 
1054
    {
 
1055
      array_stride_from_natural_sizes( m_rank , m_stride , & arg_size[0] );
 
1056
      array_stride_natural_tag( m_rank , m_tag , & arg_tag[0] );
 
1057
    }
 
1058
};
 
1059
 
 
1060
//----------------------------------------------------------------------
 
1061
// Specialization for arrays with rank and dimension tags as runtime info.
 
1062
 
 
1063
template< typename Scalar >
 
1064
class ArrayFortran< Scalar, void,void,void,void,void,void,void,void>
 
1065
{
 
1066
public:
 
1067
 
 
1068
  typedef Scalar              value_type ;
 
1069
  typedef int                 index_type ;
 
1070
  typedef size_t              size_type ;
 
1071
  typedef const ArrayDimTag * tag_type ;
 
1072
 
 
1073
  enum { Natural    = false };
 
1074
  enum { Reverse    = true };
 
1075
  enum { Contiguous = true };
 
1076
 
 
1077
private:
 
1078
 
 
1079
  Scalar  * m_ptr ;
 
1080
  size_type m_rank ;
 
1081
  size_type m_stride[8];
 
1082
  tag_type  m_tag[8] ;
 
1083
 
 
1084
  typedef ArrayHelper<8> Helper ;
 
1085
 
 
1086
  template< typename ,
 
1087
            class , class , class , class ,
 
1088
            class , class , class , class >
 
1089
  friend class ArrayNatural ;
 
1090
 
 
1091
  template< typename ,
 
1092
            class , class , class , class ,
 
1093
            class , class , class , class >
 
1094
  friend class ArrayFortran ;
 
1095
 
 
1096
public:
 
1097
 
 
1098
  //----------------------------------
 
1099
 
 
1100
  unsigned rank() const { return m_rank ; }
 
1101
  bool natural() const { return Natural ; }
 
1102
  bool reverse() const { return Reverse ; }
 
1103
  bool contiguous() const { return Contiguous ; }
 
1104
 
 
1105
  tag_type tag( unsigned ord ) const
 
1106
    {
 
1107
      array_bounds_checking( m_rank , ord );
 
1108
      return m_tag[ ord ];
 
1109
    }
 
1110
 
 
1111
  //----------------------------------
 
1112
 
 
1113
  size_type size() const { return m_stride[ m_rank - 1 ]; }
 
1114
 
 
1115
  size_type dimension( index_type ord ) const
 
1116
    {
 
1117
      array_bounds_checking( m_rank , ord );
 
1118
      return ord ? m_stride[ord] / m_stride[ord-1] : m_stride[ord] ;
 
1119
    }
 
1120
 
 
1121
  void dimensions( std::vector<size_type> & n )
 
1122
    {
 
1123
      n.resize( m_rank );
 
1124
      n[0] = m_stride[0] ;
 
1125
      for ( int i = 1 ; i < m_rank ; ++i ) {
 
1126
        n[i] = m_stride[i] / m_stride[i-1] ;
 
1127
      }
 
1128
    }
 
1129
  
 
1130
  //----------------------------------
 
1131
  /** \brief Access member data */
 
1132
  value_type * contiguous_data() const { return m_ptr ; }
 
1133
 
 
1134
  /** \brief Access member via full ordering of members. */
 
1135
  value_type & operator[]( size_type i ) const
 
1136
    {
 
1137
      ARRAY_BOUNDS_CHECKING_ORDINAL( size() , i );
 
1138
      return m_ptr[ i ];
 
1139
    }
 
1140
 
 
1141
  /** \brief Access member via Rank 8 multi-index */
 
1142
  value_type & operator()( index_type i1 , index_type i2 ,
 
1143
                           index_type i3 , index_type i4 ,
 
1144
                           index_type i5 , index_type i6 ,
 
1145
                           index_type i7 , index_type i8 ) const
 
1146
    {
 
1147
      ARRAY_BOUNDS_CHECKING_8(i1,i2,i3,i4,i5,i6,i7,i8);
 
1148
      return m_ptr[ ARRAY_INDEX_FORTRAN_8(m_stride,i1,i2,i3,i4,i5,i6,i7,i8) ];
 
1149
    }
 
1150
 
 
1151
  /** \brief Access member via Rank 7 multi-index */
 
1152
  value_type & operator()( index_type i1 , index_type i2 ,
 
1153
                           index_type i3 , index_type i4 ,
 
1154
                           index_type i5 , index_type i6 ,
 
1155
                           index_type i7 ) const
 
1156
    {
 
1157
      ARRAY_BOUNDS_CHECKING_7(i1,i2,i3,i4,i5,i6,i7);
 
1158
      return m_ptr[ ARRAY_INDEX_FORTRAN_7(m_stride,i1,i2,i3,i4,i5,i6,i7) ];
 
1159
    }
 
1160
 
 
1161
  /** \brief Access member via Rank 6 multi-index */
 
1162
  value_type & operator()( index_type i1 , index_type i2 ,
 
1163
                           index_type i3 , index_type i4 ,
 
1164
                           index_type i5 , index_type i6 ) const
 
1165
    {
 
1166
      ARRAY_BOUNDS_CHECKING_6(i1,i2,i3,i4,i5,i6);
 
1167
      return m_ptr[ ARRAY_INDEX_FORTRAN_6(m_stride,i1,i2,i3,i4,i5,i6) ];
 
1168
    }
 
1169
 
 
1170
  /** \brief Access member via Rank 5 multi-index */
 
1171
  value_type & operator()( index_type i1 , index_type i2 ,
 
1172
                           index_type i3 , index_type i4 ,
 
1173
                           index_type i5 ) const
 
1174
    {
 
1175
      ARRAY_BOUNDS_CHECKING_5(i1,i2,i3,i4,i5);
 
1176
      return m_ptr[ ARRAY_INDEX_FORTRAN_5(m_stride,i1,i2,i3,i4,i5) ];
 
1177
    }
 
1178
 
 
1179
  /** \brief Access member via Rank 4 multi-index */
 
1180
  value_type & operator()( index_type i1 , index_type i2 ,
 
1181
                           index_type i3 , index_type i4 ) const
 
1182
    {
 
1183
      ARRAY_BOUNDS_CHECKING_4(i1,i2,i3,i4);
 
1184
      return m_ptr[ ARRAY_INDEX_FORTRAN_4(m_stride,i1,i2,i3,i4) ];
 
1185
    }
 
1186
 
 
1187
  /** \brief Access member via Rank 3 multi-index */
 
1188
  value_type & operator()( index_type i1 , index_type i2 ,
 
1189
                           index_type i3 ) const
 
1190
    {
 
1191
      ARRAY_BOUNDS_CHECKING_3(i1,i2,i3);
 
1192
      return m_ptr[ ARRAY_INDEX_FORTRAN_3(m_stride,i1,i2,i3) ];
 
1193
    }
 
1194
 
 
1195
  /** \brief Access member via Rank 2 multi-index */
 
1196
  value_type & operator()( index_type i1 , index_type i2 ) const
 
1197
    {
 
1198
      ARRAY_BOUNDS_CHECKING_2(i1,i2);
 
1199
      return m_ptr[ ARRAY_INDEX_FORTRAN_2(m_stride,i1,i2) ];
 
1200
    }
 
1201
 
 
1202
  /** \brief Access member via Rank 1 multi-index */
 
1203
  value_type & operator()( index_type i1 ) const
 
1204
    {
 
1205
      ARRAY_BOUNDS_CHECKING_1(i1);
 
1206
      return m_ptr[ i1 ];
 
1207
    }
 
1208
 
 
1209
  //----------------------------------
 
1210
  /** \brief Truncate view of the array */
 
1211
 
 
1212
  typedef ArrayFortran<Scalar> TruncateViewType ;
 
1213
 
 
1214
  TruncateViewType truncated_view( index_type i ) const
 
1215
    {
 
1216
      return TruncateViewType( m_ptr + m_stride[ m_rank - 2 ] * i ,
 
1217
                               m_rank - 1 ,
 
1218
                               ArrayStride() ,
 
1219
                               m_stride ,
 
1220
                               m_tag );
 
1221
    }
 
1222
 
 
1223
  ArrayFortran( value_type      * arg_ptr ,
 
1224
                size_type         arg_rank ,
 
1225
                const ArrayStride & ,
 
1226
                const size_type * arg_stride ,
 
1227
                tag_type        * arg_tag )
 
1228
    : m_ptr( arg_ptr ), m_rank( arg_rank )
 
1229
    {
 
1230
      Helper::copy( m_stride , arg_stride );
 
1231
      Helper::copy( m_tag ,    arg_tag );
 
1232
    }
 
1233
 
 
1234
  //----------------------------------
 
1235
  // Required constructors and assignment operators:
 
1236
 
 
1237
  ArrayFortran()
 
1238
    : m_ptr(NULL), m_rank(0)
 
1239
    {
 
1240
      Helper::zero( m_stride );
 
1241
      Helper::zero( m_tag );
 
1242
    }
 
1243
 
 
1244
  ArrayFortran( const ArrayFortran & rhs )
 
1245
    : m_ptr( rhs.m_ptr ), m_rank( rhs.m_rank )
 
1246
    {
 
1247
      Helper::copy( m_stride , rhs.m_stride );
 
1248
      Helper::copy( m_tag ,    rhs.m_tag );
 
1249
    }
 
1250
 
 
1251
  ArrayFortran & operator = ( const ArrayFortran & rhs )
 
1252
    {
 
1253
      m_ptr  = rhs.m_ptr ;
 
1254
      m_rank = rhs.m_rank ;
 
1255
      Helper::copy( m_stride , rhs.m_stride );
 
1256
      Helper::copy( m_tag ,    rhs.m_tag );
 
1257
      return *this ;
 
1258
    }
 
1259
 
 
1260
  ArrayFortran( const ArrayNatural<Scalar> & rhs )
 
1261
    : m_ptr( rhs.m_ptr ), m_rank( rhs.m_rank )
 
1262
    {
 
1263
      Helper::copy( m_stride , rhs.m_stride );
 
1264
      Helper::copy( m_tag ,    rhs.m_tag );
 
1265
    }
 
1266
 
 
1267
  ArrayFortran & operator = ( const ArrayNatural<Scalar> & rhs )
 
1268
    {
 
1269
      m_ptr  = rhs.m_ptr ;
 
1270
      m_rank = rhs.m_rank ;
 
1271
      Helper::copy( m_stride , rhs.m_stride );
 
1272
      Helper::copy( m_tag ,    rhs.m_tag );
 
1273
      return *this ;
 
1274
    }
 
1275
 
 
1276
  //----------------------------------
 
1277
  // Take pointer to other compatible array types.
 
1278
 
 
1279
  template< class OtherArrayType >
 
1280
  ArrayFortran( const OtherArrayType & rhs )
 
1281
    : m_ptr( rhs.contiguous_data() ), m_rank( rhs.rank() )
 
1282
    {
 
1283
      Helper::stride( m_stride , rhs );
 
1284
      Helper::tags( m_tag , rhs );
 
1285
    }
 
1286
 
 
1287
  template< class OtherArrayType >
 
1288
  ArrayFortran & operator = ( const OtherArrayType & rhs )
 
1289
    {
 
1290
      m_ptr = rhs.contiguous_data();
 
1291
      m_rank = rhs.m_rank();
 
1292
      Helper::stride( m_stride , rhs );
 
1293
      Helper::tags( m_tag , rhs );
 
1294
      return *this ;
 
1295
    }
 
1296
 
 
1297
  //----------------------------------
 
1298
  // Class specific constructors:
 
1299
 
 
1300
  ArrayFortran( value_type * arg_ptr ,
 
1301
                const std::vector<size_type> & arg_size ,
 
1302
                const std::vector<tag_type>  & arg_tag )
 
1303
    : m_ptr( arg_ptr ),
 
1304
      m_rank( arg_size.size() == arg_tag.size() ? arg_size.size() : 0 )
 
1305
    {
 
1306
      array_stride_from_fortran_sizes( m_rank , m_stride , & arg_size[0] );
 
1307
      array_stride_fortran_tag( m_rank , m_tag , & arg_tag[0] );
 
1308
    }
 
1309
};
 
1310
 
 
1311
//----------------------------------------------------------------------
 
1312
 
 
1313
/** \endcond */
 
1314
 
 
1315
}
 
1316
 
 
1317
#endif
 
1318