/[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

temp/escript/src/Data.h revision 1387 by trankine, Fri Jan 11 07:45:26 2008 UTC trunk/escript/src/Data.h revision 2519 by jfenwick, Mon Jul 6 00:43:08 2009 UTC
# Line 1  Line 1 
1    
 /* $Id$ */  
   
2  /*******************************************************  /*******************************************************
3   *  *
4   *           Copyright 2003-2007 by ACceSS MNRF  * Copyright (c) 2003-2008 by University of Queensland
5   *       Copyright 2007 by University of Queensland  * Earth Systems Science Computational Center (ESSCC)
6   *  * http://www.uq.edu.au/esscc
7   *                http://esscc.uq.edu.au  *
8   *        Primary Business: Queensland, Australia  * Primary Business: Queensland, Australia
9   *  Licensed under the Open Software License version 3.0  * Licensed under the Open Software License version 3.0
10   *     http://www.opensource.org/licenses/osl-3.0.php  * http://www.opensource.org/licenses/osl-3.0.php
11   *  *
12   *******************************************************/  *******************************************************/
13    
14    
15  /** \file Data.h */  /** \file Data.h */
16    
# Line 19  Line 18 
18  #define DATA_H  #define DATA_H
19  #include "system_dep.h"  #include "system_dep.h"
20    
21    #include "DataTypes.h"
22  #include "DataAbstract.h"  #include "DataAbstract.h"
23  #include "DataAlgorithm.h"  #include "DataAlgorithm.h"
24  #include "FunctionSpace.h"  #include "FunctionSpace.h"
# Line 26  Line 26 
26  #include "UnaryOp.h"  #include "UnaryOp.h"
27  #include "DataException.h"  #include "DataException.h"
28    
29    
30  extern "C" {  extern "C" {
31  #include "DataC.h"  #include "DataC.h"
32  /* #include "paso/Paso.h" doesn't belong in this file...causes trouble for BruceFactory.cpp */  //#include <omp.h>
33  }  }
34    
35  #include "esysmpi.h"  #include "esysmpi.h"
36  #include <string>  #include <string>
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>
42  #include <boost/python/tuple.hpp>  #include <boost/python/tuple.hpp>
43  #include <boost/python/numeric.hpp>  
44    #include "BufferGroup.h"
45    
46  namespace escript {  namespace escript {
47    
# Line 47  namespace escript { Line 50  namespace escript {
50  class DataConstant;  class DataConstant;
51  class DataTagged;  class DataTagged;
52  class DataExpanded;  class DataExpanded;
53    class DataLazy;
54    
55  /**  /**
56     \brief     \brief
57     Data creates the appropriate Data object for the given construction     Data represents a collection of datapoints.
    arguments.  
58    
59     Description:     Description:
60     Data is essentially a factory class which creates the appropriate Data     Internally, the datapoints are actually stored by a DataAbstract object.
61     object for the given construction arguments. It retains control over     The specific instance of DataAbstract used may vary over the lifetime
62     the object created for the lifetime of the object.     of the Data object.
63     The type of Data object referred to may change during the lifetime of     Some methods on this class return references (eg getShape()).
64     the Data object.     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    
# Line 101  class Data { Line 106  class Data {
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    
      \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.  
   */  
112    ESCRIPT_DLL_API    ESCRIPT_DLL_API
113    Data(const DataArrayView& value,    Data(const DataTypes::ValueType& value,
114         const FunctionSpace& what=FunctionSpace(),           const DataTypes::ShapeType& shape,
115         bool expanded=false);                   const FunctionSpace& what=FunctionSpace(),
116                     bool expanded=false);
117    
118    /**    /**
119       \brief       \brief
120       Constructor which creates a Data from a DataArrayView shape.       Constructor which creates a Data with points having the specified shape.
121    
122       \param value - Input - Single value applied to all Data.       \param value - Input - Single value applied to all Data.
123       \param dataPointShape - Input - The shape of each data point.       \param dataPointShape - Input - The shape of each data point.
# Line 128  class Data { Line 128  class Data {
128    */    */
129    ESCRIPT_DLL_API    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 141  class Data { Line 141  class Data {
141    */    */
142    ESCRIPT_DLL_API    ESCRIPT_DLL_API
143    Data(const Data& inData,    Data(const Data& inData,
144         const DataArrayView::RegionType& region);         const DataTypes::RegionType& region);
145    
146    /**    /**
147       \brief       \brief
148       Constructor which will create Tagged data if expanded is false.       Constructor which copies data from any object that can be treated like a python array/sequence.
      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.  
     ==>*  
   */  
   ESCRIPT_DLL_API  
   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.  
   */  
   ESCRIPT_DLL_API  
   Data(const boost::python::numeric::array& value,  
        const FunctionSpace& what=FunctionSpace(),  
        bool expanded=false);  
   
   /**  
      \brief  
      Constructor which copies data from any object that can be converted into  
      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 198  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.
# Line 217  class Data { Line 180  class Data {
180         const boost::python::tuple& shape=boost::python::make_tuple(),         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       Destructor       Destructor
# Line 225  class Data { Line 204  class Data {
204    ~Data();    ~Data();
205    
206    /**    /**
207       \brief       \brief Make this object a deep copy of "other".
      Perform a deep copy.  
208    */    */
209    ESCRIPT_DLL_API    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    
# Line 247  class Data { Line 248  class Data {
248    
249    /**    /**
250       \brief       \brief
251       Returns trueif the data object is protected against update       Returns true, if the data object is protected against update
252    
253    */    */
254    ESCRIPT_DLL_API    ESCRIPT_DLL_API
255    bool    bool
256    isProtected() const;    isProtected() const;
257    
258    /**  
259       \brief  /**
260       Return the values of a data point on this process     \brief
261    */     Return the value of a data point as a python tuple.
262    */
263    ESCRIPT_DLL_API    ESCRIPT_DLL_API
264    const boost::python::numeric::array    const boost::python::object
265    getValueOfDataPoint(int dataPointNo);    getValueOfDataPointAsTuple(int dataPointNo);
266    
267    /**    /**
268       \brief       \brief
# Line 272  class Data { Line 274  class Data {
274    
275    /**    /**
276       \brief       \brief
277       sets the values of a data-point from a numarray object on this process       sets the values of a data-point from a array-like object on this process
278    */    */
279    ESCRIPT_DLL_API    ESCRIPT_DLL_API
280    void    void
281    setValueOfDataPointToArray(int dataPointNo, const boost::python::numeric::array&);    setValueOfDataPointToArray(int dataPointNo, const boost::python::object&);
282    
283    /**    /**
284       \brief       \brief
# Line 287  class Data { Line 289  class Data {
289    setValueOfDataPoint(int dataPointNo, const double);    setValueOfDataPoint(int dataPointNo, const double);
290    
291    /**    /**
292       \brief       \brief Return a data point across all processors as a python tuple.
      Return the value of the specified data-point across all processors  
293    */    */
294    ESCRIPT_DLL_API    ESCRIPT_DLL_API
295    const boost::python::numeric::array    const boost::python::object
296    getValueOfGlobalDataPoint(int procNo, int dataPointNo);    getValueOfGlobalDataPointAsTuple(int procNo, int dataPointNo);
297    
298    /**    /**
299       \brief       \brief
300       Return the tag number associated with the given data-point.       Return the tag number associated with the given data-point.
301    
      The data-point number here corresponds to the data-point number in the  
      numarray returned by convertToNumArray.  
302    */    */
303    ESCRIPT_DLL_API    ESCRIPT_DLL_API
304    int    int
# Line 313  class Data { Line 312  class Data {
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.
# Line 322  class Data { Line 323  class Data {
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    ESCRIPT_DLL_API    ESCRIPT_DLL_API
329    inline    size_t
330    std::string    getSampleBufferSize() const;
331    toString() const  
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    ESCRIPT_DLL_API    ESCRIPT_DLL_API
339    inline    std::string
340    const DataArrayView&    toString() const;
   getPointDataView() const  
   {  
      return m_data->getPointDataView();  
   }  
341    
342    /**    /**
343       \brief       \brief
# Line 360  class Data { Line 352  class Data {
352       If possible convert this Data to DataTagged. This will only allow       If possible convert this Data to DataTagged. This will only allow
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    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    ESCRIPT_DLL_API
386    bool    bool
# Line 376  class Data { Line 388  class Data {
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    ESCRIPT_DLL_API
# Line 391  class Data { Line 413  class Data {
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    ESCRIPT_DLL_API
435    bool    bool
# Line 424  class Data { Line 461  class Data {
461    */    */
462    ESCRIPT_DLL_API    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.
# Line 444  class Data { Line 497  class Data {
497    */    */
498    ESCRIPT_DLL_API    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    /**    /**
# Line 484  class Data { Line 537  class Data {
537    {    {
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       \brief
556       dumps the object into a netCDF file       dumps the object into a netCDF file
# Line 491  class Data { Line 558  class Data {
558    ESCRIPT_DLL_API    ESCRIPT_DLL_API
559    void    void
560    dump(const std::string fileName) const;    dump(const std::string fileName) const;
561    
562     /**
563      \brief returns the values of the object as a list of tuples (one for each datapoint).
564    
565      \param scalarastuple If true, scalar data will produce single valued tuples [(1,) (2,) ...]
566    If false, the result is a list of scalars [1, 2, ...]
567     */
568      ESCRIPT_DLL_API
569      const boost::python::object
570      toListOfTuples(bool scalarastuple=true);
571    
572    
573     /**
574        \brief
575        Return the sample data for the given sample no. This is not the
576        preferred interface but is provided for use by C code.
577        The bufferg parameter is only required for LazyData.
578        \param sampleNo - Input - the given sample no.
579        \param bufferg - A buffer to compute (and store) sample data in will be selected from this group.
580        \return pointer to the sample data.
581    */
582      ESCRIPT_DLL_API
583      inline
584      const DataAbstract::ValueType::value_type*
585      getSampleDataRO(DataAbstract::ValueType::size_type sampleNo, BufferGroup* bufferg=0);
586    
587    
588    /**    /**
589       \brief       \brief
590       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
591       preferred interface but is provided for use by C code.       preferred interface but is provided for use by C code.
592       \param sampleNo - Input - the given sample no.       \param sampleNo - Input - the given sample no.
593         \return pointer to the sample data.
594    */    */
595    ESCRIPT_DLL_API    ESCRIPT_DLL_API
596    inline    inline
597    DataAbstract::ValueType::value_type*    DataAbstract::ValueType::value_type*
598    getSampleData(DataAbstract::ValueType::size_type sampleNo)    getSampleDataRW(DataAbstract::ValueType::size_type sampleNo);
599    {  
     return m_data->getSampleData(sampleNo);  
   }  
600    
601    /**    /**
602       \brief       \brief
# Line 521  class Data { Line 614  class Data {
614    
615    /**    /**
616       \brief       \brief
617       Return a view into the data for the data point specified.       Return a reference into the DataVector which points to the specified data point.
      NOTE: Construction of the DataArrayView is a relatively expensive  
      operation.  
618       \param sampleNo - Input -       \param sampleNo - Input -
619       \param dataPointNo - Input -       \param dataPointNo - Input -
620    */    */
621    ESCRIPT_DLL_API    ESCRIPT_DLL_API
622      DataTypes::ValueType::const_reference
623      getDataPointRO(int sampleNo, int dataPointNo);
624    
625      /**
626         \brief
627         Return a reference into the DataVector which points to the specified data point.
628         \param sampleNo - Input -
629         \param dataPointNo - Input -
630      */
631      ESCRIPT_DLL_API
632      DataTypes::ValueType::reference
633      getDataPointRW(int sampleNo, int dataPointNo);
634    
635    
636    
637      /**
638         \brief
639         Return the offset for the given sample and point within the sample
640      */
641      ESCRIPT_DLL_API
642    inline    inline
643    DataArrayView    DataTypes::ValueType::size_type
644    getDataPoint(int sampleNo,    getDataOffset(int sampleNo,
645                 int dataPointNo)                 int dataPointNo)
646    {    {
647                  return m_data->getDataPoint(sampleNo,dataPointNo);        return m_data->getPointOffset(sampleNo,dataPointNo);
648    }    }
649    
650    /**    /**
# Line 541  class Data { Line 652  class Data {
652       Return a reference to the data point shape.       Return a reference to the data point shape.
653    */    */
654    ESCRIPT_DLL_API    ESCRIPT_DLL_API
655    const DataArrayView::ShapeType&    inline
656    getDataPointShape() const;    const DataTypes::ShapeType&
657      getDataPointShape() const
658      {
659        return m_data->getShape();
660      }
661    
662    /**    /**
663       \brief       \brief
# Line 566  class Data { Line 681  class Data {
681       Return the number of doubles stored for this Data.       Return the number of doubles stored for this Data.
682    */    */
683    ESCRIPT_DLL_API    ESCRIPT_DLL_API
684    DataArrayView::ValueType::size_type    DataTypes::ValueType::size_type
685    getLength() const;    getLength() const;
686    
687    
# Line 576  class Data { Line 691  class Data {
691       Assign the given value to the tag assocciated with name. Implicitly converts this       Assign the given value to the tag assocciated with name. Implicitly converts this
692       object to type DataTagged. Throws an exception if this object       object to type DataTagged. Throws an exception if this object
693       cannot be converted to a DataTagged object or name cannot be mapped onto a tag key.       cannot be converted to a DataTagged object or name cannot be mapped onto a tag key.
694       \param tagKey - Input - Integer key.       \param name - Input - name of tag.
695       \param value - Input - Value to associate with given key.       \param value - Input - Value to associate with given key.
     ==>*  
696    */    */
697    ESCRIPT_DLL_API    ESCRIPT_DLL_API
698    void    void
# Line 605  class Data { Line 719  class Data {
719       object to type DataTagged if it is constant.       object to type DataTagged if it is constant.
720    
721       \param tagKey - Input - Integer key.       \param tagKey - Input - Integer key.
722         \param pointshape - Input - The shape of the value parameter
723       \param value - Input - Value to associate with given key.       \param value - Input - Value to associate with given key.
724      ==>*       \param dataOffset - Input - Offset of the begining of the point within the value parameter
725    */    */
726    ESCRIPT_DLL_API    ESCRIPT_DLL_API
727    void    void
728    setTaggedValueFromCPP(int tagKey,    setTaggedValueFromCPP(int tagKey,
729                          const DataArrayView& value);              const DataTypes::ShapeType& pointshape,
730                            const DataTypes::ValueType& value,
731                int dataOffset=0);
732    
733    
734    
735    /**    /**
736      \brief      \brief
# Line 659  class Data { Line 778  class Data {
778    grad() const;    grad() const;
779    
780    /**    /**
781       \brief      \brief
782       Calculate the integral over the function space domain.       Calculate the integral over the function space domain as a python tuple.
      *  
783    */    */
784    ESCRIPT_DLL_API    ESCRIPT_DLL_API
785    boost::python::numeric::array    boost::python::object
786    integrate() const;    integrateToTuple_const() const;
787    
788    
789      /**
790        \brief
791         Calculate the integral over the function space domain as a python tuple.
792      */
793      ESCRIPT_DLL_API
794      boost::python::object
795      integrateToTuple();
796    
797    
798    
799    /**    /**
800       \brief       \brief
# Line 732  class Data { Line 861  class Data {
861    /**    /**
862       \brief       \brief
863       Return the maximum absolute value of this Data object.       Return the maximum absolute value of this Data object.
864       *  
865         The method is not const because lazy data needs to be expanded before Lsup can be computed.
866         The _const form can be used when the Data object is const, however this will only work for
867         Data which is not Lazy.
868    
869         For Data which contain no samples (or tagged Data for which no tags in use have a value)
870         zero is returned.
871    */    */
872    ESCRIPT_DLL_API    ESCRIPT_DLL_API
873    double    double
874    Lsup() const;    Lsup();
875    
876      ESCRIPT_DLL_API
877      double
878      Lsup_const() const;
879    
880    
881    /**    /**
882       \brief       \brief
883       Return the maximum value of this Data object.       Return the maximum value of this Data object.
884       *  
885         The method is not const because lazy data needs to be expanded before sup can be computed.
886         The _const form can be used when the Data object is const, however this will only work for
887         Data which is not Lazy.
888    
889         For Data which contain no samples (or tagged Data for which no tags in use have a value)
890         a large negative value is returned.
891    */    */
892    ESCRIPT_DLL_API    ESCRIPT_DLL_API
893    double    double
894    sup() const;    sup();
895    
896      ESCRIPT_DLL_API
897      double
898      sup_const() const;
899    
900    
901    /**    /**
902       \brief       \brief
903       Return the minimum value of this Data object.       Return the minimum value of this Data object.
904       *  
905         The method is not const because lazy data needs to be expanded before inf can be computed.
906         The _const form can be used when the Data object is const, however this will only work for
907         Data which is not Lazy.
908    
909         For Data which contain no samples (or tagged Data for which no tags in use have a value)
910         a large positive value is returned.
911    */    */
912    ESCRIPT_DLL_API    ESCRIPT_DLL_API
913    double    double
914    inf() const;    inf();
915    
916      ESCRIPT_DLL_API
917      double
918      inf_const() const;
919    
920    
921    
922    /**    /**
923       \brief       \brief
# Line 786  class Data { Line 949  class Data {
949    /**    /**
950       \brief       \brief
951       Return the (sample number, data-point number) of the data point with       Return the (sample number, data-point number) of the data point with
952       the minimum value in this Data object.       the minimum component value in this Data object.
953         \note If you are working in python, please consider using Locator
954    instead of manually manipulating process and point IDs.
955    */    */
956    ESCRIPT_DLL_API    ESCRIPT_DLL_API
957    const boost::python::tuple    const boost::python::tuple
958    minGlobalDataPoint() const;    minGlobalDataPoint() const;
959    
960      /**
961         \brief
962         Return the (sample number, data-point number) of the data point with
963         the minimum component value in this Data object.
964         \note If you are working in python, please consider using Locator
965    instead of manually manipulating process and point IDs.
966      */
967    ESCRIPT_DLL_API    ESCRIPT_DLL_API
968    void    const boost::python::tuple
969    calc_minGlobalDataPoint(int& ProcNo,  int& DataPointNo) const;    maxGlobalDataPoint() const;
970    
971    
972    
973    /**    /**
974       \brief       \brief
975       Return the sign of each data point of this Data object.       Return the sign of each data point of this Data object.
# Line 1212  class Data { Line 1387  class Data {
1387    */    */
1388    ESCRIPT_DLL_API    ESCRIPT_DLL_API
1389    Data    Data
1390    getSlice(const DataArrayView::RegionType& region) const;    getSlice(const DataTypes::RegionType& region) const;
1391    
1392    /**    /**
1393       \brief       \brief
# Line 1225  class Data { Line 1400  class Data {
1400    ESCRIPT_DLL_API    ESCRIPT_DLL_API
1401    void    void
1402    setSlice(const Data& value,    setSlice(const Data& value,
1403             const DataArrayView::RegionType& region);             const DataTypes::RegionType& region);
   
   /**  
      \brief  
      Archive the current Data object to the given file.  
      \param fileName - Input - file to archive to.  
   */  
   ESCRIPT_DLL_API  
   void  
   archiveData(const std::string fileName);  
   
   /**  
      \brief  
      Extract the Data object archived in the given file, overwriting  
      the current Data object.  
      Note - the current object must be of type DataEmpty.  
      \param fileName - Input - file to extract from.  
      \param fspace - Input - a suitable FunctionSpace descibing the data.  
   */  
   ESCRIPT_DLL_API  
   void  
   extractData(const std::string fileName,  
               const FunctionSpace& fspace);  
   
1404    
1405    /**    /**
1406       \brief       \brief
# Line 1290  class Data { Line 1442  class Data {
1442    /**    /**
1443       \brief       \brief
1444       return the object produced by the factory, which is a DataConstant or DataExpanded       return the object produced by the factory, which is a DataConstant or DataExpanded
1445        TODO Ownership of this object should be explained in doco.
1446    */    */
1447    ESCRIPT_DLL_API    ESCRIPT_DLL_API
1448          DataAbstract*          DataAbstract*
1449          borrowData(void) const;          borrowData(void) const;
1450    
1451      ESCRIPT_DLL_API
1452            DataAbstract_ptr
1453            borrowDataPtr(void) const;
1454    
1455      ESCRIPT_DLL_API
1456            DataReady_ptr
1457            borrowReadyPtr(void) const;
1458    
1459    
1460    
1461      /**
1462         \brief
1463         Return a pointer to the beginning of the datapoint at the specified offset.
1464         TODO Eventually these should be inlined.
1465         \param i - position(offset) in the underlying datastructure
1466      */
1467    
1468      ESCRIPT_DLL_API
1469            DataTypes::ValueType::const_reference
1470            getDataAtOffsetRO(DataTypes::ValueType::size_type i);
1471    
1472    
1473      ESCRIPT_DLL_API
1474            DataTypes::ValueType::reference
1475            getDataAtOffsetRW(DataTypes::ValueType::size_type i);
1476    
1477    
1478    
1479    /**
1480       \brief Create a buffer for use by getSample
1481       Allocates a DataVector large enough for DataLazy::resolveSample to operate on for the current Data.
1482       Do not use this buffer for other Data instances (unless you are sure they will be the same size).
1483      
1484       In multi-threaded sections, this needs to be called on each thread.
1485    
1486       \return A BufferGroup* if Data is lazy, NULL otherwise.
1487       \warning This pointer must be deallocated using freeSampleBuffer to avoid cross library memory issues.
1488    */
1489      ESCRIPT_DLL_API
1490      BufferGroup*
1491      allocSampleBuffer() const;
1492    
1493    /**
1494       \brief Free a buffer allocated with allocSampleBuffer.
1495       \param buffer Input - pointer to the buffer to deallocate.
1496    */
1497    ESCRIPT_DLL_API void freeSampleBuffer(BufferGroup* buffer);
1498    
1499   protected:   protected:
1500    
1501   private:   private:
1502    
1503      double
1504      LsupWorker() const;
1505    
1506      double
1507      supWorker() const;
1508    
1509      double
1510      infWorker() const;
1511    
1512      boost::python::object
1513      integrateWorker() const;
1514    
1515      void
1516      calc_minGlobalDataPoint(int& ProcNo,  int& DataPointNo) const;
1517    
1518      void
1519      calc_maxGlobalDataPoint(int& ProcNo,  int& DataPointNo) const;
1520    
1521    
1522    /**    /**
1523       \brief       \brief
1524       Check *this and the right operand are compatible. Throws       Check *this and the right operand are compatible. Throws
# Line 1370  class Data { Line 1590  class Data {
1590       \brief       \brief
1591       Construct a Data object of the appropriate type.       Construct a Data object of the appropriate type.
1592    */    */
1593    template <class IValueType>  
1594    void    void
1595    initialise(const IValueType& value,    initialise(const DataTypes::ValueType& value,
1596             const DataTypes::ShapeType& shape,
1597               const FunctionSpace& what,               const FunctionSpace& what,
1598               bool expanded);               bool expanded);
1599    
1600      void
1601      initialise(const WrappedArray& value,
1602                     const FunctionSpace& what,
1603                     bool expanded);
1604    
1605    //    //
1606    // flag to protect the data object against any update    // flag to protect the data object against any update
1607    bool m_protected;    bool m_protected;
1608      mutable bool m_shared;
1609      bool m_lazy;
1610    
1611    //    //
1612    // pointer to the actual data object    // pointer to the actual data object
1613    boost::shared_ptr<DataAbstract> m_data;  //   boost::shared_ptr<DataAbstract> m_data;
1614      DataAbstract_ptr m_data;
1615    
1616  };  // If possible please use getReadyPtr instead.
1617    // But see warning below.
1618      const DataReady*
1619      getReady() const;
1620    
1621  template <class IValueType>    DataReady*
1622  void    getReady();
1623  Data::initialise(const IValueType& value,  
1624                   const FunctionSpace& what,  
1625                   bool expanded)  // Be wary of using this for local operations since it (temporarily) increases reference count.
1626  {  // If you are just using this to call a method on DataReady instead of DataAbstract consider using
1627    // getReady() instead
1628      DataReady_ptr
1629      getReadyPtr();
1630    
1631      const_DataReady_ptr
1632      getReadyPtr() const;
1633    
1634    
1635      /**
1636       \brief Update the Data's shared flag
1637       This indicates that the DataAbstract used by this object is now shared (or no longer shared).
1638       For internal use only.
1639      */
1640      void updateShareStatus(bool nowshared) const
1641      {
1642        m_shared=nowshared;     // m_shared is mutable
1643      }
1644    
1645      // In the isShared() method below:
1646      // A problem would occur if m_data (the address pointed to) were being modified
1647      // while the call m_data->is_shared is being executed.
1648      //
1649      // Q: So why do I think this code can be thread safe/correct?
1650      // A: We need to make some assumptions.
1651      // 1. We assume it is acceptable to return true under some conditions when we aren't shared.
1652      // 2. We assume that no constructions or assignments which will share previously unshared
1653      //    will occur while this call is executing. This is consistent with the way Data:: and C are written.
1654    //    //
1655    // 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
1656    // 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.
1657    // 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.
1658    // within the shared_ptr constructor.    // For any threads executing before the flag switches they will assume the object is still shared.
1659    if (expanded) {    bool isShared() const
1660      DataAbstract* temp=new DataExpanded(value,what);    {
1661      boost::shared_ptr<DataAbstract> temp_data(temp);      return m_shared;
1662      m_data=temp_data;  /*  if (m_shared) return true;
1663    } else {      if (m_data->isShared())        
1664      DataAbstract* temp=new DataConstant(value,what);      {                  
1665      boost::shared_ptr<DataAbstract> temp_data(temp);          updateShareStatus(true);
1666      m_data=temp_data;          return true;
1667        }
1668        return false;*/
1669      }
1670    
1671      void forceResolve()
1672      {
1673        if (isLazy())
1674        {
1675            #ifdef _OPENMP
1676            if (omp_in_parallel())
1677            {   // Yes this is throwing an exception out of an omp thread which is forbidden.
1678            throw DataException("Please do not call forceResolve() in a parallel region.");
1679            }
1680            #endif
1681            resolve();
1682        }
1683      }
1684    
1685      /**
1686      \brief if another object is sharing out member data make a copy to work with instead.
1687      This code should only be called from single threaded sections of code.
1688      */
1689      void exclusiveWrite()
1690      {
1691    #ifdef _OPENMP
1692      if (omp_in_parallel())
1693      {
1694    // *((int*)0)=17;
1695        throw DataException("Programming error. Please do not run exclusiveWrite() in multi-threaded sections.");
1696      }
1697    #endif
1698        forceResolve();
1699        if (isShared())
1700        {
1701            DataAbstract* t=m_data->deepCopy();
1702            set_m_data(DataAbstract_ptr(t));
1703        }
1704    }    }
1705    
1706      /**
1707      \brief checks if caller can have exclusive write to the object
1708      */
1709      void checkExclusiveWrite()
1710      {
1711        if  (isLazy() || isShared())
1712        {
1713            throw DataException("Programming error. ExclusiveWrite required - please call requireWrite()");
1714        }
1715      }
1716    
1717      /**
1718      \brief Modify the data abstract hosted by this Data object
1719      For internal use only.
1720      Passing a pointer to null is permitted (do this in the destructor)
1721      \warning Only to be called in single threaded code or inside a single/critical section. This method needs to be atomic.
1722      */
1723      void set_m_data(DataAbstract_ptr p);
1724    
1725      friend class DataAbstract;        // To allow calls to updateShareStatus
1726    
1727    };
1728    
1729    }   // end namespace escript
1730    
1731    
1732    // No, this is not supposed to be at the top of the file
1733    // DataAbstact needs to be declared first, then DataReady needs to be fully declared
1734    // so that I can dynamic cast between them below.
1735    #include "DataReady.h"
1736    #include "DataLazy.h"
1737    
1738    namespace escript
1739    {
1740    
1741    inline
1742    const DataReady*
1743    Data::getReady() const
1744    {
1745       const DataReady* dr=dynamic_cast<const DataReady*>(m_data.get());
1746       EsysAssert((dr!=0), "Error - casting to DataReady.");
1747       return dr;
1748    }
1749    
1750    inline
1751    DataReady*
1752    Data::getReady()
1753    {
1754       DataReady* dr=dynamic_cast<DataReady*>(m_data.get());
1755       EsysAssert((dr!=0), "Error - casting to DataReady.");
1756       return dr;
1757    }
1758    
1759    // Be wary of using this for local operations since it (temporarily) increases reference count.
1760    // If you are just using this to call a method on DataReady instead of DataAbstract consider using
1761    // getReady() instead
1762    inline
1763    DataReady_ptr
1764    Data::getReadyPtr()
1765    {
1766       DataReady_ptr dr=boost::dynamic_pointer_cast<DataReady>(m_data);
1767       EsysAssert((dr.get()!=0), "Error - casting to DataReady.");
1768       return dr;
1769    }
1770    
1771    
1772    inline
1773    const_DataReady_ptr
1774    Data::getReadyPtr() const
1775    {
1776       const_DataReady_ptr dr=boost::dynamic_pointer_cast<const DataReady>(m_data);
1777       EsysAssert((dr.get()!=0), "Error - casting to DataReady.");
1778       return dr;
1779  }  }
1780    
1781    inline
1782    DataAbstract::ValueType::value_type*
1783    Data::getSampleDataRW(DataAbstract::ValueType::size_type sampleNo)
1784    {
1785    //    if (isLazy())
1786    //    {
1787    //  resolve();
1788    //    }
1789    //    exclusiveWrite();
1790       if (isLazy())
1791       {
1792        throw DataException("Error, attempt to acquire RW access to lazy data. Please call requireWrite() first.");
1793       }
1794       return getReady()->getSampleDataRW(sampleNo);
1795    }
1796    
1797    // inline
1798    // const DataAbstract::ValueType::value_type*
1799    // Data::getSampleDataRO(DataAbstract::ValueType::size_type sampleNo, DataTypes::ValueType* buffer)
1800    // {
1801    //    DataLazy* l=dynamic_cast<DataLazy*>(m_data.get());
1802    //    if (l!=0)
1803    //    {
1804    //  size_t offset=0;
1805    //  if (buffer==NULL)
1806    //  {
1807    //      throw DataException("Error, attempt to getSampleDataRO for lazy Data with buffer==NULL");
1808    //  }
1809    //  const DataTypes::ValueType* res=l->resolveSample(*buffer,0,sampleNo,offset);
1810    //  return &((*res)[offset]);
1811    //    }
1812    //    return getReady()->getSampleDataRO(sampleNo);
1813    // }
1814    
1815    inline
1816    const DataAbstract::ValueType::value_type*
1817    Data::getSampleDataRO(DataAbstract::ValueType::size_type sampleNo, BufferGroup* bufferg)
1818    {
1819       DataLazy* l=dynamic_cast<DataLazy*>(m_data.get());
1820       if (l!=0)
1821       {
1822        size_t offset=0;
1823        if (bufferg==NULL)
1824        {
1825            throw DataException("Error, attempt to getSampleDataRO for lazy Data with buffer==NULL");
1826        }
1827        const DataTypes::ValueType* res=l->resolveSample(*bufferg,sampleNo,offset);
1828        return &((*res)[offset]);
1829       }
1830       return getReady()->getSampleDataRO(sampleNo);
1831    }
1832    
1833    
1834    
1835    /**
1836       Modify a filename for MPI parallel output to multiple files
1837    */
1838    char *Escript_MPI_appendRankToFileName(const char *, int, int);
1839    
1840  /**  /**
1841     Binary Data object operators.     Binary Data object operators.
1842  */  */
# Line 1519  ESCRIPT_DLL_API std::ostream& operator<< Line 1948  ESCRIPT_DLL_API std::ostream& operator<<
1948  /**  /**
1949    \brief    \brief
1950    Compute a tensor product of two Data objects    Compute a tensor product of two Data objects
1951    \param arg0 - Input - Data object    \param arg_0 - Input - Data object
1952    \param arg1 - Input - Data object    \param arg_1 - Input - Data object
1953    \param axis_offset - Input - axis offset    \param axis_offset - Input - axis offset
1954    \param transpose - Input - 0: transpose neither, 1: transpose arg0, 2: transpose arg1    \param transpose - Input - 0: transpose neither, 1: transpose arg0, 2: transpose arg1
1955  */  */
1956  ESCRIPT_DLL_API  ESCRIPT_DLL_API
1957  Data  Data
1958  C_GeneralTensorProduct(Data& arg0,  C_GeneralTensorProduct(Data& arg_0,
1959                       Data& arg1,                       Data& arg_1,
1960                       int axis_offset=0,                       int axis_offset=0,
1961                       int transpose=0);                       int transpose=0);
1962    
   
 /**  
   \brief  
   Compute a tensor operation with two Data objects  
   \param arg0 - Input - Data object  
   \param arg1 - Input - Data object  
   \param operation - Input - Binary op functor  
 */  
 template <typename BinaryFunction>  
 ESCRIPT_DLL_API  
 Data  
 C_TensorBinaryOperation(Data const &arg0,  
                         Data const &arg1,  
                         BinaryFunction operation);  
   
 /**  
   \brief  
   Return true if operands are equivalent, else return false.  
   NB: this operator does very little at this point, and isn't to  
   be relied on. Requires further implementation.  
 */  
 // ESCRIPT_DLL_API bool operator==(const Data& left, const Data& right);  
   
1963  /**  /**
1964    \brief    \brief
1965    Perform the given binary operation with this and right as operands.    Perform the given binary operation with this and right as operands.
# Line 1567  Data::binaryOp(const Data& right, Line 1973  Data::binaryOp(const Data& right,
1973  {  {
1974     //     //
1975     // 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
1976     if (getPointDataView().getRank()==0 && right.getPointDataView().getRank()!=0) {     if (getDataPointRank()==0 && right.getDataPointRank()!=0) {
1977       throw DataException("Error - attempt to update rank zero object with object with rank bigger than zero.");       throw DataException("Error - attempt to update rank zero object with object with rank bigger than zero.");
1978     }     }
1979    
1980       if (isLazy() || right.isLazy())
1981       {
1982         throw DataException("Programmer error - attempt to call binaryOp with Lazy Data.");
1983       }
1984     //     //
1985     // initially make the temporary a shallow copy     // initially make the temporary a shallow copy
1986     Data tempRight(right);     Data tempRight(right);
1987    
1988     if (getFunctionSpace()!=right.getFunctionSpace()) {     if (getFunctionSpace()!=right.getFunctionSpace()) {
1989       if (right.probeInterpolation(getFunctionSpace())) {       if (right.probeInterpolation(getFunctionSpace())) {
1990         //         //
# Line 1582  Data::binaryOp(const Data& right, Line 1994  Data::binaryOp(const Data& right,
1994         //         //
1995         // interpolate onto the RHS function space         // interpolate onto the RHS function space
1996         Data tempLeft(*this,right.getFunctionSpace());         Data tempLeft(*this,right.getFunctionSpace());
1997         m_data=tempLeft.m_data;  //        m_data=tempLeft.m_data;
1998           set_m_data(tempLeft.m_data);
1999       }       }
2000     }     }
2001     operandCheck(tempRight);     operandCheck(tempRight);
# Line 1598  Data::binaryOp(const Data& right, Line 2011  Data::binaryOp(const Data& right,
2011       // of any data type       // of any data type
2012       DataExpanded* leftC=dynamic_cast<DataExpanded*>(m_data.get());       DataExpanded* leftC=dynamic_cast<DataExpanded*>(m_data.get());
2013       EsysAssert((leftC!=0), "Programming error - casting to DataExpanded.");       EsysAssert((leftC!=0), "Programming error - casting to DataExpanded.");
2014       escript::binaryOp(*leftC,*(tempRight.m_data.get()),operation);       escript::binaryOp(*leftC,*(tempRight.getReady()),operation);
2015     } else if (isTagged()) {     } else if (isTagged()) {
2016       //       //
2017       // 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 1646  Data::algorithm(BinaryFunction operation Line 2059  Data::algorithm(BinaryFunction operation
2059      DataConstant* leftC=dynamic_cast<DataConstant*>(m_data.get());      DataConstant* leftC=dynamic_cast<DataConstant*>(m_data.get());
2060      EsysAssert((leftC!=0), "Programming error - casting to DataConstant.");      EsysAssert((leftC!=0), "Programming error - casting to DataConstant.");
2061      return escript::algorithm(*leftC,operation,initial_value);      return escript::algorithm(*leftC,operation,initial_value);
2062      } else if (isEmpty()) {
2063        throw DataException("Error - Operations not permitted on instances of DataEmpty.");
2064      } else if (isLazy()) {
2065        throw DataException("Error - Operations not permitted on instances of DataLazy.");
2066      } else {
2067        throw DataException("Error - Data encapsulates an unknown type.");
2068    }    }
   return 0;  
2069  }  }
2070    
2071  /**  /**
# Line 1663  inline Line 2081  inline
2081  Data  Data
2082  Data::dp_algorithm(BinaryFunction operation, double initial_value) const  Data::dp_algorithm(BinaryFunction operation, double initial_value) const
2083  {  {
2084    if (isExpanded()) {    if (isEmpty()) {
2085      Data result(0,DataArrayView::ShapeType(),getFunctionSpace(),isExpanded());      throw DataException("Error - Operations not permitted on instances of DataEmpty.");
2086      }
2087      else if (isExpanded()) {
2088        Data result(0,DataTypes::ShapeType(),getFunctionSpace(),isExpanded());
2089      DataExpanded* dataE=dynamic_cast<DataExpanded*>(m_data.get());      DataExpanded* dataE=dynamic_cast<DataExpanded*>(m_data.get());
2090      DataExpanded* resultE=dynamic_cast<DataExpanded*>(result.m_data.get());      DataExpanded* resultE=dynamic_cast<DataExpanded*>(result.m_data.get());
2091      EsysAssert((dataE!=0), "Programming error - casting data to DataExpanded.");      EsysAssert((dataE!=0), "Programming error - casting data to DataExpanded.");
2092      EsysAssert((resultE!=0), "Programming error - casting result to DataExpanded.");      EsysAssert((resultE!=0), "Programming error - casting result to DataExpanded.");
2093      escript::dp_algorithm(*dataE,*resultE,operation,initial_value);      escript::dp_algorithm(*dataE,*resultE,operation,initial_value);
2094      return result;      return result;
2095    } else if (isTagged()) {    }
2096      else if (isTagged()) {
2097      DataTagged* dataT=dynamic_cast<DataTagged*>(m_data.get());      DataTagged* dataT=dynamic_cast<DataTagged*>(m_data.get());
     DataArrayView::ShapeType viewShape;  
     DataArrayView::ValueType viewData(1);  
     viewData[0]=0;  
     DataArrayView defaultValue(viewData,viewShape);  
     DataTagged::TagListType keys;  
     DataTagged::ValueListType values;  
     DataTagged::DataMapType::const_iterator i;  
     for (i=dataT->getTagLookup().begin();i!=dataT->getTagLookup().end();i++) {  
       keys.push_back(i->first);  
       values.push_back(defaultValue);  
     }  
     Data result(keys,values,defaultValue,getFunctionSpace());  
     DataTagged* resultT=dynamic_cast<DataTagged*>(result.m_data.get());  
2098      EsysAssert((dataT!=0), "Programming error - casting data to DataTagged.");      EsysAssert((dataT!=0), "Programming error - casting data to DataTagged.");
2099      EsysAssert((resultT!=0), "Programming error - casting result to DataTagged.");      DataTypes::ValueType defval(1);
2100        defval[0]=0;
2101        DataTagged* resultT=new DataTagged(getFunctionSpace(), DataTypes::scalarShape, defval, dataT);
2102      escript::dp_algorithm(*dataT,*resultT,operation,initial_value);      escript::dp_algorithm(*dataT,*resultT,operation,initial_value);
2103      return result;      return Data(resultT);   // note: the Data object now owns the resultT pointer
2104    } else if (isConstant()) {    }
2105      Data result(0,DataArrayView::ShapeType(),getFunctionSpace(),isExpanded());    else if (isConstant()) {
2106        Data result(0,DataTypes::ShapeType(),getFunctionSpace(),isExpanded());
2107      DataConstant* dataC=dynamic_cast<DataConstant*>(m_data.get());      DataConstant* dataC=dynamic_cast<DataConstant*>(m_data.get());
2108      DataConstant* resultC=dynamic_cast<DataConstant*>(result.m_data.get());      DataConstant* resultC=dynamic_cast<DataConstant*>(result.m_data.get());
2109      EsysAssert((dataC!=0), "Programming error - casting data to DataConstant.");      EsysAssert((dataC!=0), "Programming error - casting data to DataConstant.");
2110      EsysAssert((resultC!=0), "Programming error - casting result to DataConstant.");      EsysAssert((resultC!=0), "Programming error - casting result to DataConstant.");
2111      escript::dp_algorithm(*dataC,*resultC,operation,initial_value);      escript::dp_algorithm(*dataC,*resultC,operation,initial_value);
2112      return result;      return result;
2113      } else if (isLazy()) {
2114        throw DataException("Error - Operations not permitted on instances of DataLazy.");
2115      } else {
2116        throw DataException("Error - Data encapsulates an unknown type.");
2117    }    }
   Data falseRetVal; // to keep compiler quiet  
   return falseRetVal;  
2118  }  }
2119    
2120    /**
2121      \brief
2122      Compute a tensor operation with two Data objects
2123      \param arg_0 - Input - Data object
2124      \param arg_1 - Input - Data object
2125      \param operation - Input - Binary op functor
2126    */
2127  template <typename BinaryFunction>  template <typename BinaryFunction>
2128    inline
2129  Data  Data
2130  C_TensorBinaryOperation(Data const &arg_0,  C_TensorBinaryOperation(Data const &arg_0,
2131                          Data const &arg_1,                          Data const &arg_1,
2132                          BinaryFunction operation)                          BinaryFunction operation)
2133  {  {
2134      if (arg_0.isEmpty() || arg_1.isEmpty())
2135      {
2136         throw DataException("Error - Operations not permitted on instances of DataEmpty.");
2137      }
2138      if (arg_0.isLazy() || arg_1.isLazy())
2139      {
2140         throw DataException("Error - Operations not permitted on lazy data.");
2141      }
2142    // Interpolate if necessary and find an appropriate function space    // Interpolate if necessary and find an appropriate function space
2143    Data arg_0_Z, arg_1_Z;    Data arg_0_Z, arg_1_Z;
2144    if (arg_0.getFunctionSpace()!=arg_1.getFunctionSpace()) {    if (arg_0.getFunctionSpace()!=arg_1.getFunctionSpace()) {
# Line 1730  C_TensorBinaryOperation(Data const &arg_ Line 2160  C_TensorBinaryOperation(Data const &arg_
2160    // Get rank and shape of inputs    // Get rank and shape of inputs
2161    int rank0 = arg_0_Z.getDataPointRank();    int rank0 = arg_0_Z.getDataPointRank();
2162    int rank1 = arg_1_Z.getDataPointRank();    int rank1 = arg_1_Z.getDataPointRank();
2163    DataArrayView::ShapeType shape0 = arg_0_Z.getDataPointShape();    DataTypes::ShapeType shape0 = arg_0_Z.getDataPointShape();
2164    DataArrayView::ShapeType shape1 = arg_1_Z.getDataPointShape();    DataTypes::ShapeType shape1 = arg_1_Z.getDataPointShape();
2165    int size0 = arg_0_Z.getDataPointSize();    int size0 = arg_0_Z.getDataPointSize();
2166    int size1 = arg_1_Z.getDataPointSize();    int size1 = arg_1_Z.getDataPointSize();
   
2167    // Declare output Data object    // Declare output Data object
2168    Data res;    Data res;
2169    
2170    if (shape0 == shape1) {    if (shape0 == shape1) {
   
2171      if (arg_0_Z.isConstant()   && arg_1_Z.isConstant()) {      if (arg_0_Z.isConstant()   && arg_1_Z.isConstant()) {
2172        res = Data(0.0, shape0, arg_1_Z.getFunctionSpace());      // DataConstant output        res = Data(0.0, shape0, arg_1_Z.getFunctionSpace());      // DataConstant output
2173        double *ptr_0 = &((arg_0_Z.getPointDataView().getData())[0]);        const double *ptr_0 = &(arg_0_Z.getDataAtOffsetRO(0));
2174        double *ptr_1 = &((arg_1_Z.getPointDataView().getData())[0]);        const double *ptr_1 = &(arg_1_Z.getDataAtOffsetRO(0));
2175        double *ptr_2 = &((res.getPointDataView().getData())[0]);        double *ptr_2 = &(res.getDataAtOffsetRW(0));
2176    
2177        tensor_binary_operation(size0, ptr_0, ptr_1, ptr_2, operation);        tensor_binary_operation(size0, ptr_0, ptr_1, ptr_2, operation);
2178      }      }
2179      else if (arg_0_Z.isConstant()   && arg_1_Z.isTagged()) {      else if (arg_0_Z.isConstant()   && arg_1_Z.isTagged()) {
# Line 1762  C_TensorBinaryOperation(Data const &arg_ Line 2191  C_TensorBinaryOperation(Data const &arg_
2191    
2192        // Prepare offset into DataConstant        // Prepare offset into DataConstant
2193        int offset_0 = tmp_0->getPointOffset(0,0);        int offset_0 = tmp_0->getPointOffset(0,0);
2194        double *ptr_0 = &((arg_0_Z.getPointDataView().getData())[offset_0]);        const double *ptr_0 = &(arg_0_Z.getDataAtOffsetRO(offset_0));
2195        // Get the views  
       DataArrayView view_1 = tmp_1->getDefaultValue();  
       DataArrayView view_2 = tmp_2->getDefaultValue();  
2196        // Get the pointers to the actual data        // Get the pointers to the actual data
2197        double *ptr_1 = &((view_1.getData())[0]);        const double *ptr_1 = &(tmp_1->getDefaultValueRO(0));
2198        double *ptr_2 = &((view_2.getData())[0]);        double *ptr_2 = &(tmp_2->getDefaultValueRW(0));
2199    
2200        // Compute a result for the default        // Compute a result for the default
2201        tensor_binary_operation(size0, ptr_0, ptr_1, ptr_2, operation);        tensor_binary_operation(size0, ptr_0, ptr_1, ptr_2, operation);
2202        // Compute a result for each tag        // Compute a result for each tag
2203        const DataTagged::DataMapType& lookup_1=tmp_1->getTagLookup();        const DataTagged::DataMapType& lookup_1=tmp_1->getTagLookup();
2204        DataTagged::DataMapType::const_iterator i; // i->first is a tag, i->second is an offset into memory        DataTagged::DataMapType::const_iterator i; // i->first is a tag, i->second is an offset into memory
2205        for (i=lookup_1.begin();i!=lookup_1.end();i++) {        for (i=lookup_1.begin();i!=lookup_1.end();i++) {
2206          tmp_2->addTaggedValue(i->first,tmp_2->getDefaultValue());          tmp_2->addTag(i->first);
2207          DataArrayView view_1 = tmp_1->getDataPointByTag(i->first);          const double *ptr_1 = &(tmp_1->getDataByTagRO(i->first,0));
2208          DataArrayView view_2 = tmp_2->getDataPointByTag(i->first);          double *ptr_2 = &(tmp_2->getDataByTagRW(i->first,0));
2209          double *ptr_1 = &view_1.getData(0);  
         double *ptr_2 = &view_2.getData(0);  
2210          tensor_binary_operation(size0, ptr_0, ptr_1, ptr_2, operation);          tensor_binary_operation(size0, ptr_0, ptr_1, ptr_2, operation);
2211        }        }
2212    
2213      }      }
2214      else if (arg_0_Z.isConstant()   && arg_1_Z.isExpanded()) {      else if (arg_0_Z.isConstant()   && arg_1_Z.isExpanded()) {
   
2215        res = Data(0.0, shape0, arg_1_Z.getFunctionSpace(),true); // DataExpanded output        res = Data(0.0, shape0, arg_1_Z.getFunctionSpace(),true); // DataExpanded output
2216        DataConstant* tmp_0=dynamic_cast<DataConstant*>(arg_0_Z.borrowData());        DataConstant* tmp_0=dynamic_cast<DataConstant*>(arg_0_Z.borrowData());
2217        DataExpanded* tmp_1=dynamic_cast<DataExpanded*>(arg_1_Z.borrowData());        DataExpanded* tmp_1=dynamic_cast<DataExpanded*>(arg_1_Z.borrowData());
# Line 1795  C_TensorBinaryOperation(Data const &arg_ Line 2221  C_TensorBinaryOperation(Data const &arg_
2221        int numSamples_1 = arg_1_Z.getNumSamples();        int numSamples_1 = arg_1_Z.getNumSamples();
2222        int numDataPointsPerSample_1 = arg_1_Z.getNumDataPointsPerSample();        int numDataPointsPerSample_1 = arg_1_Z.getNumDataPointsPerSample();
2223        int offset_0 = tmp_0->getPointOffset(0,0);        int offset_0 = tmp_0->getPointOffset(0,0);
2224          res.requireWrite();
2225        #pragma omp parallel for private(sampleNo_1,dataPointNo_1) schedule(static)        #pragma omp parallel for private(sampleNo_1,dataPointNo_1) schedule(static)
2226        for (sampleNo_1 = 0; sampleNo_1 < numSamples_1; sampleNo_1++) {        for (sampleNo_1 = 0; sampleNo_1 < numSamples_1; sampleNo_1++) {
2227          for (dataPointNo_1 = 0; dataPointNo_1 < numDataPointsPerSample_1; dataPointNo_1++) {          for (dataPointNo_1 = 0; dataPointNo_1 < numDataPointsPerSample_1; dataPointNo_1++) {
2228            int offset_1 = tmp_1->getPointOffset(sampleNo_1,dataPointNo_1);            int offset_1 = tmp_1->getPointOffset(sampleNo_1,dataPointNo_1);
2229            int offset_2 = tmp_2->getPointOffset(sampleNo_1,dataPointNo_1);            int offset_2 = tmp_2->getPointOffset(sampleNo_1,dataPointNo_1);
2230            double *ptr_0 = &((arg_0_Z.getPointDataView().getData())[offset_0]);            const double *ptr_0 = &(arg_0_Z.getDataAtOffsetRO(offset_0));
2231            double *ptr_1 = &((arg_1_Z.getPointDataView().getData())[offset_1]);            const double *ptr_1 = &(arg_1_Z.getDataAtOffsetRO(offset_1));
2232            double *ptr_2 = &((res.getPointDataView().getData())[offset_2]);            double *ptr_2 = &(res.getDataAtOffsetRW(offset_2));
2233            tensor_binary_operation(size0, ptr_0, ptr_1, ptr_2, operation);            tensor_binary_operation(size0, ptr_0, ptr_1, ptr_2, operation);
2234          }          }
2235        }        }
2236    
2237      }      }
2238      else if (arg_0_Z.isTagged()     && arg_1_Z.isConstant()) {      else if (arg_0_Z.isTagged()     && arg_1_Z.isConstant()) {
   
2239        // Borrow DataTagged input from Data object        // Borrow DataTagged input from Data object
2240        DataTagged* tmp_0=dynamic_cast<DataTagged*>(arg_0_Z.borrowData());        DataTagged* tmp_0=dynamic_cast<DataTagged*>(arg_0_Z.borrowData());
2241    
# Line 1823  C_TensorBinaryOperation(Data const &arg_ Line 2249  C_TensorBinaryOperation(Data const &arg_
2249    
2250        // Prepare offset into DataConstant        // Prepare offset into DataConstant
2251        int offset_1 = tmp_1->getPointOffset(0,0);        int offset_1 = tmp_1->getPointOffset(0,0);
2252        double *ptr_1 = &((arg_1_Z.getPointDataView().getData())[offset_1]);  
2253        // Get the views        const double *ptr_1 = &(arg_1_Z.getDataAtOffsetRO(offset_1));
       DataArrayView view_0 = tmp_0->getDefaultValue();  
       DataArrayView view_2 = tmp_2->getDefaultValue();  
2254        // Get the pointers to the actual data        // Get the pointers to the actual data
2255        double *ptr_0 = &((view_0.getData())[0]);        const double *ptr_0 = &(tmp_0->getDefaultValueRO(0));
2256        double *ptr_2 = &((view_2.getData())[0]);        double *ptr_2 = &(tmp_2->getDefaultValueRW(0));
2257        // Compute a result for the default        // Compute a result for the default
2258        tensor_binary_operation(size0, ptr_0, ptr_1, ptr_2, operation);        tensor_binary_operation(size0, ptr_0, ptr_1, ptr_2, operation);
2259        // Compute a result for each tag        // Compute a result for each tag
2260        const DataTagged::DataMapType& lookup_0=tmp_0->getTagLookup();        const DataTagged::DataMapType& lookup_0=tmp_0->getTagLookup();
2261        DataTagged::DataMapType::const_iterator i; // i->first is a tag, i->second is an offset into memory        DataTagged::DataMapType::const_iterator i; // i->first is a tag, i->second is an offset into memory
2262        for (i=lookup_0.begin();i!=lookup_0.end();i++) {        for (i=lookup_0.begin();i!=lookup_0.end();i++) {
2263          tmp_2->addTaggedValue(i->first,tmp_2->getDefaultValue());          tmp_2->addTag(i->first);
2264          DataArrayView view_0 = tmp_0->getDataPointByTag(i->first);          const double *ptr_0 = &(tmp_0->getDataByTagRO(i->first,0));
2265          DataArrayView view_2 = tmp_2->getDataPointByTag(i->first);          double *ptr_2 = &(tmp_2->getDataByTagRW(i->first,0));
         double *ptr_0 = &view_0.getData(0);  
         double *ptr_2 = &view_2.getData(0);  
2266          tensor_binary_operation(size0, ptr_0, ptr_1, ptr_2, operation);          tensor_binary_operation(size0, ptr_0, ptr_1, ptr_2, operation);
2267        }        }
2268    
2269      }      }
2270      else if (arg_0_Z.isTagged()     && arg_1_Z.isTagged()) {      else if (arg_0_Z.isTagged()     && arg_1_Z.isTagged()) {
   
2271        // Borrow DataTagged input from Data object        // Borrow DataTagged input from Data object
2272        DataTagged* tmp_0=dynamic_cast<DataTagged*>(arg_0_Z.borrowData());        DataTagged* tmp_0=dynamic_cast<DataTagged*>(arg_0_Z.borrowData());
2273    
# Line 1858  C_TensorBinaryOperation(Data const &arg_ Line 2279  C_TensorBinaryOperation(Data const &arg_
2279        res.tag();        // DataTagged output        res.tag();        // DataTagged output
2280        DataTagged* tmp_2=dynamic_cast<DataTagged*>(res.borrowData());        DataTagged* tmp_2=dynamic_cast<DataTagged*>(res.borrowData());
2281    
       // Get the views  
       DataArrayView view_0 = tmp_0->getDefaultValue();  
       DataArrayView view_1 = tmp_1->getDefaultValue();  
       DataArrayView view_2 = tmp_2->getDefaultValue();  
2282        // Get the pointers to the actual data        // Get the pointers to the actual data
2283        double *ptr_0 = &((view_0.getData())[0]);        const double *ptr_0 = &(tmp_0->getDefaultValueRO(0));
2284        double *ptr_1 = &((view_1.getData())[0]);        const double *ptr_1 = &(tmp_1->getDefaultValueRO(0));
2285        double *ptr_2 = &((view_2.getData())[0]);        double *ptr_2 = &(tmp_2->getDefaultValueRW(0));
2286    
2287        // Compute a result for the default        // Compute a result for the default
2288        tensor_binary_operation(size0, ptr_0, ptr_1, ptr_2, operation);        tensor_binary_operation(size0, ptr_0, ptr_1, ptr_2, operation);
2289        // Merge the tags        // Merge the tags
# Line 1873  C_TensorBinaryOperation(Data const &arg_ Line 2291  C_TensorBinaryOperation(Data const &arg_
2291        const DataTagged::DataMapType& lookup_0=tmp_0->getTagLookup();        const DataTagged::DataMapType& lookup_0=tmp_0->getTagLookup();
2292        const DataTagged::DataMapType& lookup_1=tmp_1->getTagLookup();        const DataTagged::DataMapType& lookup_1=tmp_1->getTagLookup();
2293        for (i=lookup_0.begin();i!=lookup_0.end();i++) {        for (i=lookup_0.begin();i!=lookup_0.end();i++) {
2294          tmp_2->addTaggedValue(i->first,tmp_2->getDefaultValue()); // use tmp_2 to get correct shape          tmp_2->addTag(i->first); // use tmp_2 to get correct shape
2295        }        }
2296        for (i=lookup_1.begin();i!=lookup_1.end();i++) {        for (i=lookup_1.begin();i!=lookup_1.end();i++) {
2297          tmp_2->addTaggedValue(i->first,tmp_2->getDefaultValue());          tmp_2->addTag(i->first);
2298        }        }
2299        // Compute a result for each tag        // Compute a result for each tag
2300        const DataTagged::DataMapType& lookup_2=tmp_2->getTagLookup();        const DataTagged::DataMapType& lookup_2=tmp_2->getTagLookup();
2301        for (i=lookup_2.begin();i!=lookup_2.end();i++) {        for (i=lookup_2.begin();i!=lookup_2.end();i++) {
2302          DataArrayView view_0 = tmp_0->getDataPointByTag(i->first);  
2303          DataArrayView view_1 = tmp_1->getDataPointByTag(i->first);          const double *ptr_0 = &(tmp_0->getDataByTagRO(i->first,0));
2304          DataArrayView view_2 = tmp_2->getDataPointByTag(i->first);          const double *ptr_1 = &(tmp_1->getDataByTagRO(i->first,0));
2305          double *ptr_0 = &view_0.getData(0);          double *ptr_2 = &(tmp_2->getDataByTagRW(i->first,0));
2306          double *ptr_1 = &view_1.getData(0);  
         double *ptr_2 = &view_2.getData(0);  
2307          tensor_binary_operation(size0, ptr_0, ptr_1, ptr_2, operation);          tensor_binary_operation(size0, ptr_0, ptr_1, ptr_2, operation);
2308        }        }
2309    
2310      }      }
2311      else if (arg_0_Z.isTagged()     && arg_1_Z.isExpanded()) {      else if (arg_0_Z.isTagged()     && arg_1_Z.isExpanded()) {
   
2312        // After finding a common function space above the two inputs have the same numSamples and num DPPS        // After finding a common function space above the two inputs have the same numSamples and num DPPS
2313        res = Data(0.0, shape0, arg_1_Z.getFunctionSpace(),true); // DataExpanded output        res = Data(0.0, shape0, arg_1_Z.getFunctionSpace(),true); // DataExpanded output
2314        DataTagged*   tmp_0=dynamic_cast<DataTagged*>(arg_0_Z.borrowData());        DataTagged*   tmp_0=dynamic_cast<DataTagged*>(arg_0_Z.borrowData());
# Line 1902  C_TensorBinaryOperation(Data const &arg_ Line 2318  C_TensorBinaryOperation(Data const &arg_
2318        int sampleNo_0,dataPointNo_0;        int sampleNo_0,dataPointNo_0;
2319        int numSamples_0 = arg_0_Z.getNumSamples();        int numSamples_0 = arg_0_Z.getNumSamples();
2320        int numDataPointsPerSample_0 = arg_0_Z.getNumDataPointsPerSample();        int numDataPointsPerSample_0 = arg_0_Z.getNumDataPointsPerSample();
2321          res.requireWrite();
2322        #pragma omp parallel for private(sampleNo_0,dataPointNo_0) schedule(static)        #pragma omp parallel for private(sampleNo_0,dataPointNo_0) schedule(static)
2323        for (sampleNo_0 = 0; sampleNo_0 < numSamples_0; sampleNo_0++) {        for (sampleNo_0 = 0; sampleNo_0 < numSamples_0; sampleNo_0++) {
2324          int offset_0 = tmp_0->getPointOffset(sampleNo_0,0); // They're all the same, so just use #0          int offset_0 = tmp_0->getPointOffset(sampleNo_0,0); // They're all the same, so just use #0
2325          double *ptr_0 = &((arg_0_Z.getPointDataView().getData())[offset_0]);          const double *ptr_0 = &(arg_0_Z.getDataAtOffsetRO(offset_0));
2326          for (dataPointNo_0 = 0; dataPointNo_0 < numDataPointsPerSample_0; dataPointNo_0++) {          for (dataPointNo_0 = 0; dataPointNo_0 < numDataPointsPerSample_0; dataPointNo_0++) {
2327            int offset_1 = tmp_1->getPointOffset(sampleNo_0,dataPointNo_0);            int offset_1 = tmp_1->getPointOffset(sampleNo_0,dataPointNo_0);
2328            int offset_2 = tmp_2->getPointOffset(sampleNo_0,dataPointNo_0);            int offset_2 = tmp_2->getPointOffset(sampleNo_0,dataPointNo_0);
2329            double *ptr_1 = &((arg_1_Z.getPointDataView().getData())[offset_1]);            const double *ptr_1 = &(arg_1_Z.getDataAtOffsetRO(offset_1));
2330            double *ptr_2 = &((res.getPointDataView().getData())[offset_2]);            double *ptr_2 = &(res.getDataAtOffsetRW(offset_2));
2331            tensor_binary_operation(size0, ptr_0, ptr_1, ptr_2, operation);            tensor_binary_operation(size0, ptr_0, ptr_1, ptr_2, operation);
2332          }          }
2333        }        }
2334    
2335      }      }
2336      else if (arg_0_Z.isExpanded()   && arg_1_Z.isConstant()) {      else if (arg_0_Z.isExpanded()   && arg_1_Z.isConstant()) {
   
2337        res = Data(0.0, shape0, arg_1_Z.getFunctionSpace(),true); // DataExpanded output        res = Data(0.0, shape0, arg_1_Z.getFunctionSpace(),true); // DataExpanded output
2338        DataExpanded* tmp_0=dynamic_cast<DataExpanded*>(arg_0_Z.borrowData());        DataExpanded* tmp_0=dynamic_cast<DataExpanded*>(arg_0_Z.borrowData());
2339        DataConstant* tmp_1=dynamic_cast<DataConstant*>(arg_1_Z.borrowData());        DataConstant* tmp_1=dynamic_cast<DataConstant*>(arg_1_Z.borrowData());
# Line 1927  C_TensorBinaryOperation(Data const &arg_ Line 2343  C_TensorBinaryOperation(Data const &arg_
2343        int numSamples_0 = arg_0_Z.getNumSamples();        int numSamples_0 = arg_0_Z.getNumSamples();
2344        int numDataPointsPerSample_0 = arg_0_Z.getNumDataPointsPerSample();        int numDataPointsPerSample_0 = arg_0_Z.getNumDataPointsPerSample();
2345        int offset_1 = tmp_1->getPointOffset(0,0);        int offset_1 = tmp_1->getPointOffset(0,0);
2346          res.requireWrite();
2347        #pragma omp parallel for private(sampleNo_0,dataPointNo_0) schedule(static)        #pragma omp parallel for private(sampleNo_0,dataPointNo_0) schedule(static)
2348        for (sampleNo_0 = 0; sampleNo_0 < numSamples_0; sampleNo_0++) {        for (sampleNo_0 = 0; sampleNo_0 < numSamples_0; sampleNo_0++) {
2349          for (dataPointNo_0 = 0; dataPointNo_0 < numDataPointsPerSample_0; dataPointNo_0++) {          for (dataPointNo_0 = 0; dataPointNo_0 < numDataPointsPerSample_0; dataPointNo_0++) {
2350            int offset_0 = tmp_0->getPointOffset(sampleNo_0,dataPointNo_0);            int offset_0 = tmp_0->getPointOffset(sampleNo_0,dataPointNo_0);
2351            int offset_2 = tmp_2->getPointOffset(sampleNo_0,dataPointNo_0);            int offset_2 = tmp_2->getPointOffset(sampleNo_0,dataPointNo_0);
2352            double *ptr_0 = &((arg_0_Z.getPointDataView().getData())[offset_0]);  
2353            double *ptr_1 = &((arg_1_Z.getPointDataView().getData())[offset_1]);            const double *ptr_0 = &(arg_0_Z.getDataAtOffsetRO(offset_0));
2354            double *ptr_2 = &((res.getPointDataView().getData())[offset_2]);            const double *ptr_1 = &(arg_1_Z.getDataAtOffsetRO(offset_1));
2355              double *ptr_2 = &(res.getDataAtOffsetRW(offset_2));
2356    
2357    
2358            tensor_binary_operation(size0, ptr_0, ptr_1, ptr_2, operation);            tensor_binary_operation(size0, ptr_0, ptr_1, ptr_2, operation);
2359          }          }
2360        }        }
2361    
2362      }      }
2363      else if (arg_0_Z.isExpanded()   && arg_1_Z.isTagged()) {      else if (arg_0_Z.isExpanded()   && arg_1_Z.isTagged()) {
   
2364        // After finding a common function space above the two inputs have the same numSamples and num DPPS        // After finding a common function space above the two inputs have the same numSamples and num DPPS
2365        res = Data(0.0, shape0, arg_1_Z.getFunctionSpace(),true); // DataExpanded output        res = Data(0.0, shape0, arg_1_Z.getFunctionSpace(),true); // DataExpanded output
2366        DataExpanded* tmp_0=dynamic_cast<DataExpanded*>(arg_0_Z.borrowData());        DataExpanded* tmp_0=dynamic_cast<DataExpanded*>(arg_0_Z.borrowData());
# Line 1951  C_TensorBinaryOperation(Data const &arg_ Line 2370  C_TensorBinaryOperation(Data const &arg_
2370        int sampleNo_0,dataPointNo_0;        int sampleNo_0,dataPointNo_0;
2371        int numSamples_0 = arg_0_Z.getNumSamples();        int numSamples_0 = arg_0_Z.getNumSamples();
2372        int numDataPointsPerSample_0 = arg_0_Z.getNumDataPointsPerSample();        int numDataPointsPerSample_0 = arg_0_Z.getNumDataPointsPerSample();
2373          res.requireWrite();
2374        #pragma omp parallel for private(sampleNo_0,dataPointNo_0) schedule(static)        #pragma omp parallel for private(sampleNo_0,dataPointNo_0) schedule(static)
2375        for (sampleNo_0 = 0; sampleNo_0 < numSamples_0; sampleNo_0++) {        for (sampleNo_0 = 0; sampleNo_0 < numSamples_0; sampleNo_0++) {
2376          int offset_1 = tmp_1->getPointOffset(sampleNo_0,0);          int offset_1 = tmp_1->getPointOffset(sampleNo_0,0);
2377          double *ptr_1 = &((arg_1_Z.getPointDataView().getData())[offset_1]);          const double *ptr_1 = &(arg_1_Z.getDataAtOffsetRO(offset_1));
2378          for (dataPointNo_0 = 0; dataPointNo_0 < numDataPointsPerSample_0; dataPointNo_0++) {          for (dataPointNo_0 = 0; dataPointNo_0 < numDataPointsPerSample_0; dataPointNo_0++) {
2379            int offset_0 = tmp_0->getPointOffset(sampleNo_0,dataPointNo_0);            int offset_0 = tmp_0->getPointOffset(sampleNo_0,dataPointNo_0);
2380            int offset_2 = tmp_2->getPointOffset(sampleNo_0,dataPointNo_0);            int offset_2 = tmp_2->getPointOffset(sampleNo_0,dataPointNo_0);
2381            double *ptr_0 = &((arg_0_Z.getPointDataView().getData())[offset_0]);            const double *ptr_0 = &(arg_0_Z.getDataAtOffsetRO(offset_0));
2382            double *ptr_2 = &((res.getPointDataView().getData())[offset_2]);            double *ptr_2 = &(res.getDataAtOffsetRW(offset_2));
2383            tensor_binary_operation(size0, ptr_0, ptr_1, ptr_2, operation);            tensor_binary_operation(size0, ptr_0, ptr_1, ptr_2, operation);
2384          }          }
2385        }        }
2386    
2387      }      }
2388      else if (arg_0_Z.isExpanded()   && arg_1_Z.isExpanded()) {      else if (arg_0_Z.isExpanded()   && arg_1_Z.isExpanded()) {
   
2389        // After finding a common function space above the two inputs have the same numSamples and num DPPS        // After finding a common function space above the two inputs have the same numSamples and num DPPS
2390        res = Data(0.0, shape0, arg_1_Z.getFunctionSpace(),true); // DataExpanded output        res = Data(0.0, shape0, arg_1_Z.getFunctionSpace(),true); // DataExpanded output
2391        DataExpanded* tmp_0=dynamic_cast<DataExpanded*>(arg_0_Z.borrowData());        DataExpanded* tmp_0=dynamic_cast<DataExpanded*>(arg_0_Z.borrowData());
# Line 1976  C_TensorBinaryOperation(Data const &arg_ Line 2395  C_TensorBinaryOperation(Data const &arg_
2395        int sampleNo_0,dataPointNo_0;        int sampleNo_0,dataPointNo_0;
2396        int numSamples_0 = arg_0_Z.getNumSamples();        int numSamples_0 = arg_0_Z.getNumSamples();
2397        int numDataPointsPerSample_0 = arg_0_Z.getNumDataPointsPerSample();        int numDataPointsPerSample_0 = arg_0_Z.getNumDataPointsPerSample();
2398          res.requireWrite();
2399        #pragma omp parallel for private(sampleNo_0,dataPointNo_0) schedule(static)        #pragma omp parallel for private(sampleNo_0,dataPointNo_0) schedule(static)
2400        for (sampleNo_0 = 0; sampleNo_0 < numSamples_0; sampleNo_0++) {        for (sampleNo_0 = 0; sampleNo_0 < numSamples_0; sampleNo_0++) {
2401          for (dataPointNo_0 = 0; dataPointNo_0 < numDataPointsPerSample_0; dataPointNo_0++) {          for (dataPointNo_0 = 0; dataPointNo_0 < numDataPointsPerSample_0; dataPointNo_0++) {
2402            int offset_0 = tmp_0->getPointOffset(sampleNo_0,dataPointNo_0);            int offset_0 = tmp_0->getPointOffset(sampleNo_0,dataPointNo_0);
2403            int offset_1 = tmp_1->getPointOffset(sampleNo_0,dataPointNo_0);            int offset_1 = tmp_1->getPointOffset(sampleNo_0,dataPointNo_0);
2404            int offset_2 = tmp_2->getPointOffset(sampleNo_0,dataPointNo_0);            int offset_2 = tmp_2->getPointOffset(sampleNo_0,dataPointNo_0);
2405            double *ptr_0 = &((arg_0_Z.getPointDataView().getData())[offset_0]);            const double *ptr_0 = &(arg_0_Z.getDataAtOffsetRO(offset_0));
2406            double *ptr_1 = &((arg_1_Z.getPointDataView().getData())[offset_1]);            const double *ptr_1 = &(arg_1_Z.getDataAtOffsetRO(offset_1));
2407            double *ptr_2 = &((res.getPointDataView().getData())[offset_2]);            double *ptr_2 = &(res.getDataAtOffsetRW(offset_2));
2408            tensor_binary_operation(size0, ptr_0, ptr_1, ptr_2, operation);            tensor_binary_operation(size0, ptr_0, ptr_1, ptr_2, operation);
2409          }          }
2410        }        }
# Line 1995  C_TensorBinaryOperation(Data const &arg_ Line 2415  C_TensorBinaryOperation(Data const &arg_
2415      }      }
2416    
2417    } else if (0 == rank0) {    } else if (0 == rank0) {
   
2418      if (arg_0_Z.isConstant()   && arg_1_Z.isConstant()) {      if (arg_0_Z.isConstant()   && arg_1_Z.isConstant()) {
2419        res = Data(0.0, shape1, arg_1_Z.getFunctionSpace());      // DataConstant output        res = Data(0.0, shape1, arg_1_Z.getFunctionSpace());      // DataConstant output
2420        double *ptr_0 = &((arg_0_Z.getPointDataView().getData())[0]);        const double *ptr_0 = &(arg_0_Z.getDataAtOffsetRO(0));
2421        double *ptr_1 = &((arg_1_Z.getPointDataView().getData())[0]);        const double *ptr_1 = &(arg_1_Z.getDataAtOffsetRO(0));
2422        double *ptr_2 = &((res.getPointDataView().getData())[0]);        double *ptr_2 = &(res.getDataAtOffsetRW(0));
2423        tensor_binary_operation(size1, ptr_0[0], ptr_1, ptr_2, operation);        tensor_binary_operation(size1, ptr_0[0], ptr_1, ptr_2, operation);
2424      }      }
2425      else if (arg_0_Z.isConstant()   && arg_1_Z.isTagged()) {      else if (arg_0_Z.isConstant()   && arg_1_Z.isTagged()) {
# Line 2018  C_TensorBinaryOperation(Data const &arg_ Line 2437  C_TensorBinaryOperation(Data const &arg_
2437    
2438        // Prepare offset into DataConstant        // Prepare offset into DataConstant
2439        int offset_0 = tmp_0->getPointOffset(0,0);        int offset_0 = tmp_0->getPointOffset(0,0);
2440        double *ptr_0 = &((arg_0_Z.getPointDataView().getData())[offset_0]);        const double *ptr_0 = &(arg_0_Z.getDataAtOffsetRO(offset_0));
2441        // Get the views  
2442        DataArrayView view_1 = tmp_1->getDefaultValue();        const double *ptr_1 = &(tmp_1->getDefaultValueRO(0));
2443        DataArrayView view_2 = tmp_2->getDefaultValue();        double *ptr_2 = &(tmp_2->getDefaultValueRW(0));
2444        // Get the pointers to the actual data  
       double *ptr_1 = &((view_1.getData())[0]);  
       double *ptr_2 = &((view_2.getData())[0]);  
2445        // Compute a result for the default        // Compute a result for the default
2446        tensor_binary_operation(size1, ptr_0[0], ptr_1, ptr_2, operation);        tensor_binary_operation(size1, ptr_0[0], ptr_1, ptr_2, operation);
2447        // Compute a result for each tag        // Compute a result for each tag
2448        const DataTagged::DataMapType& lookup_1=tmp_1->getTagLookup();        const DataTagged::DataMapType& lookup_1=tmp_1->getTagLookup();
2449        DataTagged::DataMapType::const_iterator i; // i->first is a tag, i->second is an offset into memory        DataTagged::DataMapType::const_iterator i; // i->first is a tag, i->second is an offset into memory
2450        for (i=lookup_1.begin();i!=lookup_1.end();i++) {        for (i=lookup_1.begin();i!=lookup_1.end();i++) {
2451          tmp_2->addTaggedValue(i->first,tmp_2->getDefaultValue());          tmp_2->addTag(i->first);
2452          DataArrayView view_1 = tmp_1->getDataPointByTag(i->first);          const double *ptr_1 = &(tmp_1->getDataByTagRO(i->first,0));
2453          DataArrayView view_2 = tmp_2->getDataPointByTag(i->first);          double *ptr_2 = &(tmp_2->getDataByTagRW(i->first,0));
         double *ptr_1 = &view_1.getData(0);  
         double *ptr_2 = &view_2.getData(0);  
2454          tensor_binary_operation(size1, ptr_0[0], ptr_1, ptr_2, operation);          tensor_binary_operation(size1, ptr_0[0], ptr_1, ptr_2, operation);
2455        }        }
2456    
# Line 2051  C_TensorBinaryOperation(Data const &arg_ Line 2466  C_TensorBinaryOperation(Data const &arg_
2466        int numSamples_1 = arg_1_Z.getNumSamples();        int numSamples_1 = arg_1_Z.getNumSamples();
2467        int numDataPointsPerSample_1 = arg_1_Z.getNumDataPointsPerSample();        int numDataPointsPerSample_1 = arg_1_Z.getNumDataPointsPerSample();
2468        int offset_0 = tmp_0->getPointOffset(0,0);        int offset_0 = tmp_0->getPointOffset(0,0);
2469          res.requireWrite();
2470        #pragma omp parallel for private(sampleNo_1,dataPointNo_1) schedule(static)        #pragma omp parallel for private(sampleNo_1,dataPointNo_1) schedule(static)
2471        for (sampleNo_1 = 0; sampleNo_1 < numSamples_1; sampleNo_1++) {        for (sampleNo_1 = 0; sampleNo_1 < numSamples_1; sampleNo_1++) {
2472          for (dataPointNo_1 = 0; dataPointNo_1 < numDataPointsPerSample_1; dataPointNo_1++) {          for (dataPointNo_1 = 0; dataPointNo_1 < numDataPointsPerSample_1; dataPointNo_1++) {
2473            int offset_1 = tmp_1->getPointOffset(sampleNo_1,dataPointNo_1);            int offset_1 = tmp_1->getPointOffset(sampleNo_1,dataPointNo_1);
2474            int offset_2 = tmp_2->getPointOffset(sampleNo_1,dataPointNo_1);            int offset_2 = tmp_2->getPointOffset(sampleNo_1,dataPointNo_1);
2475            double *ptr_0 = &((arg_0_Z.getPointDataView().getData())[offset_0]);            const double *ptr_0 = &(arg_0_Z.getDataAtOffsetRO(offset_0));
2476            double *ptr_1 = &((arg_1_Z.getPointDataView().getData())[offset_1]);            const double *ptr_1 = &(arg_1_Z.getDataAtOffsetRO(offset_1));
2477            double *ptr_2 = &((res.getPointDataView().getData())[offset_2]);            double *ptr_2 = &(res.getDataAtOffsetRW(offset_2));
2478            tensor_binary_operation(size1, ptr_0[0], ptr_1, ptr_2, operation);            tensor_binary_operation(size1, ptr_0[0], ptr_1, ptr_2, operation);
2479    
2480          }          }
# Line 2080  C_TensorBinaryOperation(Data const &arg_ Line 2496  C_TensorBinaryOperation(Data const &arg_
2496    
2497        // Prepare offset into DataConstant        // Prepare offset into DataConstant
2498        int offset_1 = tmp_1->getPointOffset(0,0);        int offset_1 = tmp_1->getPointOffset(0,0);
2499        double *ptr_1 = &((arg_1_Z.getPointDataView().getData())[offset_1]);        const double *ptr_1 = &(arg_1_Z.getDataAtOffsetRO(offset_1));
2500        // Get the views  
       DataArrayView view_0 = tmp_0->getDefaultValue();  
       DataArrayView view_2 = tmp_2->getDefaultValue();  
2501        // Get the pointers to the actual data        // Get the pointers to the actual data
2502        double *ptr_0 = &((view_0.getData())[0]);        const double *ptr_0 = &(tmp_0->getDefaultValueRO(0));
2503        double *ptr_2 = &((view_2.getData())[0]);        double *ptr_2 = &(tmp_2->getDefaultValueRW(0));
2504    
2505    
2506        // Compute a result for the default        // Compute a result for the default
2507        tensor_binary_operation(size1, ptr_0[0], ptr_1, ptr_2, operation);        tensor_binary_operation(size1, ptr_0[0], ptr_1, ptr_2, operation);
2508        // Compute a result for each tag        // Compute a result for each tag
2509        const DataTagged::DataMapType& lookup_0=tmp_0->getTagLookup();        const DataTagged::DataMapType& lookup_0=tmp_0->getTagLookup();
2510        DataTagged::DataMapType::const_iterator i; // i->first is a tag, i->second is an offset into memory        DataTagged::DataMapType::const_iterator i; // i->first is a tag, i->second is an offset into memory
2511        for (i=lookup_0.begin();i!=lookup_0.end();i++) {        for (i=lookup_0.begin();i!=lookup_0.end();i++) {
2512          tmp_2->addTaggedValue(i->first,tmp_2->getDefaultValue());          tmp_2->addTag(i->first);
2513          DataArrayView view_0 = tmp_0->getDataPointByTag(i->first);          const double *ptr_0 = &(tmp_0->getDataByTagRO(i->first,0));
2514          DataArrayView view_2 = tmp_2->getDataPointByTag(i->first);          double *ptr_2 = &(tmp_2->getDataByTagRW(i->first,0));
2515          double *ptr_0 = &view_0.getData(0);  
         double *ptr_2 = &view_2.getData(0);  
2516          tensor_binary_operation(size1, ptr_0[0], ptr_1, ptr_2, operation);          tensor_binary_operation(size1, ptr_0[0], ptr_1, ptr_2, operation);
2517        }        }
2518    
# Line 2115  C_TensorBinaryOperation(Data const &arg_ Line 2530  C_TensorBinaryOperation(Data const &arg_
2530        res.tag();        // DataTagged output        res.tag();        // DataTagged output
2531        DataTagged* tmp_2=dynamic_cast<DataTagged*>(res.borrowData());        DataTagged* tmp_2=dynamic_cast<DataTagged*>(res.borrowData());
2532    
       // Get the views  
       DataArrayView view_0 = tmp_0->getDefaultValue();  
       DataArrayView view_1 = tmp_1->getDefaultValue();  
       DataArrayView view_2 = tmp_2->getDefaultValue();  
2533        // Get the pointers to the actual data        // Get the pointers to the actual data
2534        double *ptr_0 = &((view_0.getData())[0]);        const double *ptr_0 = &(tmp_0->getDefaultValueRO(0));
2535        double *ptr_1 = &((view_1.getData())[0]);        const double *ptr_1 = &(tmp_1->getDefaultValueRO(0));
2536        double *ptr_2 = &((view_2.getData())[0]);        double *ptr_2 = &(tmp_2->getDefaultValueRW(0));
2537    
2538        // Compute a result for the default        // Compute a result for the default
2539        tensor_binary_operation(size1, ptr_0[0], ptr_1, ptr_2, operation);        tensor_binary_operation(size1, ptr_0[0], ptr_1, ptr_2, operation);
2540        // Merge the tags        // Merge the tags
# Line 2130  C_TensorBinaryOperation(Data const &arg_ Line 2542  C_TensorBinaryOperation(Data const &arg_
2542        const DataTagged::DataMapType& lookup_0=tmp_0->getTagLookup();        const DataTagged::DataMapType& lookup_0=tmp_0->getTagLookup();
2543        const DataTagged::DataMapType& lookup_1=tmp_1->getTagLookup();        const DataTagged::DataMapType& lookup_1=tmp_1->getTagLookup();
2544        for (i=lookup_0.begin();i!=lookup_0.end();i++) {        for (i=lookup_0.begin();i!=lookup_0.end();i++) {
2545          tmp_2->addTaggedValue(i->first,tmp_2->getDefaultValue()); // use tmp_2 to get correct shape          tmp_2->addTag(i->first); // use tmp_2 to get correct shape
2546        }        }
2547        for (i=lookup_1.begin();i!=lookup_1.end();i++) {        for (i=lookup_1.begin();i!=lookup_1.end();i++) {
2548          tmp_2->addTaggedValue(i->first,tmp_2->getDefaultValue());          tmp_2->addTag(i->first);
2549        }        }
2550        // Compute a result for each tag        // Compute a result for each tag
2551        const DataTagged::DataMapType& lookup_2=tmp_2->getTagLookup();        const DataTagged::DataMapType& lookup_2=tmp_2->getTagLookup();
2552        for (i=lookup_2.begin();i!=lookup_2.end();i++) {        for (i=lookup_2.begin();i!=lookup_2.end();i++) {
2553          DataArrayView view_0 = tmp_0->getDataPointByTag(i->first);          const double *ptr_0 = &(tmp_0->getDataByTagRO(i->first,0));
2554          DataArrayView view_1 = tmp_1->getDataPointByTag(i->first);          const double *ptr_1 = &(tmp_1->getDataByTagRO(i->first,0));
2555          DataArrayView view_2 = tmp_2->getDataPointByTag(i->first);          double *ptr_2 = &(tmp_2->getDataByTagRW(i->first,0));
2556          double *ptr_0 = &view_0.getData(0);  
         double *ptr_1 = &view_1.getData(0);  
         double *ptr_2 = &view_2.getData(0);  
2557          tensor_binary_operation(size1, ptr_0[0], ptr_1, ptr_2, operation);          tensor_binary_operation(size1, ptr_0[0], ptr_1, ptr_2, operation);
2558        }        }
2559    
# Line 2159  C_TensorBinaryOperation(Data const &arg_ Line 2569  C_TensorBinaryOperation(Data const &arg_
2569        int sampleNo_0,dataPointNo_0;        int sampleNo_0,dataPointNo_0;
2570        int numSamples_0 = arg_0_Z.getNumSamples();        int numSamples_0 = arg_0_Z.getNumSamples();
2571        int numDataPointsPerSample_0 = arg_0_Z.getNumDataPointsPerSample();        int numDataPointsPerSample_0 = arg_0_Z.getNumDataPointsPerSample();
2572          res.requireWrite();
2573        #pragma omp parallel for private(sampleNo_0,dataPointNo_0) schedule(static)        #pragma omp parallel for private(sampleNo_0,dataPointNo_0) schedule(static)
2574        for (sampleNo_0 = 0; sampleNo_0 < numSamples_0; sampleNo_0++) {        for (sampleNo_0 = 0; sampleNo_0 < numSamples_0; sampleNo_0++) {
2575          int offset_0 = tmp_0->getPointOffset(sampleNo_0,0); // They're all the same, so just use #0          int offset_0 = tmp_0->getPointOffset(sampleNo_0,0); // They're all the same, so just use #0
2576          double *ptr_0 = &((arg_0_Z.getPointDataView().getData())[offset_0]);          const double *ptr_0 = &(arg_0_Z.getDataAtOffsetRO(offset_0));
2577          for (dataPointNo_0 = 0; dataPointNo_0 < numDataPointsPerSample_0; dataPointNo_0++) {          for (dataPointNo_0 = 0; dataPointNo_0 < numDataPointsPerSample_0; dataPointNo_0++) {
2578            int offset_1 = tmp_1->getPointOffset(sampleNo_0,dataPointNo_0);            int offset_1 = tmp_1->getPointOffset(sampleNo_0,dataPointNo_0);
2579            int offset_2 = tmp_2->getPointOffset(sampleNo_0,dataPointNo_0);            int offset_2 = tmp_2->getPointOffset(sampleNo_0,dataPointNo_0);
2580            double *ptr_1 = &((arg_1_Z.getPointDataView().getData())[offset_1]);            const double *ptr_1 = &(arg_1_Z.getDataAtOffsetRO(offset_1));
2581            double *ptr_2 = &((res.getPointDataView().getData())[offset_2]);            double *ptr_2 = &(res.getDataAtOffsetRW(offset_2));
2582            tensor_binary_operation(size1, ptr_0[0], ptr_1, ptr_2, operation);            tensor_binary_operation(size1, ptr_0[0], ptr_1, ptr_2, operation);
2583          }          }
2584        }        }
2585    
2586      }      }
2587      else if (arg_0_Z.isExpanded()   && arg_1_Z.isConstant()) {      else if (arg_0_Z.isExpanded()   && arg_1_Z.isConstant()) {
   
2588        res = Data(0.0, shape1, arg_1_Z.getFunctionSpace(),true); // DataExpanded output        res = Data(0.0, shape1, arg_1_Z.getFunctionSpace(),true); // DataExpanded output
2589        DataExpanded* tmp_0=dynamic_cast<DataExpanded*>(arg_0_Z.borrowData());        DataExpanded* tmp_0=dynamic_cast<DataExpanded*>(arg_0_Z.borrowData());
2590        DataConstant* tmp_1=dynamic_cast<DataConstant*>(arg_1_Z.borrowData());        DataConstant* tmp_1=dynamic_cast<DataConstant*>(arg_1_Z.borrowData());
# Line 2184  C_TensorBinaryOperation(Data const &arg_ Line 2594  C_TensorBinaryOperation(Data const &arg_
2594        int numSamples_0 = arg_0_Z.getNumSamples();        int numSamples_0 = arg_0_Z.getNumSamples();
2595        int numDataPointsPerSample_0 = arg_0_Z.getNumDataPointsPerSample();        int numDataPointsPerSample_0 = arg_0_Z.getNumDataPointsPerSample();
2596        int offset_1 = tmp_1->getPointOffset(0,0);        int offset_1 = tmp_1->getPointOffset(0,0);
2597          res.requireWrite();
2598        #pragma omp parallel for private(sampleNo_0,dataPointNo_0) schedule(static)        #pragma omp parallel for private(sampleNo_0,dataPointNo_0) schedule(static)
2599        for (sampleNo_0 = 0; sampleNo_0 < numSamples_0; sampleNo_0++) {        for (sampleNo_0 = 0; sampleNo_0 < numSamples_0; sampleNo_0++) {
2600          for (dataPointNo_0 = 0; dataPointNo_0 < numDataPointsPerSample_0; dataPointNo_0++) {          for (dataPointNo_0 = 0; dataPointNo_0 < numDataPointsPerSample_0; dataPointNo_0++) {
2601            int offset_0 = tmp_0->getPointOffset(sampleNo_0,dataPointNo_0);            int offset_0 = tmp_0->getPointOffset(sampleNo_0,dataPointNo_0);
2602            int offset_2 = tmp_2->getPointOffset(sampleNo_0,dataPointNo_0);            int offset_2 = tmp_2->getPointOffset(sampleNo_0,dataPointNo_0);
2603            double *ptr_0 = &((arg_0_Z.getPointDataView().getData())[offset_0]);            const double *ptr_0 = &(arg_0_Z.getDataAtOffsetRO(offset_0));
2604            double *ptr_1 = &((arg_1_Z.getPointDataView().getData())[offset_1]);            const double *ptr_1 = &(arg_1_Z.getDataAtOffsetRO(offset_1));
2605            double *ptr_2 = &((res.getPointDataView().getData())[offset_2]);            double *ptr_2 = &(res.getDataAtOffsetRW(offset_2));
2606            tensor_binary_operation(size1, ptr_0[0], ptr_1, ptr_2, operation);            tensor_binary_operation(size1, ptr_0[0], ptr_1, ptr_2, operation);
2607          }          }
2608        }        }
# Line 2209  C_TensorBinaryOperation(Data const &arg_ Line 2620  C_TensorBinaryOperation(Data const &arg_
2620        int sampleNo_0,dataPointNo_0;        int sampleNo_0,dataPointNo_0;
2621        int numSamples_0 = arg_0_Z.getNumSamples();        int numSamples_0 = arg_0_Z.getNumSamples();
2622        int numDataPointsPerSample_0 = arg_0_Z.getNumDataPointsPerSample();        int numDataPointsPerSample_0 = arg_0_Z.getNumDataPointsPerSample();
2623          res.requireWrite();
2624        #pragma omp parallel for private(sampleNo_0,dataPointNo_0) schedule(static)        #pragma omp parallel for private(sampleNo_0,dataPointNo_0) schedule(static)
2625        for (sampleNo_0 = 0; sampleNo_0 < numSamples_0; sampleNo_0++) {        for (sampleNo_0 = 0; sampleNo_0 < numSamples_0; sampleNo_0++) {
2626          int offset_1 = tmp_1->getPointOffset(sampleNo_0,0);          int offset_1 = tmp_1->getPointOffset(sampleNo_0,0);
2627          double *ptr_1 = &((arg_1_Z.getPointDataView().getData())[offset_1]);          const double *ptr_1 = &(arg_1_Z.getDataAtOffsetRO(offset_1));
2628          for (dataPointNo_0 = 0; dataPointNo_0 < numDataPointsPerSample_0; dataPointNo_0++) {          for (dataPointNo_0 = 0; dataPointNo_0 < numDataPointsPerSample_0; dataPointNo_0++) {
2629            int offset_0 = tmp_0->getPointOffset(sampleNo_0,dataPointNo_0);            int offset_0 = tmp_0->getPointOffset(sampleNo_0,dataPointNo_0);
2630            int offset_2 = tmp_2->getPointOffset(sampleNo_0,dataPointNo_0);            int offset_2 = tmp_2->getPointOffset(sampleNo_0,dataPointNo_0);
2631            double *ptr_0 = &((arg_0_Z.getPointDataView().getData())[offset_0]);            const double *ptr_0 = &(arg_0_Z.getDataAtOffsetRO(offset_0));
2632            double *ptr_2 = &((res.getPointDataView().getData())[offset_2]);            double *ptr_2 = &(res.getDataAtOffsetRW(offset_2));
2633            tensor_binary_operation(size1, ptr_0[0], ptr_1, ptr_2, operation);            tensor_binary_operation(size1, ptr_0[0], ptr_1, ptr_2, operation);
2634          }          }
2635        }        }
# Line 2234  C_TensorBinaryOperation(Data const &arg_ Line 2646  C_TensorBinaryOperation(Data const &arg_
2646        int sampleNo_0,dataPointNo_0;        int sampleNo_0,dataPointNo_0;
2647        int numSamples_0 = arg_0_Z.getNumSamples();        int numSamples_0 = arg_0_Z.getNumSamples();
2648        int numDataPointsPerSample_0 = arg_0_Z.getNumDataPointsPerSample();        int numDataPointsPerSample_0 = arg_0_Z.getNumDataPointsPerSample();
2649          res.requireWrite();
2650        #pragma omp parallel for private(sampleNo_0,dataPointNo_0) schedule(static)        #pragma omp parallel for private(sampleNo_0,dataPointNo_0) schedule(static)
2651        for (sampleNo_0 = 0; sampleNo_0 < numSamples_0; sampleNo_0++) {        for (sampleNo_0 = 0; sampleNo_0 < numSamples_0; sampleNo_0++) {
2652          for (dataPointNo_0 = 0; dataPointNo_0 < numDataPointsPerSample_0; dataPointNo_0++) {          for (dataPointNo_0 = 0; dataPointNo_0 < numDataPointsPerSample_0; dataPointNo_0++) {
2653            int offset_0 = tmp_0->getPointOffset(sampleNo_0,dataPointNo_0);            int offset_0 = tmp_0->getPointOffset(sampleNo_0,dataPointNo_0);
2654            int offset_1 = tmp_1->getPointOffset(sampleNo_0,dataPointNo_0);            int offset_1 = tmp_1->getPointOffset(sampleNo_0,dataPointNo_0);
2655            int offset_2 = tmp_2->getPointOffset(sampleNo_0,dataPointNo_0);            int offset_2 = tmp_2->getPointOffset(sampleNo_0,dataPointNo_0);
2656            double *ptr_0 = &((arg_0_Z.getPointDataView().getData())[offset_0]);            const double *ptr_0 = &(arg_0_Z.getDataAtOffsetRO(offset_0));
2657            double *ptr_1 = &((arg_1_Z.getPointDataView().getData())[offset_1]);            const double *ptr_1 = &(arg_1_Z.getDataAtOffsetRO(offset_1));
2658            double *ptr_2 = &((res.getPointDataView().getData())[offset_2]);            double *ptr_2 = &(res.getDataAtOffsetRW(offset_2));
2659            tensor_binary_operation(size1, ptr_0[0], ptr_1, ptr_2, operation);            tensor_binary_operation(size1, ptr_0[0], ptr_1, ptr_2, operation);
2660          }          }
2661        }        }
# Line 2253  C_TensorBinaryOperation(Data const &arg_ Line 2666  C_TensorBinaryOperation(Data const &arg_
2666      }      }
2667    
2668    } else if (0 == rank1) {    } else if (0 == rank1) {
   
2669      if (arg_0_Z.isConstant()   && arg_1_Z.isConstant()) {      if (arg_0_Z.isConstant()   && arg_1_Z.isConstant()) {
2670        res = Data(0.0, shape0, arg_1_Z.getFunctionSpace());      // DataConstant output        res = Data(0.0, shape0, arg_1_Z.getFunctionSpace());      // DataConstant output
2671        double *ptr_0 = &((arg_0_Z.getPointDataView().getData())[0]);        const double *ptr_0 = &(arg_0_Z.getDataAtOffsetRO(0));
2672        double *ptr_1 = &((arg_1_Z.getPointDataView().getData())[0]);        const double *ptr_1 = &(arg_1_Z.getDataAtOffsetRO(0));
2673        double *ptr_2 = &((res.getPointDataView().getData())[0]);        double *ptr_2 = &(res.getDataAtOffsetRW(0));
2674        tensor_binary_operation(size0, ptr_0, ptr_1[0], ptr_2, operation);        tensor_binary_operation(size0, ptr_0, ptr_1[0], ptr_2, operation);
2675      }      }
2676      else if (arg_0_Z.isConstant()   && arg_1_Z.isTagged()) {      else if (arg_0_Z.isConstant()   && arg_1_Z.isTagged()) {
# Line 2276  C_TensorBinaryOperation(Data const &arg_ Line 2688  C_TensorBinaryOperation(Data const &arg_
2688    
2689        // Prepare offset into DataConstant        // Prepare offset into DataConstant
2690        int offset_0 = tmp_0->getPointOffset(0,0);        int offset_0 = tmp_0->getPointOffset(0,0);
2691        double *ptr_0 = &((arg_0_Z.getPointDataView().getData())[offset_0]);        const double *ptr_0 = &(arg_0_Z.getDataAtOffsetRO(offset_0));
2692        // Get the views  
2693        DataArrayView view_1 = tmp_1->getDefaultValue();        //Get the pointers to the actual data
2694        DataArrayView view_2 = tmp_2->getDefaultValue();        const double *ptr_1 = &(tmp_1->getDefaultValueRO(0));
2695        // Get the pointers to the actual data        double *ptr_2 = &(tmp_2->getDefaultValueRW(0));
2696        double *ptr_1 = &((view_1.getData())[0]);  
       double *ptr_2 = &((view_2.getData())[0]);  
2697        // Compute a result for the default        // Compute a result for the default
2698        tensor_binary_operation(size0, ptr_0, ptr_1[0], ptr_2, operation);        tensor_binary_operation(size0, ptr_0, ptr_1[0], ptr_2, operation);
2699        // Compute a result for each tag        // Compute a result for each tag
2700        const DataTagged::DataMapType& lookup_1=tmp_1->getTagLookup();        const DataTagged::DataMapType& lookup_1=tmp_1->getTagLookup();
2701        DataTagged::DataMapType::const_iterator i; // i->first is a tag, i->second is an offset into memory        DataTagged::DataMapType::const_iterator i; // i->first is a tag, i->second is an offset into memory
2702        for (i=lookup_1.begin();i!=lookup_1.end();i++) {        for (i=lookup_1.begin();i!=lookup_1.end();i++) {
2703          tmp_2->addTaggedValue(i->first,tmp_2->getDefaultValue());          tmp_2->addTag(i->first);
2704          DataArrayView view_1 = tmp_1->getDataPointByTag(i->first);          const double *ptr_1 = &(tmp_1->getDataByTagRO(i->first,0));
2705          DataArrayView view_2 = tmp_2->getDataPointByTag(i->first);          double *ptr_2 = &(tmp_2->getDataByTagRW(i->first,0));
         double *ptr_1 = &view_1.getData(0);  
         double *ptr_2 = &view_2.getData(0);  
2706          tensor_binary_operation(size0, ptr_0, ptr_1[0], ptr_2, operation);          tensor_binary_operation(size0, ptr_0, ptr_1[0], ptr_2, operation);
2707        }        }
   
2708      }      }
2709      else if (arg_0_Z.isConstant()   && arg_1_Z.isExpanded()) {      else if (arg_0_Z.isConstant()   && arg_1_Z.isExpanded()) {
2710    
# Line 2309  C_TensorBinaryOperation(Data const &arg_ Line 2717  C_TensorBinaryOperation(Data const &arg_
2717        int numSamples_1 = arg_1_Z.getNumSamples();        int numSamples_1 = arg_1_Z.getNumSamples();
2718        int numDataPointsPerSample_1 = arg_1_Z.getNumDataPointsPerSample();        int numDataPointsPerSample_1 = arg_1_Z.getNumDataPointsPerSample();
2719        int offset_0 = tmp_0->getPointOffset(0,0);        int offset_0 = tmp_0->getPointOffset(0,0);
2720          res.requireWrite();
2721        #pragma omp parallel for private(sampleNo_1,dataPointNo_1) schedule(static)        #pragma omp parallel for private(sampleNo_1,dataPointNo_1) schedule(static)
2722        for (sampleNo_1 = 0; sampleNo_1 < numSamples_1; sampleNo_1++) {        for (sampleNo_1 = 0; sampleNo_1 < numSamples_1; sampleNo_1++) {
2723          for (dataPointNo_1 = 0; dataPointNo_1 < numDataPointsPerSample_1; dataPointNo_1++) {          for (dataPointNo_1 = 0; dataPointNo_1 < numDataPointsPerSample_1; dataPointNo_1++) {
2724            int offset_1 = tmp_1->getPointOffset(sampleNo_1,dataPointNo_1);            int offset_1 = tmp_1->getPointOffset(sampleNo_1,dataPointNo_1);
2725            int offset_2 = tmp_2->getPointOffset(sampleNo_1,dataPointNo_1);            int offset_2 = tmp_2->getPointOffset(sampleNo_1,dataPointNo_1);
2726            double *ptr_0 = &((arg_0_Z.getPointDataView().getData())[offset_0]);            const double *ptr_0 = &(arg_0_Z.getDataAtOffsetRO(offset_0));
2727            double *ptr_1 = &((arg_1_Z.getPointDataView().getData())[offset_1]);            const double *ptr_1 = &(arg_1_Z.getDataAtOffsetRO(offset_1));
2728            double *ptr_2 = &((res.getPointDataView().getData())[offset_2]);            double *ptr_2 = &(res.getDataAtOffsetRW(offset_2));
2729            tensor_binary_operation(size0, ptr_0, ptr_1[0], ptr_2, operation);            tensor_binary_operation(size0, ptr_0, ptr_1[0], ptr_2, operation);
2730          }          }
2731        }        }
# Line 2337  C_TensorBinaryOperation(Data const &arg_ Line 2746  C_TensorBinaryOperation(Data const &arg_
2746    
2747        // Prepare offset into DataConstant        // Prepare offset into DataConstant
2748        int offset_1 = tmp_1->getPointOffset(0,0);        int offset_1 = tmp_1->getPointOffset(0,0);
2749        double *ptr_1 = &((arg_1_Z.getPointDataView().getData())[offset_1]);        const double *ptr_1 = &(arg_1_Z.getDataAtOffsetRO(offset_1));
       // Get the views  
       DataArrayView view_0 = tmp_0->getDefaultValue();  
       DataArrayView view_2 = tmp_2->getDefaultValue();  
2750        // Get the pointers to the actual data        // Get the pointers to the actual data
2751        double *ptr_0 = &((view_0.getData())[0]);        const double *ptr_0 = &(tmp_0->getDefaultValueRO(0));
2752        double *ptr_2 = &((view_2.getData())[0]);        double *ptr_2 = &(tmp_2->getDefaultValueRW(0));
2753        // Compute a result for the default        // Compute a result for the default
2754        tensor_binary_operation(size0, ptr_0, ptr_1[0], ptr_2, operation);        tensor_binary_operation(size0, ptr_0, ptr_1[0], ptr_2, operation);
2755        // Compute a result for each tag        // Compute a result for each tag
2756        const DataTagged::DataMapType& lookup_0=tmp_0->getTagLookup();        const DataTagged::DataMapType& lookup_0=tmp_0->getTagLookup();
2757        DataTagged::DataMapType::const_iterator i; // i->first is a tag, i->second is an offset into memory        DataTagged::DataMapType::const_iterator i; // i->first is a tag, i->second is an offset into memory
2758        for (i=lookup_0.begin();i!=lookup_0.end();i++) {        for (i=lookup_0.begin();i!=lookup_0.end();i++) {
2759          tmp_2->addTaggedValue(i->first,tmp_2->getDefaultValue());          tmp_2->addTag(i->first);
2760          DataArrayView view_0 = tmp_0->getDataPointByTag(i->first);          const double *ptr_0 = &(tmp_0->getDataByTagRO(i->first,0));
2761          DataArrayView view_2 = tmp_2->getDataPointByTag(i->first);          double *ptr_2 = &(tmp_2->getDataByTagRW(i->first,0));
         double *ptr_0 = &view_0.getData(0);  
         double *ptr_2 = &view_2.getData(0);  
2762          tensor_binary_operation(size0, ptr_0, ptr_1[0], ptr_2, operation);          tensor_binary_operation(size0, ptr_0, ptr_1[0], ptr_2, operation);
2763        }        }
2764    
# Line 2372  C_TensorBinaryOperation(Data const &arg_ Line 2776  C_TensorBinaryOperation(Data const &arg_
2776        res.tag();        // DataTagged output        res.tag();        // DataTagged output
2777        DataTagged* tmp_2=dynamic_cast<DataTagged*>(res.borrowData());        DataTagged* tmp_2=dynamic_cast<DataTagged*>(res.borrowData());
2778    
       // Get the views  
       DataArrayView view_0 = tmp_0->getDefaultValue();  
       DataArrayView view_1 = tmp_1->getDefaultValue();  
       DataArrayView view_2 = tmp_2->getDefaultValue();  
2779        // Get the pointers to the actual data        // Get the pointers to the actual data
2780        double *ptr_0 = &((view_0.getData())[0]);        const double *ptr_0 = &(tmp_0->getDefaultValueRO(0));
2781        double *ptr_1 = &((view_1.getData())[0]);        const double *ptr_1 = &(tmp_1->getDefaultValueRO(0));
2782        double *ptr_2 = &((view_2.getData())[0]);        double *ptr_2 = &(tmp_2->getDefaultValueRW(0));
2783    
2784        // Compute a result for the default        // Compute a result for the default
2785        tensor_binary_operation(size0, ptr_0, ptr_1[0], ptr_2, operation);        tensor_binary_operation(size0, ptr_0, ptr_1[0], ptr_2, operation);
2786        // Merge the tags        // Merge the tags
# Line 2387  C_TensorBinaryOperation(Data const &arg_ Line 2788  C_TensorBinaryOperation(Data const &arg_
2788        const DataTagged::DataMapType& lookup_0=tmp_0->getTagLookup();        const DataTagged::DataMapType& lookup_0=tmp_0->getTagLookup();
2789        const DataTagged::DataMapType& lookup_1=tmp_1->getTagLookup();        const DataTagged::DataMapType& lookup_1=tmp_1->getTagLookup();
2790        for (i=lookup_0.begin();i!=lookup_0.end();i++) {        for (i=lookup_0.begin();i!=lookup_0.end();i++) {
2791          tmp_2->addTaggedValue(i->first,tmp_2->getDefaultValue()); // use tmp_2 to get correct shape          tmp_2->addTag(i->first); // use tmp_2 to get correct shape
2792        }        }
2793        for (i=lookup_1.begin();i!=lookup_1.end();i++) {        for (i=lookup_1.begin();i!=lookup_1.end();i++) {
2794          tmp_2->addTaggedValue(i->first,tmp_2->getDefaultValue());          tmp_2->addTag(i->first);
2795        }        }
2796        // Compute a result for each tag        // Compute a result for each tag
2797        const DataTagged::DataMapType& lookup_2=tmp_2->getTagLookup();        const DataTagged::DataMapType& lookup_2=tmp_2->getTagLookup();
2798        for (i=lookup_2.begin();i!=lookup_2.end();i++) {        for (i=lookup_2.begin();i!=lookup_2.end();i++) {
2799          DataArrayView view_0 = tmp_0->getDataPointByTag(i->first);          const double *ptr_0 = &(tmp_0->getDataByTagRO(i->first,0));
2800          DataArrayView view_1 = tmp_1->getDataPointByTag(i->first);          const double *ptr_1 = &(tmp_1->getDataByTagRO(i->first,0));
2801          DataArrayView view_2 = tmp_2->getDataPointByTag(i->first);          double *ptr_2 = &(tmp_2->getDataByTagRW(i->first,0));
         double *ptr_0 = &view_0.getData(0);  
         double *ptr_1 = &view_1.getData(0);  
         double *ptr_2 = &view_2.getData(0);  
2802          tensor_binary_operation(size0, ptr_0, ptr_1[0], ptr_2, operation);          tensor_binary_operation(size0, ptr_0, ptr_1[0], ptr_2, operation);
2803        }        }
2804    
# Line 2416  C_TensorBinaryOperation(Data const &arg_ Line 2814  C_TensorBinaryOperation(Data const &arg_
2814        int sampleNo_0,dataPointNo_0;        int sampleNo_0,dataPointNo_0;
2815        int numSamples_0 = arg_0_Z.getNumSamples();        int numSamples_0 = arg_0_Z.getNumSamples();
2816        int numDataPointsPerSample_0 = arg_0_Z.getNumDataPointsPerSample();        int numDataPointsPerSample_0 = arg_0_Z.getNumDataPointsPerSample();
2817          res.requireWrite();
2818        #pragma omp parallel for private(sampleNo_0,dataPointNo_0) schedule(static)        #pragma omp parallel for private(sampleNo_0,dataPointNo_0) schedule(static)
2819        for (sampleNo_0 = 0; sampleNo_0 < numSamples_0; sampleNo_0++) {        for (sampleNo_0 = 0; sampleNo_0 < numSamples_0; sampleNo_0++) {
2820          int offset_0 = tmp_0->getPointOffset(sampleNo_0,0); // They're all the same, so just use #0          int offset_0 = tmp_0->getPointOffset(sampleNo_0,0); // They're all the same, so just use #0
2821          double *ptr_0 = &((arg_0_Z.getPointDataView().getData())[offset_0]);          const double *ptr_0 = &(arg_0_Z.getDataAtOffsetRO(offset_0));
2822          for (dataPointNo_0 = 0; dataPointNo_0 < numDataPointsPerSample_0; dataPointNo_0++) {          for (dataPointNo_0 = 0; dataPointNo_0 < numDataPointsPerSample_0; dataPointNo_0++) {
2823            int offset_1 = tmp_1->getPointOffset(sampleNo_0,dataPointNo_0);            int offset_1 = tmp_1->getPointOffset(sampleNo_0,dataPointNo_0);
2824            int offset_2 = tmp_2->getPointOffset(sampleNo_0,dataPointNo_0);            int offset_2 = tmp_2->getPointOffset(sampleNo_0,dataPointNo_0);
2825            double *ptr_1 = &((arg_1_Z.getPointDataView().getData())[offset_1]);            const double *ptr_1 = &(arg_1_Z.getDataAtOffsetRO(offset_1));
2826            double *ptr_2 = &((res.getPointDataView().getData())[offset_2]);            double *ptr_2 = &(res.getDataAtOffsetRW(offset_2));
2827            tensor_binary_operation(size0, ptr_0, ptr_1[0], ptr_2, operation);            tensor_binary_operation(size0, ptr_0, ptr_1[0], ptr_2, operation);
2828          }          }
2829        }        }
2830    
2831      }      }
2832      else if (arg_0_Z.isExpanded()   && arg_1_Z.isConstant()) {      else if (arg_0_Z.isExpanded()   && arg_1_Z.isConstant()) {
   
2833        res = Data(0.0, shape0, arg_1_Z.getFunctionSpace(),true); // DataExpanded output        res = Data(0.0, shape0, arg_1_Z.getFunctionSpace(),true); // DataExpanded output
2834        DataExpanded* tmp_0=dynamic_cast<DataExpanded*>(arg_0_Z.borrowData());        DataExpanded* tmp_0=dynamic_cast<DataExpanded*>(arg_0_Z.borrowData());
2835        DataConstant* tmp_1=dynamic_cast<DataConstant*>(arg_1_Z.borrowData());        DataConstant* tmp_1=dynamic_cast<DataConstant*>(arg_1_Z.borrowData());
# Line 2441  C_TensorBinaryOperation(Data const &arg_ Line 2839  C_TensorBinaryOperation(Data const &arg_
2839        int numSamples_0 = arg_0_Z.getNumSamples();        int numSamples_0 = arg_0_Z.getNumSamples();
2840        int numDataPointsPerSample_0 = arg_0_Z.getNumDataPointsPerSample();        int numDataPointsPerSample_0 = arg_0_Z.getNumDataPointsPerSample();
2841        int offset_1 = tmp_1->getPointOffset(0,0);        int offset_1 = tmp_1->getPointOffset(0,0);
2842          res.requireWrite();
2843        #pragma omp parallel for private(sampleNo_0,dataPointNo_0) schedule(static)        #pragma omp parallel for private(sampleNo_0,dataPointNo_0) schedule(static)
2844        for (sampleNo_0 = 0; sampleNo_0 < numSamples_0; sampleNo_0++) {        for (sampleNo_0 = 0; sampleNo_0 < numSamples_0; sampleNo_0++) {
2845          for (dataPointNo_0 = 0; dataPointNo_0 < numDataPointsPerSample_0; dataPointNo_0++) {          for (dataPointNo_0 = 0; dataPointNo_0 < numDataPointsPerSample_0; dataPointNo_0++) {
2846            int offset_0 = tmp_0->getPointOffset(sampleNo_0,dataPointNo_0);            int offset_0 = tmp_0->getPointOffset(sampleNo_0,dataPointNo_0);
2847            int offset_2 = tmp_2->getPointOffset(sampleNo_0,dataPointNo_0);            int offset_2 = tmp_2->getPointOffset(sampleNo_0,dataPointNo_0);
2848            double *ptr_0 = &((arg_0_Z.getPointDataView().getData())[offset_0]);            const double *ptr_0 = &(arg_0_Z.getDataAtOffsetRO(offset_0));
2849            double *ptr_1 = &((arg_1_Z.getPointDataView().getData())[offset_1]);            const double *ptr_1 = &(arg_1_Z.getDataAtOffsetRO(offset_1));
2850            double *ptr_2 = &((res.getPointDataView().getData())[offset_2]);            double *ptr_2 = &(res.getDataAtOffsetRW(offset_2));
2851            tensor_binary_operation(size0, ptr_0, ptr_1[0], ptr_2, operation);            tensor_binary_operation(size0, ptr_0, ptr_1[0], ptr_2, operation);
2852          }          }
2853        }        }
# Line 2466  C_TensorBinaryOperation(Data const &arg_ Line 2865  C_TensorBinaryOperation(Data const &arg_
2865        int sampleNo_0,dataPointNo_0;        int sampleNo_0,dataPointNo_0;
2866        int numSamples_0 = arg_0_Z.getNumSamples();        int numSamples_0 = arg_0_Z.getNumSamples();
2867        int numDataPointsPerSample_0 = arg_0_Z.getNumDataPointsPerSample();        int numDataPointsPerSample_0 = arg_0_Z.getNumDataPointsPerSample();
2868          res.requireWrite();
2869        #pragma omp parallel for private(sampleNo_0,dataPointNo_0) schedule(static)        #pragma omp parallel for private(sampleNo_0,dataPointNo_0) schedule(static)
2870        for (sampleNo_0 = 0; sampleNo_0 < numSamples_0; sampleNo_0++) {        for (sampleNo_0 = 0; sampleNo_0 < numSamples_0; sampleNo_0++) {
2871          int offset_1 = tmp_1->getPointOffset(sampleNo_0,0);          int offset_1 = tmp_1->getPointOffset(sampleNo_0,0);
2872          double *ptr_1 = &((arg_1_Z.getPointDataView().getData())[offset_1]);          const double *ptr_1 = &(arg_1_Z.getDataAtOffsetRO(offset_1));
2873          for (dataPointNo_0 = 0; dataPointNo_0 < numDataPointsPerSample_0; dataPointNo_0++) {          for (dataPointNo_0 = 0; dataPointNo_0 < numDataPointsPerSample_0; dataPointNo_0++) {
2874            int offset_0 = tmp_0->getPointOffset(sampleNo_0,dataPointNo_0);            int offset_0 = tmp_0->getPointOffset(sampleNo_0,dataPointNo_0);
2875            int offset_2 = tmp_2->getPointOffset(sampleNo_0,dataPointNo_0);            int offset_2 = tmp_2->getPointOffset(sampleNo_0,dataPointNo_0);
2876            double *ptr_0 = &((arg_0_Z.getPointDataView().getData())[offset_0]);            const double *ptr_0 = &(arg_0_Z.getDataAtOffsetRO(offset_0));
2877            double *ptr_2 = &((res.getPointDataView().getData())[offset_2]);            double *ptr_2 = &(res.getDataAtOffsetRW(offset_2));
2878            tensor_binary_operation(size0, ptr_0, ptr_1[0], ptr_2, operation);            tensor_binary_operation(size0, ptr_0, ptr_1[0], ptr_2, operation);
2879          }          }
2880        }        }
# Line 2491  C_TensorBinaryOperation(Data const &arg_ Line 2891  C_TensorBinaryOperation(Data const &arg_
2891        int sampleNo_0,dataPointNo_0;        int sampleNo_0,dataPointNo_0;
2892        int numSamples_0 = arg_0_Z.getNumSamples();        int numSamples_0 = arg_0_Z.getNumSamples();
2893        int numDataPointsPerSample_0 = arg_0_Z.getNumDataPointsPerSample();        int numDataPointsPerSample_0 = arg_0_Z.getNumDataPointsPerSample();
2894          res.requireWrite();
2895        #pragma omp parallel for private(sampleNo_0,dataPointNo_0) schedule(static)        #pragma omp parallel for private(sampleNo_0,dataPointNo_0) schedule(static)
2896        for (sampleNo_0 = 0; sampleNo_0 < numSamples_0; sampleNo_0++) {        for (sampleNo_0 = 0; sampleNo_0 < numSamples_0; sampleNo_0++) {
2897          for (dataPointNo_0 = 0; dataPointNo_0 < numDataPointsPerSample_0; dataPointNo_0++) {          for (dataPointNo_0 = 0; dataPointNo_0 < numDataPointsPerSample_0; dataPointNo_0++) {
2898            int offset_0 = tmp_0->getPointOffset(sampleNo_0,dataPointNo_0);            int offset_0 = tmp_0->getPointOffset(sampleNo_0,dataPointNo_0);
2899            int offset_1 = tmp_1->getPointOffset(sampleNo_0,dataPointNo_0);            int offset_1 = tmp_1->getPointOffset(sampleNo_0,dataPointNo_0);
2900            int offset_2 = tmp_2->getPointOffset(sampleNo_0,dataPointNo_0);            int offset_2 = tmp_2->getPointOffset(sampleNo_0,dataPointNo_0);
2901            double *ptr_0 = &((arg_0_Z.getPointDataView().getData())[offset_0]);            const double *ptr_0 = &(arg_0_Z.getDataAtOffsetRO(offset_0));
2902            double *ptr_1 = &((arg_1_Z.getPointDataView().getData())[offset_1]);            const double *ptr_1 = &(arg_1_Z.getDataAtOffsetRO(offset_1));
2903            double *ptr_2 = &((res.getPointDataView().getData())[offset_2]);            double *ptr_2 = &(res.getDataAtOffsetRW(offset_2));
2904            tensor_binary_operation(size0, ptr_0, ptr_1[0], ptr_2, operation);            tensor_binary_operation(size0, ptr_0, ptr_1[0], ptr_2, operation);
2905          }          }
2906        }        }
# Line 2521  Data Line 2922  Data
2922  C_TensorUnaryOperation(Data const &arg_0,  C_TensorUnaryOperation(Data const &arg_0,
2923                         UnaryFunction operation)                         UnaryFunction operation)
2924  {  {
2925      if (arg_0.isEmpty())  // do this before we attempt to interpolate
2926      {
2927         throw DataException("Error - Operations not permitted on instances of DataEmpty.");
2928      }
2929      if (arg_0.isLazy())
2930      {
2931         throw DataException("Error - Operations not permitted on lazy data.");
2932      }
2933    // Interpolate if necessary and find an appropriate function space    // Interpolate if necessary and find an appropriate function space
2934    Data arg_0_Z = Data(arg_0);    Data arg_0_Z = Data(arg_0);
2935    
2936    // Get rank and shape of inputs    // Get rank and shape of inputs
2937    int rank0 = arg_0_Z.getDataPointRank();    const DataTypes::ShapeType& shape0 = arg_0_Z.getDataPointShape();
   DataArrayView::ShapeType shape0 = arg_0_Z.getDataPointShape();  
2938    int size0 = arg_0_Z.getDataPointSize();    int size0 = arg_0_Z.getDataPointSize();
2939    
2940    // Declare output Data object    // Declare output Data object
# Line 2534  C_TensorUnaryOperation(Data const &arg_0 Line 2942  C_TensorUnaryOperation(Data const &arg_0
2942    
2943    if (arg_0_Z.isConstant()) {    if (arg_0_Z.isConstant()) {
2944      res = Data(0.0, shape0, arg_0_Z.getFunctionSpace());      // DataConstant output      res = Data(0.0, shape0, arg_0_Z.getFunctionSpace());      // DataConstant output
2945      double *ptr_0 = &((arg_0_Z.getPointDataView().getData())[0]);      const double *ptr_0 = &(arg_0_Z.getDataAtOffsetRO(0));
2946      double *ptr_2 = &((res.getPointDataView().getData())[0]);      double *ptr_2 = &(res.getDataAtOffsetRW(0));
2947      tensor_unary_operation(size0, ptr_0, ptr_2, operation);      tensor_unary_operation(size0, ptr_0, ptr_2, operation);
2948    }    }
2949    else if (arg_0_Z.isTagged()) {    else if (arg_0_Z.isTagged()) {
# Line 2548  C_TensorUnaryOperation(Data const &arg_0 Line 2956  C_TensorUnaryOperation(Data const &arg_0
2956      res.tag();      res.tag();
2957      DataTagged* tmp_2=dynamic_cast<DataTagged*>(res.borrowData());      DataTagged* tmp_2=dynamic_cast<DataTagged*>(res.borrowData());
2958    
     // Get the views  
     DataArrayView view_0 = tmp_0->getDefaultValue();  
     DataArrayView view_2 = tmp_2->getDefaultValue();  
2959      // Get the pointers to the actual data      // Get the pointers to the actual data
2960      double *ptr_0 = &((view_0.getData())[0]);      const double *ptr_0 = &(tmp_0->getDefaultValueRO(0));
2961      double *ptr_2 = &((view_2.getData())[0]);      double *ptr_2 = &(tmp_2->getDefaultValueRW(0));
2962      // Compute a result for the default      // Compute a result for the default
2963      tensor_unary_operation(size0, ptr_0, ptr_2, operation);      tensor_unary_operation(size0, ptr_0, ptr_2, operation);
2964      // Compute a result for each tag      // Compute a result for each tag
2965      const DataTagged::DataMapType& lookup_0=tmp_0->getTagLookup();      const DataTagged::DataMapType& lookup_0=tmp_0->getTagLookup();
2966      DataTagged::DataMapType::const_iterator i; // i->first is a tag, i->second is an offset into memory      DataTagged::DataMapType::const_iterator i; // i->first is a tag, i->second is an offset into memory
2967      for (i=lookup_0.begin();i!=lookup_0.end();i++) {      for (i=lookup_0.begin();i!=lookup_0.end();i++) {
2968        tmp_2->addTaggedValue(i->first,tmp_2->getDefaultValue());        tmp_2->addTag(i->first);
2969        DataArrayView view_0 = tmp_0->getDataPointByTag(i->first);        const double *ptr_0 = &(tmp_0->getDataByTagRO(i->first,0));
2970        DataArrayView view_2 = tmp_2->getDataPointByTag(i->first);        double *ptr_2 = &(tmp_2->getDataByTagRW(i->first,0));
       double *ptr_0 = &view_0.getData(0);  
       double *ptr_2 = &view_2.getData(0);  
2971        tensor_unary_operation(size0, ptr_0, ptr_2, operation);        tensor_unary_operation(size0, ptr_0, ptr_2, operation);
2972      }      }
2973    
# Line 2583  C_TensorUnaryOperation(Data const &arg_0 Line 2986  C_TensorUnaryOperation(Data const &arg_0
2986        for (dataPointNo_0 = 0; dataPointNo_0 < numDataPointsPerSample_0; dataPointNo_0++) {        for (dataPointNo_0 = 0; dataPointNo_0 < numDataPointsPerSample_0; dataPointNo_0++) {
2987          int offset_0 = tmp_0->getPointOffset(sampleNo_0,dataPointNo_0);          int offset_0 = tmp_0->getPointOffset(sampleNo_0,dataPointNo_0);
2988          int offset_2 = tmp_2->getPointOffset(sampleNo_0,dataPointNo_0);          int offset_2 = tmp_2->getPointOffset(sampleNo_0,dataPointNo_0);
2989          double *ptr_0 = &((arg_0_Z.getPointDataView().getData())[offset_0]);          const double *ptr_0 = &(arg_0_Z.getDataAtOffsetRO(offset_0));
2990          double *ptr_2 = &((res.getPointDataView().getData())[offset_2]);          double *ptr_2 = &(res.getDataAtOffsetRW(offset_2));
2991          tensor_unary_operation(size0, ptr_0, ptr_2, operation);          tensor_unary_operation(size0, ptr_0, ptr_2, operation);
2992        }        }
2993      }      }
   
2994    }    }
2995    else {    else {
2996      throw DataException("Error - C_TensorUnaryOperation: unknown combination of inputs");      throw DataException("Error - C_TensorUnaryOperation: unknown combination of inputs");

Legend:
Removed from v.1387  
changed lines
  Added in v.2519

  ViewVC Help
Powered by ViewVC 1.1.26