/[escript]/trunk/escript/src/DataArrayView.h
ViewVC logotype

Diff of /trunk/escript/src/DataArrayView.h

Parent Directory Parent Directory | Revision Log Revision Log | View Patch Patch

trunk/esys2/escript/src/Data/DataArrayView.h revision 100 by jgs, Wed Dec 15 03:48:48 2004 UTC trunk/escript/src/DataArrayView.h revision 800 by gross, Tue Aug 8 11:23:18 2006 UTC
# Line 1  Line 1 
1    // $Id$
2    
3  /*  /*
4   ******************************************************************************   ************************************************************
5   *                                                                            *   *          Copyright 2006 by ACcESS MNRF                   *
6   *       COPYRIGHT  ACcESS 2004 -  All Rights Reserved                        *   *                                                          *
7   *                                                                            *   *              http://www.access.edu.au                    *
8   * This software is the property of ACcESS. No part of this code              *   *       Primary Business: Queensland, Australia            *
9   * may be copied in any form or by any means without the expressed written    *   *  Licensed under the Open Software License version 3.0    *
10   * consent of ACcESS.  Copying, use or modification of this software          *   *     http://www.opensource.org/licenses/osl-3.0.php       *
11   * by any unauthorised person is illegal unless that person has a software    *   *                                                          *
12   * license agreement with ACcESS.                                             *   ************************************************************
  *                                                                            *  
  ******************************************************************************  
13  */  */
14                                                                              
15  #if !defined escript_DataArrayView_20040323_H  #if !defined escript_DataArrayView_20040323_H
16  #define escript_DataArrayView_20040323_H  #define escript_DataArrayView_20040323_H
17    #include "system_dep.h"
18    
19  #include "esysUtils/EsysAssert.h"  #include "esysUtils/EsysAssert.h"
20    
21  #include <vector>  #include "DataException.h"
22    #include "DataVector.h"
23    #include "LocalOps.h"
24    
25  #include <boost/python/numeric.hpp>  #include <boost/python/numeric.hpp>
26  #include <boost/python/object.hpp>  #include <boost/python/object.hpp>
27  #include <boost/shared_ptr.hpp>  
28    #include <vector>
29    
30  namespace escript {  namespace escript {
31    
32  /**  /**
33     \brief     \brief
34     DataArrayView provides a view of data allocated externally.     DataArrayView provides a view of external data, configured according
35       to the given shape information and offset.
36    
37     Description:     Description:
38     DataArrayView provides a view of data allocated externally. The     DataArrayView provides a view of data allocated externally. The
39     external data should provide sufficient values so that the dimensions     external data should provide sufficient values so that the dimensions
40     specified for the view will be satisfied. The viewer can update     specified for the view shape will be satisfied. The viewer can update
41     values within the external data but cannot resize the external data.     values within the external data but cannot resize the external data.
42    
43       The view provided represents a single n-dimensional data-point
44       comprised of values taken from the given data array, starting at the
45       specified offset and extending for as many values as are necessary to
46       satisfy the given shape. The default offset can be changed, or different
47       offsets specified, in order to provide views of other data-points in
48       the underlying data array.
49  */  */
50    
51  class DataArrayView {  class DataArrayView {
52    
53    friend bool operator==(const DataArrayView& left, const DataArrayView& right);    ESCRIPT_DLL_API friend bool operator==(const DataArrayView& left, const DataArrayView& right);
54    friend bool operator!=(const DataArrayView& left, const DataArrayView& right);    ESCRIPT_DLL_API friend bool operator!=(const DataArrayView& left, const DataArrayView& right);
55    
56   public:   public:
57    
58    typedef std::vector<double> ValueType;    //
59    typedef std::vector<int> ShapeType;    // Some basic types which define the data values and view shapes.
60      typedef DataVector                        ValueType;
61      //typedef std::vector<double>               ValueType;
62      typedef std::vector<int>                  ShapeType;
63    typedef std::vector<std::pair<int, int> > RegionType;    typedef std::vector<std::pair<int, int> > RegionType;
64      typedef std::vector<std::pair<int, int> > RegionLoopRangeType;
65    
66    /**    /**
67       \brief       \brief
68       Default constructor for DataArrayView.       Default constructor for DataArrayView.
69    
70       Description:       Description:
71       Default constructor for DataArrayView. Creates an       Default constructor for DataArrayView.
72       empty view with no associated data.       Creates an empty view with no associated data array.
73    
74         This is essentially useless but here for completeness.
75    */    */
76      ESCRIPT_DLL_API
77    DataArrayView();    DataArrayView();
78    
79    /**    /**
# Line 63  class DataArrayView { Line 83  class DataArrayView {
83       Description:       Description:
84       Constructor for DataArrayView.       Constructor for DataArrayView.
85    
86       \param data - Input - Container holding data to be viewed. This must       \param data - Input -
87                             be large enough for the specified shape.                  Array holding data to be viewed. This must be large enough
88       \param viewShape - Input - The shape of the view.                  to supply sufficient values for the specified shape starting at
89       \param offset - Input - Offset into the data at which view should start.                  the given offset.
90         \param viewShape - Input -
91                    The shape of the view to create.
92         \param offset - Input -
93                    Offset into the data at which view should start.
94    */    */
95      ESCRIPT_DLL_API
96    DataArrayView(ValueType& data,    DataArrayView(ValueType& data,
97                  const ShapeType& viewShape,                  const ShapeType& viewShape,
98                  int offset=0);                  int offset=0);
# Line 79  class DataArrayView { Line 104  class DataArrayView {
104       Description:       Description:
105       Copy constructor for DataArrayView.       Copy constructor for DataArrayView.
106    
107       \param other - Input - DataArrayView to copy.       \param other - Input -
108                    DataArrayView to copy.
109    
110       NOTE: The copy references the same data container.       NOTE: The copy references the same data array.
111    */    */
112      ESCRIPT_DLL_API
113    DataArrayView(const DataArrayView& other);    DataArrayView(const DataArrayView& other);
114    
115    /**    /**
116       \brief       \brief
117       Check if view is empty.       Copy from a numarray into the data array viewed by this.
118    */       This must have same shape as the given value - will throw if not.
   bool  
   isEmpty() const;  
   
   /**  
      \brief  
      Copy from a numarray into the data managed by DataArrayView.  
119    */    */
120      ESCRIPT_DLL_API
121    void    void
122    copy(const boost::python::numeric::array& value);    copy(const boost::python::numeric::array& value);
123    
124    /**    /**
125       \brief       \brief
126       Copy from another DataArrayViewinto the data managed by this       Copy data from another DataArrayView into the data array viewed by this
127       DataArrayView at the default offset.       starting at the default offset in both views.
128         The shapes of both views must be the same - will throw if not.
129         NB: views may or may not reference same underlying data array!
130    */    */
131      ESCRIPT_DLL_API
132    void    void
133    copy(const DataArrayView& other);    copy(const DataArrayView& other);
134    
135    /**    /**
136       \brief       \brief
137       Copy from another DataArrayView into this view at the       Copy data from another DataArrayView into this view starting at the
138       given offset.       given offset in this view and the default offset in the other view.
139         The shapes of both views must be the same - will throw if not.
140         NB: views may or may not reference same underlying data array!
141    */    */
142      ESCRIPT_DLL_API
143    void    void
144    copy(ValueType::size_type offset,    copy(ValueType::size_type offset,
145         const DataArrayView& other);         const DataArrayView& other);
146    
147    /**    /**
148       \brief       \brief
149       Copy from another DataArrayView into this view at the       Copy data from another DataArrayView into this view starting at the
150       given offsets.       given offsets in each view.
151         The shapes of both views must be compatible - will throw if not.
152       \param thisOffset - Input - Offset into this view's data array to copy to.       NB: views may or may not reference same underlying data array!
153       \param other - Input - View for the data to copy.  
154       \param otherOffset - Input - Offset into the other view to copy data from.       \param thisOffset - Input -
155                       Offset into this view's data array to copy to.
156         \param other - Input -
157                       View to copy data from.
158         \param otherOffset - Input -
159                       Offset into the other view's data array to copy from.
160    */    */
161      ESCRIPT_DLL_API
162    void    void
163    copy(ValueType::size_type thisOffset,    copy(ValueType::size_type thisOffset,
164         const DataArrayView& other,         const DataArrayView& other,
# Line 132  class DataArrayView { Line 166  class DataArrayView {
166    
167    /**    /**
168       \brief       \brief
169       Copy the given single value over the entire view.       Copy the given single value over the entire view starting at the given
170         offset in the view's data array.
171    
172       \param offset - Input - Offset into this view's data array to copy to.       \param offset - Input -
173       \param value - Input - Value to copy.                     Offset into this view's data array to copy to.
174         \param value - Input -
175                       Value to copy.
176    */    */
177      ESCRIPT_DLL_API
178    void    void
179    copy(ValueType::size_type offset,    copy(ValueType::size_type offset,
180         ValueType::value_type value);         ValueType::value_type value);
181    
182    /**    /**
183         \brief
184         Check if view is empty. ie: does not point to any actual data.
185      */
186      ESCRIPT_DLL_API
187      bool
188      isEmpty() const;
189    
190      /**
191      \brief      \brief
192      Return this view's offset to the start of data.      Return this view's offset into the viewed data array.
193    */    */
194      ESCRIPT_DLL_API
195    ValueType::size_type    ValueType::size_type
196    getOffset() const;    getOffset() const;
197    
198    /**    /**
199      \brief       \brief
200       Return the rank of the data.       Set the offset into the data array at which this view is to start.
201         This could be used to step through the underlying data array by incrementing
202         the offset by noValues successively. Thus this view would provide a moving
203         window on the underlying data with the given shape.
204    */    */
205    int    ESCRIPT_DLL_API
206    getRank() const;    void
207      setOffset(ValueType::size_type offset);
208    
209    /**    /**
210       \brief       \brief
211       Return the shape of the data.       Increment the offset by the number of values in the shape, if possible. Thus
212         moving the view onto the next data point of the given shape in the underlying
213         data array.
214    */    */
215    const ShapeType&    ESCRIPT_DLL_API
216    getShape() const;    void
217      incrOffset();
218    
219    /**    /**
220       \brief       \brief
221       Calculate the number of values for the given shape.       Check the (given) offset will not result in two few values being available in
222         the underlying data array for this view's shape.
223    */    */
224    static int    ESCRIPT_DLL_API
225    noValues(const ShapeType& shape);    bool
226      checkOffset() const;
227    
228      ESCRIPT_DLL_API
229      bool
230      checkOffset(ValueType::size_type offset) const;
231    
232      /**
233        \brief
234         Return the rank of the shape of this view.
235      */
236      ESCRIPT_DLL_API
237      int
238      getRank() const;
239    
240    /**    /**
241       \brief       \brief
242       Return the number of values for the current view.       Return the number of values for the shape of this view.
243    */    */
244      ESCRIPT_DLL_API
245    int    int
246    noValues() const;    noValues() const;
247    
248    /**    /**
249       \brief       \brief
250       Return true if the shapes are the same.       Calculate the number of values for the given shape or region.
251         This is purely a utility method and has no bearing on this view.
252    */    */
253    bool    ESCRIPT_DLL_API
254    checkShape(const DataArrayView::ShapeType& other) const;    static
255      int
256      noValues(const ShapeType& shape);
257    
258      ESCRIPT_DLL_API
259      static
260      int
261      noValues(const RegionLoopRangeType& region);
262    
263    /**    /**
264       \brief       \brief
265       Return a reference to the underlying data.       Return a reference to the underlying data array.
266    */    */
267      ESCRIPT_DLL_API
268    ValueType&    ValueType&
269    getData() const;    getData() const;
270    
271    /**    /**
272       \brief       \brief
273       Return the value with the given index for the view. This takes into account       Return a reference to the data value with the given
274       the offset. Effectively this returns each value of the view in sequential order.       index in this view. This takes into account the offset.
275    */    */
276      ESCRIPT_DLL_API
277    ValueType::reference    ValueType::reference
278    getData(ValueType::size_type i) const;    getData(ValueType::size_type i) const;
279    
280    /**    /**
281       \brief       \brief
282       Create a shape error message. Normally used when there is a shape       Return the shape of this view.
      mismatch.  
   
      \param messagePrefix - Input - First part of the error message.  
      \param other - Input - The other shape.  
283    */    */
284    std::string    ESCRIPT_DLL_API
285    createShapeErrorMessage(const std::string& messagePrefix,    const
286                            const DataArrayView::ShapeType& other) const;    ShapeType&
287      getShape() const;
288    
289    /**    /**
290       \brief       \brief
291       Set the offset into the array.       Return true if the given shape is the same as this view's shape.
292    */    */
293    void    ESCRIPT_DLL_API
294    setOffset(ValueType::size_type offset);    bool
295      checkShape(const ShapeType& other) const;
296    
297    /**    /**
298       \brief       \brief
299       Return the shape of a data view following the slice specified.       Create a shape error message. Normally used when there is a shape
300         mismatch between this shape and the other shape.
301    
302       \param region - Input - Slice region.       \param messagePrefix - Input -
303                           First part of the error message.
304         \param other - Input -
305                           The other shape.
306    */    */
307    static DataArrayView::ShapeType    ESCRIPT_DLL_API
308    getResultSliceShape(const RegionType& region);    std::string
309      createShapeErrorMessage(const std::string& messagePrefix,
310                              const ShapeType& other) const;
311    
312    /**    /**
313       \brief       \brief
314       Copy a slice from the given view. This view must be the right shape for the slice.       Return the viewed data as a formatted string.
315         Not recommended for very large objects!
316    
317       \param other - Input - Data to copy from.       \param suffix - Input -
318       \param region - Input - Slice region.                         Suffix appended to index display.
319    */    */
320    void    ESCRIPT_DLL_API
321    copySlice(const DataArrayView& other,    std::string
322              const RegionType& region);    toString(const std::string& suffix=std::string("")) const;
323    
324    /**    /**
325       \brief       \brief
326       Copy a slice from the given view. This view must be the right shape for the slice.       Return the given shape as a string.
327         This is purely a utility method and has no bearing on this view.
328    
329       \param thisOffset - Input - use this view offset instead of the default.       \param shape - Input.
      \param other - Input - Data to copy from.  
      \param otherOffset - Input - use this slice offset instead of the default.  
      \param region - Input - Slice region.  
330    */    */
331    void    ESCRIPT_DLL_API
332    copySlice(ValueType::size_type thisOffset,    static
333              const DataArrayView& other,    std::string
334              ValueType::size_type otherOffset,    shapeToString(const ShapeType& shape);
             const RegionType& region);  
   
  /**  
      \brief  
      Copy into a slice from the given view. This view must have the same rank as the slice region.  
335    
336       \param other - Input - Data to copy from.    /**
337       \param region - Input - Slice region.      \brief
338        Return the 1 dimensional index into the data array of the only element
339        in the view, *ignoring the offset*.
340        Assumes a rank 0 view.
341    */    */
342    void    ESCRIPT_DLL_API
343    copySliceFrom(const DataArrayView& other,    ValueType::size_type
344                  const RegionType& region);    relIndex() const;
345    
346    /**    /**
347       \brief      \brief
348       Copy into a slice from the given value. This view must have the same rank as the slice region.      Return the 1 dimensional index into the data array of the element at
349        position i in the view, *ignoring the offset*.
350       \param thisOffset - Input - use this view offset instead of the default.      Assumes a rank 1 view.
      \param other - Input - Data to copy from.  
      \param otherOffset - Input - use this slice offset instead of the default.  
      \param region - Input - Slice region.  
351    */    */
352    void    ESCRIPT_DLL_API
353    copySliceFrom(ValueType::size_type thisOffset,    ValueType::size_type
354                  const DataArrayView& other,    relIndex(ValueType::size_type i) const;
                 ValueType::size_type otherOffset,  
                 const RegionType& region);  
355    
356    /**    /**
357       \brief      \brief
358       Determine the shape of the result array for a matrix multiplication.      Return the 1 dimensional index into the data array of the element at
359        position i,j in the view, *ignoring the offset*.
360        Assumes a rank 2 view.
361    */    */
362    static ShapeType    ESCRIPT_DLL_API
363    determineResultShape(const DataArrayView& left,    ValueType::size_type
364                         const DataArrayView& right);    relIndex(ValueType::size_type i,
365               ValueType::size_type j) const;
366    
367    /**    /**
368       \brief      \brief
369       Returns the range for the slices defined by the key which is a Python      Return the 1 dimensional index into the data array of the element at
370       slice object or a tuple of Python slice object.      position i,j,k in the view, *ignoring the offset*.
371    */      Assumes a rank 3 view.
   DataArrayView::RegionType  
   getSliceRegion(const boost::python::object& key) const;  
   /*  
   DataArrayView::RegionType  
   getSliceRegion2(const boost::python::object& key) const;  
372    */    */
373      ESCRIPT_DLL_API
374    // *******************************************************************    ValueType::size_type
375    // NOTE: The following relIndex functions are a hack. The indexing    relIndex(ValueType::size_type i,
376    // mechanism should be split out from DataArrayView to get the flexability             ValueType::size_type j,
377    // needed.             ValueType::size_type k) const;
378    
379    /**    /**
380      \brief      \brief
381      Return the 1 dimensional index of the element at position i,j,k,m      Return the 1 dimensional index into the data array of the element at
382      ignoring the offset.      position i,j,k,m in the view, *ignoring the offset*.
383        Assumes a rank 4 view.
384    */    */
385      ESCRIPT_DLL_API
386    ValueType::size_type    ValueType::size_type
387    relIndex(ValueType::size_type i,    relIndex(ValueType::size_type i,
388             ValueType::size_type j,             ValueType::size_type j,
# Line 316  class DataArrayView { Line 391  class DataArrayView {
391    
392    /**    /**
393      \brief      \brief
394      Return the 1 dimensional index of the element at position i,j,k      Return the 1 dimensional index into the data array of the only element
395      ignoring the offset.      in the view.
396        Assumes a rank 0 view.
397    */    */
398      ESCRIPT_DLL_API
399    ValueType::size_type    ValueType::size_type
400    relIndex(ValueType::size_type i,    index() const;
            ValueType::size_type j,  
            ValueType::size_type k) const;  
401    
402    /**    /**
403      \brief      \brief
404      Return the 1 dimensional index of the element at position i,j      Return the 1 dimensional index into the data array of the element at
405      ignoring the offset.      position i in the view.
406        Assumes a rank 1 view.
407    */    */
408      ESCRIPT_DLL_API
409    ValueType::size_type    ValueType::size_type
410    relIndex(ValueType::size_type i,    index(ValueType::size_type i) const;
            ValueType::size_type j) const;  
   
   // ********************************************************************  
411    
412    /**    /**
413      \brief      \brief
414      Return the 1 dimensional index of the element at position i,j,k,m.      Return the 1 dimensional index into the data array of the element at
415        position i,j in the view.
416        Assumes a rank 2 view.
417    */    */
418      ESCRIPT_DLL_API
419    ValueType::size_type    ValueType::size_type
420    index(ValueType::size_type i,    index(ValueType::size_type i,
421          ValueType::size_type j,          ValueType::size_type j) const;
         ValueType::size_type k,  
         ValueType::size_type m) const;  
422    
423    /**    /**
424      \brief      \brief
425      Return the 1 dimensional index of the element at position i,j,k.      Return the 1 dimensional index into the data array of the element at
426        position i,j,k in the view.
427        Assumes a rank 3 view.
428    */    */
429      ESCRIPT_DLL_API
430    ValueType::size_type    ValueType::size_type
431    index(ValueType::size_type i,    index(ValueType::size_type i,
432          ValueType::size_type j,          ValueType::size_type j,
# Line 356  class DataArrayView { Line 434  class DataArrayView {
434    
435    /**    /**
436      \brief      \brief
437      Return the 1 dimensional index of the element at position i,j.      Return the 1 dimensional index into the data array of the element at
438        position i,j,k,m in the view.
439        Assumes a rank 4 view.
440    */    */
441      ESCRIPT_DLL_API
442    ValueType::size_type    ValueType::size_type
443    index(ValueType::size_type i,    index(ValueType::size_type i,
444          ValueType::size_type j) const;          ValueType::size_type j,
445            ValueType::size_type k,
446            ValueType::size_type m) const;
447    
448    /**    /**
449      \brief      \brief
450      Return the 1 dimensional index of the element at position i.      Return a reference for the only element in the view.
451        Assumes a rank 0 view.
452    */    */
453    ValueType::size_type    ESCRIPT_DLL_API
454    index(ValueType::size_type i) const;    ValueType::reference
455      operator()();
456    
457      ESCRIPT_DLL_API
458      ValueType::const_reference
459      operator()() const;
460    
461    /**    /**
462      \brief      \brief
463      Return a reference to the element at position i.      Return a reference to the element at position i in the view.
464        Assumes a rank 1 view.
465    */    */
466      ESCRIPT_DLL_API
467    ValueType::reference    ValueType::reference
468    operator()(ValueType::size_type i);    operator()(ValueType::size_type i);
469    
470      ESCRIPT_DLL_API
471    ValueType::const_reference    ValueType::const_reference
472    operator()(ValueType::size_type i) const;    operator()(ValueType::size_type i) const;
473    
474    /**    /**
475      \brief      \brief
476      Return a reference to the element at position i,j.      Return a reference to the element at position i,j in the view.
477        Assumes a rank 2 view.
478    */    */
479      ESCRIPT_DLL_API
480    ValueType::reference    ValueType::reference
481    operator()(ValueType::size_type i,    operator()(ValueType::size_type i,
482               ValueType::size_type j);               ValueType::size_type j);
483    
484      ESCRIPT_DLL_API
485    ValueType::const_reference    ValueType::const_reference
486    operator()(ValueType::size_type i,    operator()(ValueType::size_type i,
487               ValueType::size_type j) const;               ValueType::size_type j) const;
488    
489    /**    /**
490      \brief      \brief
491      Return a reference to the element at position i,j,k.      Return a reference to the element at position i,j,k in the view.
492        Assumes a rank 3 view.
493    */    */
494      ESCRIPT_DLL_API
495    ValueType::reference    ValueType::reference
496    operator()(ValueType::size_type i,    operator()(ValueType::size_type i,
497               ValueType::size_type j,               ValueType::size_type j,
498               ValueType::size_type k);               ValueType::size_type k);
499    
500      ESCRIPT_DLL_API
501    ValueType::const_reference    ValueType::const_reference
502    operator()(ValueType::size_type i,    operator()(ValueType::size_type i,
503               ValueType::size_type j,               ValueType::size_type j,
# Line 407  class DataArrayView { Line 505  class DataArrayView {
505    
506   /**   /**
507      \brief      \brief
508      Return a reference to the element at position i,j,k,m.      Return a reference to the element at position i,j,k,m in the view.
509        Assumes a rank 4 view.
510    */    */
511      ESCRIPT_DLL_API
512    ValueType::reference    ValueType::reference
513    operator()(ValueType::size_type i,    operator()(ValueType::size_type i,
514               ValueType::size_type j,               ValueType::size_type j,
515               ValueType::size_type k,               ValueType::size_type k,
516               ValueType::size_type m);               ValueType::size_type m);
517    
518      ESCRIPT_DLL_API
519    ValueType::const_reference    ValueType::const_reference
520    operator()(ValueType::size_type i,    operator()(ValueType::size_type i,
521               ValueType::size_type j,               ValueType::size_type j,
# Line 422  class DataArrayView { Line 523  class DataArrayView {
523               ValueType::size_type m) const;               ValueType::size_type m) const;
524    
525    /**    /**
526      \brief       \brief
527      Return a reference for the only element, assuming rank 0.       Determine the shape of the specified slice region.
528         This is purely a utility method and has no bearing on this view.
529    
530         \param region - Input -
531                           Slice region.
532    */    */
533    ValueType::reference    ESCRIPT_DLL_API
534    operator()();    static
535      ShapeType
536      getResultSliceShape(const RegionType& region);
537    
538    ValueType::const_reference    /**
539    operator()() const;       \brief
540         Determine the region specified by the given python slice object.
541    
542         \param key - Input -
543                        python slice object specifying region to be returned.
544    
545         The slice object is a tuple of n python slice specifiers, where
546         n <= the rank of this Data object. Each slice specifier specifies the
547         range of indexes to be sliced from the corresponding dimension. The
548         first specifier corresponds to the first dimension, the second to the
549         second and so on. Where n < the rank, the remaining dimensions are
550         sliced across the full range of their indicies.
551    
552         Each slice specifier is of the form "a:b", which specifies a slice
553         from index a, up to but not including index b. Where index a is ommitted
554         a is assumed to be 0. Where index b is ommitted, b is assumed to be the
555         length of this dimension. Where both are ommitted (eg: ":") the slice is
556         assumed to encompass that entire dimension.
557    
558         Where one of the slice specifiers is a single integer, eg: [1], we
559         want to generate a rank-1 dimension object, as opposed to eg: [1,2]
560         which implies we want to take a rank dimensional object with one
561         dimension of size 1.
562    
563         The return value is a vector of pairs with length equal to the rank of
564         this object. Each pair corresponds to the range of indicies from the
565         corresponding dimension to be sliced from, as specified in the input
566         slice object.
567    
568         Examples:
569    
570           For a rank 1 object of shape(5):
571    
572             getSliceRegion(:)   => < <0,5> >
573             getSliceRegion(2:3) => < <2,3> >
574             getSliceRegion(:3)  => < <0,3> >
575             getSliceRegion(2:)  => < <2,5> >
576    
577           For a rank 2 object of shape(4,5):
578    
579             getSliceRegion(2:3) => < <2,3> <0,5> >
580             getSliceRegion(2)   => < <2,3> <0,5> >
581               NB: but return object requested will have rank 1, shape(5), with
582                   values taken from index 2 of this object's first dimension.
583    
584           For a rank 3 object of shape (2,4,6):
585    
586             getSliceRegion(0:2,0:4,0:6) => < <0,2> <0,4> <0,6> >
587             getSliceRegion(:,:,:)       => < <0,2> <0,4> <0,6> >
588             getSliceRegion(0:1)         => < <0,1> <0,4> <0,6> >
589             getSliceRegion(:1,0:2)      => < <0,1> <0,2> <0,6> >
590    
591      */
592      ESCRIPT_DLL_API
593      RegionType
594      getSliceRegion(const boost::python::object& key) const;
595    
596    /**    /**
597       \brief       \brief
598       Perform the unary operation using the given offset       Copy a data slice specified by the given region from the given view
599       instead of the offset defined within the view.       into this view, using the default offsets in both views.
600    
601         \param other - Input -
602                          View to copy data from.
603         \param region - Input -
604                          Region in other view to copy data from.
605    */    */
606    template <class UnaryFunction>    ESCRIPT_DLL_API
607    void    void
608    unaryOp(ValueType::size_type leftOffset,    copySlice(const DataArrayView& other,
609            UnaryFunction operation);              const RegionLoopRangeType& region);
610    
611    /**    /**
612       \brief       \brief
613       Perform the unary operation using the view's offset.       Copy a data slice specified by the given region and offset from the
614         given view into this view at the given offset.
615    
616         \param thisOffset - Input -
617                          Copy the slice to this offset in this view.
618         \param other - Input -
619                          View to copy data from.
620         \param otherOffset - Input -
621                          Copy the slice from this offset in the given view.
622         \param region - Input -
623                          Region in other view to copy data from.
624    */    */
625    template <class UnaryFunction>    ESCRIPT_DLL_API
626    void    void
627    unaryOp(UnaryFunction operation);    copySlice(ValueType::size_type thisOffset,
628                const DataArrayView& other,
629                ValueType::size_type otherOffset,
630                const RegionLoopRangeType& region);
631    
632    /**    /**
633       \brief       \brief
634       Perform the given operation and return a result.       Copy data into a slice specified by the given region in this view from
635         the given view, using the default offsets in both views.
636    
637         \param other - Input -
638                      View to copy data from.
639         \param region - Input -
640                      Region in this view to copy data to.
641    */    */
642    template <class UnaryFunction>    ESCRIPT_DLL_API
643    double    void
644    algorithm(ValueType::size_type leftOffset,    copySliceFrom(const DataArrayView& other,
645              UnaryFunction operation) const;                  const RegionLoopRangeType& region);
646    
647    /**    /**
648       \brief       \brief
649       Perform the given operation and return a result. Use the default offset.       Copy data into a slice specified by the given region and offset in
650         this view from the given view at the given offset.
651    
652         \param thisOffset - Input -
653                        Copy the slice to this offset in this view.
654         \param other - Input -
655                        View to copy data from.
656         \param otherOffset - Input -
657                        Copy the slice from this offset in the given view.
658         \param region - Input -
659                        Region in this view to copy data to.
660      */
661      ESCRIPT_DLL_API
662      void
663      copySliceFrom(ValueType::size_type thisOffset,
664                    const DataArrayView& other,
665                    ValueType::size_type otherOffset,
666                    const RegionLoopRangeType& region);
667    
668      /**
669         \brief
670         Perform the unary operation on the data point specified by the view's
671         default offset. Applies the specified operation to each value in the data
672         point.
673    
674         Called by escript::unaryOp.
675    
676         \param operation - Input -
677                      Operation to apply. Operation must be a pointer to a function.
678    */    */
679    template <class UnaryFunction>    template <class UnaryFunction>
680    double    ESCRIPT_DLL_API
681    algorithm(UnaryFunction operation) const;    void
682      unaryOp(UnaryFunction operation);
683    
684    /**    /**
685       \brief       \brief
686       Perform the binary operation. Version which applies a double value       Perform the unary operation on the data point specified by the given
687       to all values within the view. The given offset is used instead of       offset. Applies the specified operation to each value in the data
688       the default offset specified within the view.       point. Operation must be a pointer to a function.
689    
690         Called by escript::unaryOp.
691    
692         \param offset - Input -
693                      Apply the operation to data starting at this offset in this view.
694         \param operation - Input -
695                      Operation to apply. Must be a pointer to a function.
696    */    */
697    template <class BinaryFunction>    template <class UnaryFunction>
698      ESCRIPT_DLL_API
699    void    void
700    binaryOp(ValueType::size_type leftOffset,    unaryOp(ValueType::size_type offset,
701             double right,            UnaryFunction operation);
            BinaryFunction operation);  
702    
703    /**    /**
704       \brief       \brief
705       Perform the binary operation. Version which applies a double value       Perform the binary operation on the data points specified by the default
706       to all values within the view.       offsets in this view and in view "right". Applies the specified operation
707         to corresponding values in both data points. Operation must be a pointer
708         to a function.
709    
710         Called by escript::binaryOp.
711    
712         \param right - Input -
713                      View to act as RHS in given binary operation.
714         \param operation - Input -
715                      Operation to apply. Must be a pointer to a function.
716    */    */
717    template <class BinaryFunction>    template <class BinaryFunction>
718      ESCRIPT_DLL_API
719    void    void
720    binaryOp(double right,    binaryOp(const DataArrayView& right,
721             BinaryFunction operation);             BinaryFunction operation);
722    
723    /**    /**
724       \brief       \brief
725       Perform the binary operation. The given offsets override the default       Perform the binary operation on the data points specified by the given
726       offsets specified within the views.       offsets in this view and in view "right". Applies the specified operation
727         to corresponding values in both data points. Operation must be a pointer
728         to a function.
729    
730         Called by escript::binaryOp.
731    
732         \param leftOffset - Input -
733                      Apply the operation to data starting at this offset in this view.
734         \param right - Input -
735                      View to act as RHS in given binary operation.
736         \param rightOffset - Input -
737                      Apply the operation to data starting at this offset in right.
738         \param operation - Input -
739                      Operation to apply. Must be a pointer to a function.
740    */    */
741    template <class BinaryFunction>    template <class BinaryFunction>
742      ESCRIPT_DLL_API
743    void    void
744    binaryOp(ValueType::size_type leftOffset,    binaryOp(ValueType::size_type leftOffset,
745             const DataArrayView& right,             const DataArrayView& right,
# Line 502  class DataArrayView { Line 748  class DataArrayView {
748    
749    /**    /**
750       \brief       \brief
751       Perform the binary operation.       Perform the binary operation on the data point specified by the default
752         offset in this view using the scalar value "right". Applies the specified
753         operation to values in the data point. Operation must be a pointer to
754         a function.
755    
756         Called by escript::binaryOp.
757    
758         \param right - Input -
759                      Value to act as RHS in given binary operation.
760         \param operation - Input -
761                      Operation to apply. Must be a pointer to a function.
762    */    */
763    template <class BinaryFunction>    template <class BinaryFunction>
764      ESCRIPT_DLL_API
765    void    void
766    binaryOp(const DataArrayView& right,    binaryOp(double right,
767             BinaryFunction operation);             BinaryFunction operation);
768    
769    /**    /**
770       \brief       \brief
771       Return the data as a string. Not recommended for very large objects.       Perform the binary operation on the data point specified by the given
772       \param suffix - Input - Suffix appended to index display.       offset in this view using the scalar value "right". Applies the specified
773         operation to values in the data point. Operation must be a pointer
774         to a function.
775    
776         Called by escript::binaryOp.
777    
778         \param offset - Input -
779                      Apply the operation to data starting at this offset in this view.
780         \param right - Input -
781                      Value to act as RHS in given binary operation.
782         \param operation - Input -
783                      Operation to apply. Must be a pointer to a function.
784    */    */
785    std::string    template <class BinaryFunction>
786    toString(const std::string& suffix=std::string("")) const;    ESCRIPT_DLL_API
787      void
788      binaryOp(ValueType::size_type offset,
789               double right,
790               BinaryFunction operation);
791    
792    /**    /**
793       \brief       \brief
794       Return the given shape as a string.       Perform the given data point reduction operation on the data point
795         specified by the default offset into the view. Reduces all elements of
796         the data point using the given operation, returning the result as a
797         scalar. Operation must be a pointer to a function.
798    
799       \param shape - Input.       Called by escript::algorithm.
800    
801         \param operation - Input -
802                      Operation to apply. Must be a pointer to a function.
803    */    */
804    static std::string    template <class BinaryFunction>
805    shapeToString(const DataArrayView::ShapeType& shape);    ESCRIPT_DLL_API
806      double
807      reductionOp(BinaryFunction operation,
808                  double initial_value) const;
809    
810    /**    /**
811       \brief       \brief
812       Perform matrix multiply.       Perform the given data point reduction operation on the data point
813         specified by the given offset into the view. Reduces all elements of
814         the data point using the given operation, returning the result as a
815         scalar. Operation must be a pointer to a function.
816    
817       Description:       Called by escript::algorithm.
818       Perform matrix multiply.  
819         \param offset - Input -
820                      Apply the operation to data starting at this offset in this view.
821         \param operation - Input -
822                      Operation to apply. Must be a pointer to a function.
823      */
824      template <class BinaryFunction>
825      ESCRIPT_DLL_API
826      double
827      reductionOp(ValueType::size_type offset,
828                  BinaryFunction operation,
829                  double initial_value) const;
830    
831     /**
832         \brief
833         Perform a matrix multiply of the given views.
834         This is purely a utility method and has no bearing on this view.
835    
836         NB: Only multiplies together the two given views, ie: two data-points,
837         would need to call this over all data-points to multiply the entire
838         Data objects involved.
839    
840       \param left - Input - The left hand side.       \param left - Input - The left hand side.
841       \param right - Input - The right hand side.       \param right - Input - The right hand side.
842       \param result - Output - The result of the operation.       \param result - Output - The result of the operation.
843    */    */
844    static void    ESCRIPT_DLL_API
845      static
846      void
847    matMult(const DataArrayView& left,    matMult(const DataArrayView& left,
848            const DataArrayView& right,            const DataArrayView& right,
849            DataArrayView& result);            DataArrayView& result);
850    
851      /**
852         \brief
853         Determine the shape of the result array for a matrix multiplication
854         of the given views.
855         This is purely a utility method and has no bearing on this view.
856      */
857      ESCRIPT_DLL_API
858      static
859      ShapeType
860      determineResultShape(const DataArrayView& left,
861                           const DataArrayView& right);
862    
863      /**
864         \brief
865         computes a symmetric matrix from your square matrix A: (A + transpose(A)) / 2
866    
867         \param in - Input - matrix
868         \param inOffset - Input - offset into in
869         \param ev - Output - The symmetric matrix
870         \param inOffset - Input - offset into ev
871      */
872      ESCRIPT_DLL_API
873      static
874      inline
875      void
876      symmetric(DataArrayView& in,
877                ValueType::size_type inOffset,
878                DataArrayView& ev,
879                ValueType::size_type evOffset)
880      {
881       if (in.getRank() == 2) {
882         int i0, i1;
883         int s0=in.getShape()[0];
884         int s1=in.getShape()[1];
885         for (i0=0; i0<s0; i0++) {
886           for (i1=0; i1<s1; i1++) {
887             (*(ev.m_data))[evOffset+ev.index(i0,i1)] = ((*(in.m_data))[inOffset+in.index(i0,i1)] + (*(in.m_data))[inOffset+in.index(i1,i0)]) / 2.0;
888           }
889         }
890       }
891       if (in.getRank() == 4) {
892         int i0, i1, i2, i3;
893         int s0=in.getShape()[0];
894         int s1=in.getShape()[1];
895         int s2=in.getShape()[2];
896         int s3=in.getShape()[3];
897         for (i0=0; i0<s0; i0++) {
898           for (i1=0; i1<s1; i1++) {
899             for (i2=0; i2<s2; i2++) {
900               for (i3=0; i3<s3; i3++) {
901                 (*(ev.m_data))[evOffset+ev.index(i0,i1,i2,i3)] = ((*(in.m_data))[inOffset+in.index(i0,i1,i2,i3)] + (*(in.m_data))[inOffset+in.index(i2,i3,i0,i1)]) / 2.0;
902               }
903             }
904           }
905         }
906       }
907      }
908    
909      /**
910         \brief
911         computes a nonsymmetric matrix from your square matrix A: (A - transpose(A)) / 2
912    
913         \param in - Input - matrix
914         \param inOffset - Input - offset into in
915         \param ev - Output - The nonsymmetric matrix
916         \param inOffset - Input - offset into ev
917      */
918      ESCRIPT_DLL_API
919      static
920      inline
921      void
922      nonsymmetric(DataArrayView& in,
923                ValueType::size_type inOffset,
924                DataArrayView& ev,
925                ValueType::size_type evOffset)
926      {
927       if (in.getRank() == 2) {
928         int i0, i1;
929         int s0=in.getShape()[0];
930         int s1=in.getShape()[1];
931         for (i0=0; i0<s0; i0++) {
932           for (i1=0; i1<s1; i1++) {
933             (*(ev.m_data))[evOffset+ev.index(i0,i1)] = ((*(in.m_data))[inOffset+in.index(i0,i1)] - (*(in.m_data))[inOffset+in.index(i1,i0)]) / 2.0;
934           }
935         }
936       }
937       if (in.getRank() == 4) {
938         int i0, i1, i2, i3;
939         int s0=in.getShape()[0];
940         int s1=in.getShape()[1];
941         int s2=in.getShape()[2];
942         int s3=in.getShape()[3];
943         for (i0=0; i0<s0; i0++) {
944           for (i1=0; i1<s1; i1++) {
945             for (i2=0; i2<s2; i2++) {
946               for (i3=0; i3<s3; i3++) {
947                 (*(ev.m_data))[evOffset+ev.index(i0,i1,i2,i3)] = ((*(in.m_data))[inOffset+in.index(i0,i1,i2,i3)] - (*(in.m_data))[inOffset+in.index(i2,i3,i0,i1)]) / 2.0;
948               }
949             }
950           }
951         }
952       }
953      }
954    
955      /**
956         \brief
957         computes the trace of a matrix
958    
959         \param in - Input - matrix
960         \param inOffset - Input - offset into in
961         \param ev - Output - The nonsymmetric matrix
962         \param inOffset - Input - offset into ev
963      */
964      static
965      inline
966      void
967      trace(DataArrayView& in,
968                ValueType::size_type inOffset,
969                DataArrayView& ev,
970                ValueType::size_type evOffset,
971            int axis_offset)
972      {
973       if (in.getRank() == 2) {
974         int s0=in.getShape()[0]; // Python wrapper limits to square matrix
975         int i;
976         for (i=0; i<s0; i++) {
977           (*(ev.m_data))[evOffset+ev.index()] += (*(in.m_data))[inOffset+in.index(i,i)];
978         }
979       }
980       else if (in.getRank() == 3) {
981         if (axis_offset==0) {
982           int s0=in.getShape()[0];
983           int s2=in.getShape()[2];
984           int i0, i2;
985           for (i0=0; i0<s0; i0++) {
986             for (i2=0; i2<s2; i2++) {
987               (*(ev.m_data))[evOffset+ev.index(i2)] += (*(in.m_data))[inOffset+in.index(i0,i0,i2)];
988             }
989           }
990         }
991         else if (axis_offset==1) {
992           int s0=in.getShape()[0];
993           int s1=in.getShape()[1];
994           int i0, i1;
995           for (i0=0; i0<s0; i0++) {
996             for (i1=0; i1<s1; i1++) {
997               (*(ev.m_data))[evOffset+ev.index(i0)] += (*(in.m_data))[inOffset+in.index(i0,i1,i1)];
998             }
999           }
1000         }
1001       }
1002       else if (in.getRank() == 4) {
1003         if (axis_offset==0) {
1004           int s0=in.getShape()[0];
1005           int s2=in.getShape()[2];
1006           int s3=in.getShape()[3];
1007           int i0, i2, i3;
1008           for (i0=0; i0<s0; i0++) {
1009             for (i2=0; i2<s2; i2++) {
1010               for (i3=0; i3<s3; i3++) {
1011                 (*(ev.m_data))[evOffset+ev.index(i2,i3)] += (*(in.m_data))[inOffset+in.index(i0,i0,i2,i3)];
1012               }
1013             }
1014           }
1015         }
1016         else if (axis_offset==1) {
1017           int s0=in.getShape()[0];
1018           int s1=in.getShape()[1];
1019           int s3=in.getShape()[3];
1020           int i0, i1, i3;
1021           for (i0=0; i0<s0; i0++) {
1022             for (i1=0; i1<s1; i1++) {
1023               for (i3=0; i3<s3; i3++) {
1024                 (*(ev.m_data))[evOffset+ev.index(i0,i3)] += (*(in.m_data))[inOffset+in.index(i0,i1,i1,i3)];
1025               }
1026             }
1027           }
1028         }
1029         else if (axis_offset==2) {
1030           int s0=in.getShape()[0];
1031           int s1=in.getShape()[1];
1032           int s2=in.getShape()[2];
1033           int i0, i1, i2;
1034           for (i0=0; i0<s0; i0++) {
1035             for (i1=0; i1<s1; i1++) {
1036               for (i2=0; i2<s2; i2++) {
1037                 (*(ev.m_data))[evOffset+ev.index(i0,i1)] += (*(in.m_data))[inOffset+in.index(i0,i1,i2,i2)];
1038               }
1039             }
1040           }
1041         }
1042       }
1043      }
1044    
1045      /**
1046         \brief
1047         Transpose each data point of this Data object around the given axis.
1048    
1049         \param in - Input - matrix
1050         \param inOffset - Input - offset into in
1051         \param ev - Output - The nonsymmetric matrix
1052         \param inOffset - Input - offset into ev
1053      */
1054      ESCRIPT_DLL_API
1055      static
1056      inline
1057      void
1058      transpose(DataArrayView& in,
1059                ValueType::size_type inOffset,
1060                DataArrayView& ev,
1061                ValueType::size_type evOffset,
1062            int axis_offset)
1063      {
1064       if (in.getRank() == 4) {
1065         int s0=ev.getShape()[0];
1066         int s1=ev.getShape()[1];
1067         int s2=ev.getShape()[2];
1068         int s3=ev.getShape()[3];
1069         int i0, i1, i2, i3;
1070         if (axis_offset==1) {
1071           for (i0=0; i0<s0; i0++) {
1072             for (i1=0; i1<s1; i1++) {
1073               for (i2=0; i2<s2; i2++) {
1074                 for (i3=0; i3<s3; i3++) {
1075                   (*(ev.m_data))[evOffset+ev.index(i0,i1,i2,i3)] = (*(in.m_data))[inOffset+in.index(i3,i0,i1,i2)];
1076                 }
1077               }
1078             }
1079           }
1080         }
1081         else if (axis_offset==2) {
1082           for (i0=0; i0<s0; i0++) {
1083             for (i1=0; i1<s1; i1++) {
1084               for (i2=0; i2<s2; i2++) {
1085                 for (i3=0; i3<s3; i3++) {
1086                   (*(ev.m_data))[evOffset+ev.index(i0,i1,i2,i3)] = (*(in.m_data))[inOffset+in.index(i2,i3,i0,i1)];
1087                 }
1088               }
1089             }
1090           }
1091         }
1092         else if (axis_offset==3) {
1093           for (i0=0; i0<s0; i0++) {
1094             for (i1=0; i1<s1; i1++) {
1095               for (i2=0; i2<s2; i2++) {
1096                 for (i3=0; i3<s3; i3++) {
1097                   (*(ev.m_data))[evOffset+ev.index(i0,i1,i2,i3)] = (*(in.m_data))[inOffset+in.index(i1,i2,i3,i0)];
1098                 }
1099               }
1100             }
1101           }
1102         }
1103         else {
1104           for (i0=0; i0<s0; i0++) {
1105             for (i1=0; i1<s1; i1++) {
1106               for (i2=0; i2<s2; i2++) {
1107                 for (i3=0; i3<s3; i3++) {
1108                   (*(ev.m_data))[evOffset+ev.index(i0,i1,i2,i3)] = (*(in.m_data))[inOffset+in.index(i0,i1,i2,i3)];
1109                 }
1110               }
1111             }
1112           }
1113         }
1114       }
1115       else if (in.getRank() == 3) {
1116         int s0=ev.getShape()[0];
1117         int s1=ev.getShape()[1];
1118         int s2=ev.getShape()[2];
1119         int i0, i1, i2;
1120         if (axis_offset==1) {
1121           for (i0=0; i0<s0; i0++) {
1122             for (i1=0; i1<s1; i1++) {
1123               for (i2=0; i2<s2; i2++) {
1124                 (*(ev.m_data))[evOffset+ev.index(i0,i1,i2)] = (*(in.m_data))[inOffset+in.index(i2,i0,i1)];
1125               }
1126             }
1127           }
1128         }
1129         else if (axis_offset==2) {
1130           for (i0=0; i0<s0; i0++) {
1131             for (i1=0; i1<s1; i1++) {
1132               for (i2=0; i2<s2; i2++) {
1133                 (*(ev.m_data))[evOffset+ev.index(i0,i1,i2)] = (*(in.m_data))[inOffset+in.index(i1,i2,i0)];
1134               }
1135             }
1136           }
1137         }
1138         else {
1139           // Copy the matrix unchanged
1140           for (i0=0; i0<s0; i0++) {
1141             for (i1=0; i1<s1; i1++) {
1142               for (i2=0; i2<s2; i2++) {
1143                 (*(ev.m_data))[evOffset+ev.index(i0,i1,i2)] = (*(in.m_data))[inOffset+in.index(i0,i1,i2)];
1144               }
1145             }
1146           }
1147         }
1148       }
1149       else if (in.getRank() == 2) {
1150         int s0=ev.getShape()[0];
1151         int s1=ev.getShape()[1];
1152         int i0, i1;
1153         if (axis_offset==1) {
1154           for (i0=0; i0<s0; i0++) {
1155             for (i1=0; i1<s1; i1++) {
1156               (*(ev.m_data))[evOffset+ev.index(i0,i1)] = (*(in.m_data))[inOffset+in.index(i1,i0)];
1157             }
1158           }
1159         }
1160         else {
1161           for (i0=0; i0<s0; i0++) {
1162             for (i1=0; i1<s1; i1++) {
1163               (*(ev.m_data))[evOffset+ev.index(i0,i1)] = (*(in.m_data))[inOffset+in.index(i0,i1)];
1164             }
1165           }
1166         }
1167       }
1168       else if (in.getRank() == 1) {
1169         int s0=ev.getShape()[0];
1170         int i0;
1171         for (i0=0; i0<s0; i0++) {
1172           (*(ev.m_data))[evOffset+ev.index(i0)] = (*(in.m_data))[inOffset+in.index(i0)];
1173         }
1174       }
1175       else if (in.getRank() == 0) {
1176         (*(ev.m_data))[evOffset+ev.index()] = (*(in.m_data))[inOffset+in.index()];
1177       }
1178       else {
1179          throw DataException("Error - DataArrayView::transpose can only be calculated for rank 0, 1, 2, 3 or 4 objects.");
1180       }
1181      }
1182    
1183      /**
1184         \brief
1185         swaps the components axis_offset and axis_offset+1
1186    
1187         \param in - Input - matrix
1188         \param inOffset - Input - offset into in
1189         \param ev - Output - The nonsymmetric matrix
1190         \param inOffset - Input - offset into ev
1191      */
1192      ESCRIPT_DLL_API
1193      static
1194      inline
1195      void
1196      swap(DataArrayView& in,
1197           ValueType::size_type inOffset,
1198           DataArrayView& ev,
1199           ValueType::size_type evOffset,
1200           int axis_offset)
1201      {
1202       if (in.getRank() == 4) {
1203         int s0=ev.getShape()[0];
1204         int s1=ev.getShape()[1];
1205         int s2=ev.getShape()[2];
1206         int s3=ev.getShape()[3];
1207         int i0, i1, i2, i3;
1208         if (axis_offset==0) {
1209           for (i0=0; i0<s0; i0++) {
1210             for (i1=0; i1<s1; i1++) {
1211               for (i2=0; i2<s2; i2++) {
1212                 for (i3=0; i3<s3; i3++) {
1213                   (*(ev.m_data))[evOffset+ev.index(i0,i1,i2,i3)] = (*(in.m_data))[inOffset+in.index(i1,i0,i2,i3)];
1214                 }
1215               }
1216             }
1217           }
1218         }
1219         else if (axis_offset==1) {
1220           for (i0=0; i0<s0; i0++) {
1221             for (i1=0; i1<s1; i1++) {
1222               for (i2=0; i2<s2; i2++) {
1223                 for (i3=0; i3<s3; i3++) {
1224                   (*(ev.m_data))[evOffset+ev.index(i0,i1,i2,i3)] = (*(in.m_data))[inOffset+in.index(i0,i2,i1,i3)];
1225                 }
1226               }
1227             }
1228           }
1229         }
1230         else {
1231           for (i0=0; i0<s0; i0++) {
1232             for (i1=0; i1<s1; i1++) {
1233               for (i2=0; i2<s2; i2++) {
1234                 for (i3=0; i3<s3; i3++) {
1235                   (*(ev.m_data))[evOffset+ev.index(i0,i1,i2,i3)] = (*(in.m_data))[inOffset+in.index(i0,i1,i3,i2)];
1236                 }
1237               }
1238             }
1239           }
1240         }
1241       }
1242       else if (in.getRank() == 3) {
1243         int s0=ev.getShape()[0];
1244         int s1=ev.getShape()[1];
1245         int s2=ev.getShape()[2];
1246         int i0, i1, i2;
1247         if (axis_offset==0) {
1248           for (i0=0; i0<s0; i0++) {
1249             for (i1=0; i1<s1; i1++) {
1250               for (i2=0; i2<s2; i2++) {
1251                 (*(ev.m_data))[evOffset+ev.index(i0,i1,i2)] = (*(in.m_data))[inOffset+in.index(i1,i0,i2)];
1252               }
1253             }
1254           }
1255         }
1256         else {
1257           for (i0=0; i0<s0; i0++) {
1258             for (i1=0; i1<s1; i1++) {
1259               for (i2=0; i2<s2; i2++) {
1260                 (*(ev.m_data))[evOffset+ev.index(i0,i1,i2)] = (*(in.m_data))[inOffset+in.index(i0,i2,i1)];
1261               }
1262             }
1263           }
1264         }
1265       }
1266       else if (in.getRank() == 2) {
1267         int s0=ev.getShape()[0];
1268         int s1=ev.getShape()[1];
1269         int i0, i1;
1270         for (i0=0; i0<s0; i0++) {
1271             for (i1=0; i1<s1; i1++) {
1272               (*(ev.m_data))[evOffset+ev.index(i0,i1)] = (*(in.m_data))[inOffset+in.index(i1,i0)];
1273             }
1274          }
1275       }
1276       else if (in.getRank() == 1) {
1277         int s0=ev.getShape()[0];
1278         int i0;
1279         for (i0=0; i0<s0; i0++) {
1280           (*(ev.m_data))[evOffset+ev.index(i0)] = (*(in.m_data))[inOffset+in.index(i0)];
1281         }
1282       }
1283       else if (in.getRank() == 0) {
1284         (*(ev.m_data))[evOffset+ev.index()] = (*(in.m_data))[inOffset+in.index()];
1285       }
1286       else {
1287          throw DataException("Error - DataArrayView::transpose can only be calculated for rank 0, 1, 2, 3 or 4 objects.");
1288       }
1289      }
1290    
1291      /**
1292         \brief
1293         solves a local eigenvalue problem
1294    
1295         \param in - Input - matrix
1296         \param inOffset - Input - offset into in
1297         \param ev - Output - The eigenvalues
1298         \param inOffset - Input - offset into ev
1299      */
1300      ESCRIPT_DLL_API
1301      static
1302      inline
1303      void
1304      eigenvalues(DataArrayView& in,
1305                  ValueType::size_type inOffset,
1306                  DataArrayView& ev,
1307                  ValueType::size_type evOffset)
1308      {
1309       double in00,in10,in20,in01,in11,in21,in02,in12,in22;
1310       double ev0,ev1,ev2;
1311       int s=in.getShape()[0];
1312       if (s==1) {
1313          in00=(*(in.m_data))[inOffset+in.index(0,0)];
1314          eigenvalues1(in00,&ev0);
1315          (*(ev.m_data))[evOffset+ev.index(0)]=ev0;
1316    
1317       } else  if (s==2) {
1318          in00=(*(in.m_data))[inOffset+in.index(0,0)];
1319          in10=(*(in.m_data))[inOffset+in.index(1,0)];
1320          in01=(*(in.m_data))[inOffset+in.index(0,1)];
1321          in11=(*(in.m_data))[inOffset+in.index(1,1)];
1322          eigenvalues2(in00,(in01+in10)/2.,in11,&ev0,&ev1);
1323          (*(ev.m_data))[evOffset+ev.index(0)]=ev0;
1324          (*(ev.m_data))[evOffset+ev.index(1)]=ev1;
1325    
1326       } else  if (s==3) {
1327          in00=(*(in.m_data))[inOffset+in.index(0,0)];
1328          in10=(*(in.m_data))[inOffset+in.index(1,0)];
1329          in20=(*(in.m_data))[inOffset+in.index(2,0)];
1330          in01=(*(in.m_data))[inOffset+in.index(0,1)];
1331          in11=(*(in.m_data))[inOffset+in.index(1,1)];
1332          in21=(*(in.m_data))[inOffset+in.index(2,1)];
1333          in02=(*(in.m_data))[inOffset+in.index(0,2)];
1334          in12=(*(in.m_data))[inOffset+in.index(1,2)];
1335          in22=(*(in.m_data))[inOffset+in.index(2,2)];
1336          eigenvalues3(in00,(in01+in10)/2.,(in02+in20)/2.,in11,(in21+in12)/2.,in22,
1337                     &ev0,&ev1,&ev2);
1338          (*(ev.m_data))[evOffset+ev.index(0)]=ev0;
1339          (*(ev.m_data))[evOffset+ev.index(1)]=ev1;
1340          (*(ev.m_data))[evOffset+ev.index(2)]=ev2;
1341    
1342       }
1343      }
1344    
1345      /**
1346         \brief
1347         solves a local eigenvalue problem
1348    
1349         \param in - Input - matrix
1350         \param inOffset - Input - offset into in
1351         \param ev - Output - The eigenvalues
1352         \param evOffset - Input - offset into ev
1353         \param V - Output - The eigenvectors
1354         \param VOffset - Input - offset into V
1355         \param tol - Input - eigenvalues with relative difference tol are treated as equal
1356      */
1357      ESCRIPT_DLL_API
1358      static
1359      inline
1360      void
1361      eigenvalues_and_eigenvectors(DataArrayView& in,
1362                                   ValueType::size_type inOffset,
1363                                   DataArrayView& ev,
1364                                   ValueType::size_type evOffset,
1365                                   DataArrayView& V,
1366                                   ValueType::size_type VOffset,
1367                                   const double tol=1.e-13)
1368      {
1369       double in00,in10,in20,in01,in11,in21,in02,in12,in22;
1370       double V00,V10,V20,V01,V11,V21,V02,V12,V22;
1371       double ev0,ev1,ev2;
1372       int s=in.getShape()[0];
1373       if (s==1) {
1374          in00=(*(in.m_data))[inOffset+in.index(0,0)];
1375          eigenvalues_and_eigenvectors1(in00,&ev0,&V00,tol);
1376          (*(ev.m_data))[evOffset+ev.index(0)]=ev0;
1377          (*(V.m_data))[inOffset+V.index(0,0)]=V00;
1378       } else  if (s==2) {
1379          in00=(*(in.m_data))[inOffset+in.index(0,0)];
1380          in10=(*(in.m_data))[inOffset+in.index(1,0)];
1381          in01=(*(in.m_data))[inOffset+in.index(0,1)];
1382          in11=(*(in.m_data))[inOffset+in.index(1,1)];
1383          eigenvalues_and_eigenvectors2(in00,(in01+in10)/2.,in11,
1384                       &ev0,&ev1,&V00,&V10,&V01,&V11,tol);
1385          (*(ev.m_data))[evOffset+ev.index(0)]=ev0;
1386          (*(ev.m_data))[evOffset+ev.index(1)]=ev1;
1387          (*(V.m_data))[inOffset+V.index(0,0)]=V00;
1388          (*(V.m_data))[inOffset+V.index(1,0)]=V10;
1389          (*(V.m_data))[inOffset+V.index(0,1)]=V01;
1390          (*(V.m_data))[inOffset+V.index(1,1)]=V11;
1391       } else  if (s==3) {
1392          in00=(*(in.m_data))[inOffset+in.index(0,0)];
1393          in10=(*(in.m_data))[inOffset+in.index(1,0)];
1394          in20=(*(in.m_data))[inOffset+in.index(2,0)];
1395          in01=(*(in.m_data))[inOffset+in.index(0,1)];
1396          in11=(*(in.m_data))[inOffset+in.index(1,1)];
1397          in21=(*(in.m_data))[inOffset+in.index(2,1)];
1398          in02=(*(in.m_data))[inOffset+in.index(0,2)];
1399          in12=(*(in.m_data))[inOffset+in.index(1,2)];
1400          in22=(*(in.m_data))[inOffset+in.index(2,2)];
1401          eigenvalues_and_eigenvectors3(in00,(in01+in10)/2.,(in02+in20)/2.,in11,(in21+in12)/2.,in22,
1402                     &ev0,&ev1,&ev2,
1403                     &V00,&V10,&V20,&V01,&V11,&V21,&V02,&V12,&V22,tol);
1404          (*(ev.m_data))[evOffset+ev.index(0)]=ev0;
1405          (*(ev.m_data))[evOffset+ev.index(1)]=ev1;
1406          (*(ev.m_data))[evOffset+ev.index(2)]=ev2;
1407          (*(V.m_data))[inOffset+V.index(0,0)]=V00;
1408          (*(V.m_data))[inOffset+V.index(1,0)]=V10;
1409          (*(V.m_data))[inOffset+V.index(2,0)]=V20;
1410          (*(V.m_data))[inOffset+V.index(0,1)]=V01;
1411          (*(V.m_data))[inOffset+V.index(1,1)]=V11;
1412          (*(V.m_data))[inOffset+V.index(2,1)]=V21;
1413          (*(V.m_data))[inOffset+V.index(0,2)]=V02;
1414          (*(V.m_data))[inOffset+V.index(1,2)]=V12;
1415          (*(V.m_data))[inOffset+V.index(2,2)]=V22;
1416    
1417       }
1418     }
1419   protected:   protected:
1420    
1421   private:   private:
1422    
1423    //    //
1424      // The maximum rank allowed for the shape of any view.
1425    static const int m_maxRank=4;    static const int m_maxRank=4;
1426    
1427    //    //
1428    // The data values for the view. NOTE: This points to data external to the view.    // The data values for the view.
1429      // NOTE: This points to data external to the view.
1430      // This is just a pointer to an array of ValueType.
1431    ValueType* m_data;    ValueType* m_data;
1432    
1433    //    //
1434    // The offset into the data array used by different views.    // The offset into the data array used by different views.
1435      // This is simply an integer specifying a position in the data array
1436      // pointed to by m_data.
1437    ValueType::size_type m_offset;    ValueType::size_type m_offset;
1438    
1439    //    //
1440    // The shape of the data.    // The shape of the data.
1441      // This is simply an STL vector specifying the lengths of each dimension
1442      // of the shape as ints.
1443    ShapeType m_shape;    ShapeType m_shape;
1444    
1445    //    //
1446    // The number of values needed for the array.    // The number of values needed for the array.
1447      // This can be derived from m_shape by multiplying the size of each dimension, but
1448      // is stored here for convenience.
1449    int m_noValues;    int m_noValues;
1450    
1451  };  };
1452    
1453  inline  ESCRIPT_DLL_API bool operator==(const DataArrayView& left, const DataArrayView& right);
1454  DataArrayView::ValueType::size_type  ESCRIPT_DLL_API bool operator!=(const DataArrayView& left, const DataArrayView& right);
 DataArrayView::getOffset() const  
 {  
   return m_offset;  
 }  
1455    
1456  inline  /**
1457  DataArrayView::ValueType&    \brief
1458  DataArrayView::getData() const     Modify region to copy from in order to
1459  {     deal with the case where one range in the region contains identical indexes,
1460    EsysAssert(!isEmpty(),"Error - View is empty");     eg: <<1,1><0,3><0,3>>
1461    return *m_data;     This situation implies we want to copy from an object with rank greater than that of this
1462  }     object. eg: we want to copy the values from a two dimensional slice out of a three
1463       dimensional object into a two dimensional object.
1464       We do this by taking a slice from the other object where one dimension of
1465       the slice region is of size 1. So in the above example, we modify the above
1466       region like so: <<1,2><0,3><0,3>> and take this slice.
1467    */
1468    DataArrayView::RegionLoopRangeType
1469    getSliceRegionLoopRange(const DataArrayView::RegionType& region);
1470    
1471  inline  /**
1472  DataArrayView::ValueType::reference    \brief
1473  DataArrayView::getData(ValueType::size_type i) const    Calculate the slice range from the given python key object
1474  {    Used by DataArrayView::getSliceRegion.
1475    //    Returns the python slice object key as a pair of ints where the first
1476    // don't do any checking to allow one past the end of the vector providing    member is the start and the second member is the end. the presence of a possible
1477    // the equivalent of end()    step attribute with value different from one will throw an exception
   return (*m_data)[i+m_offset];  
 }  
1478    
1479  template <class UnaryFunction>    /param key - Input - key object specifying slice range.
1480  inline  */
1481  void  std::pair<int,int>
1482  DataArrayView::unaryOp(ValueType::size_type leftOffset,  getSliceRange(const boost::python::object& key,
1483                         UnaryFunction operation)                const int shape);
1484  {  
1485    for (ValueType::size_type i=0;i<(noValues(m_shape));++i) {  /**
1486      (*m_data)[i+leftOffset]=operation((*m_data)[i+leftOffset]);     Inline function definitions.
1487    }  */
 }  
1488    
1489  template <class UnaryFunction>  template <class UnaryFunction>
1490  inline  inline
# Line 613  DataArrayView::unaryOp(UnaryFunction ope Line 1496  DataArrayView::unaryOp(UnaryFunction ope
1496    
1497  template <class UnaryFunction>  template <class UnaryFunction>
1498  inline  inline
1499  double  void
1500  DataArrayView::algorithm(ValueType::size_type leftOffset,  DataArrayView::unaryOp(ValueType::size_type offset,
1501                           UnaryFunction operation) const                         UnaryFunction operation)
1502  {  {
1503    for (ValueType::size_type i=0;i<(noValues(m_shape));++i) {    EsysAssert((!isEmpty()&&checkOffset(offset)),
1504      operation((*m_data)[i+leftOffset]);                 "Error - Couldn't perform unaryOp due to insufficient storage.");
1505      for (ValueType::size_type i=0;i<noValues();i++) {
1506        (*m_data)[offset+i]=operation((*m_data)[offset+i]);
1507    }    }
   return operation.getResult();  
1508  }  }
1509    
1510  template <class UnaryFunction>  template <class BinaryFunction>
1511  inline  inline
1512  double  void
1513  DataArrayView::algorithm(UnaryFunction operation) const  DataArrayView::binaryOp(const DataArrayView& right,
1514                            BinaryFunction operation)
1515  {  {
1516    return algorithm(m_offset,operation);    binaryOp(m_offset,right,right.getOffset(),operation);
1517  }  }
1518    
1519  template <class BinaryFunction>  template <class BinaryFunction>
1520  inline  inline
1521  void  void
1522  DataArrayView::binaryOp(ValueType::size_type leftOffset,  DataArrayView::binaryOp(ValueType::size_type leftOffset,
1523                          const DataArrayView& right,                          const DataArrayView& right,
1524                          ValueType::size_type rightOffset,                          ValueType::size_type rightOffset,
1525                          BinaryFunction operation)                          BinaryFunction operation)
1526  {  {
1527    EsysAssert(getShape()==right.getShape(),    EsysAssert(getShape()==right.getShape(),
1528           "Error - Right hand shape: " << shapeToString(right.getShape()) << " doesn't match the left: " << shapeToString(getShape()));           "Error - Couldn't perform binaryOp due to shape mismatch,");
1529    for (ValueType::size_type i=0;i<noValues();++i) {    EsysAssert((!isEmpty()&&checkOffset(leftOffset)),
1530      (*m_data)[i+leftOffset]=operation((*m_data)[i+leftOffset],(*right.m_data)[i+rightOffset]);               "Error - Couldn't perform binaryOp due to insufficient storage in left object.");
1531      EsysAssert((!right.isEmpty()&&right.checkOffset(rightOffset)),
1532                 "Error - Couldn't perform binaryOp due to insufficient storage in right object.");
1533      for (ValueType::size_type i=0;i<noValues();i++) {
1534        (*m_data)[leftOffset+i]=operation((*m_data)[leftOffset+i],(*right.m_data)[rightOffset+i]);
1535    }    }
1536  }  }
1537    
1538  template <class BinaryFunction>  template <class BinaryFunction>
1539  inline  inline
1540  void  void
1541  DataArrayView::binaryOp(const DataArrayView& right,  DataArrayView::binaryOp(double right,
1542                          BinaryFunction operation)                          BinaryFunction operation)
1543  {  {
1544    //    binaryOp(m_offset,right,operation);
   // version using the offsets specified within the view  
   binaryOp(m_offset,right,right.getOffset(),operation);  
1545  }  }
1546    
1547  template <class BinaryFunction>  template <class BinaryFunction>
1548  inline  inline
1549  void  void
1550  DataArrayView::binaryOp(ValueType::size_type leftOffset,  DataArrayView::binaryOp(ValueType::size_type offset,
1551                          double right,                          double right,
1552                          BinaryFunction operation)                          BinaryFunction operation)
1553  {  {
1554    //    EsysAssert((!isEmpty()&&checkOffset(offset)),
1555    // If a scalar is to be applied to the entire array force the caller to               "Error - Couldn't perform binaryOp due to insufficient storage in left object.");
1556    // explicitly specify a single value    for (ValueType::size_type i=0;i<noValues();i++) {
1557    for (ValueType::size_type i=0;i<(noValues(m_shape));++i) {      (*m_data)[offset+i]=operation((*m_data)[offset+i],right);
     (*m_data)[i+leftOffset]=operation((*m_data)[i+leftOffset],right);  
1558    }    }
1559  }  }
1560    
1561  template <class BinaryFunction>  template <class BinaryFunction>
1562  inline  inline
1563  void  double
1564  DataArrayView::binaryOp(double right,  DataArrayView::reductionOp(BinaryFunction operation,
1565                          BinaryFunction operation)                             double initial_value) const
1566  {  {
1567    //    return reductionOp(m_offset,operation,initial_value);
1568    // use the default offset  }
1569    binaryOp(m_offset,right,operation);  
1570    template <class BinaryFunction>
1571    inline
1572    double
1573    DataArrayView::reductionOp(ValueType::size_type offset,
1574                               BinaryFunction operation,
1575                               double initial_value) const
1576    {
1577      EsysAssert((!isEmpty()&&checkOffset(offset)),
1578                   "Error - Couldn't perform reductionOp due to insufficient storage.");
1579      double current_value=initial_value;
1580      for (ValueType::size_type i=0;i<noValues();i++) {
1581        current_value=operation(current_value,(*m_data)[offset+i]);
1582      }
1583      return current_value;
1584    }
1585    
1586    inline
1587    DataArrayView::ValueType::size_type
1588    DataArrayView::relIndex() const
1589    {
1590      EsysAssert((getRank()==0),"Incorrect number of indices for the rank of this object.");
1591      return 0;
1592    }
1593    
1594    inline
1595    DataArrayView::ValueType::size_type
1596    DataArrayView::index() const
1597    {
1598      EsysAssert((getRank()==0),"Incorrect number of indices for the rank of this object.");
1599      return (m_offset);
1600    }
1601    
1602    inline
1603    DataArrayView::ValueType::size_type
1604    DataArrayView::relIndex(ValueType::size_type i) const
1605    {
1606      EsysAssert((getRank()==1),"Incorrect number of indices for the rank of this object.");
1607      EsysAssert((i < noValues(m_shape)), "Error - Invalid index.");
1608      return i;
1609  }  }
1610    
1611  inline  inline
1612  DataArrayView::ValueType::size_type  DataArrayView::ValueType::size_type
1613  DataArrayView::index(ValueType::size_type i) const  DataArrayView::index(ValueType::size_type i) const
1614  {  {
1615      //EsysAssert((i >= 0) && (i < noValues(m_shape)), "Invalid index.");    EsysAssert((getRank()==1),"Incorrect number of indices for the rank of this object.");
1616      EsysAssert((i < noValues(m_shape)), "Invalid index.");    EsysAssert((i < noValues(m_shape)), "Error - Invalid index.");
1617      return (i+m_offset);    return (m_offset+i);
1618  }  }
1619    
1620  inline  inline
# Line 697  DataArrayView::ValueType::size_type Line 1622  DataArrayView::ValueType::size_type
1622  DataArrayView::relIndex(ValueType::size_type i,  DataArrayView::relIndex(ValueType::size_type i,
1623                          ValueType::size_type j) const                          ValueType::size_type j) const
1624  {  {
1625      EsysAssert((getRank()==2),"Incorrect number of indices for the rank of this object.");
1626    ValueType::size_type temp=i+j*m_shape[0];    ValueType::size_type temp=i+j*m_shape[0];
1627    //EsysAssert((temp >= 0 || temp < noValues(m_shape)), "Invalid index.");    EsysAssert((temp < noValues(m_shape)), "Error - Invalid index.");
   EsysAssert((temp < noValues(m_shape)), "Invalid index.");  
   //  
   // no offset  
1628    return temp;    return temp;
1629  }  }
1630    
# Line 710  DataArrayView::ValueType::size_type Line 1633  DataArrayView::ValueType::size_type
1633  DataArrayView::index(ValueType::size_type i,  DataArrayView::index(ValueType::size_type i,
1634               ValueType::size_type j) const               ValueType::size_type j) const
1635  {  {
1636      EsysAssert((getRank()==2),"Incorrect number of indices for the rank of this object.");
1637    ValueType::size_type temp=i+j*m_shape[0];    ValueType::size_type temp=i+j*m_shape[0];
1638    //EsysAssert((temp >= 0 || temp < noValues(m_shape)), "Invalid index.");    EsysAssert((temp < noValues(m_shape)), "Error - Invalid index.");
1639    EsysAssert((temp < noValues(m_shape)), "Invalid index.");    return (m_offset+temp);
   return (temp+m_offset);  
1640  }  }
1641    
1642  inline  inline
# Line 722  DataArrayView::relIndex(ValueType::size_ Line 1645  DataArrayView::relIndex(ValueType::size_
1645              ValueType::size_type j,              ValueType::size_type j,
1646              ValueType::size_type k) const              ValueType::size_type k) const
1647  {  {
1648      ValueType::size_type temp=i+j*m_shape[0]+k*m_shape[1]*m_shape[0];    EsysAssert((getRank()==3),"Incorrect number of indices for the rank of this object.");
1649      //EsysAssert((temp >= 0 || temp < noValues(m_shape)), "Invalid index.");    ValueType::size_type temp=i+j*m_shape[0]+k*m_shape[1]*m_shape[0];
1650      EsysAssert((temp < noValues(m_shape)), "Invalid index.");    EsysAssert((temp < noValues(m_shape)), "Error - Invalid index.");
1651      //    return temp;
     // no offset  
     return temp;  
1652  }  }
1653    
1654  inline  inline
# Line 736  DataArrayView::index(ValueType::size_typ Line 1657  DataArrayView::index(ValueType::size_typ
1657               ValueType::size_type j,               ValueType::size_type j,
1658               ValueType::size_type k) const               ValueType::size_type k) const
1659  {  {
1660      ValueType::size_type temp=i+j*m_shape[0]+k*m_shape[1]*m_shape[0];    EsysAssert((getRank()==3),"Incorrect number of indices for the rank of this object.");
1661      //EsysAssert((temp >= 0 || temp < noValues(m_shape)), "Invalid index.");    ValueType::size_type temp=i+j*m_shape[0]+k*m_shape[1]*m_shape[0];
1662      EsysAssert((temp < noValues(m_shape)), "Invalid index.");    EsysAssert((temp < noValues(m_shape)), "Error - Invalid index.");
1663      return (temp+m_offset);    return (m_offset+temp);
1664  }  }
1665    
1666  inline  inline
# Line 749  DataArrayView::relIndex(ValueType::size_ Line 1670  DataArrayView::relIndex(ValueType::size_
1670                          ValueType::size_type k,                          ValueType::size_type k,
1671                          ValueType::size_type m) const                          ValueType::size_type m) const
1672  {  {
1673      EsysAssert((getRank()==4),"Incorrect number of indices for the rank of this object.");
1674    ValueType::size_type temp=i+j*m_shape[0]+k*m_shape[1]*m_shape[0]+m*m_shape[2]*m_shape[1]*m_shape[0];    ValueType::size_type temp=i+j*m_shape[0]+k*m_shape[1]*m_shape[0]+m*m_shape[2]*m_shape[1]*m_shape[0];
1675    //EsysAssert((temp >= 0 || temp < noValues(m_shape)), "Invalid index.");    EsysAssert((temp < noValues(m_shape)), "Error - Invalid index.");
   EsysAssert((temp < noValues(m_shape)), "Invalid index.");  
   //  
   // no offset  
1676    return temp;    return temp;
1677  }  }
1678    
# Line 764  DataArrayView::index(ValueType::size_typ Line 1683  DataArrayView::index(ValueType::size_typ
1683               ValueType::size_type k,               ValueType::size_type k,
1684               ValueType::size_type m) const               ValueType::size_type m) const
1685  {  {
1686      EsysAssert((getRank()==4),"Incorrect number of indices for the rank of this object.");
1687    ValueType::size_type temp=i+j*m_shape[0]+k*m_shape[1]*m_shape[0]+m*m_shape[2]*m_shape[1]*m_shape[0];    ValueType::size_type temp=i+j*m_shape[0]+k*m_shape[1]*m_shape[0]+m*m_shape[2]*m_shape[1]*m_shape[0];
1688    //EsysAssert((temp >= 0 || temp < noValues(m_shape)), "Invalid index.");    EsysAssert((temp < noValues(m_shape)), "Error - Invalid index.");
1689    EsysAssert((temp < noValues(m_shape)), "Invalid index.");    return (m_offset+temp);
   return (temp+m_offset);  
1690  }  }
1691    
1692  inline  inline
1693  DataArrayView::ValueType::reference  DataArrayView::ValueType::reference
1694  DataArrayView::operator()(ValueType::size_type i,  DataArrayView::operator()()
                           ValueType::size_type j,  
                           ValueType::size_type k,  
                           ValueType::size_type m)  
1695  {  {
1696      EsysAssert((getRank()==4),    return (*m_data)[index()];
            "Incorrect number of indices for the rank of this object.");  
     return (*m_data)[index(i,j,k,m)];  
1697  }  }
1698    
1699  inline  inline
1700  DataArrayView::ValueType::const_reference  DataArrayView::ValueType::const_reference
1701  DataArrayView::operator()(ValueType::size_type i,  DataArrayView::operator()() const
                           ValueType::size_type j,  
                           ValueType::size_type k,  
                           ValueType::size_type m) const  
1702  {  {
1703      EsysAssert((getRank()==4),    return (*m_data)[index()];
            "Incorrect number of indices for the rank of this object.");  
     return (*m_data)[index(i,j,k,m)];  
1704  }  }
1705    
1706  inline  inline
1707  DataArrayView::ValueType::reference  DataArrayView::ValueType::reference
1708  DataArrayView::operator()(ValueType::size_type i,  DataArrayView::operator()(ValueType::size_type i)
                           ValueType::size_type j,  
                           ValueType::size_type k)  
1709  {  {
1710      EsysAssert((getRank()==3),    return (*m_data)[index(i)];
            "Incorrect number of indices for the rank of this object.");  
     return (*m_data)[index(i,j,k)];  
1711  }  }
1712    
1713  inline  inline
1714  DataArrayView::ValueType::const_reference  DataArrayView::ValueType::const_reference
1715  DataArrayView::operator()(ValueType::size_type i,  DataArrayView::operator()(ValueType::size_type i) const
                           ValueType::size_type j,  
                           ValueType::size_type k) const  
1716  {  {
1717      EsysAssert((getRank()==3),    return (*m_data)[index(i)];
            "Incorrect number of indices for the rank of this object.");  
     return (*m_data)[index(i,j,k)];  
1718  }  }
1719    
1720  inline  inline
# Line 821  DataArrayView::ValueType::reference Line 1722  DataArrayView::ValueType::reference
1722  DataArrayView::operator()(ValueType::size_type i,  DataArrayView::operator()(ValueType::size_type i,
1723                            ValueType::size_type j)                            ValueType::size_type j)
1724  {  {
1725      EsysAssert((getRank()==2),    return (*m_data)[index(i,j)];
            "Incorrect number of indices for the rank of this object.");  
     return (*m_data)[index(i,j)];  
1726  }  }
1727    
1728  inline  inline
# Line 831  DataArrayView::ValueType::const_referenc Line 1730  DataArrayView::ValueType::const_referenc
1730  DataArrayView::operator()(ValueType::size_type i,  DataArrayView::operator()(ValueType::size_type i,
1731                            ValueType::size_type j) const                            ValueType::size_type j) const
1732  {  {
1733      EsysAssert((getRank()==2),    return (*m_data)[index(i,j)];
            "Incorrect number of indices for the rank of this object.");  
     return (*m_data)[index(i,j)];  
1734  }  }
1735    
1736  inline  inline
1737  DataArrayView::ValueType::reference  DataArrayView::ValueType::reference
1738  DataArrayView::operator()(ValueType::size_type i)  DataArrayView::operator()(ValueType::size_type i,
1739                              ValueType::size_type j,
1740                              ValueType::size_type k)
1741  {  {
1742      EsysAssert((getRank()==1),    return (*m_data)[index(i,j,k)];
            "Incorrect number of indices for the rank of this object.");  
     return (*m_data)[index(i)];  
1743  }  }
1744    
1745  inline  inline
1746  DataArrayView::ValueType::const_reference  DataArrayView::ValueType::const_reference
1747  DataArrayView::operator()(ValueType::size_type i) const  DataArrayView::operator()(ValueType::size_type i,
1748                              ValueType::size_type j,
1749                              ValueType::size_type k) const
1750  {  {
1751      EsysAssert((getRank()==1),    return (*m_data)[index(i,j,k)];
            "Incorrect number of indices for the rank of this object.");  
     return (*m_data)[index(i)];  
1752  }  }
1753    
1754  inline  inline
1755  DataArrayView::ValueType::reference  DataArrayView::ValueType::reference
1756  DataArrayView::operator()()  DataArrayView::operator()(ValueType::size_type i,
1757                              ValueType::size_type j,
1758                              ValueType::size_type k,
1759                              ValueType::size_type m)
1760  {  {
1761      EsysAssert((getRank()==0),    return (*m_data)[index(i,j,k,m)];
            "Incorrect number of indices for the rank of this object.");  
     return (*m_data)[m_offset];  
1762  }  }
1763    
1764  inline  inline
1765  DataArrayView::ValueType::const_reference  DataArrayView::ValueType::const_reference
1766  DataArrayView::operator()() const  DataArrayView::operator()(ValueType::size_type i,
1767                              ValueType::size_type j,
1768                              ValueType::size_type k,
1769                              ValueType::size_type m) const
1770  {  {
1771      EsysAssert((getRank()==0),    return (*m_data)[index(i,j,k,m)];
            "Incorrect number of indices for the rank of this object.");  
     return (*m_data)[m_offset];  
1772  }  }
1773    
 bool operator==(const DataArrayView& left, const DataArrayView& right);  
 bool operator!=(const DataArrayView& left, const DataArrayView& right);  
   
 /* returns the python slice object key as a pair of ints where the first  
     member is the start and the second member is the end. the presence of a possible  
     step attribute with value different from one will truow an exception */  
 std::pair<int,int>  
 getSliceRange(const int s,  
               const boost::python::object&  key);  
   
1774  } // end of namespace  } // end of namespace
1775    
1776  #endif  #endif

Legend:
Removed from v.100  
changed lines
  Added in v.800

  ViewVC Help
Powered by ViewVC 1.1.26