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

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

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

trunk/esys2/escript/src/Data/Data.h revision 97 by jgs, Tue Dec 14 05:39:33 2004 UTC trunk/escript/src/Data.h revision 2458 by jfenwick, Wed Jun 3 06:18:21 2009 UTC
# Line 1  Line 1 
 // $Id$  
 /*=============================================================================  
1    
2   ******************************************************************************  /*******************************************************
3   *                                                                            *  *
4   *       COPYRIGHT ACcESS 2004 -  All Rights Reserved                         *  * Copyright (c) 2003-2008 by University of Queensland
5   *                                                                            *  * Earth Systems Science Computational Center (ESSCC)
6   * This software is the property of ACcESS.  No part of this code             *  * http://www.uq.edu.au/esscc
7   * may be copied in any form or by any means without the expressed written    *  *
8   * consent of ACcESS.  Copying, use or modification of this software          *  * Primary Business: Queensland, Australia
9   * by any unauthorised person is illegal unless that                          *  * Licensed under the Open Software License version 3.0
10   * person has a software license agreement with ACcESS.                       *  * http://www.opensource.org/licenses/osl-3.0.php
11   *                                                                            *  *
12   ******************************************************************************  *******************************************************/
13    
14  ******************************************************************************/  
15    /** \file Data.h */
16    
17  #ifndef DATA_H  #ifndef DATA_H
18  #define DATA_H  #define DATA_H
19    #include "system_dep.h"
20    
21    #include "DataTypes.h"
22    #include "DataAbstract.h"
23    #include "DataAlgorithm.h"
24    #include "FunctionSpace.h"
25    #include "BinaryOp.h"
26    #include "UnaryOp.h"
27    #include "DataException.h"
28    
 #include "escript/Data/DataAbstract.h"  
 #include "escript/Data/DataTagged.h"  
 #include "escript/Data/FunctionSpace.h"  
 #include "escript/Data/BinaryOp.h"  
 #include "escript/Data/UnaryOp.h"  
29    
30  extern "C" {  extern "C" {
31  #include "escript/Data/DataC.h"  #include "DataC.h"
32    //#include <omp.h>
33  }  }
34    
35  #include <iostream>  #include "esysmpi.h"
36  #include <string>  #include <string>
 #include <memory>  
37  #include <algorithm>  #include <algorithm>
38    #include <sstream>
39    
40  #include <boost/shared_ptr.hpp>  #include <boost/shared_ptr.hpp>
41  #include <boost/python/object.hpp>  #include <boost/python/object.hpp>
 #include <boost/python/list.hpp>  
42  #include <boost/python/tuple.hpp>  #include <boost/python/tuple.hpp>
 #include <boost/python/numeric.hpp>  
   
 /**  
    \brief  
    Data is essentially a factory class which creates the appropriate Data  
    object for the given construction arguments. It retains control over  
    the object created for the lifetime of the object.  
    The type of Data object referred to may change during the lifetime of  
    the Data object.  
43    
44     Description:  #include "BufferGroup.h"
    Data is essentially a factory class which creates the appropriate Data  
    object for the given construction arguments. It retains control over  
    the object created for the lifetime of the object.  
    The type of Data object referred to may change during the lifetime of  
    the Data object.  
 */  
45    
46  namespace escript {  namespace escript {
47    
48  //  //
49  // Forward declaration for various implimentations of Data.  // Forward declaration for various implementations of Data.
 class DataEmpty;  
50  class DataConstant;  class DataConstant;
51  class DataTagged;  class DataTagged;
52  class DataExpanded;  class DataExpanded;
53    class DataLazy;
54    
55    /**
56       \brief
57       Data represents a collection of datapoints.
58    
59       Description:
60       Internally, the datapoints are actually stored by a DataAbstract object.
61       The specific instance of DataAbstract used may vary over the lifetime
62       of the Data object.
63       Some methods on this class return references (eg getShape()).
64       These references should not be used after an operation which changes the underlying DataAbstract object.
65       Doing so will lead to invalid memory access.
66       This should not affect any methods exposed via boost::python.
67    */
68  class Data {  class Data {
69    
70    public:    public:
71    
72      // These typedefs allow function names to be cast to pointers
73      // to functions of the appropriate type when calling unaryOp etc.
74    typedef double (*UnaryDFunPtr)(double);    typedef double (*UnaryDFunPtr)(double);
75    typedef double (*BinaryDFunPtr)(double,double);    typedef double (*BinaryDFunPtr)(double,double);
76    
77    
78    /**    /**
79       Constructors.       Constructors.
80    */    */
# Line 80  class Data { Line 84  class Data {
84       Default constructor.       Default constructor.
85       Creates a DataEmpty object.       Creates a DataEmpty object.
86    */    */
87      ESCRIPT_DLL_API
88    Data();    Data();
89    
90    /**    /**
# Line 87  class Data { Line 92  class Data {
92       Copy constructor.       Copy constructor.
93       WARNING: Only performs a shallow copy.       WARNING: Only performs a shallow copy.
94    */    */
95      ESCRIPT_DLL_API
96    Data(const Data& inData);    Data(const Data& inData);
97    
98    /**    /**
# Line 95  class Data { Line 101  class Data {
101       function space of inData the inData are tried to be interpolated to what,       function space of inData the inData are tried to be interpolated to what,
102       otherwise a shallow copy of inData is returned.       otherwise a shallow copy of inData is returned.
103    */    */
104      ESCRIPT_DLL_API
105    Data(const Data& inData,    Data(const Data& inData,
106         const FunctionSpace& what);         const FunctionSpace& what);
107    
108    /**    /**
109       \brief      \brief Copy Data from an existing vector
110       Constructor which copies data from a DataArrayView.    */
111    
112       \param value - Input - Data value for a single point.    ESCRIPT_DLL_API
113       \param what - Input - A description of what this data represents.    Data(const DataTypes::ValueType& value,
114       \param expanded - Input - Flag, if true fill the entire container with           const DataTypes::ShapeType& shape,
115                         the value. Otherwise a more efficient storage                   const FunctionSpace& what=FunctionSpace(),
116                         mechanism will be used.                   bool expanded=false);
   */  
   Data(const DataArrayView& value,  
        const FunctionSpace& what=FunctionSpace(),  
        bool expanded=false);  
117    
118    /**    /**
119       \brief       \brief
# Line 123  class Data { Line 126  class Data {
126                         the given value. Otherwise a more efficient storage                         the given value. Otherwise a more efficient storage
127                         mechanism will be used.                         mechanism will be used.
128    */    */
129      ESCRIPT_DLL_API
130    Data(double value,    Data(double value,
131         const DataArrayView::ShapeType& dataPointShape=DataArrayView::ShapeType(),         const DataTypes::ShapeType& dataPointShape=DataTypes::ShapeType(),
132         const FunctionSpace& what=FunctionSpace(),         const FunctionSpace& what=FunctionSpace(),
133         bool expanded=false);         bool expanded=false);
134    
# Line 135  class Data { Line 139  class Data {
139       \param inData - Input - Input Data object.       \param inData - Input - Input Data object.
140       \param region - Input - Region to copy.       \param region - Input - Region to copy.
141    */    */
142      ESCRIPT_DLL_API
143    Data(const Data& inData,    Data(const Data& inData,
144         const DataArrayView::RegionType& region);         const DataTypes::RegionType& region);
   
   /**  
      \brief  
      Constructor which will create Tagged data if expanded is false.  
      No attempt is made to ensure the tag keys match the tag keys  
      within the function space.  
   
      \param tagKeys - Input - List of tag values.  
      \param values - Input - List of values, one for each tag.  
      \param defaultValue - Input - A default value, used if tag doesn't exist.  
      \param what - Input - A description of what this data represents.  
      \param expanded - Input - Flag, if true fill the entire container with  
                        the appropriate values.  
   */  
   Data(const DataTagged::TagListType& tagKeys,  
        const DataTagged::ValueListType& values,  
        const DataArrayView& defaultValue,  
        const FunctionSpace& what=FunctionSpace(),  
        bool expanded=false);  
   
   /**  
      \brief  
      Constructor which copies data from a python numarray.  
   
      \param value - Input - Data value for a single point.  
      \param what - Input - A description of what this data represents.  
      \param expanded - Input - Flag, if true fill the entire container with  
                        the value. Otherwise a more efficient storage  
                        mechanism will be used.  
   */  
   Data(const boost::python::numeric::array& value,  
        const FunctionSpace& what=FunctionSpace(),  
        bool expanded=false);  
145    
146    /**    /**
147       \brief       \brief
148       Constructor which copies data from any object that can be converted into       Constructor which copies data from any object that can be treated like a python array/sequence.
      a python numarray.  
149    
150       \param value - Input - Input data.       \param value - Input - Input data.
151       \param what - Input - A description of what this data represents.       \param what - Input - A description of what this data represents.
# Line 182  class Data { Line 153  class Data {
153                         the value. Otherwise a more efficient storage                         the value. Otherwise a more efficient storage
154                         mechanism will be used.                         mechanism will be used.
155    */    */
156      ESCRIPT_DLL_API
157    Data(const boost::python::object& value,    Data(const boost::python::object& value,
158         const FunctionSpace& what=FunctionSpace(),         const FunctionSpace& what=FunctionSpace(),
159         bool expanded=false);         bool expanded=false);
# Line 189  class Data { Line 161  class Data {
161    /**    /**
162       \brief       \brief
163       Constructor which creates a DataConstant.       Constructor which creates a DataConstant.
164       Copies data from any object that can be converted       Copies data from any object that can be treated like a python array/sequence.
165       into a numarray. All other parameters are copied from other.       All other parameters are copied from other.
166    
167       \param value - Input - Input data.       \param value - Input - Input data.
168       \param other - Input - contains all other parameters.       \param other - Input - contains all other parameters.
169    */    */
170      ESCRIPT_DLL_API
171    Data(const boost::python::object& value,    Data(const boost::python::object& value,
172         const Data& other);         const Data& other);
173    
# Line 202  class Data { Line 175  class Data {
175       \brief       \brief
176       Constructor which creates a DataConstant of "shape" with constant value.       Constructor which creates a DataConstant of "shape" with constant value.
177    */    */
178    Data(double value,    ESCRIPT_DLL_API
179         const boost::python::tuple& shape=boost::python::make_tuple(),    Data(double value,
180           const boost::python::tuple& shape=boost::python::make_tuple(),
181         const FunctionSpace& what=FunctionSpace(),         const FunctionSpace& what=FunctionSpace(),
182         bool expanded=false);         bool expanded=false);
183    
184    
185    
186      /**
187        \brief Create a Data using an existing DataAbstract. Warning: The new object assumes ownership of the pointer!
188        Once you have passed the pointer, do not delete it.
189      */
190      ESCRIPT_DLL_API
191      explicit Data(DataAbstract* underlyingdata);
192    
193      /**
194        \brief Create a Data based on the supplied DataAbstract
195      */
196      ESCRIPT_DLL_API
197      explicit Data(DataAbstract_ptr underlyingdata);
198    
199    /**    /**
200       \brief       \brief
201       Perform a deep copy.       Destructor
202    */    */
203      ESCRIPT_DLL_API
204      ~Data();
205    
206      /**
207         \brief Make this object a deep copy of "other".
208      */
209      ESCRIPT_DLL_API
210    void    void
211    copy(const Data& other);    copy(const Data& other);
212    
213    /**    /**
214         \brief Return a pointer to a deep copy of this object.
215      */
216      ESCRIPT_DLL_API
217      Data
218      copySelf();
219    
220    
221      /**
222         \brief produce a delayed evaluation version of this Data.
223      */
224      ESCRIPT_DLL_API
225      Data
226      delay();
227    
228      /**
229         \brief convert the current data into lazy data.
230      */
231      ESCRIPT_DLL_API
232      void
233      delaySelf();
234    
235    
236      /**
237       Member access methods.       Member access methods.
238    */    */
239    
240    /**    /**
241       \brief       \brief
242         switches on update protection
243    
244      */
245      ESCRIPT_DLL_API
246      void
247      setProtection();
248    
249      /**
250         \brief
251         Returns trueif the data object is protected against update
252    
253      */
254      ESCRIPT_DLL_API
255      bool
256      isProtected() const;
257    
258    
259    /**
260       \brief
261       Return teh value of a data point as a python tuple.
262    */
263      ESCRIPT_DLL_API
264      const boost::python::object
265      getValueOfDataPointAsTuple(int dataPointNo);
266    
267      /**
268         \brief
269         sets the values of a data-point from a python object on this process
270      */
271      ESCRIPT_DLL_API
272      void
273      setValueOfDataPointToPyObject(int dataPointNo, const boost::python::object& py_object);
274    
275      /**
276         \brief
277         sets the values of a data-point from a array-like object on this process
278      */
279      ESCRIPT_DLL_API
280      void
281      setValueOfDataPointToArray(int dataPointNo, const boost::python::object&);
282    
283      /**
284         \brief
285         sets the values of a data-point on this process
286      */
287      ESCRIPT_DLL_API
288      void
289      setValueOfDataPoint(int dataPointNo, const double);
290    
291      /**
292         \brief Return a data point across all processors as a python tuple.
293      */
294      ESCRIPT_DLL_API
295      const boost::python::object
296      getValueOfGlobalDataPointAsTuple(int procNo, int dataPointNo);
297    
298      /**
299         \brief
300         Return the tag number associated with the given data-point.
301    
302      */
303      ESCRIPT_DLL_API
304      int
305      getTagNumber(int dpno);
306    
307      /**
308         \brief
309       Return the C wrapper for the Data object.       Return the C wrapper for the Data object.
310    */    */
311      ESCRIPT_DLL_API
312    escriptDataC    escriptDataC
313    getDataC();    getDataC();
314    
315    
316    
317    /**    /**
318       \brief       \brief
319       Return the C wrapper for the Data object - const version.       Return the C wrapper for the Data object - const version.
320    */    */
321      ESCRIPT_DLL_API
322    escriptDataC    escriptDataC
323    getDataC() const;    getDataC() const;
324    
325    /**    /**
326       \brief      \brief How much space is required to evaulate a sample of the Data.
      Write the data as a string.  
327    */    */
328    inline    ESCRIPT_DLL_API
329    std::string    size_t
330    toString() const    getSampleBufferSize() const;
331    {  
332      return m_data->toString();  
   }  
333    
334    /**    /**
335       \brief       \brief
336       Return the DataArrayView of the point data. This essentially contains       Write the data as a string. For large amounts of data, a summary is printed.
      the shape information for each data point although it also may be used  
      to manipulate the point data.  
337    */    */
338    inline    ESCRIPT_DLL_API
339    const DataArrayView&    std::string
340    getPointDataView() const    toString() const;
   {  
      return m_data->getPointDataView();  
   }  
341    
342    /**    /**
343       \brief       \brief
344       Whatever the current Data type make this into a DataExpanded.       Whatever the current Data type make this into a DataExpanded.
345    */    */
346      ESCRIPT_DLL_API
347    void    void
348    expand();    expand();
349    
# Line 269  class Data { Line 353  class Data {
353       Constant data to be converted to tagged. An attempt to convert       Constant data to be converted to tagged. An attempt to convert
354       Expanded data to tagged will throw an exception.       Expanded data to tagged will throw an exception.
355    */    */
356      ESCRIPT_DLL_API
357    void    void
358    tag();    tag();
359    
360    /**    /**
361        \brief If this data is lazy, then convert it to ready data.
362        What type of ready data depends on the expression. For example, Constant+Tagged==Tagged.
363      */
364      ESCRIPT_DLL_API
365      void
366      resolve();
367    
368    
369      /**
370       \brief Ensures data is ready for write access.
371      This means that the data will be resolved if lazy and will be copied if shared with another Data object.
372      \warning This method should only be called in single threaded sections of code. (It modifies m_data).
373      Do not create any Data objects from this one between calling requireWrite and getSampleDataRW.
374      Doing so might introduce additional sharing.
375      */
376      ESCRIPT_DLL_API
377      void
378      requireWrite();
379    
380      /**
381       \brief       \brief
382       Return true if this Data is expanded.       Return true if this Data is expanded.
383         \note To determine if a sample will contain separate values for each datapoint. Use actsExpanded instead.
384    */    */
385      ESCRIPT_DLL_API
386    bool    bool
387    isExpanded() const;    isExpanded() const;
388    
389    /**    /**
390       \brief       \brief
391         Return true if this Data is expanded or resolves to expanded.
392         That is, if it has a separate value for each datapoint in the sample.
393      */
394      ESCRIPT_DLL_API
395      bool
396      actsExpanded() const;
397      
398    
399      /**
400         \brief
401       Return true if this Data is tagged.       Return true if this Data is tagged.
402    */    */
403      ESCRIPT_DLL_API
404    bool    bool
405    isTagged() const;    isTagged() const;
406    
# Line 290  class Data { Line 408  class Data {
408       \brief       \brief
409       Return true if this Data is constant.       Return true if this Data is constant.
410    */    */
411      ESCRIPT_DLL_API
412    bool    bool
413    isConstant() const;    isConstant() const;
414    
415    /**    /**
416         \brief Return true if this Data is lazy.
417      */
418      ESCRIPT_DLL_API
419      bool
420      isLazy() const;
421    
422      /**
423         \brief Return true if this data is ready.
424      */
425      ESCRIPT_DLL_API
426      bool
427      isReady() const;
428    
429      /**
430       \brief       \brief
431       Return true if this Data is empty.       Return true if this Data holds an instance of DataEmpty. This is _not_ the same as asking if the object
432    contains datapoints.
433    */    */
434      ESCRIPT_DLL_API
435    bool    bool
436    isEmpty() const;    isEmpty() const;
437    
# Line 304  class Data { Line 439  class Data {
439       \brief       \brief
440       Return the function space.       Return the function space.
441    */    */
442      ESCRIPT_DLL_API
443    inline    inline
444    const FunctionSpace&    const FunctionSpace&
445    getFunctionSpace() const    getFunctionSpace() const
# Line 315  class Data { Line 451  class Data {
451       \brief       \brief
452       Return a copy of the function space.       Return a copy of the function space.
453    */    */
454      ESCRIPT_DLL_API
455    const FunctionSpace    const FunctionSpace
456    getCopyOfFunctionSpace() const;    getCopyOfFunctionSpace() const;
457    
# Line 322  class Data { Line 459  class Data {
459       \brief       \brief
460       Return the domain.       Return the domain.
461    */    */
462      ESCRIPT_DLL_API
463    inline    inline
464    const AbstractDomain&  //   const AbstractDomain&
465      const_Domain_ptr
466    getDomain() const    getDomain() const
467    {    {
468       return getFunctionSpace().getDomain();       return getFunctionSpace().getDomain();
469    }    }
470    
471    
472      /**
473         \brief
474         Return the domain.
475         TODO: For internal use only.   This should be removed.
476      */
477      ESCRIPT_DLL_API
478      inline
479    //   const AbstractDomain&
480      Domain_ptr
481      getDomainPython() const
482      {
483         return getFunctionSpace().getDomainPython();
484      }
485    
486    /**    /**
487       \brief       \brief
488       Return a copy of the domain.       Return a copy of the domain.
489    */    */
490      ESCRIPT_DLL_API
491    const AbstractDomain    const AbstractDomain
492    getCopyOfDomain() const;    getCopyOfDomain() const;
493    
# Line 340  class Data { Line 495  class Data {
495       \brief       \brief
496       Return the rank of the point data.       Return the rank of the point data.
497    */    */
498      ESCRIPT_DLL_API
499    inline    inline
500    int    unsigned int
501    getDataPointRank() const    getDataPointRank() const
502    {    {
503      return m_data->getPointDataView().getRank();      return m_data->getRank();
504    }    }
505    
506    /**    /**
507       \brief       \brief
508         Return the number of data points
509      */
510      ESCRIPT_DLL_API
511      inline
512      int
513      getNumDataPoints() const
514      {
515        return getNumSamples() * getNumDataPointsPerSample();
516      }
517      /**
518         \brief
519       Return the number of samples.       Return the number of samples.
520    */    */
521      ESCRIPT_DLL_API
522    inline    inline
523    int    int
524    getNumSamples() const    getNumSamples() const
# Line 362  class Data { Line 530  class Data {
530       \brief       \brief
531       Return the number of data points per sample.       Return the number of data points per sample.
532    */    */
533      ESCRIPT_DLL_API
534    inline    inline
535    int    int
536    getNumDataPointsPerSample() const    getNumDataPointsPerSample() const
# Line 369  class Data { Line 538  class Data {
538      return m_data->getNumDPPSample();      return m_data->getNumDPPSample();
539    }    }
540    
541    
542      /**
543        \brief
544        Return the number of values in the shape for this object.
545      */
546      ESCRIPT_DLL_API
547      int
548      getNoValues() const
549      {
550        return m_data->getNoValues();
551      }
552    
553    
554      /**
555         \brief
556         dumps the object into a netCDF file
557      */
558      ESCRIPT_DLL_API
559      void
560      dump(const std::string fileName) const;
561    
562    //  /**
563    //     \brief
564    //     Return the sample data for the given sample no. This is not the
565    //     preferred interface but is provided for use by C code.
566    //     The buffer parameter is only required for LazyData.
567    //     \param sampleNo - Input - the given sample no.
568    //     \param buffer - Vector to compute (and store) sample data in.
569    //     \return pointer to the sample data.
570    //*/
571    //   ESCRIPT_DLL_API
572    //   inline
573    //   const DataAbstract::ValueType::value_type*
574    //   getSampleDataRO(DataAbstract::ValueType::size_type sampleNo, DataTypes::ValueType* buffer=0);
575    
576    
577     /**
578        \brief
579        Return the sample data for the given sample no. This is not the
580        preferred interface but is provided for use by C code.
581        The bufferg parameter is only required for LazyData.
582        \param sampleNo - Input - the given sample no.
583        \param buffer - A buffer to to compute (and store) sample data in will be selected from this group.
584        \return pointer to the sample data.
585    */
586      ESCRIPT_DLL_API
587      inline
588      const DataAbstract::ValueType::value_type*
589      getSampleDataRO(DataAbstract::ValueType::size_type sampleNo, BufferGroup* bufferg=0);
590    
591    
592    /**    /**
593       \brief       \brief
594       Return the sample data for the given sample no. This is not the       Return the sample data for the given sample no. This is not the
595       preferred interface but is provided for use by C code.       preferred interface but is provided for use by C code.
596       \param sampleNo - Input - the given sample no.       \param sampleNo - Input - the given sample no.
597         \return pointer to the sample data.
598    */    */
599      ESCRIPT_DLL_API
600    inline    inline
601    DataAbstract::ValueType::value_type*    DataAbstract::ValueType::value_type*
602    getSampleData(DataAbstract::ValueType::size_type sampleNo)    getSampleDataRW(DataAbstract::ValueType::size_type sampleNo);
603    {  
     return m_data->getSampleData(sampleNo);  
   }  
604    
605    /**    /**
606       \brief       \brief
# Line 388  class Data { Line 608  class Data {
608       access data that isn't tagged an exception will be thrown.       access data that isn't tagged an exception will be thrown.
609       \param tag - Input - the tag key.       \param tag - Input - the tag key.
610    */    */
611      ESCRIPT_DLL_API
612    inline    inline
613    DataAbstract::ValueType::value_type*    DataAbstract::ValueType::value_type*
614    getSampleDataByTag(int tag)    getSampleDataByTag(int tag)
# Line 403  class Data { Line 624  class Data {
624       \param sampleNo - Input -       \param sampleNo - Input -
625       \param dataPointNo - Input -       \param dataPointNo - Input -
626    */    */
627      ESCRIPT_DLL_API
628      DataTypes::ValueType::const_reference
629      getDataPointRO(int sampleNo, int dataPointNo);
630    
631    
632      ESCRIPT_DLL_API
633      DataTypes::ValueType::reference
634      getDataPointRW(int sampleNo, int dataPointNo);
635    
636    
637    
638      /**
639         \brief
640         Return the offset for the given sample and point within the sample
641      */
642      ESCRIPT_DLL_API
643    inline    inline
644    DataArrayView    DataTypes::ValueType::size_type
645    getDataPoint(int sampleNo,    getDataOffset(int sampleNo,
646                 int dataPointNo)                 int dataPointNo)
647    {    {
648      return m_data->getDataPoint(sampleNo,dataPointNo);        return m_data->getPointOffset(sampleNo,dataPointNo);
649    }    }
650    
651    /**    /**
652       \brief       \brief
653       Return a reference to the data point shape.       Return a reference to the data point shape.
654    */    */
655    const DataArrayView::ShapeType&    ESCRIPT_DLL_API
656    getDataPointShape() const;    inline
657      const DataTypes::ShapeType&
658      getDataPointShape() const
659      {
660        return m_data->getShape();
661      }
662    
663    /**    /**
664       \brief       \brief
665       Return the data point shape as a tuple of integers.       Return the data point shape as a tuple of integers.
666    */    */
667    boost::python::tuple    ESCRIPT_DLL_API
668      const boost::python::tuple
669    getShapeTuple() const;    getShapeTuple() const;
670    
671    /**    /**
# Line 430  class Data { Line 673  class Data {
673       Return the size of the data point. It is the product of the       Return the size of the data point. It is the product of the
674       data point shape dimensions.       data point shape dimensions.
675    */    */
676      ESCRIPT_DLL_API
677    int    int
678    getDataPointSize() const;    getDataPointSize() const;
679    
# Line 437  class Data { Line 681  class Data {
681       \brief       \brief
682       Return the number of doubles stored for this Data.       Return the number of doubles stored for this Data.
683    */    */
684    DataArrayView::ValueType::size_type    ESCRIPT_DLL_API
685      DataTypes::ValueType::size_type
686    getLength() const;    getLength() const;
687    
688    
689    
690    /**    /**
691       \brief       \brief
692       Assign the given value to the tag. Implicitly converts this       Assign the given value to the tag assocciated with name. Implicitly converts this
693       object to type DataTagged. Throws an exception if this object       object to type DataTagged. Throws an exception if this object
694       cannot be converted to a DataTagged object.       cannot be converted to a DataTagged object or name cannot be mapped onto a tag key.
695       \param tagKey - Input - Integer key.       \param tagKey - Input - Integer key.
696       \param value - Input - Value to associate with given key.       \param value - Input - Value to associate with given key.
697        ==>*
698    */    */
699      ESCRIPT_DLL_API
700    void    void
701    setTaggedValue(int tagKey,    setTaggedValueByName(std::string name,
702                   const boost::python::object& value);                         const boost::python::object& value);
703    
704    /**    /**
705       \brief       \brief
706       Assign the given value to the tag. Implicitly converts this       Assign the given value to the tag. Implicitly converts this
707       object to type DataTagged. Throws an exception if this object       object to type DataTagged if it is constant.
708       cannot be converted to a DataTagged object.  
709       \param tagKey - Input - Integer key.       \param tagKey - Input - Integer key.
710       \param value - Input - Value to associate with given key.       \param value - Input - Value to associate with given key.
711       Note: removed for now - this version not needed, and breaks escript.cpp      ==>*
712    */    */
713    /*    ESCRIPT_DLL_API
714    void    void
715    setTaggedValue(int tagKey,    setTaggedValue(int tagKey,
716                   const DataArrayView& value);                   const boost::python::object& value);
717    
718      /**
719         \brief
720         Assign the given value to the tag. Implicitly converts this
721         object to type DataTagged if it is constant.
722    
723         \param tagKey - Input - Integer key.
724         \param pointshape - Input - The shape of the value parameter
725         \param value - Input - Value to associate with given key.
726         \param dataOffset - Input - Offset of the begining of the point within the value parameter
727    */    */
728      ESCRIPT_DLL_API
729      void
730      setTaggedValueFromCPP(int tagKey,
731                const DataTypes::ShapeType& pointshape,
732                            const DataTypes::ValueType& value,
733                int dataOffset=0);
734    
735    
736    
737    /**    /**
738      \brief      \brief
739      Copy other Data object into this Data object where mask is positive.      Copy other Data object into this Data object where mask is positive.
740    */    */
741      ESCRIPT_DLL_API
742    void    void
743    copyWithMask(const Data& other,    copyWithMask(const Data& other,
744                 const Data& mask);                 const Data& mask);
# Line 481  class Data { Line 749  class Data {
749    
750    /**    /**
751       \brief       \brief
752         set all values to zero
753         *
754      */
755      ESCRIPT_DLL_API
756      void
757      setToZero();
758    
759      /**
760         \brief
761       Interpolates this onto the given functionspace and returns       Interpolates this onto the given functionspace and returns
762       the result as a Data object.       the result as a Data object.
763         *
764    */    */
765      ESCRIPT_DLL_API
766    Data    Data
767    interpolate(const FunctionSpace& functionspace) const;    interpolate(const FunctionSpace& functionspace) const;
   
768    /**    /**
769       \brief       \brief
770       Calculates the gradient of the data at the data points of functionspace.       Calculates the gradient of the data at the data points of functionspace.
771       If functionspace is not present the function space of Function(getDomain()) is used.       If functionspace is not present the function space of Function(getDomain()) is used.
772         *
773    */    */
774      ESCRIPT_DLL_API
775    Data    Data
776    gradOn(const FunctionSpace& functionspace) const;    gradOn(const FunctionSpace& functionspace) const;
777    
778      ESCRIPT_DLL_API
779    Data    Data
780    grad() const;    grad() const;
781    
782    /**    /**
783       \brief      \brief
784       Calculate the integral over the function space domain.       Calculate the integral over the function space domain as a python tuple.
785      */
786      ESCRIPT_DLL_API
787      boost::python::object
788      integrateToTuple_const() const;
789    
790    
791      /**
792        \brief
793         Calculate the integral over the function space domain as a python tuple.
794    */    */
795    boost::python::numeric::array    ESCRIPT_DLL_API
796    integrate() const;    boost::python::object
797      integrateToTuple();
798    
799    
800    
801    /**    /**
802       \brief       \brief
803         Returns 1./ Data object
804         *
805      */
806      ESCRIPT_DLL_API
807      Data
808      oneOver() const;
809      /**
810         \brief
811       Return a Data with a 1 for +ive values and a 0 for 0 or -ive values.       Return a Data with a 1 for +ive values and a 0 for 0 or -ive values.
812         *
813    */    */
814      ESCRIPT_DLL_API
815    Data    Data
816    wherePositive() const;    wherePositive() const;
817    
818    /**    /**
819       \brief       \brief
820       Return a Data with a 1 for -ive values and a 0 for +ive or 0 values.       Return a Data with a 1 for -ive values and a 0 for +ive or 0 values.
821         *
822    */    */
823      ESCRIPT_DLL_API
824    Data    Data
825    whereNegative() const;    whereNegative() const;
826    
827    /**    /**
828       \brief       \brief
829       Return a Data with a 1 for +ive or 0 values and a 0 for -ive values.       Return a Data with a 1 for +ive or 0 values and a 0 for -ive values.
830         *
831    */    */
832      ESCRIPT_DLL_API
833    Data    Data
834    whereNonNegative() const;    whereNonNegative() const;
835    
836    /**    /**
837       \brief       \brief
838       Return a Data with a 1 for -ive or 0 values and a 0 for +ive values.       Return a Data with a 1 for -ive or 0 values and a 0 for +ive values.
839         *
840    */    */
841      ESCRIPT_DLL_API
842    Data    Data
843    whereNonPositive() const;    whereNonPositive() const;
844    
845    /**    /**
846       \brief       \brief
847       Return a Data with a 1 for 0 values and a 0 for +ive or -ive values.       Return a Data with a 1 for 0 values and a 0 for +ive or -ive values.
848         *
849    */    */
850      ESCRIPT_DLL_API
851    Data    Data
852    whereZero() const;    whereZero(double tol=0.0) const;
853    
854    /**    /**
855       \brief       \brief
856       Return a Data with a 0 for 0 values and a 1 for +ive or -ive values.       Return a Data with a 0 for 0 values and a 1 for +ive or -ive values.
857         *
858      */
859      ESCRIPT_DLL_API
860      Data
861      whereNonZero(double tol=0.0) const;
862    
863      /**
864         \brief
865         Return the maximum absolute value of this Data object.
866    
867         The method is not const because lazy data needs to be expanded before Lsup can be computed.
868         The _const form can be used when the Data object is const, however this will only work for
869         Data which is not Lazy.
870    
871         For Data which contain no samples (or tagged Data for which no tags in use have a value)
872         zero is returned.
873      */
874      ESCRIPT_DLL_API
875      double
876      Lsup();
877    
878      ESCRIPT_DLL_API
879      double
880      Lsup_const() const;
881    
882    
883      /**
884         \brief
885         Return the maximum value of this Data object.
886    
887         The method is not const because lazy data needs to be expanded before sup can be computed.
888         The _const form can be used when the Data object is const, however this will only work for
889         Data which is not Lazy.
890    
891         For Data which contain no samples (or tagged Data for which no tags in use have a value)
892         a large negative value is returned.
893      */
894      ESCRIPT_DLL_API
895      double
896      sup();
897    
898      ESCRIPT_DLL_API
899      double
900      sup_const() const;
901    
902    
903      /**
904         \brief
905         Return the minimum value of this Data object.
906    
907         The method is not const because lazy data needs to be expanded before inf can be computed.
908         The _const form can be used when the Data object is const, however this will only work for
909         Data which is not Lazy.
910    
911         For Data which contain no samples (or tagged Data for which no tags in use have a value)
912         a large positive value is returned.
913      */
914      ESCRIPT_DLL_API
915      double
916      inf();
917    
918      ESCRIPT_DLL_API
919      double
920      inf_const() const;
921    
922    
923    
924      /**
925         \brief
926         Return the absolute value of each data point of this Data object.
927         *
928      */
929      ESCRIPT_DLL_API
930      Data
931      abs() const;
932    
933      /**
934         \brief
935         Return the maximum value of each data point of this Data object.
936         *
937      */
938      ESCRIPT_DLL_API
939      Data
940      maxval() const;
941    
942      /**
943         \brief
944         Return the minimum value of each data point of this Data object.
945         *
946      */
947      ESCRIPT_DLL_API
948      Data
949      minval() const;
950    
951      /**
952         \brief
953         Return the (sample number, data-point number) of the data point with
954         the minimum value in this Data object.
955      */
956      ESCRIPT_DLL_API
957      const boost::python::tuple
958      minGlobalDataPoint() const;
959    
960      ESCRIPT_DLL_API
961      void
962      calc_minGlobalDataPoint(int& ProcNo,  int& DataPointNo) const;
963      /**
964         \brief
965         Return the sign of each data point of this Data object.
966         -1 for negative values, zero for zero values, 1 for positive values.
967         *
968      */
969      ESCRIPT_DLL_API
970      Data
971      sign() const;
972    
973      /**
974         \brief
975         Return the symmetric part of a matrix which is half the matrix plus its transpose.
976         *
977      */
978      ESCRIPT_DLL_API
979      Data
980      symmetric() const;
981    
982      /**
983         \brief
984         Return the nonsymmetric part of a matrix which is half the matrix minus its transpose.
985         *
986      */
987      ESCRIPT_DLL_API
988      Data
989      nonsymmetric() const;
990    
991      /**
992         \brief
993         Return the trace of a matrix
994         *
995      */
996      ESCRIPT_DLL_API
997      Data
998      trace(int axis_offset) const;
999    
1000      /**
1001         \brief
1002         Transpose each data point of this Data object around the given axis.
1003         *
1004      */
1005      ESCRIPT_DLL_API
1006      Data
1007      transpose(int axis_offset) const;
1008    
1009      /**
1010         \brief
1011         Return the eigenvalues of the symmetric part at each data point of this Data object in increasing values.
1012         Currently this function is restricted to rank 2, square shape, and dimension 3.
1013         *
1014      */
1015      ESCRIPT_DLL_API
1016      Data
1017      eigenvalues() const;
1018    
1019      /**
1020         \brief
1021         Return the eigenvalues and corresponding eigenvcetors of the symmetric part at each data point of this Data object.
1022         the eigenvalues are ordered in increasing size where eigenvalues with relative difference less than
1023         tol are treated as equal. The eigenvectors are orthogonal, normalized and the sclaed such that the
1024         first non-zero entry is positive.
1025         Currently this function is restricted to rank 2, square shape, and dimension 3
1026         *
1027      */
1028      ESCRIPT_DLL_API
1029      const boost::python::tuple
1030      eigenvalues_and_eigenvectors(const double tol=1.e-12) const;
1031    
1032      /**
1033         \brief
1034         swaps the components axis0 and axis1
1035         *
1036      */
1037      ESCRIPT_DLL_API
1038      Data
1039      swapaxes(const int axis0, const int axis1) const;
1040    
1041      /**
1042         \brief
1043         Return the error function erf of each data point of this Data object.
1044         *
1045    */    */
1046      ESCRIPT_DLL_API
1047    Data    Data
1048    whereNonZero() const;    erf() const;
1049    
1050    /**    /**
1051       \brief       \brief
1052       Return the sin of each data point of this Data object.       Return the sin of each data point of this Data object.
1053         *
1054    */    */
1055      ESCRIPT_DLL_API
1056    Data    Data
1057    sin() const;    sin() const;
1058    
1059    /**    /**
1060       \brief       \brief
1061       Return the cos of each data point of this Data object.       Return the cos of each data point of this Data object.
1062         *
1063    */    */
1064      ESCRIPT_DLL_API
1065    Data    Data
1066    cos() const;    cos() const;
1067    
1068    /**    /**
1069       \brief       \brief
1070       Return the tan of each data point of this Data object.       Return the tan of each data point of this Data object.
1071         *
1072    */    */
1073      ESCRIPT_DLL_API
1074    Data    Data
1075    tan() const;    tan() const;
1076    
1077    /**    /**
1078       \brief       \brief
1079       Return the log to base 10 of each data point of this Data object.       Return the asin of each data point of this Data object.
1080         *
1081    */    */
1082      ESCRIPT_DLL_API
1083    Data    Data
1084    log() const;    asin() const;
1085    
1086    /**    /**
1087       \brief       \brief
1088       Return the natural log of each data point of this Data object.       Return the acos of each data point of this Data object.
1089         *
1090    */    */
1091      ESCRIPT_DLL_API
1092    Data    Data
1093    ln() const;    acos() const;
1094    
1095    /**    /**
1096       \brief       \brief
1097       Return the maximum absolute value of this Data object.       Return the atan of each data point of this Data object.
1098         *
1099    */    */
1100    double    ESCRIPT_DLL_API
1101    Lsup() const;    Data
1102      atan() const;
1103    
1104    /**    /**
1105       \brief       \brief
1106       Return the maximum value of this Data object.       Return the sinh of each data point of this Data object.
1107         *
1108    */    */
1109    double    ESCRIPT_DLL_API
1110    sup() const;    Data
1111      sinh() const;
1112    
1113    /**    /**
1114       \brief       \brief
1115       Return the minimum value of this Data object.       Return the cosh of each data point of this Data object.
1116         *
1117    */    */
1118    double    ESCRIPT_DLL_API
1119    inf() const;    Data
1120      cosh() const;
1121    
1122    /**    /**
1123       \brief       \brief
1124       Return the absolute value of each data point of this Data object.       Return the tanh of each data point of this Data object.
1125         *
1126    */    */
1127      ESCRIPT_DLL_API
1128    Data    Data
1129    abs() const;    tanh() const;
1130    
1131    /**    /**
1132       \brief       \brief
1133       Return the maximum value of each data point of this Data object.       Return the asinh of each data point of this Data object.
1134         *
1135    */    */
1136      ESCRIPT_DLL_API
1137    Data    Data
1138    maxval() const;    asinh() const;
1139    
1140    /**    /**
1141       \brief       \brief
1142       Return the minimum value of each data point of this Data object.       Return the acosh of each data point of this Data object.
1143         *
1144    */    */
1145      ESCRIPT_DLL_API
1146    Data    Data
1147    minval() const;    acosh() const;
1148    
1149    /**    /**
1150       \brief       \brief
1151       Return the length of each data point of this Data object.       Return the atanh of each data point of this Data object.
1152       sqrt(sum(A[i,j,k,l]^2))       *
1153    */    */
1154      ESCRIPT_DLL_API
1155    Data    Data
1156    length() const;    atanh() const;
1157    
1158    /**    /**
1159       \brief       \brief
1160       Return the sign of each data point of this Data object.       Return the log to base 10 of each data point of this Data object.
1161       -1 for negative values, zero for zero values, 1 for positive values.       *
1162    */    */
1163      ESCRIPT_DLL_API
1164    Data    Data
1165    sign() const;    log10() const;
1166    
1167      /**
1168         \brief
1169         Return the natural log of each data point of this Data object.
1170         *
1171      */
1172      ESCRIPT_DLL_API
1173      Data
1174      log() const;
1175    
1176    /**    /**
1177      \transpose       \brief
1178      Transpose each data point of this Data object around the given axis.       Return the exponential function of each data point of this Data object.
1179         *
1180    */    */
1181      ESCRIPT_DLL_API
1182    Data    Data
1183    transpose(int axis) const;    exp() const;
1184    
1185    /**    /**
1186      \trace       \brief
1187      Calculate the trace of each data point of this Data object.       Return the square root of each data point of this Data object.
1188      sum(A[i,i,i,i])       *
1189    */    */
1190      ESCRIPT_DLL_API
1191    Data    Data
1192    trace() const;    sqrt() const;
1193    
1194    /**    /**
1195      \exp       \brief
1196      Return the exponential function of each data point of this Data object.       Return the negation of each data point of this Data object.
1197         *
1198    */    */
1199      ESCRIPT_DLL_API
1200    Data    Data
1201    exp() const;    neg() const;
1202    
1203    /**    /**
1204      \sqrt       \brief
1205      Return the square root of each data point of this Data object.       Return the identity of each data point of this Data object.
1206         Simply returns this object unmodified.
1207         *
1208    */    */
1209      ESCRIPT_DLL_API
1210    Data    Data
1211    sqrt() const;    pos() const;
1212    
1213    /**    /**
1214       \brief       \brief
1215       Return the given power of each data point of this Data object.       Return the given power of each data point of this Data object.
1216    
1217         \param right Input - the power to raise the object to.
1218         *
1219    */    */
1220      ESCRIPT_DLL_API
1221    Data    Data
1222    powD(const Data& right) const;    powD(const Data& right) const;
1223    
1224      /**
1225         \brief
1226         Return the given power of each data point of this boost python object.
1227    
1228         \param right Input - the power to raise the object to.
1229         *
1230       */
1231      ESCRIPT_DLL_API
1232    Data    Data
1233    powO(const boost::python::object& right) const;    powO(const boost::python::object& right) const;
1234    
1235    /**    /**
1236      \brief       \brief
1237      Return the negation of each data point of this Data object.       Return the given power of each data point of this boost python object.
1238    */  
1239         \param left Input - the bases
1240         *
1241       */
1242    
1243      ESCRIPT_DLL_API
1244    Data    Data
1245    neg() const;    rpowO(const boost::python::object& left) const;
1246    
1247    /**    /**
1248      \brief       \brief
1249      Return the identity of each data point of this Data object.       writes the object to a file in the DX file format
     Simply returns this object unmodified.  
1250    */    */
1251    Data    ESCRIPT_DLL_API
1252    pos() const;    void
1253      saveDX(std::string fileName) const;
1254    
1255      /**
1256         \brief
1257         writes the object to a file in the VTK file format
1258      */
1259      ESCRIPT_DLL_API
1260      void
1261      saveVTK(std::string fileName) const;
1262    
1263    /**    /**
1264       \brief       \brief
1265       Overloaded operator +=       Overloaded operator +=
1266       \param right - Input - The right hand side.       \param right - Input - The right hand side.
1267         *
1268    */    */
1269      ESCRIPT_DLL_API
1270    Data& operator+=(const Data& right);    Data& operator+=(const Data& right);
1271      ESCRIPT_DLL_API
1272    Data& operator+=(const boost::python::object& right);    Data& operator+=(const boost::python::object& right);
1273    
1274      ESCRIPT_DLL_API
1275      Data& operator=(const Data& other);
1276    
1277    /**    /**
1278       \brief       \brief
1279       Overloaded operator -=       Overloaded operator -=
1280       \param right - Input - The right hand side.       \param right - Input - The right hand side.
1281         *
1282    */    */
1283      ESCRIPT_DLL_API
1284    Data& operator-=(const Data& right);    Data& operator-=(const Data& right);
1285      ESCRIPT_DLL_API
1286    Data& operator-=(const boost::python::object& right);    Data& operator-=(const boost::python::object& right);
1287    
1288   /**   /**
1289       \brief       \brief
1290       Overloaded operator *=       Overloaded operator *=
1291       \param right - Input - The right hand side.       \param right - Input - The right hand side.
1292         *
1293    */    */
1294      ESCRIPT_DLL_API
1295    Data& operator*=(const Data& right);    Data& operator*=(const Data& right);
1296      ESCRIPT_DLL_API
1297    Data& operator*=(const boost::python::object& right);    Data& operator*=(const boost::python::object& right);
1298    
1299   /**   /**
1300       \brief       \brief
1301       Overloaded operator /=       Overloaded operator /=
1302       \param right - Input - The right hand side.       \param right - Input - The right hand side.
1303         *
1304    */    */
1305      ESCRIPT_DLL_API
1306    Data& operator/=(const Data& right);    Data& operator/=(const Data& right);
1307      ESCRIPT_DLL_API
1308    Data& operator/=(const boost::python::object& right);    Data& operator/=(const boost::python::object& right);
1309    
1310    /**    /**
1311       \brief       \brief
1312       Returns true if this can be interpolated to functionspace.       Returns true if this can be interpolated to functionspace.
1313    */    */
1314      ESCRIPT_DLL_API
1315    bool    bool
1316    probeInterpolation(const FunctionSpace& functionspace) const;    probeInterpolation(const FunctionSpace& functionspace) const;
1317    
# Line 748  class Data { Line 1330  class Data {
1330       \param key - Input - python slice tuple specifying       \param key - Input - python slice tuple specifying
1331       slice to return.       slice to return.
1332    */    */
1333      ESCRIPT_DLL_API
1334    Data    Data
1335    getItem(const boost::python::object& key) const;    getItem(const boost::python::object& key) const;
1336    
# Line 755  class Data { Line 1338  class Data {
1338       \brief       \brief
1339       Copies slice from value into this Data object.       Copies slice from value into this Data object.
1340    
      \description  
1341       Implements the [] set operator in python.       Implements the [] set operator in python.
1342       Calls setSlice.       Calls setSlice.
1343    
# Line 763  class Data { Line 1345  class Data {
1345       slice to copy from value.       slice to copy from value.
1346       \param value - Input - Data object to copy from.       \param value - Input - Data object to copy from.
1347    */    */
1348      ESCRIPT_DLL_API
1349    void    void
1350    setItemD(const boost::python::object& key,    setItemD(const boost::python::object& key,
1351             const Data& value);             const Data& value);
1352    
1353      ESCRIPT_DLL_API
1354    void    void
1355    setItemO(const boost::python::object& key,    setItemO(const boost::python::object& key,
1356             const boost::python::object& value);             const boost::python::object& value);
# Line 779  class Data { Line 1363  class Data {
1363       this Data object.       this Data object.
1364    */    */
1365    template <class UnaryFunction>    template <class UnaryFunction>
1366      ESCRIPT_DLL_API
1367    inline    inline
1368    void    void
1369    unaryOp(UnaryFunction operation);    unaryOp2(UnaryFunction operation);
1370    
1371    /**    /**
1372       \brief       \brief
1373       Return a Data object containing the specified slice of       Return a Data object containing the specified slice of
1374       this Data object.       this Data object.
1375       \param region - Input - Region to copy.       \param region - Input - Region to copy.
1376         *
1377    */    */
1378      ESCRIPT_DLL_API
1379    Data    Data
1380    getSlice(const DataArrayView::RegionType& region) const;    getSlice(const DataTypes::RegionType& region) const;
1381    
1382    /**    /**
1383       \brief       \brief
# Line 798  class Data { Line 1385  class Data {
1385       Data object.       Data object.
1386       \param value - Input - Data to copy from.       \param value - Input - Data to copy from.
1387       \param region - Input - Region to copy.       \param region - Input - Region to copy.
1388         *
1389    */    */
1390      ESCRIPT_DLL_API
1391    void    void
1392    setSlice(const Data& value,    setSlice(const Data& value,
1393             const DataArrayView::RegionType& region);             const DataTypes::RegionType& region);
1394    
1395      /**
1396         \brief
1397         print the data values to stdout. Used for debugging
1398      */
1399      ESCRIPT_DLL_API
1400      void
1401            print(void);
1402    
1403      /**
1404         \brief
1405         return the MPI rank number of the local data
1406                     MPI_COMM_WORLD is assumed and the result of MPI_Comm_size()
1407                     is returned
1408      */
1409      ESCRIPT_DLL_API
1410            int
1411            get_MPIRank(void) const;
1412    
1413      /**
1414         \brief
1415         return the MPI rank number of the local data
1416                     MPI_COMM_WORLD is assumed and the result of MPI_Comm_rank()
1417                     is returned
1418      */
1419      ESCRIPT_DLL_API
1420            int
1421            get_MPISize(void) const;
1422    
1423      /**
1424         \brief
1425         return the MPI rank number of the local data
1426                     MPI_COMM_WORLD is assumed and returned.
1427      */
1428      ESCRIPT_DLL_API
1429            MPI_Comm
1430            get_MPIComm(void) const;
1431    
1432      /**
1433         \brief
1434         return the object produced by the factory, which is a DataConstant or DataExpanded
1435        TODO Ownership of this object should be explained in doco.
1436      */
1437      ESCRIPT_DLL_API
1438            DataAbstract*
1439            borrowData(void) const;
1440    
1441      ESCRIPT_DLL_API
1442            DataAbstract_ptr
1443            borrowDataPtr(void) const;
1444    
1445      ESCRIPT_DLL_API
1446            DataReady_ptr
1447            borrowReadyPtr(void) const;
1448    
1449    
1450    
1451      /**
1452         \brief
1453         Return a pointer to the beginning of the datapoint at the specified offset.
1454         TODO Eventually these should be inlined.
1455         \param i - position(offset) in the underlying datastructure
1456      */
1457    
1458      ESCRIPT_DLL_API
1459            DataTypes::ValueType::const_reference
1460            getDataAtOffsetRO(DataTypes::ValueType::size_type i);
1461    
1462    
1463      ESCRIPT_DLL_API
1464            DataTypes::ValueType::reference
1465            getDataAtOffsetRW(DataTypes::ValueType::size_type i);
1466    
1467    
1468    
1469    /**
1470       \brief Create a buffer for use by getSample
1471       Allocates a DataVector large enough for DataLazy::resolveSample to operate on for the current Data.
1472       Do not use this buffer for other Data instances (unless you are sure they will be the same size).
1473      
1474       In multi-threaded sections, this needs to be called on each thread.
1475    
1476       \return A BufferGroup* if Data is lazy, NULL otherwise.
1477       \warning This pointer must be deallocated using freeSampleBuffer to avoid cross library memory issues.
1478    */
1479      ESCRIPT_DLL_API
1480      BufferGroup*
1481      allocSampleBuffer() const;
1482    
1483    /**
1484       \brief Free a buffer allocated with allocSampleBuffer.
1485       \param buffer Input - pointer to the buffer to deallocate.
1486    */
1487    ESCRIPT_DLL_API void freeSampleBuffer(BufferGroup* buffer);
1488    
1489   protected:   protected:
1490    
1491   private:   private:
1492    
1493      double
1494      LsupWorker() const;
1495    
1496      double
1497      supWorker() const;
1498    
1499      double
1500      infWorker() const;
1501    
1502      boost::python::object
1503      integrateWorker() const;
1504    
1505    /**    /**
1506       \brief       \brief
1507       Check *this and the right operand are compatible. Throws       Check *this and the right operand are compatible. Throws
# Line 823  class Data { Line 1518  class Data {
1518    /**    /**
1519       \brief       \brief
1520       Perform the specified reduction algorithm on every element of every data point in       Perform the specified reduction algorithm on every element of every data point in
1521       this Data object and return the single double value result.       this Data object according to the given function and return the single value result.
1522    */    */
1523    template <class UnaryFunction>    template <class BinaryFunction>
1524    inline    inline
1525    double    double
1526    algorithm(UnaryFunction operation) const;    algorithm(BinaryFunction operation,
1527                double initial_value) const;
1528    
1529    /**    /**
1530       \brief       \brief
1531       Perform the given binary operation on all of the data's elements.       Reduce each data-point in this Data object using the given operation. Return a Data
1532       The underlying type of the right hand side (right) determines the final       object with the same number of data-points, but with each data-point containing only
1533       type of *this after the operation. For example if the right hand side       one value - the result of the reduction operation on the corresponding data-point in
1534       is expanded *this will be expanded if necessary.       this Data object
      RHS is a Data object.  
1535    */    */
1536    template <class BinaryFunction>    template <class BinaryFunction>
1537    inline    inline
1538    void    Data
1539    binaryOp(const Data& right,    dp_algorithm(BinaryFunction operation,
1540             BinaryFunction operation);                 double initial_value) const;
1541    
1542    /**    /**
1543       \brief       \brief
1544       Perform the given binary operation on all of the data's elements.       Perform the given binary operation on all of the data's elements.
1545       RHS is a boost::python object.       The underlying type of the right hand side (right) determines the final
1546         type of *this after the operation. For example if the right hand side
1547         is expanded *this will be expanded if necessary.
1548         RHS is a Data object.
1549    */    */
1550    template <class BinaryFunction>    template <class BinaryFunction>
1551    inline    inline
1552    void    void
1553    binaryOp(const boost::python::object& right,    binaryOp(const Data& right,
1554             BinaryFunction operation);             BinaryFunction operation);
1555    
1556    /**    /**
# Line 875  class Data { Line 1573  class Data {
1573       \brief       \brief
1574       Construct a Data object of the appropriate type.       Construct a Data object of the appropriate type.
1575    */    */
1576    template <class IValueType>  
1577    void    void
1578    initialise(const IValueType& value,    initialise(const DataTypes::ValueType& value,
1579             const DataTypes::ShapeType& shape,
1580               const FunctionSpace& what,               const FunctionSpace& what,
1581               bool expanded);               bool expanded);
1582    
   /**  
      \brief  
      Reshape the data point if the data point is currently rank 0.  
      Will throw an exception if the data points are not rank 0.  
      The original data point value is used for all values of the new  
      data point.  
   */  
1583    void    void
1584    reshapeDataPoint(const DataArrayView::ShapeType& shape);    initialise(const WrappedArray& value,
1585                     const FunctionSpace& what,
1586                     bool expanded);
1587    
1588      //
1589      // flag to protect the data object against any update
1590      bool m_protected;
1591      mutable bool m_shared;
1592      bool m_lazy;
1593    
1594    //    //
1595    // pointer to the actual data object    // pointer to the actual data object
1596    boost::shared_ptr<DataAbstract> m_data;  //   boost::shared_ptr<DataAbstract> m_data;
1597      DataAbstract_ptr m_data;
1598    
1599  };  // If possible please use getReadyPtr instead.
1600    // But see warning below.
1601      const DataReady*
1602      getReady() const;
1603    
1604  template <class IValueType>    DataReady*
1605  void    getReady();
1606  Data::initialise(const IValueType& value,  
1607                   const FunctionSpace& what,  
1608                   bool expanded)  // Be wary of using this for local operations since it (temporarily) increases reference count.
1609  {  // If you are just using this to call a method on DataReady instead of DataAbstract consider using
1610    // getReady() instead
1611      DataReady_ptr
1612      getReadyPtr();
1613    
1614      const_DataReady_ptr
1615      getReadyPtr() const;
1616    
1617    
1618      /**
1619       \brief Update the Data's shared flag
1620       This indicates that the DataAbstract used by this object is now shared (or no longer shared).
1621       For internal use only.
1622      */
1623      void updateShareStatus(bool nowshared) const
1624      {
1625        m_shared=nowshared;     // m_shared is mutable
1626      }
1627    
1628      // In the isShared() method below:
1629      // A problem would occur if m_data (the address pointed to) were being modified
1630      // while the call m_data->is_shared is being executed.
1631      //
1632      // Q: So why do I think this code can be thread safe/correct?
1633      // A: We need to make some assumptions.
1634      // 1. We assume it is acceptable to return true under some conditions when we aren't shared.
1635      // 2. We assume that no constructions or assignments which will share previously unshared
1636      //    will occur while this call is executing. This is consistent with the way Data:: and C are written.
1637    //    //
1638    // Construct a Data object of the appropriate type.    // This means that the only transition we need to consider, is when a previously shared object is
1639    // Construct the object first as there seems to be a bug which causes    // not shared anymore. ie. the other objects have been destroyed or a deep copy has been made.
1640    // undefined behaviour if an exception is thrown during construction    // In those cases the m_shared flag changes to false after m_data has completed changing.
1641    // within the shared_ptr constructor.    // For any threads executing before the flag switches they will assume the object is still shared.
1642    if (expanded) {    bool isShared() const
1643      DataAbstract* temp=new DataExpanded(value,what);    {
1644      boost::shared_ptr<DataAbstract> temp_data(temp);      return m_shared;
1645      m_data=temp_data;  /*  if (m_shared) return true;
1646    } else {      if (m_data->isShared())        
1647      DataAbstract* temp=new DataConstant(value,what);      {                  
1648      boost::shared_ptr<DataAbstract> temp_data(temp);          updateShareStatus(true);
1649      m_data=temp_data;          return true;
1650        }
1651        return false;*/
1652      }
1653    
1654      void forceResolve()
1655      {
1656        if (isLazy())
1657        {
1658            #ifdef _OPENMP
1659            if (omp_in_parallel())
1660            {   // Yes this is throwing an exception out of an omp thread which is forbidden.
1661            throw DataException("Please do not call forceResolve() in a parallel region.");
1662            }
1663            #endif
1664            resolve();
1665        }
1666      }
1667    
1668      /**
1669      \brief if another object is sharing out member data make a copy to work with instead.
1670      This code should only be called from single threaded sections of code.
1671      */
1672      void exclusiveWrite()
1673      {
1674    #ifdef _OPENMP
1675      if (omp_in_parallel())
1676      {
1677    // *((int*)0)=17;
1678        throw DataException("Programming error. Please do not run exclusiveWrite() in multi-threaded sections.");
1679      }
1680    #endif
1681        forceResolve();
1682        if (isShared())
1683        {
1684            DataAbstract* t=m_data->deepCopy();
1685            set_m_data(DataAbstract_ptr(t));
1686        }
1687      }
1688    
1689      /**
1690      \brief checks if caller can have exclusive write to the object
1691      */
1692      void checkExclusiveWrite()
1693      {
1694        if  (isLazy() || isShared())
1695        {
1696            throw DataException("Programming error. ExclusiveWrite required - please call requireWrite()");
1697        }
1698    }    }
1699    
1700      /**
1701      \brief Modify the data abstract hosted by this Data object
1702      For internal use only.
1703      Passing a pointer to null is permitted (do this in the destructor)
1704      \warning Only to be called in single threaded code or inside a single/critical section. This method needs to be atomic.
1705      */
1706      void set_m_data(DataAbstract_ptr p);
1707    
1708      friend class DataAbstract;        // To allow calls to updateShareStatus
1709    
1710    };
1711    
1712    }   // end namespace escript
1713    
1714    
1715    // No, this is not supposed to be at the top of the file
1716    // DataAbstact needs to be declared first, then DataReady needs to be fully declared
1717    // so that I can dynamic cast between them below.
1718    #include "DataReady.h"
1719    #include "DataLazy.h"
1720    
1721    namespace escript
1722    {
1723    
1724    inline
1725    const DataReady*
1726    Data::getReady() const
1727    {
1728       const DataReady* dr=dynamic_cast<const DataReady*>(m_data.get());
1729       EsysAssert((dr!=0), "Error - casting to DataReady.");
1730       return dr;
1731    }
1732    
1733    inline
1734    DataReady*
1735    Data::getReady()
1736    {
1737       DataReady* dr=dynamic_cast<DataReady*>(m_data.get());
1738       EsysAssert((dr!=0), "Error - casting to DataReady.");
1739       return dr;
1740    }
1741    
1742    // Be wary of using this for local operations since it (temporarily) increases reference count.
1743    // If you are just using this to call a method on DataReady instead of DataAbstract consider using
1744    // getReady() instead
1745    inline
1746    DataReady_ptr
1747    Data::getReadyPtr()
1748    {
1749       DataReady_ptr dr=boost::dynamic_pointer_cast<DataReady>(m_data);
1750       EsysAssert((dr.get()!=0), "Error - casting to DataReady.");
1751       return dr;
1752  }  }
1753    
1754    
1755    inline
1756    const_DataReady_ptr
1757    Data::getReadyPtr() const
1758    {
1759       const_DataReady_ptr dr=boost::dynamic_pointer_cast<const DataReady>(m_data);
1760       EsysAssert((dr.get()!=0), "Error - casting to DataReady.");
1761       return dr;
1762    }
1763    
1764    inline
1765    DataAbstract::ValueType::value_type*
1766    Data::getSampleDataRW(DataAbstract::ValueType::size_type sampleNo)
1767    {
1768    //    if (isLazy())
1769    //    {
1770    //  resolve();
1771    //    }
1772    //    exclusiveWrite();
1773       if (isLazy())
1774       {
1775        throw DataException("Error, attempt to acquire RW access to lazy data. Please call requireWrite() first.");
1776       }
1777       return getReady()->getSampleDataRW(sampleNo);
1778    }
1779    
1780    // inline
1781    // const DataAbstract::ValueType::value_type*
1782    // Data::getSampleDataRO(DataAbstract::ValueType::size_type sampleNo, DataTypes::ValueType* buffer)
1783    // {
1784    //    DataLazy* l=dynamic_cast<DataLazy*>(m_data.get());
1785    //    if (l!=0)
1786    //    {
1787    //  size_t offset=0;
1788    //  if (buffer==NULL)
1789    //  {
1790    //      throw DataException("Error, attempt to getSampleDataRO for lazy Data with buffer==NULL");
1791    //  }
1792    //  const DataTypes::ValueType* res=l->resolveSample(*buffer,0,sampleNo,offset);
1793    //  return &((*res)[offset]);
1794    //    }
1795    //    return getReady()->getSampleDataRO(sampleNo);
1796    // }
1797    
1798    inline
1799    const DataAbstract::ValueType::value_type*
1800    Data::getSampleDataRO(DataAbstract::ValueType::size_type sampleNo, BufferGroup* bufferg)
1801    {
1802       DataLazy* l=dynamic_cast<DataLazy*>(m_data.get());
1803       if (l!=0)
1804       {
1805        size_t offset=0;
1806        if (bufferg==NULL)
1807        {
1808            throw DataException("Error, attempt to getSampleDataRO for lazy Data with buffer==NULL");
1809        }
1810        const DataTypes::ValueType* res=l->resolveSample(*bufferg,sampleNo,offset);
1811        return &((*res)[offset]);
1812       }
1813       return getReady()->getSampleDataRO(sampleNo);
1814    }
1815    
1816    
1817    
1818    /**
1819       Modify a filename for MPI parallel output to multiple files
1820    */
1821    char *Escript_MPI_appendRankToFileName(const char *, int, int);
1822    
1823  /**  /**
1824     Binary Data object operators.     Binary Data object operators.
1825  */  */
1826    inline double rpow(double x,double y)
1827    {
1828        return pow(y,x);
1829    }
1830    
1831  /**  /**
1832    \brief    \brief
1833    Operator+    Operator+
1834    Takes two Data objects.    Takes two Data objects.
1835  */  */
1836  Data operator+(const Data& left, const Data& right);  ESCRIPT_DLL_API Data operator+(const Data& left, const Data& right);
1837    
1838  /**  /**
1839    \brief    \brief
1840    Operator-    Operator-
1841    Takes two Data objects.    Takes two Data objects.
1842  */  */
1843  Data operator-(const Data& left, const Data& right);  ESCRIPT_DLL_API Data operator-(const Data& left, const Data& right);
1844    
1845  /**  /**
1846    \brief    \brief
1847    Operator*    Operator*
1848    Takes two Data objects.    Takes two Data objects.
1849  */  */
1850  Data operator*(const Data& left, const Data& right);  ESCRIPT_DLL_API Data operator*(const Data& left, const Data& right);
1851    
1852  /**  /**
1853    \brief    \brief
1854    Operator/    Operator/
1855    Takes two Data objects.    Takes two Data objects.
1856  */  */
1857  Data operator/(const Data& left, const Data& right);  ESCRIPT_DLL_API Data operator/(const Data& left, const Data& right);
1858    
1859  /**  /**
1860    \brief    \brief
# Line 957  Data operator/(const Data& left, const D Line 1862  Data operator/(const Data& left, const D
1862    Takes LHS Data object and RHS python::object.    Takes LHS Data object and RHS python::object.
1863    python::object must be convertable to Data type.    python::object must be convertable to Data type.
1864  */  */
1865  Data operator+(const Data& left, const boost::python::object& right);  ESCRIPT_DLL_API Data operator+(const Data& left, const boost::python::object& right);
1866    
1867  /**  /**
1868    \brief    \brief
# Line 965  Data operator+(const Data& left, const b Line 1870  Data operator+(const Data& left, const b
1870    Takes LHS Data object and RHS python::object.    Takes LHS Data object and RHS python::object.
1871    python::object must be convertable to Data type.    python::object must be convertable to Data type.
1872  */  */
1873  Data operator-(const Data& left, const boost::python::object& right);  ESCRIPT_DLL_API Data operator-(const Data& left, const boost::python::object& right);
1874    
1875  /**  /**
1876    \brief    \brief
# Line 973  Data operator-(const Data& left, const b Line 1878  Data operator-(const Data& left, const b
1878    Takes LHS Data object and RHS python::object.    Takes LHS Data object and RHS python::object.
1879    python::object must be convertable to Data type.    python::object must be convertable to Data type.
1880  */  */
1881  Data operator*(const Data& left, const boost::python::object& right);  ESCRIPT_DLL_API Data operator*(const Data& left, const boost::python::object& right);
1882    
1883  /**  /**
1884    \brief    \brief
# Line 981  Data operator*(const Data& left, const b Line 1886  Data operator*(const Data& left, const b
1886    Takes LHS Data object and RHS python::object.    Takes LHS Data object and RHS python::object.
1887    python::object must be convertable to Data type.    python::object must be convertable to Data type.
1888  */  */
1889  Data operator/(const Data& left, const boost::python::object& right);  ESCRIPT_DLL_API Data operator/(const Data& left, const boost::python::object& right);
1890    
1891  /**  /**
1892    \brief    \brief
# Line 989  Data operator/(const Data& left, const b Line 1894  Data operator/(const Data& left, const b
1894    Takes LHS python::object and RHS Data object.    Takes LHS python::object and RHS Data object.
1895    python::object must be convertable to Data type.    python::object must be convertable to Data type.
1896  */  */
1897  Data operator+(const boost::python::object& left, const Data& right);  ESCRIPT_DLL_API Data operator+(const boost::python::object& left, const Data& right);
1898    
1899  /**  /**
1900    \brief    \brief
# Line 997  Data operator+(const boost::python::obje Line 1902  Data operator+(const boost::python::obje
1902    Takes LHS python::object and RHS Data object.    Takes LHS python::object and RHS Data object.
1903    python::object must be convertable to Data type.    python::object must be convertable to Data type.
1904  */  */
1905  Data operator-(const boost::python::object& left, const Data& right);  ESCRIPT_DLL_API Data operator-(const boost::python::object& left, const Data& right);
1906    
1907  /**  /**
1908    \brief    \brief
# Line 1005  Data operator-(const boost::python::obje Line 1910  Data operator-(const boost::python::obje
1910    Takes LHS python::object and RHS Data object.    Takes LHS python::object and RHS Data object.
1911    python::object must be convertable to Data type.    python::object must be convertable to Data type.
1912  */  */
1913  Data operator*(const boost::python::object& left, const Data& right);  ESCRIPT_DLL_API Data operator*(const boost::python::object& left, const Data& right);
1914    
1915  /**  /**
1916    \brief    \brief
# Line 1013  Data operator*(const boost::python::obje Line 1918  Data operator*(const boost::python::obje
1918    Takes LHS python::object and RHS Data object.    Takes LHS python::object and RHS Data object.
1919    python::object must be convertable to Data type.    python::object must be convertable to Data type.
1920  */  */
1921  Data operator/(const boost::python::object& left, const Data& right);  ESCRIPT_DLL_API Data operator/(const boost::python::object& left, const Data& right);
1922    
1923    
1924    
1925  /**  /**
1926    \brief    \brief
1927    Output operator    Output operator
1928  */  */
1929  std::ostream& operator<<(std::ostream& o, const Data& data);  ESCRIPT_DLL_API std::ostream& operator<<(std::ostream& o, const Data& data);
1930    
1931  /**  /**
1932    \brief    \brief
1933    Return true if operands are equivalent, else return false.    Compute a tensor product of two Data objects
1934    NB: this operator does very little at this point, and isn't to    \param arg0 - Input - Data object
1935    be relied on. Requires further implementation.    \param arg1 - Input - Data object
1936      \param axis_offset - Input - axis offset
1937      \param transpose - Input - 0: transpose neither, 1: transpose arg0, 2: transpose arg1
1938  */  */
1939  //bool operator==(const Data& left, const Data& right);  ESCRIPT_DLL_API
1940    Data
1941    C_GeneralTensorProduct(Data& arg0,
1942                         Data& arg1,
1943                         int axis_offset=0,
1944                         int transpose=0);
1945    
1946  /**  /**
1947    \brief    \brief
# Line 1042  Data::binaryOp(const Data& right, Line 1956  Data::binaryOp(const Data& right,
1956  {  {
1957     //     //
1958     // if this has a rank of zero promote it to the rank of the RHS     // if this has a rank of zero promote it to the rank of the RHS
1959     if (getPointDataView().getRank()==0 && right.getPointDataView().getRank()!=0) {     if (getDataPointRank()==0 && right.getDataPointRank()!=0) {
1960       reshapeDataPoint(right.getPointDataView().getShape());       throw DataException("Error - attempt to update rank zero object with object with rank bigger than zero.");
1961       }
1962    
1963       if (isLazy() || right.isLazy())
1964       {
1965         throw DataException("Programmer error - attempt to call binaryOp with Lazy Data.");
1966     }     }
1967     //     //
1968     // initially make the temporary a shallow copy     // initially make the temporary a shallow copy
1969     Data tempRight(right);     Data tempRight(right);
1970    
1971     if (getFunctionSpace()!=right.getFunctionSpace()) {     if (getFunctionSpace()!=right.getFunctionSpace()) {
1972       if (right.probeInterpolation(getFunctionSpace())) {       if (right.probeInterpolation(getFunctionSpace())) {
1973         //         //
1974         // an interpolation is required so create a new Data         // an interpolation is required so create a new Data
1975         tempRight=Data(right,this->getFunctionSpace());         tempRight=Data(right,this->getFunctionSpace());
1976       } else if (probeInterpolation(right.getFunctionSpace())) {       } else if (probeInterpolation(right.getFunctionSpace())) {
1977         //         //
1978         // interpolate onto the RHS function space         // interpolate onto the RHS function space
1979         Data tempLeft(*this,right.getFunctionSpace());         Data tempLeft(*this,right.getFunctionSpace());
1980         m_data=tempLeft.m_data;  //        m_data=tempLeft.m_data;
1981           set_m_data(tempLeft.m_data);
1982       }       }
1983     }     }
1984     operandCheck(tempRight);     operandCheck(tempRight);
# Line 1073  Data::binaryOp(const Data& right, Line 1994  Data::binaryOp(const Data& right,
1994       // of any data type       // of any data type
1995       DataExpanded* leftC=dynamic_cast<DataExpanded*>(m_data.get());       DataExpanded* leftC=dynamic_cast<DataExpanded*>(m_data.get());
1996       EsysAssert((leftC!=0), "Programming error - casting to DataExpanded.");       EsysAssert((leftC!=0), "Programming error - casting to DataExpanded.");
1997       escript::binaryOp(*leftC,*(tempRight.m_data.get()),operation);       escript::binaryOp(*leftC,*(tempRight.getReady()),operation);
1998     } else if (isTagged()) {     } else if (isTagged()) {
1999       //       //
2000       // Tagged data is operated on serially, the right hand side can be       // Tagged data is operated on serially, the right hand side can be
# Line 1099  Data::binaryOp(const Data& right, Line 2020  Data::binaryOp(const Data& right,
2020    
2021  /**  /**
2022    \brief    \brief
   Perform the given binary operation with this and right as operands.  
   Right is a boost::python object.  
 */  
 template <class BinaryFunction>  
 inline  
 void  
 Data::binaryOp(const boost::python::object& right,  
                BinaryFunction operation)  
 {  
    DataArray temp(right);  
    //  
    // if this has a rank of zero promote it to the rank of the RHS.  
    if (getPointDataView().getRank()==0 && temp.getView().getRank()!=0) {  
       reshapeDataPoint(temp.getView().getShape());  
    }  
    //  
    // Always allow scalar values for the RHS but check other shapes  
    if (temp.getView().getRank()!=0) {  
      if (!getPointDataView().checkShape(temp.getView().getShape())) {  
        throw DataException(getPointDataView().createShapeErrorMessage(  
                   "Error - RHS shape doesn't match LHS shape.",temp.getView().getShape()));  
      }  
    }  
    if (isExpanded()) {  
      DataExpanded* leftC=dynamic_cast<DataExpanded*>(m_data.get());  
      EsysAssert((leftC!=0),"Programming error - casting to DataExpanded.");  
      escript::binaryOp(*leftC,temp.getView(),operation);  
    } else if (isTagged()) {  
      DataTagged* leftC=dynamic_cast<DataTagged*>(m_data.get());  
      EsysAssert((leftC!=0), "Programming error - casting to DataTagged.");  
      escript::binaryOp(*leftC,temp.getView(),operation);  
    } else if (isConstant()) {  
      DataConstant* leftC=dynamic_cast<DataConstant*>(m_data.get());  
      EsysAssert((leftC!=0),"Programming error - casting to DataConstant.");  
      escript::binaryOp(*leftC,temp.getView(),operation);  
    }  
 }  
   
 /**  
   \brief  
   Perform the given unary operation on other and return the result.  
   Given operation is performed on each element of each data point, thus  
   argument object is a rank n Data object, and returned object is a rank n  
   Data object.  
   Calls Data::unaryOp.  
 */  
 template <class UnaryFunction>  
 inline  
 Data  
 unaryOp(const Data& other,  
         UnaryFunction operation)  
 {  
   Data result;  
   result.copy(other);  
   result.unaryOp(operation);  
   return result;  
 }  
   
 /**  
   \brief  
   Perform the given unary operation on this.  
   Given operation is performed on each element of each data point.  
   Calls escript::unaryOp.  
 */  
 template <class UnaryFunction>  
 inline  
 void  
 Data::unaryOp(UnaryFunction operation)  
 {  
   if (isExpanded()) {  
     DataExpanded* leftC=dynamic_cast<DataExpanded*>(m_data.get());  
     EsysAssert((leftC!=0), "Programming error - casting to DataExpanded.");  
     escript::unaryOp(*leftC,operation);  
   } else if (isTagged()) {  
     DataTagged* leftC=dynamic_cast<DataTagged*>(m_data.get());  
     EsysAssert((leftC!=0), "Programming error - casting to DataTagged.");  
     escript::unaryOp(*leftC,operation);  
   } else if (isConstant()) {  
     DataConstant* leftC=dynamic_cast<DataConstant*>(m_data.get());  
     EsysAssert((leftC!=0), "Programming error - casting to DataConstant.");  
     escript::unaryOp(*leftC,operation);  
   }  
 }  
   
 /**  
   \brief  
2023    Perform the given Data object reduction algorithm on this and return the result.    Perform the given Data object reduction algorithm on this and return the result.
2024    Given operation combines each element of each data point, thus argument    Given operation combines each element of each data point, thus argument
2025    object (*this) is a rank n Data object, and returned object is a scalar.    object (*this) is a rank n Data object, and returned object is a scalar.
2026    Calls escript::algorithm.    Calls escript::algorithm.
2027  */  */
2028  template <class UnaryFunction>  template <class BinaryFunction>
2029  inline  inline
2030  double  double
2031  Data::algorithm(UnaryFunction operation) const  Data::algorithm(BinaryFunction operation, double initial_value) const
2032  {  {
2033    if (isExpanded()) {    if (isExpanded()) {
2034      DataExpanded* leftC=dynamic_cast<DataExpanded*>(m_data.get());      DataExpanded* leftC=dynamic_cast<DataExpanded*>(m_data.get());
2035      EsysAssert((leftC!=0), "Programming error - casting to DataExpanded.");      EsysAssert((leftC!=0), "Programming error - casting to DataExpanded.");
2036      return escript::algorithm(*leftC,operation);      return escript::algorithm(*leftC,operation,initial_value);
2037    } else if (isTagged()) {    } else if (isTagged()) {
2038      DataTagged* leftC=dynamic_cast<DataTagged*>(m_data.get());      DataTagged* leftC=dynamic_cast<DataTagged*>(m_data.get());
2039      EsysAssert((leftC!=0), "Programming error - casting to DataTagged.");      EsysAssert((leftC!=0), "Programming error - casting to DataTagged.");
2040      return escript::algorithm(*leftC,operation);      return escript::algorithm(*leftC,operation,initial_value);
2041    } else if (isConstant()) {    } else if (isConstant()) {
2042      DataConstant* leftC=dynamic_cast<DataConstant*>(m_data.get());      DataConstant* leftC=dynamic_cast<DataConstant*>(m_data.get());
2043      EsysAssert((leftC!=0), "Programming error - casting to DataConstant.");      EsysAssert((leftC!=0), "Programming error - casting to DataConstant.");
2044      return escript::algorithm(*leftC,operation);      return escript::algorithm(*leftC,operation,initial_value);
2045      } else if (isEmpty()) {
2046        throw DataException("Error - Operations not permitted on instances of DataEmpty.");
2047      } else if (isLazy()) {
2048        throw DataException("Error - Operations not permitted on instances of DataLazy.");
2049      } else {
2050        throw DataException("Error - Data encapsulates an unknown type.");
2051    }    }
   return 0;  
2052  }  }
2053    
2054  /**  /**
2055    \brief    \brief
2056    Perform the given data point reduction algorithm on data and return the result.    Perform the given data point reduction algorithm on data and return the result.
2057    Given operation combines each element within each data point into a scalar,    Given operation combines each element within each data point into a scalar,
2058    thus argument object is a rank n Data object, and returned object is a    thus argument object is a rank n Data object, and returned object is a
2059    rank 0 Data object.    rank 0 Data object.
2060    Calls escript::dp_algorithm.    Calls escript::dp_algorithm.
2061  */  */
2062  template <class UnaryFunction>  template <class BinaryFunction>
2063  inline  inline
2064  Data  Data
2065  dp_algorithm(const Data& data,  Data::dp_algorithm(BinaryFunction operation, double initial_value) const
              UnaryFunction operation)  
2066  {  {
2067    Data result(0,DataArrayView::ShapeType(),data.getFunctionSpace(),data.isExpanded());    if (isEmpty()) {
2068    if (data.isExpanded()) {      throw DataException("Error - Operations not permitted on instances of DataEmpty.");
2069      DataExpanded* dataE=dynamic_cast<DataExpanded*>(data.m_data.get());    }
2070      else if (isExpanded()) {
2071        Data result(0,DataTypes::ShapeType(),getFunctionSpace(),isExpanded());
2072        DataExpanded* dataE=dynamic_cast<DataExpanded*>(m_data.get());
2073      DataExpanded* resultE=dynamic_cast<DataExpanded*>(result.m_data.get());      DataExpanded* resultE=dynamic_cast<DataExpanded*>(result.m_data.get());
2074      EsysAssert((dataE!=0), "Programming error - casting data to DataExpanded.");      EsysAssert((dataE!=0), "Programming error - casting data to DataExpanded.");
2075      EsysAssert((resultE!=0), "Programming error - casting result to DataExpanded.");      EsysAssert((resultE!=0), "Programming error - casting result to DataExpanded.");
2076      escript::dp_algorithm(*dataE,*resultE,operation);      escript::dp_algorithm(*dataE,*resultE,operation,initial_value);
2077    } else if (data.isTagged()) {      return result;
2078      DataTagged* dataT=dynamic_cast<DataTagged*>(data.m_data.get());    }
2079      DataTagged* resultT=dynamic_cast<DataTagged*>(result.m_data.get());    else if (isTagged()) {
2080        DataTagged* dataT=dynamic_cast<DataTagged*>(m_data.get());
2081      EsysAssert((dataT!=0), "Programming error - casting data to DataTagged.");      EsysAssert((dataT!=0), "Programming error - casting data to DataTagged.");
2082      EsysAssert((resultT!=0), "Programming error - casting result to DataTagged.");      DataTypes::ValueType defval(1);
2083      escript::dp_algorithm(*dataT,*resultT,operation);      defval[0]=0;
2084    } else if (data.isConstant()) {      DataTagged* resultT=new DataTagged(getFunctionSpace(), DataTypes::scalarShape, defval, dataT);
2085      DataConstant* dataC=dynamic_cast<DataConstant*>(data.m_data.get());      escript::dp_algorithm(*dataT,*resultT,operation,initial_value);
2086        return Data(resultT);   // note: the Data object now owns the resultT pointer
2087      }
2088      else if (isConstant()) {
2089        Data result(0,DataTypes::ShapeType(),getFunctionSpace(),isExpanded());
2090        DataConstant* dataC=dynamic_cast<DataConstant*>(m_data.get());
2091      DataConstant* resultC=dynamic_cast<DataConstant*>(result.m_data.get());      DataConstant* resultC=dynamic_cast<DataConstant*>(result.m_data.get());
2092      EsysAssert((dataC!=0), "Programming error - casting data to DataConstant.");      EsysAssert((dataC!=0), "Programming error - casting data to DataConstant.");
2093      EsysAssert((resultC!=0), "Programming error - casting result to DataConstant.");      EsysAssert((resultC!=0), "Programming error - casting result to DataConstant.");
2094      escript::dp_algorithm(*dataC,*resultC,operation);      escript::dp_algorithm(*dataC,*resultC,operation,initial_value);
2095        return result;
2096      } else if (isLazy()) {
2097        throw DataException("Error - Operations not permitted on instances of DataLazy.");
2098      } else {
2099        throw DataException("Error - Data encapsulates an unknown type.");
2100      }
2101    }
2102    
2103    /**
2104      \brief
2105      Compute a tensor operation with two Data objects
2106      \param arg0 - Input - Data object
2107      \param arg1 - Input - Data object
2108      \param operation - Input - Binary op functor
2109    */
2110    template <typename BinaryFunction>
2111    inline
2112    Data
2113    C_TensorBinaryOperation(Data const &arg_0,
2114                            Data const &arg_1,
2115                            BinaryFunction operation)
2116    {
2117      if (arg_0.isEmpty() || arg_1.isEmpty())
2118      {
2119         throw DataException("Error - Operations not permitted on instances of DataEmpty.");
2120      }
2121      if (arg_0.isLazy() || arg_1.isLazy())
2122      {
2123         throw DataException("Error - Operations not permitted on lazy data.");
2124      }
2125      // Interpolate if necessary and find an appropriate function space
2126      Data arg_0_Z, arg_1_Z;
2127      if (arg_0.getFunctionSpace()!=arg_1.getFunctionSpace()) {
2128        if (arg_0.probeInterpolation(arg_1.getFunctionSpace())) {
2129          arg_0_Z = arg_0.interpolate(arg_1.getFunctionSpace());
2130          arg_1_Z = Data(arg_1);
2131        }
2132        else if (arg_1.probeInterpolation(arg_0.getFunctionSpace())) {
2133          arg_1_Z=arg_1.interpolate(arg_0.getFunctionSpace());
2134          arg_0_Z =Data(arg_0);
2135        }
2136        else {
2137          throw DataException("Error - C_TensorBinaryOperation: arguments have incompatible function spaces.");
2138        }
2139      } else {
2140          arg_0_Z = Data(arg_0);
2141          arg_1_Z = Data(arg_1);
2142      }
2143      // Get rank and shape of inputs
2144      int rank0 = arg_0_Z.getDataPointRank();
2145      int rank1 = arg_1_Z.getDataPointRank();
2146      DataTypes::ShapeType shape0 = arg_0_Z.getDataPointShape();
2147      DataTypes::ShapeType shape1 = arg_1_Z.getDataPointShape();
2148      int size0 = arg_0_Z.getDataPointSize();
2149      int size1 = arg_1_Z.getDataPointSize();
2150      // Declare output Data object
2151      Data res;
2152    
2153      if (shape0 == shape1) {
2154        if (arg_0_Z.isConstant()   && arg_1_Z.isConstant()) {
2155          res = Data(0.0, shape0, arg_1_Z.getFunctionSpace());      // DataConstant output
2156          const double *ptr_0 = &(arg_0_Z.getDataAtOffsetRO(0));
2157          const double *ptr_1 = &(arg_1_Z.getDataAtOffsetRO(0));
2158          double *ptr_2 = &(res.getDataAtOffsetRW(0));
2159    
2160          tensor_binary_operation(size0, ptr_0, ptr_1, ptr_2, operation);
2161        }
2162        else if (arg_0_Z.isConstant()   && arg_1_Z.isTagged()) {
2163    
2164          // Prepare the DataConstant input
2165          DataConstant* tmp_0=dynamic_cast<DataConstant*>(arg_0_Z.borrowData());
2166    
2167          // Borrow DataTagged input from Data object
2168          DataTagged* tmp_1=dynamic_cast<DataTagged*>(arg_1_Z.borrowData());
2169    
2170          // Prepare a DataTagged output 2
2171          res = Data(0.0, shape0, arg_1_Z.getFunctionSpace());      // DataTagged output
2172          res.tag();
2173          DataTagged* tmp_2=dynamic_cast<DataTagged*>(res.borrowData());
2174    
2175          // Prepare offset into DataConstant
2176          int offset_0 = tmp_0->getPointOffset(0,0);
2177          const double *ptr_0 = &(arg_0_Z.getDataAtOffsetRO(offset_0));
2178    
2179          // Get the pointers to the actual data
2180          const double *ptr_1 = &(tmp_1->getDefaultValueRO(0));
2181          double *ptr_2 = &(tmp_2->getDefaultValueRW(0));
2182    
2183          // Compute a result for the default
2184          tensor_binary_operation(size0, ptr_0, ptr_1, ptr_2, operation);
2185          // Compute a result for each tag
2186          const DataTagged::DataMapType& lookup_1=tmp_1->getTagLookup();
2187          DataTagged::DataMapType::const_iterator i; // i->first is a tag, i->second is an offset into memory
2188          for (i=lookup_1.begin();i!=lookup_1.end();i++) {
2189            tmp_2->addTag(i->first);
2190            const double *ptr_1 = &(tmp_1->getDataByTagRO(i->first,0));
2191            double *ptr_2 = &(tmp_2->getDataByTagRW(i->first,0));
2192    
2193            tensor_binary_operation(size0, ptr_0, ptr_1, ptr_2, operation);
2194          }
2195    
2196        }
2197        else if (arg_0_Z.isConstant()   && arg_1_Z.isExpanded()) {
2198          res = Data(0.0, shape0, arg_1_Z.getFunctionSpace(),true); // DataExpanded output
2199          DataConstant* tmp_0=dynamic_cast<DataConstant*>(arg_0_Z.borrowData());
2200          DataExpanded* tmp_1=dynamic_cast<DataExpanded*>(arg_1_Z.borrowData());
2201          DataExpanded* tmp_2=dynamic_cast<DataExpanded*>(res.borrowData());
2202    
2203          int sampleNo_1,dataPointNo_1;
2204          int numSamples_1 = arg_1_Z.getNumSamples();
2205          int numDataPointsPerSample_1 = arg_1_Z.getNumDataPointsPerSample();
2206          int offset_0 = tmp_0->getPointOffset(0,0);
2207          res.requireWrite();
2208          #pragma omp parallel for private(sampleNo_1,dataPointNo_1) schedule(static)
2209          for (sampleNo_1 = 0; sampleNo_1 < numSamples_1; sampleNo_1++) {
2210            for (dataPointNo_1 = 0; dataPointNo_1 < numDataPointsPerSample_1; dataPointNo_1++) {
2211              int offset_1 = tmp_1->getPointOffset(sampleNo_1,dataPointNo_1);
2212              int offset_2 = tmp_2->getPointOffset(sampleNo_1,dataPointNo_1);
2213              const double *ptr_0 = &(arg_0_Z.getDataAtOffsetRO(offset_0));
2214              const double *ptr_1 = &(arg_1_Z.getDataAtOffsetRO(offset_1));
2215              double *ptr_2 = &(res.getDataAtOffsetRW(offset_2));
2216              tensor_binary_operation(size0, ptr_0, ptr_1, ptr_2, operation);
2217            }
2218          }
2219    
2220        }
2221        else if (arg_0_Z.isTagged()     && arg_1_Z.isConstant()) {
2222          // Borrow DataTagged input from Data object
2223          DataTagged* tmp_0=dynamic_cast<DataTagged*>(arg_0_Z.borrowData());
2224    
2225          // Prepare the DataConstant input
2226          DataConstant* tmp_1=dynamic_cast<DataConstant*>(arg_1_Z.borrowData());
2227    
2228          // Prepare a DataTagged output 2
2229          res = Data(0.0, shape0, arg_0_Z.getFunctionSpace());      // DataTagged output
2230          res.tag();
2231          DataTagged* tmp_2=dynamic_cast<DataTagged*>(res.borrowData());
2232    
2233          // Prepare offset into DataConstant
2234          int offset_1 = tmp_1->getPointOffset(0,0);
2235    
2236          const double *ptr_1 = &(arg_1_Z.getDataAtOffsetRO(offset_1));
2237          // Get the pointers to the actual data
2238          const double *ptr_0 = &(tmp_0->getDefaultValueRO(0));
2239          double *ptr_2 = &(tmp_2->getDefaultValueRW(0));
2240          // Compute a result for the default
2241          tensor_binary_operation(size0, ptr_0, ptr_1, ptr_2, operation);
2242          // Compute a result for each tag
2243          const DataTagged::DataMapType& lookup_0=tmp_0->getTagLookup();
2244          DataTagged::DataMapType::const_iterator i; // i->first is a tag, i->second is an offset into memory
2245          for (i=lookup_0.begin();i!=lookup_0.end();i++) {
2246            tmp_2->addTag(i->first);
2247            const double *ptr_0 = &(tmp_0->getDataByTagRO(i->first,0));
2248            double *ptr_2 = &(tmp_2->getDataByTagRW(i->first,0));
2249            tensor_binary_operation(size0, ptr_0, ptr_1, ptr_2, operation);
2250          }
2251    
2252        }
2253        else if (arg_0_Z.isTagged()     && arg_1_Z.isTagged()) {
2254          // Borrow DataTagged input from Data object
2255          DataTagged* tmp_0=dynamic_cast<DataTagged*>(arg_0_Z.borrowData());
2256    
2257          // Borrow DataTagged input from Data object
2258          DataTagged* tmp_1=dynamic_cast<DataTagged*>(arg_1_Z.borrowData());
2259    
2260          // Prepare a DataTagged output 2
2261          res = Data(0.0, shape0, arg_1_Z.getFunctionSpace());
2262          res.tag();        // DataTagged output
2263          DataTagged* tmp_2=dynamic_cast<DataTagged*>(res.borrowData());
2264    
2265          // Get the pointers to the actual data
2266          const double *ptr_0 = &(tmp_0->getDefaultValueRO(0));
2267          const double *ptr_1 = &(tmp_1->getDefaultValueRO(0));
2268          double *ptr_2 = &(tmp_2->getDefaultValueRW(0));
2269    
2270          // Compute a result for the default
2271          tensor_binary_operation(size0, ptr_0, ptr_1, ptr_2, operation);
2272          // Merge the tags
2273          DataTagged::DataMapType::const_iterator i; // i->first is a tag, i->second is an offset into memory
2274          const DataTagged::DataMapType& lookup_0=tmp_0->getTagLookup();
2275          const DataTagged::DataMapType& lookup_1=tmp_1->getTagLookup();
2276          for (i=lookup_0.begin();i!=lookup_0.end();i++) {
2277            tmp_2->addTag(i->first); // use tmp_2 to get correct shape
2278          }
2279          for (i=lookup_1.begin();i!=lookup_1.end();i++) {
2280            tmp_2->addTag(i->first);
2281          }
2282          // Compute a result for each tag
2283          const DataTagged::DataMapType& lookup_2=tmp_2->getTagLookup();
2284          for (i=lookup_2.begin();i!=lookup_2.end();i++) {
2285    
2286            const double *ptr_0 = &(tmp_0->getDataByTagRO(i->first,0));
2287            const double *ptr_1 = &(tmp_1->getDataByTagRO(i->first,0));
2288            double *ptr_2 = &(tmp_2->getDataByTagRW(i->first,0));
2289    
2290            tensor_binary_operation(size0, ptr_0, ptr_1, ptr_2, operation);
2291          }
2292    
2293        }
2294        else if (arg_0_Z.isTagged()     && arg_1_Z.isExpanded()) {
2295          // After finding a common function space above the two inputs have the same numSamples and num DPPS
2296          res = Data(0.0, shape0, arg_1_Z.getFunctionSpace(),true); // DataExpanded output
2297          DataTagged*   tmp_0=dynamic_cast<DataTagged*>(arg_0_Z.borrowData());
2298          DataExpanded* tmp_1=dynamic_cast<DataExpanded*>(arg_1_Z.borrowData());
2299          DataExpanded* tmp_2=dynamic_cast<DataExpanded*>(res.borrowData());
2300    
2301          int sampleNo_0,dataPointNo_0;
2302          int numSamples_0 = arg_0_Z.getNumSamples();
2303          int numDataPointsPerSample_0 = arg_0_Z.getNumDataPointsPerSample();
2304          res.requireWrite();
2305          #pragma omp parallel for private(sampleNo_0,dataPointNo_0) schedule(static)
2306          for (sampleNo_0 = 0; sampleNo_0 < numSamples_0; sampleNo_0++) {
2307            int offset_0 = tmp_0->getPointOffset(sampleNo_0,0); // They're all the same, so just use #0
2308            const double *ptr_0 = &(arg_0_Z.getDataAtOffsetRO(offset_0));
2309            for (dataPointNo_0 = 0; dataPointNo_0 < numDataPointsPerSample_0; dataPointNo_0++) {
2310              int offset_1 = tmp_1->getPointOffset(sampleNo_0,dataPointNo_0);
2311              int offset_2 = tmp_2->getPointOffset(sampleNo_0,dataPointNo_0);
2312              const double *ptr_1 = &(arg_1_Z.getDataAtOffsetRO(offset_1));
2313              double *ptr_2 = &(res.getDataAtOffsetRW(offset_2));
2314              tensor_binary_operation(size0, ptr_0, ptr_1, ptr_2, operation);
2315            }
2316          }
2317    
2318        }
2319        else if (arg_0_Z.isExpanded()   && arg_1_Z.isConstant()) {
2320          res = Data(0.0, shape0, arg_1_Z.getFunctionSpace(),true); // DataExpanded output
2321          DataExpanded* tmp_0=dynamic_cast<DataExpanded*>(arg_0_Z.borrowData());
2322          DataConstant* tmp_1=dynamic_cast<DataConstant*>(arg_1_Z.borrowData());
2323          DataExpanded* tmp_2=dynamic_cast<DataExpanded*>(res.borrowData());
2324    
2325          int sampleNo_0,dataPointNo_0;
2326          int numSamples_0 = arg_0_Z.getNumSamples();
2327          int numDataPointsPerSample_0 = arg_0_Z.getNumDataPointsPerSample();
2328          int offset_1 = tmp_1->getPointOffset(0,0);
2329          res.requireWrite();
2330          #pragma omp parallel for private(sampleNo_0,dataPointNo_0) schedule(static)
2331          for (sampleNo_0 = 0; sampleNo_0 < numSamples_0; sampleNo_0++) {
2332            for (dataPointNo_0 = 0; dataPointNo_0 < numDataPointsPerSample_0; dataPointNo_0++) {
2333              int offset_0 = tmp_0->getPointOffset(sampleNo_0,dataPointNo_0);
2334              int offset_2 = tmp_2->getPointOffset(sampleNo_0,dataPointNo_0);
2335    
2336              const double *ptr_0 = &(arg_0_Z.getDataAtOffsetRO(offset_0));
2337              const double *ptr_1 = &(arg_1_Z.getDataAtOffsetRO(offset_1));
2338              double *ptr_2 = &(res.getDataAtOffsetRW(offset_2));
2339    
2340    
2341              tensor_binary_operation(size0, ptr_0, ptr_1, ptr_2, operation);
2342            }
2343          }
2344    
2345        }
2346        else if (arg_0_Z.isExpanded()   && arg_1_Z.isTagged()) {
2347          // After finding a common function space above the two inputs have the same numSamples and num DPPS
2348          res = Data(0.0, shape0, arg_1_Z.getFunctionSpace(),true); // DataExpanded output
2349          DataExpanded* tmp_0=dynamic_cast<DataExpanded*>(arg_0_Z.borrowData());
2350          DataTagged*   tmp_1=dynamic_cast<DataTagged*>(arg_1_Z.borrowData());
2351          DataExpanded* tmp_2=dynamic_cast<DataExpanded*>(res.borrowData());
2352    
2353          int sampleNo_0,dataPointNo_0;
2354          int numSamples_0 = arg_0_Z.getNumSamples();
2355          int numDataPointsPerSample_0 = arg_0_Z.getNumDataPointsPerSample();
2356          res.requireWrite();
2357          #pragma omp parallel for private(sampleNo_0,dataPointNo_0) schedule(static)
2358          for (sampleNo_0 = 0; sampleNo_0 < numSamples_0; sampleNo_0++) {
2359            int offset_1 = tmp_1->getPointOffset(sampleNo_0,0);
2360            const double *ptr_1 = &(arg_1_Z.getDataAtOffsetRO(offset_1));
2361            for (dataPointNo_0 = 0; dataPointNo_0 < numDataPointsPerSample_0; dataPointNo_0++) {
2362              int offset_0 = tmp_0->getPointOffset(sampleNo_0,dataPointNo_0);
2363              int offset_2 = tmp_2->getPointOffset(sampleNo_0,dataPointNo_0);
2364              const double *ptr_0 = &(arg_0_Z.getDataAtOffsetRO(offset_0));
2365              double *ptr_2 = &(res.getDataAtOffsetRW(offset_2));
2366              tensor_binary_operation(size0, ptr_0, ptr_1, ptr_2, operation);
2367            }
2368          }
2369    
2370        }
2371        else if (arg_0_Z.isExpanded()   && arg_1_Z.isExpanded()) {
2372          // After finding a common function space above the two inputs have the same numSamples and num DPPS
2373          res = Data(0.0, shape0, arg_1_Z.getFunctionSpace(),true); // DataExpanded output
2374          DataExpanded* tmp_0=dynamic_cast<DataExpanded*>(arg_0_Z.borrowData());
2375          DataExpanded* tmp_1=dynamic_cast<DataExpanded*>(arg_1_Z.borrowData());
2376          DataExpanded* tmp_2=dynamic_cast<DataExpanded*>(res.borrowData());
2377    
2378          int sampleNo_0,dataPointNo_0;
2379          int numSamples_0 = arg_0_Z.getNumSamples();
2380          int numDataPointsPerSample_0 = arg_0_Z.getNumDataPointsPerSample();
2381          res.requireWrite();
2382          #pragma omp parallel for private(sampleNo_0,dataPointNo_0) schedule(static)
2383          for (sampleNo_0 = 0; sampleNo_0 < numSamples_0; sampleNo_0++) {
2384            for (dataPointNo_0 = 0; dataPointNo_0 < numDataPointsPerSample_0; dataPointNo_0++) {
2385              int offset_0 = tmp_0->getPointOffset(sampleNo_0,dataPointNo_0);
2386              int offset_1 = tmp_1->getPointOffset(sampleNo_0,dataPointNo_0);
2387              int offset_2 = tmp_2->getPointOffset(sampleNo_0,dataPointNo_0);
2388              const double *ptr_0 = &(arg_0_Z.getDataAtOffsetRO(offset_0));
2389              const double *ptr_1 = &(arg_1_Z.getDataAtOffsetRO(offset_1));
2390              double *ptr_2 = &(res.getDataAtOffsetRW(offset_2));
2391              tensor_binary_operation(size0, ptr_0, ptr_1, ptr_2, operation);
2392            }
2393          }
2394    
2395        }
2396        else {
2397          throw DataException("Error - C_TensorBinaryOperation: unknown combination of inputs");
2398        }
2399    
2400      } else if (0 == rank0) {
2401        if (arg_0_Z.isConstant()   && arg_1_Z.isConstant()) {
2402          res = Data(0.0, shape1, arg_1_Z.getFunctionSpace());      // DataConstant output
2403          const double *ptr_0 = &(arg_0_Z.getDataAtOffsetRO(0));
2404          const double *ptr_1 = &(arg_1_Z.getDataAtOffsetRO(0));
2405          double *ptr_2 = &(res.getDataAtOffsetRW(0));
2406          tensor_binary_operation(size1, ptr_0[0], ptr_1, ptr_2, operation);
2407        }
2408        else if (arg_0_Z.isConstant()   && arg_1_Z.isTagged()) {
2409    
2410          // Prepare the DataConstant input
2411          DataConstant* tmp_0=dynamic_cast<DataConstant*>(arg_0_Z.borrowData());
2412    
2413          // Borrow DataTagged input from Data object
2414          DataTagged* tmp_1=dynamic_cast<DataTagged*>(arg_1_Z.borrowData());
2415    
2416          // Prepare a DataTagged output 2
2417          res = Data(0.0, shape1, arg_1_Z.getFunctionSpace());      // DataTagged output
2418          res.tag();
2419          DataTagged* tmp_2=dynamic_cast<DataTagged*>(res.borrowData());
2420    
2421          // Prepare offset into DataConstant
2422          int offset_0 = tmp_0->getPointOffset(0,0);
2423          const double *ptr_0 = &(arg_0_Z.getDataAtOffsetRO(offset_0));
2424    
2425          const double *ptr_1 = &(tmp_1->getDefaultValueRO(0));
2426          double *ptr_2 = &(tmp_2->getDefaultValueRW(0));
2427    
2428          // Compute a result for the default
2429          tensor_binary_operation(size1, ptr_0[0], ptr_1, ptr_2, operation);
2430          // Compute a result for each tag
2431          const DataTagged::DataMapType& lookup_1=tmp_1->getTagLookup();
2432          DataTagged::DataMapType::const_iterator i; // i->first is a tag, i->second is an offset into memory
2433          for (i=lookup_1.begin();i!=lookup_1.end();i++) {
2434            tmp_2->addTag(i->first);
2435            const double *ptr_1 = &(tmp_1->getDataByTagRO(i->first,0));
2436            double *ptr_2 = &(tmp_2->getDataByTagRW(i->first,0));
2437            tensor_binary_operation(size1, ptr_0[0], ptr_1, ptr_2, operation);
2438          }
2439    
2440        }
2441        else if (arg_0_Z.isConstant()   && arg_1_Z.isExpanded()) {
2442    
2443          res = Data(0.0, shape1, arg_1_Z.getFunctionSpace(),true); // DataExpanded output
2444          DataConstant* tmp_0=dynamic_cast<DataConstant*>(arg_0_Z.borrowData());
2445          DataExpanded* tmp_1=dynamic_cast<DataExpanded*>(arg_1_Z.borrowData());
2446          DataExpanded* tmp_2=dynamic_cast<DataExpanded*>(res.borrowData());
2447    
2448          int sampleNo_1,dataPointNo_1;
2449          int numSamples_1 = arg_1_Z.getNumSamples();
2450          int numDataPointsPerSample_1 = arg_1_Z.getNumDataPointsPerSample();
2451          int offset_0 = tmp_0->getPointOffset(0,0);
2452          res.requireWrite();
2453          #pragma omp parallel for private(sampleNo_1,dataPointNo_1) schedule(static)
2454          for (sampleNo_1 = 0; sampleNo_1 < numSamples_1; sampleNo_1++) {
2455            for (dataPointNo_1 = 0; dataPointNo_1 < numDataPointsPerSample_1; dataPointNo_1++) {
2456              int offset_1 = tmp_1->getPointOffset(sampleNo_1,dataPointNo_1);
2457              int offset_2 = tmp_2->getPointOffset(sampleNo_1,dataPointNo_1);
2458              const double *ptr_0 = &(arg_0_Z.getDataAtOffsetRO(offset_0));
2459              const double *ptr_1 = &(arg_1_Z.getDataAtOffsetRO(offset_1));
2460              double *ptr_2 = &(res.getDataAtOffsetRW(offset_2));
2461              tensor_binary_operation(size1, ptr_0[0], ptr_1, ptr_2, operation);
2462    
2463            }
2464          }
2465    
2466        }
2467        else if (arg_0_Z.isTagged()     && arg_1_Z.isConstant()) {
2468    
2469          // Borrow DataTagged input from Data object
2470          DataTagged* tmp_0=dynamic_cast<DataTagged*>(arg_0_Z.borrowData());
2471    
2472          // Prepare the DataConstant input
2473          DataConstant* tmp_1=dynamic_cast<DataConstant*>(arg_1_Z.borrowData());
2474    
2475          // Prepare a DataTagged output 2
2476          res = Data(0.0, shape1, arg_0_Z.getFunctionSpace());      // DataTagged output
2477          res.tag();
2478          DataTagged* tmp_2=dynamic_cast<DataTagged*>(res.borrowData());
2479    
2480          // Prepare offset into DataConstant
2481          int offset_1 = tmp_1->getPointOffset(0,0);
2482          const double *ptr_1 = &(arg_1_Z.getDataAtOffsetRO(offset_1));
2483    
2484          // Get the pointers to the actual data
2485          const double *ptr_0 = &(tmp_0->getDefaultValueRO(0));
2486          double *ptr_2 = &(tmp_2->getDefaultValueRW(0));
2487    
2488    
2489          // Compute a result for the default
2490          tensor_binary_operation(size1, ptr_0[0], ptr_1, ptr_2, operation);
2491          // Compute a result for each tag
2492          const DataTagged::DataMapType& lookup_0=tmp_0->getTagLookup();
2493          DataTagged::DataMapType::const_iterator i; // i->first is a tag, i->second is an offset into memory
2494          for (i=lookup_0.begin();i!=lookup_0.end();i++) {
2495            tmp_2->addTag(i->first);
2496            const double *ptr_0 = &(tmp_0->getDataByTagRO(i->first,0));
2497            double *ptr_2 = &(tmp_2->getDataByTagRW(i->first,0));
2498    
2499            tensor_binary_operation(size1, ptr_0[0], ptr_1, ptr_2, operation);
2500          }
2501    
2502        }
2503        else if (arg_0_Z.isTagged()     && arg_1_Z.isTagged()) {
2504    
2505          // Borrow DataTagged input from Data object
2506          DataTagged* tmp_0=dynamic_cast<DataTagged*>(arg_0_Z.borrowData());
2507    
2508          // Borrow DataTagged input from Data object
2509          DataTagged* tmp_1=dynamic_cast<DataTagged*>(arg_1_Z.borrowData());
2510    
2511          // Prepare a DataTagged output 2
2512          res = Data(0.0, shape1, arg_1_Z.getFunctionSpace());
2513          res.tag();        // DataTagged output
2514          DataTagged* tmp_2=dynamic_cast<DataTagged*>(res.borrowData());
2515    
2516          // Get the pointers to the actual data
2517          const double *ptr_0 = &(tmp_0->getDefaultValueRO(0));
2518          const double *ptr_1 = &(tmp_1->getDefaultValueRO(0));
2519          double *ptr_2 = &(tmp_2->getDefaultValueRW(0));
2520    
2521          // Compute a result for the default
2522          tensor_binary_operation(size1, ptr_0[0], ptr_1, ptr_2, operation);
2523          // Merge the tags
2524          DataTagged::DataMapType::const_iterator i; // i->first is a tag, i->second is an offset into memory
2525          const DataTagged::DataMapType& lookup_0=tmp_0->getTagLookup();
2526          const DataTagged::DataMapType& lookup_1=tmp_1->getTagLookup();
2527          for (i=lookup_0.begin();i!=lookup_0.end();i++) {
2528            tmp_2->addTag(i->first); // use tmp_2 to get correct shape
2529          }
2530          for (i=lookup_1.begin();i!=lookup_1.end();i++) {
2531            tmp_2->addTag(i->first);
2532          }
2533          // Compute a result for each tag
2534          const DataTagged::DataMapType& lookup_2=tmp_2->getTagLookup();
2535          for (i=lookup_2.begin();i!=lookup_2.end();i++) {
2536    
2537    /*        DataArrayView view_0 = tmp_0->getDataPointByTag(i->first);
2538            DataArrayView view_1 = tmp_1->getDataPointByTag(i->first);
2539            DataArrayView view_2 = tmp_2->getDataPointByTag(i->first);
2540            double *ptr_0 = &view_0.getData(0);
2541            double *ptr_1 = &view_1.getData(0);
2542            double *ptr_2 = &view_2.getData(0);*/
2543    
2544            const double *ptr_0 = &(tmp_0->getDataByTagRO(i->first,0));
2545            const double *ptr_1 = &(tmp_1->getDataByTagRO(i->first,0));
2546            double *ptr_2 = &(tmp_2->getDataByTagRW(i->first,0));
2547    
2548            tensor_binary_operation(size1, ptr_0[0], ptr_1, ptr_2, operation);
2549          }
2550    
2551        }
2552        else if (arg_0_Z.isTagged()     && arg_1_Z.isExpanded()) {
2553    
2554          // After finding a common function space above the two inputs have the same numSamples and num DPPS
2555          res = Data(0.0, shape1, arg_1_Z.getFunctionSpace(),true); // DataExpanded output
2556          DataTagged*   tmp_0=dynamic_cast<DataTagged*>(arg_0_Z.borrowData());
2557          DataExpanded* tmp_1=dynamic_cast<DataExpanded*>(arg_1_Z.borrowData());
2558          DataExpanded* tmp_2=dynamic_cast<DataExpanded*>(res.borrowData());
2559    
2560          int sampleNo_0,dataPointNo_0;
2561          int numSamples_0 = arg_0_Z.getNumSamples();
2562          int numDataPointsPerSample_0 = arg_0_Z.getNumDataPointsPerSample();
2563          res.requireWrite();
2564          #pragma omp parallel for private(sampleNo_0,dataPointNo_0) schedule(static)
2565          for (sampleNo_0 = 0; sampleNo_0 < numSamples_0; sampleNo_0++) {
2566            int offset_0 = tmp_0->getPointOffset(sampleNo_0,0); // They're all the same, so just use #0
2567            const double *ptr_0 = &(arg_0_Z.getDataAtOffsetRO(offset_0));
2568            for (dataPointNo_0 = 0; dataPointNo_0 < numDataPointsPerSample_0; dataPointNo_0++) {
2569              int offset_1 = tmp_1->getPointOffset(sampleNo_0,dataPointNo_0);
2570              int offset_2 = tmp_2->getPointOffset(sampleNo_0,dataPointNo_0);
2571              const double *ptr_1 = &(arg_1_Z.getDataAtOffsetRO(offset_1));
2572              double *ptr_2 = &(res.getDataAtOffsetRW(offset_2));
2573              tensor_binary_operation(size1, ptr_0[0], ptr_1, ptr_2, operation);
2574            }
2575          }
2576    
2577        }
2578        else if (arg_0_Z.isExpanded()   && arg_1_Z.isConstant()) {
2579          res = Data(0.0, shape1, arg_1_Z.getFunctionSpace(),true); // DataExpanded output
2580          DataExpanded* tmp_0=dynamic_cast<DataExpanded*>(arg_0_Z.borrowData());
2581          DataConstant* tmp_1=dynamic_cast<DataConstant*>(arg_1_Z.borrowData());
2582          DataExpanded* tmp_2=dynamic_cast<DataExpanded*>(res.borrowData());
2583    
2584          int sampleNo_0,dataPointNo_0;
2585          int numSamples_0 = arg_0_Z.getNumSamples();
2586          int numDataPointsPerSample_0 = arg_0_Z.getNumDataPointsPerSample();
2587          int offset_1 = tmp_1->getPointOffset(0,0);
2588          res.requireWrite();
2589          #pragma omp parallel for private(sampleNo_0,dataPointNo_0) schedule(static)
2590          for (sampleNo_0 = 0; sampleNo_0 < numSamples_0; sampleNo_0++) {
2591            for (dataPointNo_0 = 0; dataPointNo_0 < numDataPointsPerSample_0; dataPointNo_0++) {
2592              int offset_0 = tmp_0->getPointOffset(sampleNo_0,dataPointNo_0);
2593              int offset_2 = tmp_2->getPointOffset(sampleNo_0,dataPointNo_0);
2594              const double *ptr_0 = &(arg_0_Z.getDataAtOffsetRO(offset_0));
2595              const double *ptr_1 = &(arg_1_Z.getDataAtOffsetRO(offset_1));
2596              double *ptr_2 = &(res.getDataAtOffsetRW(offset_2));
2597              tensor_binary_operation(size1, ptr_0[0], ptr_1, ptr_2, operation);
2598            }
2599          }
2600    
2601    
2602        }
2603        else if (arg_0_Z.isExpanded()   && arg_1_Z.isTagged()) {
2604    
2605          // After finding a common function space above the two inputs have the same numSamples and num DPPS
2606          res = Data(0.0, shape1, arg_1_Z.getFunctionSpace(),true); // DataExpanded output
2607          DataExpanded* tmp_0=dynamic_cast<DataExpanded*>(arg_0_Z.borrowData());
2608          DataTagged*   tmp_1=dynamic_cast<DataTagged*>(arg_1_Z.borrowData());
2609          DataExpanded* tmp_2=dynamic_cast<DataExpanded*>(res.borrowData());
2610    
2611          int sampleNo_0,dataPointNo_0;
2612          int numSamples_0 = arg_0_Z.getNumSamples();
2613          int numDataPointsPerSample_0 = arg_0_Z.getNumDataPointsPerSample();
2614          res.requireWrite();
2615          #pragma omp parallel for private(sampleNo_0,dataPointNo_0) schedule(static)
2616          for (sampleNo_0 = 0; sampleNo_0 < numSamples_0; sampleNo_0++) {
2617            int offset_1 = tmp_1->getPointOffset(sampleNo_0,0);
2618            const double *ptr_1 = &(arg_1_Z.getDataAtOffsetRO(offset_1));
2619            for (dataPointNo_0 = 0; dataPointNo_0 < numDataPointsPerSample_0; dataPointNo_0++) {
2620              int offset_0 = tmp_0->getPointOffset(sampleNo_0,dataPointNo_0);
2621              int offset_2 = tmp_2->getPointOffset(sampleNo_0,dataPointNo_0);
2622              const double *ptr_0 = &(arg_0_Z.getDataAtOffsetRO(offset_0));
2623              double *ptr_2 = &(res.getDataAtOffsetRW(offset_2));
2624              tensor_binary_operation(size1, ptr_0[0], ptr_1, ptr_2, operation);
2625            }
2626          }
2627    
2628        }
2629        else if (arg_0_Z.isExpanded()   && arg_1_Z.isExpanded()) {
2630    
2631          // After finding a common function space above the two inputs have the same numSamples and num DPPS
2632          res = Data(0.0, shape1, arg_1_Z.getFunctionSpace(),true); // DataExpanded output
2633          DataExpanded* tmp_0=dynamic_cast<DataExpanded*>(arg_0_Z.borrowData());
2634          DataExpanded* tmp_1=dynamic_cast<DataExpanded*>(arg_1_Z.borrowData());
2635          DataExpanded* tmp_2=dynamic_cast<DataExpanded*>(res.borrowData());
2636    
2637          int sampleNo_0,dataPointNo_0;
2638          int numSamples_0 = arg_0_Z.getNumSamples();
2639          int numDataPointsPerSample_0 = arg_0_Z.getNumDataPointsPerSample();
2640          res.requireWrite();
2641          #pragma omp parallel for private(sampleNo_0,dataPointNo_0) schedule(static)
2642          for (sampleNo_0 = 0; sampleNo_0 < numSamples_0; sampleNo_0++) {
2643            for (dataPointNo_0 = 0; dataPointNo_0 < numDataPointsPerSample_0; dataPointNo_0++) {
2644              int offset_0 = tmp_0->getPointOffset(sampleNo_0,dataPointNo_0);
2645              int offset_1 = tmp_1->getPointOffset(sampleNo_0,dataPointNo_0);
2646              int offset_2 = tmp_2->getPointOffset(sampleNo_0,dataPointNo_0);
2647              const double *ptr_0 = &(arg_0_Z.getDataAtOffsetRO(offset_0));
2648              const double *ptr_1 = &(arg_1_Z.getDataAtOffsetRO(offset_1));
2649              double *ptr_2 = &(res.getDataAtOffsetRW(offset_2));
2650              tensor_binary_operation(size1, ptr_0[0], ptr_1, ptr_2, operation);
2651            }
2652          }
2653    
2654        }
2655        else {
2656          throw DataException("Error - C_TensorBinaryOperation: unknown combination of inputs");
2657        }
2658    
2659      } else if (0 == rank1) {
2660        if (arg_0_Z.isConstant()   && arg_1_Z.isConstant()) {
2661          res = Data(0.0, shape0, arg_1_Z.getFunctionSpace());      // DataConstant output
2662          const double *ptr_0 = &(arg_0_Z.getDataAtOffsetRO(0));
2663          const double *ptr_1 = &(arg_1_Z.getDataAtOffsetRO(0));
2664          double *ptr_2 = &(res.getDataAtOffsetRW(0));
2665          tensor_binary_operation(size0, ptr_0, ptr_1[0], ptr_2, operation);
2666        }
2667        else if (arg_0_Z.isConstant()   && arg_1_Z.isTagged()) {
2668    
2669          // Prepare the DataConstant input
2670          DataConstant* tmp_0=dynamic_cast<DataConstant*>(arg_0_Z.borrowData());
2671    
2672          // Borrow DataTagged input from Data object
2673          DataTagged* tmp_1=dynamic_cast<DataTagged*>(arg_1_Z.borrowData());
2674    
2675          // Prepare a DataTagged output 2
2676          res = Data(0.0, shape0, arg_1_Z.getFunctionSpace());      // DataTagged output
2677          res.tag();
2678          DataTagged* tmp_2=dynamic_cast<DataTagged*>(res.borrowData());
2679    
2680          // Prepare offset into DataConstant
2681          int offset_0 = tmp_0->getPointOffset(0,0);
2682          const double *ptr_0 = &(arg_0_Z.getDataAtOffsetRO(offset_0));
2683    
2684          //Get the pointers to the actual data
2685          const double *ptr_1 = &(tmp_1->getDefaultValueRO(0));
2686          double *ptr_2 = &(tmp_2->getDefaultValueRW(0));
2687    
2688          // Compute a result for the default
2689          tensor_binary_operation(size0, ptr_0, ptr_1[0], ptr_2, operation);
2690          // Compute a result for each tag
2691          const DataTagged::DataMapType& lookup_1=tmp_1->getTagLookup();
2692          DataTagged::DataMapType::const_iterator i; // i->first is a tag, i->second is an offset into memory
2693          for (i=lookup_1.begin();i!=lookup_1.end();i++) {
2694            tmp_2->addTag(i->first);
2695            const double *ptr_1 = &(tmp_1->getDataByTagRO(i->first,0));
2696            double *ptr_2 = &(tmp_2->getDataByTagRW(i->first,0));
2697            tensor_binary_operation(size0, ptr_0, ptr_1[0], ptr_2, operation);
2698          }
2699        }
2700        else if (arg_0_Z.isConstant()   && arg_1_Z.isExpanded()) {
2701    
2702          res = Data(0.0, shape0, arg_1_Z.getFunctionSpace(),true); // DataExpanded output
2703          DataConstant* tmp_0=dynamic_cast<DataConstant*>(arg_0_Z.borrowData());
2704          DataExpanded* tmp_1=dynamic_cast<DataExpanded*>(arg_1_Z.borrowData());
2705          DataExpanded* tmp_2=dynamic_cast<DataExpanded*>(res.borrowData());
2706    
2707          int sampleNo_1,dataPointNo_1;
2708          int numSamples_1 = arg_1_Z.getNumSamples();
2709          int numDataPointsPerSample_1 = arg_1_Z.getNumDataPointsPerSample();
2710          int offset_0 = tmp_0->getPointOffset(0,0);
2711          res.requireWrite();
2712          #pragma omp parallel for private(sampleNo_1,dataPointNo_1) schedule(static)
2713          for (sampleNo_1 = 0; sampleNo_1 < numSamples_1; sampleNo_1++) {
2714            for (dataPointNo_1 = 0; dataPointNo_1 < numDataPointsPerSample_1; dataPointNo_1++) {
2715              int offset_1 = tmp_1->getPointOffset(sampleNo_1,dataPointNo_1);
2716              int offset_2 = tmp_2->getPointOffset(sampleNo_1,dataPointNo_1);
2717              const double *ptr_0 = &(arg_0_Z.getDataAtOffsetRO(offset_0));
2718              const double *ptr_1 = &(arg_1_Z.getDataAtOffsetRO(offset_1));
2719              double *ptr_2 = &(res.getDataAtOffsetRW(offset_2));
2720              tensor_binary_operation(size0, ptr_0, ptr_1[0], ptr_2, operation);
2721            }
2722          }
2723    
2724        }
2725        else if (arg_0_Z.isTagged()     && arg_1_Z.isConstant()) {
2726    
2727          // Borrow DataTagged input from Data object
2728          DataTagged* tmp_0=dynamic_cast<DataTagged*>(arg_0_Z.borrowData());
2729    
2730          // Prepare the DataConstant input
2731          DataConstant* tmp_1=dynamic_cast<DataConstant*>(arg_1_Z.borrowData());
2732    
2733          // Prepare a DataTagged output 2
2734          res = Data(0.0, shape0, arg_0_Z.getFunctionSpace());      // DataTagged output
2735          res.tag();
2736          DataTagged* tmp_2=dynamic_cast<DataTagged*>(res.borrowData());
2737    
2738          // Prepare offset into DataConstant
2739          int offset_1 = tmp_1->getPointOffset(0,0);
2740          const double *ptr_1 = &(arg_1_Z.getDataAtOffsetRO(offset_1));
2741          // Get the pointers to the actual data
2742          const double *ptr_0 = &(tmp_0->getDefaultValueRO(0));
2743          double *ptr_2 = &(tmp_2->getDefaultValueRW(0));
2744          // Compute a result for the default
2745          tensor_binary_operation(size0, ptr_0, ptr_1[0], ptr_2, operation);
2746          // Compute a result for each tag
2747          const DataTagged::DataMapType& lookup_0=tmp_0->getTagLookup();
2748          DataTagged::DataMapType::const_iterator i; // i->first is a tag, i->second is an offset into memory
2749          for (i=lookup_0.begin();i!=lookup_0.end();i++) {
2750            tmp_2->addTag(i->first);
2751    /*        DataArrayView view_0 = tmp_0->getDataPointByTag(i->first);
2752            DataArrayView view_2 = tmp_2->getDataPointByTag(i->first);
2753            double *ptr_0 = &view_0.getData(0);
2754            double *ptr_2 = &view_2.getData(0);*/
2755            const double *ptr_0 = &(tmp_0->getDataByTagRO(i->first,0));
2756            double *ptr_2 = &(tmp_2->getDataByTagRW(i->first,0));
2757            tensor_binary_operation(size0, ptr_0, ptr_1[0], ptr_2, operation);
2758          }
2759    
2760        }
2761        else if (arg_0_Z.isTagged()     && arg_1_Z.isTagged()) {
2762    
2763          // Borrow DataTagged input from Data object
2764          DataTagged* tmp_0=dynamic_cast<DataTagged*>(arg_0_Z.borrowData());
2765    
2766          // Borrow DataTagged input from Data object
2767          DataTagged* tmp_1=dynamic_cast<DataTagged*>(arg_1_Z.borrowData());
2768    
2769          // Prepare a DataTagged output 2
2770          res = Data(0.0, shape0, arg_1_Z.getFunctionSpace());
2771          res.tag();        // DataTagged output
2772          DataTagged* tmp_2=dynamic_cast<DataTagged*>(res.borrowData());
2773    
2774          // Get the pointers to the actual data
2775          const double *ptr_0 = &(tmp_0->getDefaultValueRO(0));
2776          const double *ptr_1 = &(tmp_1->getDefaultValueRO(0));
2777          double *ptr_2 = &(tmp_2->getDefaultValueRW(0));
2778    
2779          // Compute a result for the default
2780          tensor_binary_operation(size0, ptr_0, ptr_1[0], ptr_2, operation);
2781          // Merge the tags
2782          DataTagged::DataMapType::const_iterator i; // i->first is a tag, i->second is an offset into memory
2783          const DataTagged::DataMapType& lookup_0=tmp_0->getTagLookup();
2784          const DataTagged::DataMapType& lookup_1=tmp_1->getTagLookup();
2785          for (i=lookup_0.begin();i!=lookup_0.end();i++) {
2786            tmp_2->addTag(i->first); // use tmp_2 to get correct shape
2787          }
2788          for (i=lookup_1.begin();i!=lookup_1.end();i++) {
2789            tmp_2->addTag(i->first);
2790          }
2791          // Compute a result for each tag
2792          const DataTagged::DataMapType& lookup_2=tmp_2->getTagLookup();
2793          for (i=lookup_2.begin();i!=lookup_2.end();i++) {
2794    //         DataArrayView view_0 = tmp_0->getDataPointByTag(i->first);
2795    //         DataArrayView view_1 = tmp_1->getDataPointByTag(i->first);
2796    //         DataArrayView view_2 = tmp_2->getDataPointByTag(i->first);
2797    //         double *ptr_0 = &view_0.getData(0);
2798    //         double *ptr_1 = &view_1.getData(0);
2799    //         double *ptr_2 = &view_2.getData(0);
2800    
2801            const double *ptr_0 = &(tmp_0->getDataByTagRO(i->first,0));
2802            const double *ptr_1 = &(tmp_1->getDataByTagRO(i->first,0));
2803            double *ptr_2 = &(tmp_2->getDataByTagRW(i->first,0));
2804            tensor_binary_operation(size0, ptr_0, ptr_1[0], ptr_2, operation);
2805          }
2806    
2807        }
2808        else if (arg_0_Z.isTagged()     && arg_1_Z.isExpanded()) {
2809    
2810          // After finding a common function space above the two inputs have the same numSamples and num DPPS
2811          res = Data(0.0, shape0, arg_1_Z.getFunctionSpace(),true); // DataExpanded output
2812          DataTagged*   tmp_0=dynamic_cast<DataTagged*>(arg_0_Z.borrowData());
2813          DataExpanded* tmp_1=dynamic_cast<DataExpanded*>(arg_1_Z.borrowData());
2814          DataExpanded* tmp_2=dynamic_cast<DataExpanded*>(res.borrowData());
2815    
2816          int sampleNo_0,dataPointNo_0;
2817          int numSamples_0 = arg_0_Z.getNumSamples();
2818          int numDataPointsPerSample_0 = arg_0_Z.getNumDataPointsPerSample();
2819          res.requireWrite();
2820          #pragma omp parallel for private(sampleNo_0,dataPointNo_0) schedule(static)
2821          for (sampleNo_0 = 0; sampleNo_0 < numSamples_0; sampleNo_0++) {
2822            int offset_0 = tmp_0->getPointOffset(sampleNo_0,0); // They're all the same, so just use #0
2823            const double *ptr_0 = &(arg_0_Z.getDataAtOffsetRO(offset_0));
2824            for (dataPointNo_0 = 0; dataPointNo_0 < numDataPointsPerSample_0; dataPointNo_0++) {
2825              int offset_1 = tmp_1->getPointOffset(sampleNo_0,dataPointNo_0);
2826              int offset_2 = tmp_2->getPointOffset(sampleNo_0,dataPointNo_0);
2827              const double *ptr_1 = &(arg_1_Z.getDataAtOffsetRO(offset_1));
2828              double *ptr_2 = &(res.getDataAtOffsetRW(offset_2));
2829              tensor_binary_operation(size0, ptr_0, ptr_1[0], ptr_2, operation);
2830            }
2831          }
2832    
2833        }
2834        else if (arg_0_Z.isExpanded()   && arg_1_Z.isConstant()) {
2835          res = Data(0.0, shape0, arg_1_Z.getFunctionSpace(),true); // DataExpanded output
2836          DataExpanded* tmp_0=dynamic_cast<DataExpanded*>(arg_0_Z.borrowData());
2837          DataConstant* tmp_1=dynamic_cast<DataConstant*>(arg_1_Z.borrowData());
2838          DataExpanded* tmp_2=dynamic_cast<DataExpanded*>(res.borrowData());
2839    
2840          int sampleNo_0,dataPointNo_0;
2841          int numSamples_0 = arg_0_Z.getNumSamples();
2842          int numDataPointsPerSample_0 = arg_0_Z.getNumDataPointsPerSample();
2843          int offset_1 = tmp_1->getPointOffset(0,0);
2844          res.requireWrite();
2845          #pragma omp parallel for private(sampleNo_0,dataPointNo_0) schedule(static)
2846          for (sampleNo_0 = 0; sampleNo_0 < numSamples_0; sampleNo_0++) {
2847            for (dataPointNo_0 = 0; dataPointNo_0 < numDataPointsPerSample_0; dataPointNo_0++) {
2848              int offset_0 = tmp_0->getPointOffset(sampleNo_0,dataPointNo_0);
2849              int offset_2 = tmp_2->getPointOffset(sampleNo_0,dataPointNo_0);
2850              const double *ptr_0 = &(arg_0_Z.getDataAtOffsetRO(offset_0));
2851              const double *ptr_1 = &(arg_1_Z.getDataAtOffsetRO(offset_1));
2852              double *ptr_2 = &(res.getDataAtOffsetRW(offset_2));
2853              tensor_binary_operation(size0, ptr_0, ptr_1[0], ptr_2, operation);
2854            }
2855          }
2856    
2857    
2858        }
2859        else if (arg_0_Z.isExpanded()   && arg_1_Z.isTagged()) {
2860    
2861          // After finding a common function space above the two inputs have the same numSamples and num DPPS
2862          res = Data(0.0, shape0, arg_1_Z.getFunctionSpace(),true); // DataExpanded output
2863          DataExpanded* tmp_0=dynamic_cast<DataExpanded*>(arg_0_Z.borrowData());
2864          DataTagged*   tmp_1=dynamic_cast<DataTagged*>(arg_1_Z.borrowData());
2865          DataExpanded* tmp_2=dynamic_cast<DataExpanded*>(res.borrowData());
2866    
2867          int sampleNo_0,dataPointNo_0;
2868          int numSamples_0 = arg_0_Z.getNumSamples();
2869          int numDataPointsPerSample_0 = arg_0_Z.getNumDataPointsPerSample();
2870          res.requireWrite();
2871          #pragma omp parallel for private(sampleNo_0,dataPointNo_0) schedule(static)
2872          for (sampleNo_0 = 0; sampleNo_0 < numSamples_0; sampleNo_0++) {
2873            int offset_1 = tmp_1->getPointOffset(sampleNo_0,0);
2874            const double *ptr_1 = &(arg_1_Z.getDataAtOffsetRO(offset_1));
2875            for (dataPointNo_0 = 0; dataPointNo_0 < numDataPointsPerSample_0; dataPointNo_0++) {
2876              int offset_0 = tmp_0->getPointOffset(sampleNo_0,dataPointNo_0);
2877              int offset_2 = tmp_2->getPointOffset(sampleNo_0,dataPointNo_0);
2878              const double *ptr_0 = &(arg_0_Z.getDataAtOffsetRO(offset_0));
2879              double *ptr_2 = &(res.getDataAtOffsetRW(offset_2));
2880              tensor_binary_operation(size0, ptr_0, ptr_1[0], ptr_2, operation);
2881            }
2882          }
2883    
2884        }
2885        else if (arg_0_Z.isExpanded()   && arg_1_Z.isExpanded()) {
2886    
2887          // After finding a common function space above the two inputs have the same numSamples and num DPPS
2888          res = Data(0.0, shape0, arg_1_Z.getFunctionSpace(),true); // DataExpanded output
2889          DataExpanded* tmp_0=dynamic_cast<DataExpanded*>(arg_0_Z.borrowData());
2890          DataExpanded* tmp_1=dynamic_cast<DataExpanded*>(arg_1_Z.borrowData());
2891          DataExpanded* tmp_2=dynamic_cast<DataExpanded*>(res.borrowData());
2892    
2893          int sampleNo_0,dataPointNo_0;
2894          int numSamples_0 = arg_0_Z.getNumSamples();
2895          int numDataPointsPerSample_0 = arg_0_Z.getNumDataPointsPerSample();
2896          res.requireWrite();
2897          #pragma omp parallel for private(sampleNo_0,dataPointNo_0) schedule(static)
2898          for (sampleNo_0 = 0; sampleNo_0 < numSamples_0; sampleNo_0++) {
2899            for (dataPointNo_0 = 0; dataPointNo_0 < numDataPointsPerSample_0; dataPointNo_0++) {
2900              int offset_0 = tmp_0->getPointOffset(sampleNo_0,dataPointNo_0);
2901              int offset_1 = tmp_1->getPointOffset(sampleNo_0,dataPointNo_0);
2902              int offset_2 = tmp_2->getPointOffset(sampleNo_0,dataPointNo_0);
2903              const double *ptr_0 = &(arg_0_Z.getDataAtOffsetRO(offset_0));
2904              const double *ptr_1 = &(arg_1_Z.getDataAtOffsetRO(offset_1));
2905              double *ptr_2 = &(res.getDataAtOffsetRW(offset_2));
2906              tensor_binary_operation(size0, ptr_0, ptr_1[0], ptr_2, operation);
2907            }
2908          }
2909    
2910        }
2911        else {
2912          throw DataException("Error - C_TensorBinaryOperation: unknown combination of inputs");
2913        }
2914    
2915      } else {
2916        throw DataException("Error - C_TensorBinaryOperation: arguments have incompatible shapes");
2917      }
2918    
2919      return res;
2920    }
2921    
2922    template <typename UnaryFunction>
2923    Data
2924    C_TensorUnaryOperation(Data const &arg_0,
2925                           UnaryFunction operation)
2926    {
2927      if (arg_0.isEmpty())  // do this before we attempt to interpolate
2928      {
2929         throw DataException("Error - Operations not permitted on instances of DataEmpty.");
2930      }
2931      if (arg_0.isLazy())
2932      {
2933         throw DataException("Error - Operations not permitted on lazy data.");
2934    }    }
2935    return result;    // Interpolate if necessary and find an appropriate function space
2936      Data arg_0_Z = Data(arg_0);
2937    
2938      // Get rank and shape of inputs
2939    //  int rank0 = arg_0_Z.getDataPointRank();
2940      const DataTypes::ShapeType& shape0 = arg_0_Z.getDataPointShape();
2941      int size0 = arg_0_Z.getDataPointSize();
2942    
2943      // Declare output Data object
2944      Data res;
2945    
2946      if (arg_0_Z.isConstant()) {
2947        res = Data(0.0, shape0, arg_0_Z.getFunctionSpace());      // DataConstant output
2948    //     double *ptr_0 = &((arg_0_Z.getPointDataView().getData())[0]);
2949    //     double *ptr_2 = &((res.getPointDataView().getData())[0]);
2950        const double *ptr_0 = &(arg_0_Z.getDataAtOffsetRO(0));
2951        double *ptr_2 = &(res.getDataAtOffsetRW(0));
2952        tensor_unary_operation(size0, ptr_0, ptr_2, operation);
2953      }
2954      else if (arg_0_Z.isTagged()) {
2955    
2956        // Borrow DataTagged input from Data object
2957        DataTagged* tmp_0=dynamic_cast<DataTagged*>(arg_0_Z.borrowData());
2958    
2959        // Prepare a DataTagged output 2
2960        res = Data(0.0, shape0, arg_0_Z.getFunctionSpace());   // DataTagged output
2961        res.tag();
2962        DataTagged* tmp_2=dynamic_cast<DataTagged*>(res.borrowData());
2963    
2964        // Get the pointers to the actual data
2965        const double *ptr_0 = &(tmp_0->getDefaultValueRO(0));
2966        double *ptr_2 = &(tmp_2->getDefaultValueRW(0));
2967        // Compute a result for the default
2968        tensor_unary_operation(size0, ptr_0, ptr_2, operation);
2969        // Compute a result for each tag
2970        const DataTagged::DataMapType& lookup_0=tmp_0->getTagLookup();
2971        DataTagged::DataMapType::const_iterator i; // i->first is a tag, i->second is an offset into memory
2972        for (i=lookup_0.begin();i!=lookup_0.end();i++) {
2973          tmp_2->addTag(i->first);
2974          const double *ptr_0 = &(tmp_0->getDataByTagRO(i->first,0));
2975          double *ptr_2 = &(tmp_2->getDataByTagRW(i->first,0));
2976          tensor_unary_operation(size0, ptr_0, ptr_2, operation);
2977        }
2978    
2979      }
2980      else if (arg_0_Z.isExpanded()) {
2981    
2982        res = Data(0.0, shape0, arg_0_Z.getFunctionSpace(),true); // DataExpanded output
2983        DataExpanded* tmp_0=dynamic_cast<DataExpanded*>(arg_0_Z.borrowData());
2984        DataExpanded* tmp_2=dynamic_cast<DataExpanded*>(res.borrowData());
2985    
2986        int sampleNo_0,dataPointNo_0;
2987        int numSamples_0 = arg_0_Z.getNumSamples();
2988        int numDataPointsPerSample_0 = arg_0_Z.getNumDataPointsPerSample();
2989        #pragma omp parallel for private(sampleNo_0,dataPointNo_0) schedule(static)
2990        for (sampleNo_0 = 0; sampleNo_0 < numSamples_0; sampleNo_0++) {
2991          for (dataPointNo_0 = 0; dataPointNo_0 < numDataPointsPerSample_0; dataPointNo_0++) {
2992            int offset_0 = tmp_0->getPointOffset(sampleNo_0,dataPointNo_0);
2993            int offset_2 = tmp_2->getPointOffset(sampleNo_0,dataPointNo_0);
2994            const double *ptr_0 = &(arg_0_Z.getDataAtOffsetRO(offset_0));
2995            double *ptr_2 = &(res.getDataAtOffsetRW(offset_2));
2996            tensor_unary_operation(size0, ptr_0, ptr_2, operation);
2997          }
2998        }
2999      }
3000      else {
3001        throw DataException("Error - C_TensorUnaryOperation: unknown combination of inputs");
3002      }
3003    
3004      return res;
3005  }  }
3006    
3007  }  }

Legend:
Removed from v.97  
changed lines
  Added in v.2458

  ViewVC Help
Powered by ViewVC 1.1.26