/[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 2628 by jfenwick, Tue Aug 25 03:50:00 2009 UTC
# Line 1  Line 1 
1    
 /* $Id$ */  
   
2  /*******************************************************  /*******************************************************
3   *  *
4   *           Copyright 2003-2007 by ACceSS MNRF  * Copyright (c) 2003-2009 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);
   
   /**  
      \brief  
      Constructor which will create Tagged data if expanded is false.  
      No attempt is made to ensure the tag keys match the tag keys  
      within the function space.  
   
      \param tagKeys - Input - List of tag values.  
      \param values - Input - List of values, one for each tag.  
      \param defaultValue - Input - A default value, used if tag doesn't exist.  
      \param what - Input - A description of what this data represents.  
      \param expanded - Input - Flag, if true fill the entire container with  
                        the appropriate values.  
     ==>*  
   */  
   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);  
145    
146    /**    /**
147       \brief       \brief
148       Constructor which copies data from any object that can be converted into       Constructor which copies data from any object that can be treated like a python array/sequence.
      a python numarray.  
149    
150       \param value - Input - Input data.       \param value - Input - Input data.
151       \param what - Input - A description of what this data represents.       \param what - Input - A description of what this data represents.
# Line 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 644  class Data { Line 763  class Data {
763    ESCRIPT_DLL_API    ESCRIPT_DLL_API
764    Data    Data
765    interpolate(const FunctionSpace& functionspace) const;    interpolate(const FunctionSpace& functionspace) const;
766    
767    
768      ESCRIPT_DLL_API
769      Data
770      interpolateFromTable2D(const WrappedArray& table, double Amin, double Astep,
771                           double undef, Data& B, double Bmin, double Bstep);
772    
773    
774      // This signature needs to be tuned
775      ESCRIPT_DLL_API
776      Data
777      interpolateFromTable(boost::python::object table, double Amin, double Astep,
778                           double undef, Data& B, double Bmin, double Bstep);
779    
780    /**    /**
781       \brief       \brief
782       Calculates the gradient of the data at the data points of functionspace.       Calculates the gradient of the data at the data points of functionspace.
# Line 659  class Data { Line 792  class Data {
792    grad() const;    grad() const;
793    
794    /**    /**
795       \brief      \brief
796       Calculate the integral over the function space domain.       Calculate the integral over the function space domain as a python tuple.
      *  
797    */    */
798    ESCRIPT_DLL_API    ESCRIPT_DLL_API
799    boost::python::numeric::array    boost::python::object
800    integrate() const;    integrateToTuple_const() const;
801    
802    
803      /**
804        \brief
805         Calculate the integral over the function space domain as a python tuple.
806      */
807      ESCRIPT_DLL_API
808      boost::python::object
809      integrateToTuple();
810    
811    
812    
813    /**    /**
814       \brief       \brief
# Line 732  class Data { Line 875  class Data {
875    /**    /**
876       \brief       \brief
877       Return the maximum absolute value of this Data object.       Return the maximum absolute value of this Data object.
878       *  
879         The method is not const because lazy data needs to be expanded before Lsup can be computed.
880         The _const form can be used when the Data object is const, however this will only work for
881         Data which is not Lazy.
882    
883         For Data which contain no samples (or tagged Data for which no tags in use have a value)
884         zero is returned.
885    */    */
886    ESCRIPT_DLL_API    ESCRIPT_DLL_API
887    double    double
888    Lsup() const;    Lsup();
889    
890      ESCRIPT_DLL_API
891      double
892      Lsup_const() const;
893    
894    
895    /**    /**
896       \brief       \brief
897       Return the maximum value of this Data object.       Return the maximum value of this Data object.
898       *  
899         The method is not const because lazy data needs to be expanded before sup can be computed.
900         The _const form can be used when the Data object is const, however this will only work for
901         Data which is not Lazy.
902    
903         For Data which contain no samples (or tagged Data for which no tags in use have a value)
904         a large negative value is returned.
905    */    */
906    ESCRIPT_DLL_API    ESCRIPT_DLL_API
907    double    double
908    sup() const;    sup();
909    
910      ESCRIPT_DLL_API
911      double
912      sup_const() const;
913    
914    
915    /**    /**
916       \brief       \brief
917       Return the minimum value of this Data object.       Return the minimum value of this Data object.
918       *  
919         The method is not const because lazy data needs to be expanded before inf can be computed.
920         The _const form can be used when the Data object is const, however this will only work for
921         Data which is not Lazy.
922    
923         For Data which contain no samples (or tagged Data for which no tags in use have a value)
924         a large positive value is returned.
925    */    */
926    ESCRIPT_DLL_API    ESCRIPT_DLL_API
927    double    double
928    inf() const;    inf();
929    
930      ESCRIPT_DLL_API
931      double
932      inf_const() const;
933    
934    
935    
936    /**    /**
937       \brief       \brief
# Line 786  class Data { Line 963  class Data {
963    /**    /**
964       \brief       \brief
965       Return the (sample number, data-point number) of the data point with       Return the (sample number, data-point number) of the data point with
966       the minimum value in this Data object.       the minimum component value in this Data object.
967         \note If you are working in python, please consider using Locator
968    instead of manually manipulating process and point IDs.
969    */    */
970    ESCRIPT_DLL_API    ESCRIPT_DLL_API
971    const boost::python::tuple    const boost::python::tuple
972    minGlobalDataPoint() const;    minGlobalDataPoint() const;
973    
974      /**
975         \brief
976         Return the (sample number, data-point number) of the data point with
977         the minimum component value in this Data object.
978         \note If you are working in python, please consider using Locator
979    instead of manually manipulating process and point IDs.
980      */
981    ESCRIPT_DLL_API    ESCRIPT_DLL_API
982    void    const boost::python::tuple
983    calc_minGlobalDataPoint(int& ProcNo,  int& DataPointNo) const;    maxGlobalDataPoint() const;
984    
985    
986    
987    /**    /**
988       \brief       \brief
989       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 1401  class Data {
1401    */    */
1402    ESCRIPT_DLL_API    ESCRIPT_DLL_API
1403    Data    Data
1404    getSlice(const DataArrayView::RegionType& region) const;    getSlice(const DataTypes::RegionType& region) const;
1405    
1406    /**    /**
1407       \brief       \brief
# Line 1225  class Data { Line 1414  class Data {
1414    ESCRIPT_DLL_API    ESCRIPT_DLL_API
1415    void    void
1416    setSlice(const Data& value,    setSlice(const Data& value,
1417             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);  
   
1418    
1419    /**    /**
1420       \brief       \brief
# Line 1290  class Data { Line 1456  class Data {
1456    /**    /**
1457       \brief       \brief
1458       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
1459        TODO Ownership of this object should be explained in doco.
1460    */    */
1461    ESCRIPT_DLL_API    ESCRIPT_DLL_API
1462          DataAbstract*          DataAbstract*
1463          borrowData(void) const;          borrowData(void) const;
1464    
1465      ESCRIPT_DLL_API
1466            DataAbstract_ptr
1467            borrowDataPtr(void) const;
1468    
1469      ESCRIPT_DLL_API
1470            DataReady_ptr
1471            borrowReadyPtr(void) const;
1472    
1473    
1474    
1475      /**
1476         \brief
1477         Return a pointer to the beginning of the datapoint at the specified offset.
1478         TODO Eventually these should be inlined.
1479         \param i - position(offset) in the underlying datastructure
1480      */
1481    
1482      ESCRIPT_DLL_API
1483            DataTypes::ValueType::const_reference
1484            getDataAtOffsetRO(DataTypes::ValueType::size_type i);
1485    
1486    
1487      ESCRIPT_DLL_API
1488            DataTypes::ValueType::reference
1489            getDataAtOffsetRW(DataTypes::ValueType::size_type i);
1490    
1491    
1492    
1493    /**
1494       \brief Create a buffer for use by getSample
1495       Allocates a DataVector large enough for DataLazy::resolveSample to operate on for the current Data.
1496       Do not use this buffer for other Data instances (unless you are sure they will be the same size).
1497      
1498       In multi-threaded sections, this needs to be called on each thread.
1499    
1500       \return A BufferGroup* if Data is lazy, NULL otherwise.
1501       \warning This pointer must be deallocated using freeSampleBuffer to avoid cross library memory issues.
1502    */
1503      ESCRIPT_DLL_API
1504      BufferGroup*
1505      allocSampleBuffer() const;
1506    
1507    /**
1508       \brief Free a buffer allocated with allocSampleBuffer.
1509       \param buffer Input - pointer to the buffer to deallocate.
1510    */
1511    ESCRIPT_DLL_API void freeSampleBuffer(BufferGroup* buffer);
1512    
1513   protected:   protected:
1514    
1515   private:   private:
1516    
1517      double
1518      LsupWorker() const;
1519    
1520      double
1521      supWorker() const;
1522    
1523      double
1524      infWorker() const;
1525    
1526      boost::python::object
1527      integrateWorker() const;
1528    
1529      void
1530      calc_minGlobalDataPoint(int& ProcNo,  int& DataPointNo) const;
1531    
1532      void
1533      calc_maxGlobalDataPoint(int& ProcNo,  int& DataPointNo) const;
1534    
1535    
1536    /**    /**
1537       \brief       \brief
1538       Check *this and the right operand are compatible. Throws       Check *this and the right operand are compatible. Throws
# Line 1370  class Data { Line 1604  class Data {
1604       \brief       \brief
1605       Construct a Data object of the appropriate type.       Construct a Data object of the appropriate type.
1606    */    */
1607    template <class IValueType>  
1608    void    void
1609    initialise(const IValueType& value,    initialise(const DataTypes::ValueType& value,
1610             const DataTypes::ShapeType& shape,
1611               const FunctionSpace& what,               const FunctionSpace& what,
1612               bool expanded);               bool expanded);
1613    
1614      void
1615      initialise(const WrappedArray& value,
1616                     const FunctionSpace& what,
1617                     bool expanded);
1618    
1619    //    //
1620    // flag to protect the data object against any update    // flag to protect the data object against any update
1621    bool m_protected;    bool m_protected;
1622      mutable bool m_shared;
1623      bool m_lazy;
1624    
1625    //    //
1626    // pointer to the actual data object    // pointer to the actual data object
1627    boost::shared_ptr<DataAbstract> m_data;  //   boost::shared_ptr<DataAbstract> m_data;
1628      DataAbstract_ptr m_data;
1629    
1630  };  // If possible please use getReadyPtr instead.
1631    // But see warning below.
1632      const DataReady*
1633      getReady() const;
1634    
1635  template <class IValueType>    DataReady*
1636  void    getReady();
1637  Data::initialise(const IValueType& value,  
1638                   const FunctionSpace& what,  
1639                   bool expanded)  // Be wary of using this for local operations since it (temporarily) increases reference count.
1640  {  // If you are just using this to call a method on DataReady instead of DataAbstract consider using
1641    // getReady() instead
1642      DataReady_ptr
1643      getReadyPtr();
1644    
1645      const_DataReady_ptr
1646      getReadyPtr() const;
1647    
1648    
1649      /**
1650       \brief Update the Data's shared flag
1651       This indicates that the DataAbstract used by this object is now shared (or no longer shared).
1652       For internal use only.
1653      */
1654      void updateShareStatus(bool nowshared) const
1655      {
1656        m_shared=nowshared;     // m_shared is mutable
1657      }
1658    
1659      // In the isShared() method below:
1660      // A problem would occur if m_data (the address pointed to) were being modified
1661      // while the call m_data->is_shared is being executed.
1662      //
1663      // Q: So why do I think this code can be thread safe/correct?
1664      // A: We need to make some assumptions.
1665      // 1. We assume it is acceptable to return true under some conditions when we aren't shared.
1666      // 2. We assume that no constructions or assignments which will share previously unshared
1667      //    will occur while this call is executing. This is consistent with the way Data:: and C are written.
1668    //    //
1669    // 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
1670    // 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.
1671    // 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.
1672    // within the shared_ptr constructor.    // For any threads executing before the flag switches they will assume the object is still shared.
1673    if (expanded) {    bool isShared() const
1674      DataAbstract* temp=new DataExpanded(value,what);    {
1675      boost::shared_ptr<DataAbstract> temp_data(temp);      return m_shared;
1676      m_data=temp_data;  /*  if (m_shared) return true;
1677    } else {      if (m_data->isShared())        
1678      DataAbstract* temp=new DataConstant(value,what);      {                  
1679      boost::shared_ptr<DataAbstract> temp_data(temp);          updateShareStatus(true);
1680      m_data=temp_data;          return true;
1681        }
1682        return false;*/
1683      }
1684    
1685      void forceResolve()
1686      {
1687        if (isLazy())
1688        {
1689            #ifdef _OPENMP
1690            if (omp_in_parallel())
1691            {   // Yes this is throwing an exception out of an omp thread which is forbidden.
1692            throw DataException("Please do not call forceResolve() in a parallel region.");
1693            }
1694            #endif
1695            resolve();
1696        }
1697    }    }
1698    
1699      /**
1700      \brief if another object is sharing out member data make a copy to work with instead.
1701      This code should only be called from single threaded sections of code.
1702      */
1703      void exclusiveWrite()
1704      {
1705    #ifdef _OPENMP
1706      if (omp_in_parallel())
1707      {
1708    // *((int*)0)=17;
1709        throw DataException("Programming error. Please do not run exclusiveWrite() in multi-threaded sections.");
1710      }
1711    #endif
1712        forceResolve();
1713        if (isShared())
1714        {
1715            DataAbstract* t=m_data->deepCopy();
1716            set_m_data(DataAbstract_ptr(t));
1717        }
1718      }
1719    
1720      /**
1721      \brief checks if caller can have exclusive write to the object
1722      */
1723      void checkExclusiveWrite()
1724      {
1725        if  (isLazy() || isShared())
1726        {
1727            throw DataException("Programming error. ExclusiveWrite required - please call requireWrite()");
1728        }
1729      }
1730    
1731      /**
1732      \brief Modify the data abstract hosted by this Data object
1733      For internal use only.
1734      Passing a pointer to null is permitted (do this in the destructor)
1735      \warning Only to be called in single threaded code or inside a single/critical section. This method needs to be atomic.
1736      */
1737      void set_m_data(DataAbstract_ptr p);
1738    
1739      friend class DataAbstract;        // To allow calls to updateShareStatus
1740    
1741    };
1742    
1743    }   // end namespace escript
1744    
1745    
1746    // No, this is not supposed to be at the top of the file
1747    // DataAbstact needs to be declared first, then DataReady needs to be fully declared
1748    // so that I can dynamic cast between them below.
1749    #include "DataReady.h"
1750    #include "DataLazy.h"
1751    
1752    namespace escript
1753    {
1754    
1755    inline
1756    const DataReady*
1757    Data::getReady() const
1758    {
1759       const DataReady* dr=dynamic_cast<const DataReady*>(m_data.get());
1760       EsysAssert((dr!=0), "Error - casting to DataReady.");
1761       return dr;
1762  }  }
1763    
1764    inline
1765    DataReady*
1766    Data::getReady()
1767    {
1768       DataReady* dr=dynamic_cast<DataReady*>(m_data.get());
1769       EsysAssert((dr!=0), "Error - casting to DataReady.");
1770       return dr;
1771    }
1772    
1773    // Be wary of using this for local operations since it (temporarily) increases reference count.
1774    // If you are just using this to call a method on DataReady instead of DataAbstract consider using
1775    // getReady() instead
1776    inline
1777    DataReady_ptr
1778    Data::getReadyPtr()
1779    {
1780       DataReady_ptr dr=boost::dynamic_pointer_cast<DataReady>(m_data);
1781       EsysAssert((dr.get()!=0), "Error - casting to DataReady.");
1782       return dr;
1783    }
1784    
1785    
1786    inline
1787    const_DataReady_ptr
1788    Data::getReadyPtr() const
1789    {
1790       const_DataReady_ptr dr=boost::dynamic_pointer_cast<const DataReady>(m_data);
1791       EsysAssert((dr.get()!=0), "Error - casting to DataReady.");
1792       return dr;
1793    }
1794    
1795    inline
1796    DataAbstract::ValueType::value_type*
1797    Data::getSampleDataRW(DataAbstract::ValueType::size_type sampleNo)
1798    {
1799       if (isLazy())
1800       {
1801        throw DataException("Error, attempt to acquire RW access to lazy data. Please call requireWrite() first.");
1802       }
1803       return getReady()->getSampleDataRW(sampleNo);
1804    }
1805    
1806    inline
1807    const DataAbstract::ValueType::value_type*
1808    Data::getSampleDataRO(DataAbstract::ValueType::size_type sampleNo, BufferGroup* bufferg)
1809    {
1810       DataLazy* l=dynamic_cast<DataLazy*>(m_data.get());
1811       if (l!=0)
1812       {
1813        size_t offset=0;
1814        if (bufferg==NULL)
1815        {
1816            throw DataException("Error, attempt to getSampleDataRO for lazy Data with buffer==NULL");
1817        }
1818        const DataTypes::ValueType* res=l->resolveSample(*bufferg,sampleNo,offset);
1819        return &((*res)[offset]);
1820       }
1821       return getReady()->getSampleDataRO(sampleNo);
1822    }
1823    
1824    
1825    
1826    /**
1827       Modify a filename for MPI parallel output to multiple files
1828    */
1829    char *Escript_MPI_appendRankToFileName(const char *, int, int);
1830    
1831  /**  /**
1832     Binary Data object operators.     Binary Data object operators.
1833  */  */
# Line 1519  ESCRIPT_DLL_API std::ostream& operator<< Line 1939  ESCRIPT_DLL_API std::ostream& operator<<
1939  /**  /**
1940    \brief    \brief
1941    Compute a tensor product of two Data objects    Compute a tensor product of two Data objects
1942    \param arg0 - Input - Data object    \param arg_0 - Input - Data object
1943    \param arg1 - Input - Data object    \param arg_1 - Input - Data object
1944    \param axis_offset - Input - axis offset    \param axis_offset - Input - axis offset
1945    \param transpose - Input - 0: transpose neither, 1: transpose arg0, 2: transpose arg1    \param transpose - Input - 0: transpose neither, 1: transpose arg0, 2: transpose arg1
1946  */  */
1947  ESCRIPT_DLL_API  ESCRIPT_DLL_API
1948  Data  Data
1949  C_GeneralTensorProduct(Data& arg0,  C_GeneralTensorProduct(Data& arg_0,
1950                       Data& arg1,                       Data& arg_1,
1951                       int axis_offset=0,                       int axis_offset=0,
1952                       int transpose=0);                       int transpose=0);
1953    
   
 /**  
   \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);  
   
1954  /**  /**
1955    \brief    \brief
1956    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 1964  Data::binaryOp(const Data& right,
1964  {  {
1965     //     //
1966     // 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
1967     if (getPointDataView().getRank()==0 && right.getPointDataView().getRank()!=0) {     if (getDataPointRank()==0 && right.getDataPointRank()!=0) {
1968       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.");
1969     }     }
1970    
1971       if (isLazy() || right.isLazy())
1972       {
1973         throw DataException("Programmer error - attempt to call binaryOp with Lazy Data.");
1974       }
1975     //     //
1976     // initially make the temporary a shallow copy     // initially make the temporary a shallow copy
1977     Data tempRight(right);     Data tempRight(right);
1978    
1979     if (getFunctionSpace()!=right.getFunctionSpace()) {     if (getFunctionSpace()!=right.getFunctionSpace()) {
1980       if (right.probeInterpolation(getFunctionSpace())) {       if (right.probeInterpolation(getFunctionSpace())) {
1981         //         //
# Line 1582  Data::binaryOp(const Data& right, Line 1985  Data::binaryOp(const Data& right,
1985         //         //
1986         // interpolate onto the RHS function space         // interpolate onto the RHS function space
1987         Data tempLeft(*this,right.getFunctionSpace());         Data tempLeft(*this,right.getFunctionSpace());
1988         m_data=tempLeft.m_data;  //        m_data=tempLeft.m_data;
1989           set_m_data(tempLeft.m_data);
1990       }       }
1991     }     }
1992     operandCheck(tempRight);     operandCheck(tempRight);
# Line 1598  Data::binaryOp(const Data& right, Line 2002  Data::binaryOp(const Data& right,
2002       // of any data type       // of any data type
2003       DataExpanded* leftC=dynamic_cast<DataExpanded*>(m_data.get());       DataExpanded* leftC=dynamic_cast<DataExpanded*>(m_data.get());
2004       EsysAssert((leftC!=0), "Programming error - casting to DataExpanded.");       EsysAssert((leftC!=0), "Programming error - casting to DataExpanded.");
2005       escript::binaryOp(*leftC,*(tempRight.m_data.get()),operation);       escript::binaryOp(*leftC,*(tempRight.getReady()),operation);
2006     } else if (isTagged()) {     } else if (isTagged()) {
2007       //       //
2008       // 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 2050  Data::algorithm(BinaryFunction operation
2050      DataConstant* leftC=dynamic_cast<DataConstant*>(m_data.get());      DataConstant* leftC=dynamic_cast<DataConstant*>(m_data.get());
2051      EsysAssert((leftC!=0), "Programming error - casting to DataConstant.");      EsysAssert((leftC!=0), "Programming error - casting to DataConstant.");
2052      return escript::algorithm(*leftC,operation,initial_value);      return escript::algorithm(*leftC,operation,initial_value);
2053      } else if (isEmpty()) {
2054        throw DataException("Error - Operations not permitted on instances of DataEmpty.");
2055      } else if (isLazy()) {
2056        throw DataException("Error - Operations not permitted on instances of DataLazy.");
2057      } else {
2058        throw DataException("Error - Data encapsulates an unknown type.");
2059    }    }
   return 0;  
2060  }  }
2061    
2062  /**  /**
# Line 1663  inline Line 2072  inline
2072  Data  Data
2073  Data::dp_algorithm(BinaryFunction operation, double initial_value) const  Data::dp_algorithm(BinaryFunction operation, double initial_value) const
2074  {  {
2075    if (isExpanded()) {    if (isEmpty()) {
2076      Data result(0,DataArrayView::ShapeType(),getFunctionSpace(),isExpanded());      throw DataException("Error - Operations not permitted on instances of DataEmpty.");
2077      }
2078      else if (isExpanded()) {
2079        Data result(0,DataTypes::ShapeType(),getFunctionSpace(),isExpanded());
2080      DataExpanded* dataE=dynamic_cast<DataExpanded*>(m_data.get());      DataExpanded* dataE=dynamic_cast<DataExpanded*>(m_data.get());
2081      DataExpanded* resultE=dynamic_cast<DataExpanded*>(result.m_data.get());      DataExpanded* resultE=dynamic_cast<DataExpanded*>(result.m_data.get());
2082      EsysAssert((dataE!=0), "Programming error - casting data to DataExpanded.");      EsysAssert((dataE!=0), "Programming error - casting data to DataExpanded.");
2083      EsysAssert((resultE!=0), "Programming error - casting result to DataExpanded.");      EsysAssert((resultE!=0), "Programming error - casting result to DataExpanded.");
2084      escript::dp_algorithm(*dataE,*resultE,operation,initial_value);      escript::dp_algorithm(*dataE,*resultE,operation,initial_value);
2085      return result;      return result;
2086    } else if (isTagged()) {    }
2087      else if (isTagged()) {
2088      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());  
2089      EsysAssert((dataT!=0), "Programming error - casting data to DataTagged.");      EsysAssert((dataT!=0), "Programming error - casting data to DataTagged.");
2090      EsysAssert((resultT!=0), "Programming error - casting result to DataTagged.");      DataTypes::ValueType defval(1);
2091        defval[0]=0;
2092        DataTagged* resultT=new DataTagged(getFunctionSpace(), DataTypes::scalarShape, defval, dataT);
2093      escript::dp_algorithm(*dataT,*resultT,operation,initial_value);      escript::dp_algorithm(*dataT,*resultT,operation,initial_value);
2094      return result;      return Data(resultT);   // note: the Data object now owns the resultT pointer
2095    } else if (isConstant()) {    }
2096      Data result(0,DataArrayView::ShapeType(),getFunctionSpace(),isExpanded());    else if (isConstant()) {
2097        Data result(0,DataTypes::ShapeType(),getFunctionSpace(),isExpanded());
2098      DataConstant* dataC=dynamic_cast<DataConstant*>(m_data.get());      DataConstant* dataC=dynamic_cast<DataConstant*>(m_data.get());
2099      DataConstant* resultC=dynamic_cast<DataConstant*>(result.m_data.get());      DataConstant* resultC=dynamic_cast<DataConstant*>(result.m_data.get());
2100      EsysAssert((dataC!=0), "Programming error - casting data to DataConstant.");      EsysAssert((dataC!=0), "Programming error - casting data to DataConstant.");
2101      EsysAssert((resultC!=0), "Programming error - casting result to DataConstant.");      EsysAssert((resultC!=0), "Programming error - casting result to DataConstant.");
2102      escript::dp_algorithm(*dataC,*resultC,operation,initial_value);      escript::dp_algorithm(*dataC,*resultC,operation,initial_value);
2103      return result;      return result;
2104      } else if (isLazy()) {
2105        throw DataException("Error - Operations not permitted on instances of DataLazy.");
2106      } else {
2107        throw DataException("Error - Data encapsulates an unknown type.");
2108    }    }
   Data falseRetVal; // to keep compiler quiet  
   return falseRetVal;  
2109  }  }
2110    
2111    /**
2112      \brief
2113      Compute a tensor operation with two Data objects
2114      \param arg_0 - Input - Data object
2115      \param arg_1 - Input - Data object
2116      \param operation - Input - Binary op functor
2117    */
2118  template <typename BinaryFunction>  template <typename BinaryFunction>
2119    inline
2120  Data  Data
2121  C_TensorBinaryOperation(Data const &arg_0,  C_TensorBinaryOperation(Data const &arg_0,
2122                          Data const &arg_1,                          Data const &arg_1,
2123                          BinaryFunction operation)                          BinaryFunction operation)
2124  {  {
2125      if (arg_0.isEmpty() || arg_1.isEmpty())
2126      {
2127         throw DataException("Error - Operations not permitted on instances of DataEmpty.");
2128      }
2129      if (arg_0.isLazy() || arg_1.isLazy())
2130      {
2131         throw DataException("Error - Operations not permitted on lazy data.");
2132      }
2133    // Interpolate if necessary and find an appropriate function space    // Interpolate if necessary and find an appropriate function space
2134    Data arg_0_Z, arg_1_Z;    Data arg_0_Z, arg_1_Z;
2135    if (arg_0.getFunctionSpace()!=arg_1.getFunctionSpace()) {    if (arg_0.getFunctionSpace()!=arg_1.getFunctionSpace()) {
# Line 1730  C_TensorBinaryOperation(Data const &arg_ Line 2151  C_TensorBinaryOperation(Data const &arg_
2151    // Get rank and shape of inputs    // Get rank and shape of inputs
2152    int rank0 = arg_0_Z.getDataPointRank();    int rank0 = arg_0_Z.getDataPointRank();
2153    int rank1 = arg_1_Z.getDataPointRank();    int rank1 = arg_1_Z.getDataPointRank();
2154    DataArrayView::ShapeType shape0 = arg_0_Z.getDataPointShape();    DataTypes::ShapeType shape0 = arg_0_Z.getDataPointShape();
2155    DataArrayView::ShapeType shape1 = arg_1_Z.getDataPointShape();    DataTypes::ShapeType shape1 = arg_1_Z.getDataPointShape();
2156    int size0 = arg_0_Z.getDataPointSize();    int size0 = arg_0_Z.getDataPointSize();
2157    int size1 = arg_1_Z.getDataPointSize();    int size1 = arg_1_Z.getDataPointSize();
   
2158    // Declare output Data object    // Declare output Data object
2159    Data res;    Data res;
2160    
2161    if (shape0 == shape1) {    if (shape0 == shape1) {
   
2162      if (arg_0_Z.isConstant()   && arg_1_Z.isConstant()) {      if (arg_0_Z.isConstant()   && arg_1_Z.isConstant()) {
2163        res = Data(0.0, shape0, arg_1_Z.getFunctionSpace());      // DataConstant output        res = Data(0.0, shape0, arg_1_Z.getFunctionSpace());      // DataConstant output
2164        double *ptr_0 = &((arg_0_Z.getPointDataView().getData())[0]);        const double *ptr_0 = &(arg_0_Z.getDataAtOffsetRO(0));
2165        double *ptr_1 = &((arg_1_Z.getPointDataView().getData())[0]);        const double *ptr_1 = &(arg_1_Z.getDataAtOffsetRO(0));
2166        double *ptr_2 = &((res.getPointDataView().getData())[0]);        double *ptr_2 = &(res.getDataAtOffsetRW(0));
2167    
2168        tensor_binary_operation(size0, ptr_0, ptr_1, ptr_2, operation);        tensor_binary_operation(size0, ptr_0, ptr_1, ptr_2, operation);
2169      }      }
2170      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 2182  C_TensorBinaryOperation(Data const &arg_
2182    
2183        // Prepare offset into DataConstant        // Prepare offset into DataConstant
2184        int offset_0 = tmp_0->getPointOffset(0,0);        int offset_0 = tmp_0->getPointOffset(0,0);
2185        double *ptr_0 = &((arg_0_Z.getPointDataView().getData())[offset_0]);        const double *ptr_0 = &(arg_0_Z.getDataAtOffsetRO(offset_0));
2186        // Get the views  
       DataArrayView view_1 = tmp_1->getDefaultValue();  
       DataArrayView view_2 = tmp_2->getDefaultValue();  
2187        // Get the pointers to the actual data        // Get the pointers to the actual data
2188        double *ptr_1 = &((view_1.getData())[0]);        const double *ptr_1 = &(tmp_1->getDefaultValueRO(0));
2189        double *ptr_2 = &((view_2.getData())[0]);        double *ptr_2 = &(tmp_2->getDefaultValueRW(0));
2190    
2191        // Compute a result for the default        // Compute a result for the default
2192        tensor_binary_operation(size0, ptr_0, ptr_1, ptr_2, operation);        tensor_binary_operation(size0, ptr_0, ptr_1, ptr_2, operation);
2193        // Compute a result for each tag        // Compute a result for each tag
2194        const DataTagged::DataMapType& lookup_1=tmp_1->getTagLookup();        const DataTagged::DataMapType& lookup_1=tmp_1->getTagLookup();
2195        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
2196        for (i=lookup_1.begin();i!=lookup_1.end();i++) {        for (i=lookup_1.begin();i!=lookup_1.end();i++) {
2197          tmp_2->addTaggedValue(i->first,tmp_2->getDefaultValue());          tmp_2->addTag(i->first);
2198          DataArrayView view_1 = tmp_1->getDataPointByTag(i->first);          const double *ptr_1 = &(tmp_1->getDataByTagRO(i->first,0));
2199          DataArrayView view_2 = tmp_2->getDataPointByTag(i->first);          double *ptr_2 = &(tmp_2->getDataByTagRW(i->first,0));
2200          double *ptr_1 = &view_1.getData(0);  
         double *ptr_2 = &view_2.getData(0);  
2201          tensor_binary_operation(size0, ptr_0, ptr_1, ptr_2, operation);          tensor_binary_operation(size0, ptr_0, ptr_1, ptr_2, operation);
2202        }        }
2203    
2204      }      }
2205      else if (arg_0_Z.isConstant()   && arg_1_Z.isExpanded()) {      else if (arg_0_Z.isConstant()   && arg_1_Z.isExpanded()) {
   
2206        res = Data(0.0, shape0, arg_1_Z.getFunctionSpace(),true); // DataExpanded output        res = Data(0.0, shape0, arg_1_Z.getFunctionSpace(),true); // DataExpanded output
2207        DataConstant* tmp_0=dynamic_cast<DataConstant*>(arg_0_Z.borrowData());        DataConstant* tmp_0=dynamic_cast<DataConstant*>(arg_0_Z.borrowData());
2208        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 2212  C_TensorBinaryOperation(Data const &arg_
2212        int numSamples_1 = arg_1_Z.getNumSamples();        int numSamples_1 = arg_1_Z.getNumSamples();
2213        int numDataPointsPerSample_1 = arg_1_Z.getNumDataPointsPerSample();        int numDataPointsPerSample_1 = arg_1_Z.getNumDataPointsPerSample();
2214        int offset_0 = tmp_0->getPointOffset(0,0);        int offset_0 = tmp_0->getPointOffset(0,0);
2215          res.requireWrite();
2216        #pragma omp parallel for private(sampleNo_1,dataPointNo_1) schedule(static)        #pragma omp parallel for private(sampleNo_1,dataPointNo_1) schedule(static)
2217        for (sampleNo_1 = 0; sampleNo_1 < numSamples_1; sampleNo_1++) {        for (sampleNo_1 = 0; sampleNo_1 < numSamples_1; sampleNo_1++) {
2218          for (dataPointNo_1 = 0; dataPointNo_1 < numDataPointsPerSample_1; dataPointNo_1++) {          for (dataPointNo_1 = 0; dataPointNo_1 < numDataPointsPerSample_1; dataPointNo_1++) {
2219            int offset_1 = tmp_1->getPointOffset(sampleNo_1,dataPointNo_1);            int offset_1 = tmp_1->getPointOffset(sampleNo_1,dataPointNo_1);
2220            int offset_2 = tmp_2->getPointOffset(sampleNo_1,dataPointNo_1);            int offset_2 = tmp_2->getPointOffset(sampleNo_1,dataPointNo_1);
2221            double *ptr_0 = &((arg_0_Z.getPointDataView().getData())[offset_0]);            const double *ptr_0 = &(arg_0_Z.getDataAtOffsetRO(offset_0));
2222            double *ptr_1 = &((arg_1_Z.getPointDataView().getData())[offset_1]);            const double *ptr_1 = &(arg_1_Z.getDataAtOffsetRO(offset_1));
2223            double *ptr_2 = &((res.getPointDataView().getData())[offset_2]);            double *ptr_2 = &(res.getDataAtOffsetRW(offset_2));
2224            tensor_binary_operation(size0, ptr_0, ptr_1, ptr_2, operation);            tensor_binary_operation(size0, ptr_0, ptr_1, ptr_2, operation);
2225          }          }
2226        }        }
2227    
2228      }      }
2229      else if (arg_0_Z.isTagged()     && arg_1_Z.isConstant()) {      else if (arg_0_Z.isTagged()     && arg_1_Z.isConstant()) {
   
2230        // Borrow DataTagged input from Data object        // Borrow DataTagged input from Data object
2231        DataTagged* tmp_0=dynamic_cast<DataTagged*>(arg_0_Z.borrowData());        DataTagged* tmp_0=dynamic_cast<DataTagged*>(arg_0_Z.borrowData());
2232    
# Line 1823  C_TensorBinaryOperation(Data const &arg_ Line 2240  C_TensorBinaryOperation(Data const &arg_
2240    
2241        // Prepare offset into DataConstant        // Prepare offset into DataConstant
2242        int offset_1 = tmp_1->getPointOffset(0,0);        int offset_1 = tmp_1->getPointOffset(0,0);
2243        double *ptr_1 = &((arg_1_Z.getPointDataView().getData())[offset_1]);  
2244        // 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();  
2245        // Get the pointers to the actual data        // Get the pointers to the actual data
2246        double *ptr_0 = &((view_0.getData())[0]);        const double *ptr_0 = &(tmp_0->getDefaultValueRO(0));
2247        double *ptr_2 = &((view_2.getData())[0]);        double *ptr_2 = &(tmp_2->getDefaultValueRW(0));
2248        // Compute a result for the default        // Compute a result for the default
2249        tensor_binary_operation(size0, ptr_0, ptr_1, ptr_2, operation);        tensor_binary_operation(size0, ptr_0, ptr_1, ptr_2, operation);
2250        // Compute a result for each tag        // Compute a result for each tag
2251        const DataTagged::DataMapType& lookup_0=tmp_0->getTagLookup();        const DataTagged::DataMapType& lookup_0=tmp_0->getTagLookup();
2252        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
2253        for (i=lookup_0.begin();i!=lookup_0.end();i++) {        for (i=lookup_0.begin();i!=lookup_0.end();i++) {
2254          tmp_2->addTaggedValue(i->first,tmp_2->getDefaultValue());          tmp_2->addTag(i->first);
2255          DataArrayView view_0 = tmp_0->getDataPointByTag(i->first);          const double *ptr_0 = &(tmp_0->getDataByTagRO(i->first,0));
2256          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);  
2257          tensor_binary_operation(size0, ptr_0, ptr_1, ptr_2, operation);          tensor_binary_operation(size0, ptr_0, ptr_1, ptr_2, operation);
2258        }        }
2259    
2260      }      }
2261      else if (arg_0_Z.isTagged()     && arg_1_Z.isTagged()) {      else if (arg_0_Z.isTagged()     && arg_1_Z.isTagged()) {
   
2262        // Borrow DataTagged input from Data object        // Borrow DataTagged input from Data object
2263        DataTagged* tmp_0=dynamic_cast<DataTagged*>(arg_0_Z.borrowData());        DataTagged* tmp_0=dynamic_cast<DataTagged*>(arg_0_Z.borrowData());
2264    
# Line 1858  C_TensorBinaryOperation(Data const &arg_ Line 2270  C_TensorBinaryOperation(Data const &arg_
2270        res.tag();        // DataTagged output        res.tag();        // DataTagged output
2271        DataTagged* tmp_2=dynamic_cast<DataTagged*>(res.borrowData());        DataTagged* tmp_2=dynamic_cast<DataTagged*>(res.borrowData());
2272    
       // Get the views  
       DataArrayView view_0 = tmp_0->getDefaultValue();  
       DataArrayView view_1 = tmp_1->getDefaultValue();  
       DataArrayView view_2 = tmp_2->getDefaultValue();  
2273        // Get the pointers to the actual data        // Get the pointers to the actual data
2274        double *ptr_0 = &((view_0.getData())[0]);        const double *ptr_0 = &(tmp_0->getDefaultValueRO(0));
2275        double *ptr_1 = &((view_1.getData())[0]);        const double *ptr_1 = &(tmp_1->getDefaultValueRO(0));
2276        double *ptr_2 = &((view_2.getData())[0]);        double *ptr_2 = &(tmp_2->getDefaultValueRW(0));
2277    
2278        // Compute a result for the default        // Compute a result for the default
2279        tensor_binary_operation(size0, ptr_0, ptr_1, ptr_2, operation);        tensor_binary_operation(size0, ptr_0, ptr_1, ptr_2, operation);
2280        // Merge the tags        // Merge the tags
# Line 1873  C_TensorBinaryOperation(Data const &arg_ Line 2282  C_TensorBinaryOperation(Data const &arg_
2282        const DataTagged::DataMapType& lookup_0=tmp_0->getTagLookup();        const DataTagged::DataMapType& lookup_0=tmp_0->getTagLookup();
2283        const DataTagged::DataMapType& lookup_1=tmp_1->getTagLookup();        const DataTagged::DataMapType& lookup_1=tmp_1->getTagLookup();
2284        for (i=lookup_0.begin();i!=lookup_0.end();i++) {        for (i=lookup_0.begin();i!=lookup_0.end();i++) {
2285          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
2286        }        }
2287        for (i=lookup_1.begin();i!=lookup_1.end();i++) {        for (i=lookup_1.begin();i!=lookup_1.end();i++) {
2288          tmp_2->addTaggedValue(i->first,tmp_2->getDefaultValue());          tmp_2->addTag(i->first);
2289        }        }
2290        // Compute a result for each tag        // Compute a result for each tag
2291        const DataTagged::DataMapType& lookup_2=tmp_2->getTagLookup();        const DataTagged::DataMapType& lookup_2=tmp_2->getTagLookup();
2292        for (i=lookup_2.begin();i!=lookup_2.end();i++) {        for (i=lookup_2.begin();i!=lookup_2.end();i++) {
2293          DataArrayView view_0 = tmp_0->getDataPointByTag(i->first);  
2294          DataArrayView view_1 = tmp_1->getDataPointByTag(i->first);          const double *ptr_0 = &(tmp_0->getDataByTagRO(i->first,0));
2295          DataArrayView view_2 = tmp_2->getDataPointByTag(i->first);          const double *ptr_1 = &(tmp_1->getDataByTagRO(i->first,0));
2296          double *ptr_0 = &view_0.getData(0);          double *ptr_2 = &(tmp_2->getDataByTagRW(i->first,0));
2297          double *ptr_1 = &view_1.getData(0);  
         double *ptr_2 = &view_2.getData(0);  
2298          tensor_binary_operation(size0, ptr_0, ptr_1, ptr_2, operation);          tensor_binary_operation(size0, ptr_0, ptr_1, ptr_2, operation);
2299        }        }
2300    
2301      }      }
2302      else if (arg_0_Z.isTagged()     && arg_1_Z.isExpanded()) {      else if (arg_0_Z.isTagged()     && arg_1_Z.isExpanded()) {
   
2303        // 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
2304        res = Data(0.0, shape0, arg_1_Z.getFunctionSpace(),true); // DataExpanded output        res = Data(0.0, shape0, arg_1_Z.getFunctionSpace(),true); // DataExpanded output
2305        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 2309  C_TensorBinaryOperation(Data const &arg_
2309        int sampleNo_0,dataPointNo_0;        int sampleNo_0,dataPointNo_0;
2310        int numSamples_0 = arg_0_Z.getNumSamples();        int numSamples_0 = arg_0_Z.getNumSamples();
2311        int numDataPointsPerSample_0 = arg_0_Z.getNumDataPointsPerSample();        int numDataPointsPerSample_0 = arg_0_Z.getNumDataPointsPerSample();
2312          res.requireWrite();
2313        #pragma omp parallel for private(sampleNo_0,dataPointNo_0) schedule(static)        #pragma omp parallel for private(sampleNo_0,dataPointNo_0) schedule(static)
2314        for (sampleNo_0 = 0; sampleNo_0 < numSamples_0; sampleNo_0++) {        for (sampleNo_0 = 0; sampleNo_0 < numSamples_0; sampleNo_0++) {
2315          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
2316          double *ptr_0 = &((arg_0_Z.getPointDataView().getData())[offset_0]);          const double *ptr_0 = &(arg_0_Z.getDataAtOffsetRO(offset_0));
2317          for (dataPointNo_0 = 0; dataPointNo_0 < numDataPointsPerSample_0; dataPointNo_0++) {          for (dataPointNo_0 = 0; dataPointNo_0 < numDataPointsPerSample_0; dataPointNo_0++) {
2318            int offset_1 = tmp_1->getPointOffset(sampleNo_0,dataPointNo_0);            int offset_1 = tmp_1->getPointOffset(sampleNo_0,dataPointNo_0);
2319            int offset_2 = tmp_2->getPointOffset(sampleNo_0,dataPointNo_0);            int offset_2 = tmp_2->getPointOffset(sampleNo_0,dataPointNo_0);
2320            double *ptr_1 = &((arg_1_Z.getPointDataView().getData())[offset_1]);            const double *ptr_1 = &(arg_1_Z.getDataAtOffsetRO(offset_1));
2321            double *ptr_2 = &((res.getPointDataView().getData())[offset_2]);            double *ptr_2 = &(res.getDataAtOffsetRW(offset_2));
2322            tensor_binary_operation(size0, ptr_0, ptr_1, ptr_2, operation);            tensor_binary_operation(size0, ptr_0, ptr_1, ptr_2, operation);
2323          }          }
2324        }        }
2325    
2326      }      }
2327      else if (arg_0_Z.isExpanded()   && arg_1_Z.isConstant()) {      else if (arg_0_Z.isExpanded()   && arg_1_Z.isConstant()) {
   
2328        res = Data(0.0, shape0, arg_1_Z.getFunctionSpace(),true); // DataExpanded output        res = Data(0.0, shape0, arg_1_Z.getFunctionSpace(),true); // DataExpanded output
2329        DataExpanded* tmp_0=dynamic_cast<DataExpanded*>(arg_0_Z.borrowData());        DataExpanded* tmp_0=dynamic_cast<DataExpanded*>(arg_0_Z.borrowData());
2330        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 2334  C_TensorBinaryOperation(Data const &arg_
2334        int numSamples_0 = arg_0_Z.getNumSamples();        int numSamples_0 = arg_0_Z.getNumSamples();
2335        int numDataPointsPerSample_0 = arg_0_Z.getNumDataPointsPerSample();        int numDataPointsPerSample_0 = arg_0_Z.getNumDataPointsPerSample();
2336        int offset_1 = tmp_1->getPointOffset(0,0);        int offset_1 = tmp_1->getPointOffset(0,0);
2337          res.requireWrite();
2338        #pragma omp parallel for private(sampleNo_0,dataPointNo_0) schedule(static)        #pragma omp parallel for private(sampleNo_0,dataPointNo_0) schedule(static)
2339        for (sampleNo_0 = 0; sampleNo_0 < numSamples_0; sampleNo_0++) {        for (sampleNo_0 = 0; sampleNo_0 < numSamples_0; sampleNo_0++) {
2340          for (dataPointNo_0 = 0; dataPointNo_0 < numDataPointsPerSample_0; dataPointNo_0++) {          for (dataPointNo_0 = 0; dataPointNo_0 < numDataPointsPerSample_0; dataPointNo_0++) {
2341            int offset_0 = tmp_0->getPointOffset(sampleNo_0,dataPointNo_0);            int offset_0 = tmp_0->getPointOffset(sampleNo_0,dataPointNo_0);
2342            int offset_2 = tmp_2->getPointOffset(sampleNo_0,dataPointNo_0);            int offset_2 = tmp_2->getPointOffset(sampleNo_0,dataPointNo_0);
2343            double *ptr_0 = &((arg_0_Z.getPointDataView().getData())[offset_0]);  
2344            double *ptr_1 = &((arg_1_Z.getPointDataView().getData())[offset_1]);            const double *ptr_0 = &(arg_0_Z.getDataAtOffsetRO(offset_0));
2345            double *ptr_2 = &((res.getPointDataView().getData())[offset_2]);            const double *ptr_1 = &(arg_1_Z.getDataAtOffsetRO(offset_1));
2346              double *ptr_2 = &(res.getDataAtOffsetRW(offset_2));
2347    
2348    
2349            tensor_binary_operation(size0, ptr_0, ptr_1, ptr_2, operation);            tensor_binary_operation(size0, ptr_0, ptr_1, ptr_2, operation);
2350          }          }
2351        }        }
2352    
2353      }      }
2354      else if (arg_0_Z.isExpanded()   && arg_1_Z.isTagged()) {      else if (arg_0_Z.isExpanded()   && arg_1_Z.isTagged()) {
   
2355        // 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
2356        res = Data(0.0, shape0, arg_1_Z.getFunctionSpace(),true); // DataExpanded output        res = Data(0.0, shape0, arg_1_Z.getFunctionSpace(),true); // DataExpanded output
2357        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 2361  C_TensorBinaryOperation(Data const &arg_
2361        int sampleNo_0,dataPointNo_0;        int sampleNo_0,dataPointNo_0;
2362        int numSamples_0 = arg_0_Z.getNumSamples();        int numSamples_0 = arg_0_Z.getNumSamples();
2363        int numDataPointsPerSample_0 = arg_0_Z.getNumDataPointsPerSample();        int numDataPointsPerSample_0 = arg_0_Z.getNumDataPointsPerSample();
2364          res.requireWrite();
2365        #pragma omp parallel for private(sampleNo_0,dataPointNo_0) schedule(static)        #pragma omp parallel for private(sampleNo_0,dataPointNo_0) schedule(static)
2366        for (sampleNo_0 = 0; sampleNo_0 < numSamples_0; sampleNo_0++) {        for (sampleNo_0 = 0; sampleNo_0 < numSamples_0; sampleNo_0++) {
2367          int offset_1 = tmp_1->getPointOffset(sampleNo_0,0);          int offset_1 = tmp_1->getPointOffset(sampleNo_0,0);
2368          double *ptr_1 = &((arg_1_Z.getPointDataView().getData())[offset_1]);          const double *ptr_1 = &(arg_1_Z.getDataAtOffsetRO(offset_1));
2369          for (dataPointNo_0 = 0; dataPointNo_0 < numDataPointsPerSample_0; dataPointNo_0++) {          for (dataPointNo_0 = 0; dataPointNo_0 < numDataPointsPerSample_0; dataPointNo_0++) {
2370            int offset_0 = tmp_0->getPointOffset(sampleNo_0,dataPointNo_0);            int offset_0 = tmp_0->getPointOffset(sampleNo_0,dataPointNo_0);
2371            int offset_2 = tmp_2->getPointOffset(sampleNo_0,dataPointNo_0);            int offset_2 = tmp_2->getPointOffset(sampleNo_0,dataPointNo_0);
2372            double *ptr_0 = &((arg_0_Z.getPointDataView().getData())[offset_0]);            const double *ptr_0 = &(arg_0_Z.getDataAtOffsetRO(offset_0));
2373            double *ptr_2 = &((res.getPointDataView().getData())[offset_2]);            double *ptr_2 = &(res.getDataAtOffsetRW(offset_2));
2374            tensor_binary_operation(size0, ptr_0, ptr_1, ptr_2, operation);            tensor_binary_operation(size0, ptr_0, ptr_1, ptr_2, operation);
2375          }          }
2376        }        }
2377    
2378      }      }
2379      else if (arg_0_Z.isExpanded()   && arg_1_Z.isExpanded()) {      else if (arg_0_Z.isExpanded()   && arg_1_Z.isExpanded()) {
   
2380        // 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
2381        res = Data(0.0, shape0, arg_1_Z.getFunctionSpace(),true); // DataExpanded output        res = Data(0.0, shape0, arg_1_Z.getFunctionSpace(),true); // DataExpanded output
2382        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 2386  C_TensorBinaryOperation(Data const &arg_
2386        int sampleNo_0,dataPointNo_0;        int sampleNo_0,dataPointNo_0;
2387        int numSamples_0 = arg_0_Z.getNumSamples();        int numSamples_0 = arg_0_Z.getNumSamples();
2388        int numDataPointsPerSample_0 = arg_0_Z.getNumDataPointsPerSample();        int numDataPointsPerSample_0 = arg_0_Z.getNumDataPointsPerSample();
2389          res.requireWrite();
2390        #pragma omp parallel for private(sampleNo_0,dataPointNo_0) schedule(static)        #pragma omp parallel for private(sampleNo_0,dataPointNo_0) schedule(static)
2391        for (sampleNo_0 = 0; sampleNo_0 < numSamples_0; sampleNo_0++) {        for (sampleNo_0 = 0; sampleNo_0 < numSamples_0; sampleNo_0++) {
2392          for (dataPointNo_0 = 0; dataPointNo_0 < numDataPointsPerSample_0; dataPointNo_0++) {          for (dataPointNo_0 = 0; dataPointNo_0 < numDataPointsPerSample_0; dataPointNo_0++) {
2393            int offset_0 = tmp_0->getPointOffset(sampleNo_0,dataPointNo_0);            int offset_0 = tmp_0->getPointOffset(sampleNo_0,dataPointNo_0);
2394            int offset_1 = tmp_1->getPointOffset(sampleNo_0,dataPointNo_0);            int offset_1 = tmp_1->getPointOffset(sampleNo_0,dataPointNo_0);
2395            int offset_2 = tmp_2->getPointOffset(sampleNo_0,dataPointNo_0);            int offset_2 = tmp_2->getPointOffset(sampleNo_0,dataPointNo_0);
2396            double *ptr_0 = &((arg_0_Z.getPointDataView().getData())[offset_0]);            const double *ptr_0 = &(arg_0_Z.getDataAtOffsetRO(offset_0));
2397            double *ptr_1 = &((arg_1_Z.getPointDataView().getData())[offset_1]);            const double *ptr_1 = &(arg_1_Z.getDataAtOffsetRO(offset_1));
2398            double *ptr_2 = &((res.getPointDataView().getData())[offset_2]);            double *ptr_2 = &(res.getDataAtOffsetRW(offset_2));
2399            tensor_binary_operation(size0, ptr_0, ptr_1, ptr_2, operation);            tensor_binary_operation(size0, ptr_0, ptr_1, ptr_2, operation);
2400          }          }
2401        }        }
# Line 1995  C_TensorBinaryOperation(Data const &arg_ Line 2406  C_TensorBinaryOperation(Data const &arg_
2406      }      }
2407    
2408    } else if (0 == rank0) {    } else if (0 == rank0) {
   
2409      if (arg_0_Z.isConstant()   && arg_1_Z.isConstant()) {      if (arg_0_Z.isConstant()   && arg_1_Z.isConstant()) {
2410        res = Data(0.0, shape1, arg_1_Z.getFunctionSpace());      // DataConstant output        res = Data(0.0, shape1, arg_1_Z.getFunctionSpace());      // DataConstant output
2411        double *ptr_0 = &((arg_0_Z.getPointDataView().getData())[0]);        const double *ptr_0 = &(arg_0_Z.getDataAtOffsetRO(0));
2412        double *ptr_1 = &((arg_1_Z.getPointDataView().getData())[0]);        const double *ptr_1 = &(arg_1_Z.getDataAtOffsetRO(0));
2413        double *ptr_2 = &((res.getPointDataView().getData())[0]);        double *ptr_2 = &(res.getDataAtOffsetRW(0));
2414        tensor_binary_operation(size1, ptr_0[0], ptr_1, ptr_2, operation);        tensor_binary_operation(size1, ptr_0[0], ptr_1, ptr_2, operation);
2415      }      }
2416      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 2428  C_TensorBinaryOperation(Data const &arg_
2428    
2429        // Prepare offset into DataConstant        // Prepare offset into DataConstant
2430        int offset_0 = tmp_0->getPointOffset(0,0);        int offset_0 = tmp_0->getPointOffset(0,0);
2431        double *ptr_0 = &((arg_0_Z.getPointDataView().getData())[offset_0]);        const double *ptr_0 = &(arg_0_Z.getDataAtOffsetRO(offset_0));
2432        // Get the views  
2433        DataArrayView view_1 = tmp_1->getDefaultValue();        const double *ptr_1 = &(tmp_1->getDefaultValueRO(0));
2434        DataArrayView view_2 = tmp_2->getDefaultValue();        double *ptr_2 = &(tmp_2->getDefaultValueRW(0));
2435        // Get the pointers to the actual data  
       double *ptr_1 = &((view_1.getData())[0]);  
       double *ptr_2 = &((view_2.getData())[0]);  
2436        // Compute a result for the default        // Compute a result for the default
2437        tensor_binary_operation(size1, ptr_0[0], ptr_1, ptr_2, operation);        tensor_binary_operation(size1, ptr_0[0], ptr_1, ptr_2, operation);
2438        // Compute a result for each tag        // Compute a result for each tag
2439        const DataTagged::DataMapType& lookup_1=tmp_1->getTagLookup();        const DataTagged::DataMapType& lookup_1=tmp_1->getTagLookup();
2440        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
2441        for (i=lookup_1.begin();i!=lookup_1.end();i++) {        for (i=lookup_1.begin();i!=lookup_1.end();i++) {
2442          tmp_2->addTaggedValue(i->first,tmp_2->getDefaultValue());          tmp_2->addTag(i->first);
2443          DataArrayView view_1 = tmp_1->getDataPointByTag(i->first);          const double *ptr_1 = &(tmp_1->getDataByTagRO(i->first,0));
2444          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);  
2445          tensor_binary_operation(size1, ptr_0[0], ptr_1, ptr_2, operation);          tensor_binary_operation(size1, ptr_0[0], ptr_1, ptr_2, operation);
2446        }        }
2447    
# Line 2051  C_TensorBinaryOperation(Data const &arg_ Line 2457  C_TensorBinaryOperation(Data const &arg_
2457        int numSamples_1 = arg_1_Z.getNumSamples();        int numSamples_1 = arg_1_Z.getNumSamples();
2458        int numDataPointsPerSample_1 = arg_1_Z.getNumDataPointsPerSample();        int numDataPointsPerSample_1 = arg_1_Z.getNumDataPointsPerSample();
2459        int offset_0 = tmp_0->getPointOffset(0,0);        int offset_0 = tmp_0->getPointOffset(0,0);
2460          res.requireWrite();
2461        #pragma omp parallel for private(sampleNo_1,dataPointNo_1) schedule(static)        #pragma omp parallel for private(sampleNo_1,dataPointNo_1) schedule(static)
2462        for (sampleNo_1 = 0; sampleNo_1 < numSamples_1; sampleNo_1++) {        for (sampleNo_1 = 0; sampleNo_1 < numSamples_1; sampleNo_1++) {
2463          for (dataPointNo_1 = 0; dataPointNo_1 < numDataPointsPerSample_1; dataPointNo_1++) {          for (dataPointNo_1 = 0; dataPointNo_1 < numDataPointsPerSample_1; dataPointNo_1++) {
2464            int offset_1 = tmp_1->getPointOffset(sampleNo_1,dataPointNo_1);            int offset_1 = tmp_1->getPointOffset(sampleNo_1,dataPointNo_1);
2465            int offset_2 = tmp_2->getPointOffset(sampleNo_1,dataPointNo_1);            int offset_2 = tmp_2->getPointOffset(sampleNo_1,dataPointNo_1);
2466            double *ptr_0 = &((arg_0_Z.getPointDataView().getData())[offset_0]);            const double *ptr_0 = &(arg_0_Z.getDataAtOffsetRO(offset_0));
2467            double *ptr_1 = &((arg_1_Z.getPointDataView().getData())[offset_1]);            const double *ptr_1 = &(arg_1_Z.getDataAtOffsetRO(offset_1));
2468            double *ptr_2 = &((res.getPointDataView().getData())[offset_2]);            double *ptr_2 = &(res.getDataAtOffsetRW(offset_2));
2469            tensor_binary_operation(size1, ptr_0[0], ptr_1, ptr_2, operation);            tensor_binary_operation(size1, ptr_0[0], ptr_1, ptr_2, operation);
2470    
2471          }          }
# Line 2080  C_TensorBinaryOperation(Data const &arg_ Line 2487  C_TensorBinaryOperation(Data const &arg_
2487    
2488        // Prepare offset into DataConstant        // Prepare offset into DataConstant
2489        int offset_1 = tmp_1->getPointOffset(0,0);        int offset_1 = tmp_1->getPointOffset(0,0);
2490        double *ptr_1 = &((arg_1_Z.getPointDataView().getData())[offset_1]);        const double *ptr_1 = &(arg_1_Z.getDataAtOffsetRO(offset_1));
2491        // Get the views  
       DataArrayView view_0 = tmp_0->getDefaultValue();  
       DataArrayView view_2 = tmp_2->getDefaultValue();  
2492        // Get the pointers to the actual data        // Get the pointers to the actual data
2493        double *ptr_0 = &((view_0.getData())[0]);        const double *ptr_0 = &(tmp_0->getDefaultValueRO(0));
2494        double *ptr_2 = &((view_2.getData())[0]);        double *ptr_2 = &(tmp_2->getDefaultValueRW(0));
2495    
2496    
2497        // Compute a result for the default        // Compute a result for the default
2498        tensor_binary_operation(size1, ptr_0[0], ptr_1, ptr_2, operation);        tensor_binary_operation(size1, ptr_0[0], ptr_1, ptr_2, operation);
2499        // Compute a result for each tag        // Compute a result for each tag
2500        const DataTagged::DataMapType& lookup_0=tmp_0->getTagLookup();        const DataTagged::DataMapType& lookup_0=tmp_0->getTagLookup();
2501        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
2502        for (i=lookup_0.begin();i!=lookup_0.end();i++) {        for (i=lookup_0.begin();i!=lookup_0.end();i++) {
2503          tmp_2->addTaggedValue(i->first,tmp_2->getDefaultValue());          tmp_2->addTag(i->first);
2504          DataArrayView view_0 = tmp_0->getDataPointByTag(i->first);          const double *ptr_0 = &(tmp_0->getDataByTagRO(i->first,0));
2505          DataArrayView view_2 = tmp_2->getDataPointByTag(i->first);          double *ptr_2 = &(tmp_2->getDataByTagRW(i->first,0));
2506          double *ptr_0 = &view_0.getData(0);  
         double *ptr_2 = &view_2.getData(0);  
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        }        }
2509    
# Line 2115  C_TensorBinaryOperation(Data const &arg_ Line 2521  C_TensorBinaryOperation(Data const &arg_
2521        res.tag();        // DataTagged output        res.tag();        // DataTagged output
2522        DataTagged* tmp_2=dynamic_cast<DataTagged*>(res.borrowData());        DataTagged* tmp_2=dynamic_cast<DataTagged*>(res.borrowData());
2523    
       // Get the views  
       DataArrayView view_0 = tmp_0->getDefaultValue();  
       DataArrayView view_1 = tmp_1->getDefaultValue();  
       DataArrayView view_2 = tmp_2->getDefaultValue();  
2524        // Get the pointers to the actual data        // Get the pointers to the actual data
2525        double *ptr_0 = &((view_0.getData())[0]);        const double *ptr_0 = &(tmp_0->getDefaultValueRO(0));
2526        double *ptr_1 = &((view_1.getData())[0]);        const double *ptr_1 = &(tmp_1->getDefaultValueRO(0));
2527        double *ptr_2 = &((view_2.getData())[0]);        double *ptr_2 = &(tmp_2->getDefaultValueRW(0));
2528    
2529        // Compute a result for the default        // Compute a result for the default
2530        tensor_binary_operation(size1, ptr_0[0], ptr_1, ptr_2, operation);        tensor_binary_operation(size1, ptr_0[0], ptr_1, ptr_2, operation);
2531        // Merge the tags        // Merge the tags
# Line 2130  C_TensorBinaryOperation(Data const &arg_ Line 2533  C_TensorBinaryOperation(Data const &arg_
2533        const DataTagged::DataMapType& lookup_0=tmp_0->getTagLookup();        const DataTagged::DataMapType& lookup_0=tmp_0->getTagLookup();
2534        const DataTagged::DataMapType& lookup_1=tmp_1->getTagLookup();        const DataTagged::DataMapType& lookup_1=tmp_1->getTagLookup();
2535        for (i=lookup_0.begin();i!=lookup_0.end();i++) {        for (i=lookup_0.begin();i!=lookup_0.end();i++) {
2536          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
2537        }        }
2538        for (i=lookup_1.begin();i!=lookup_1.end();i++) {        for (i=lookup_1.begin();i!=lookup_1.end();i++) {
2539          tmp_2->addTaggedValue(i->first,tmp_2->getDefaultValue());          tmp_2->addTag(i->first);
2540        }        }
2541        // Compute a result for each tag        // Compute a result for each tag
2542        const DataTagged::DataMapType& lookup_2=tmp_2->getTagLookup();        const DataTagged::DataMapType& lookup_2=tmp_2->getTagLookup();
2543        for (i=lookup_2.begin();i!=lookup_2.end();i++) {        for (i=lookup_2.begin();i!=lookup_2.end();i++) {
2544          DataArrayView view_0 = tmp_0->getDataPointByTag(i->first);          const double *ptr_0 = &(tmp_0->getDataByTagRO(i->first,0));
2545          DataArrayView view_1 = tmp_1->getDataPointByTag(i->first);          const double *ptr_1 = &(tmp_1->getDataByTagRO(i->first,0));
2546          DataArrayView view_2 = tmp_2->getDataPointByTag(i->first);          double *ptr_2 = &(tmp_2->getDataByTagRW(i->first,0));
2547          double *ptr_0 = &view_0.getData(0);  
         double *ptr_1 = &view_1.getData(0);  
         double *ptr_2 = &view_2.getData(0);  
2548          tensor_binary_operation(size1, ptr_0[0], ptr_1, ptr_2, operation);          tensor_binary_operation(size1, ptr_0[0], ptr_1, ptr_2, operation);
2549        }        }
2550    
# Line 2159  C_TensorBinaryOperation(Data const &arg_ Line 2560  C_TensorBinaryOperation(Data const &arg_
2560        int sampleNo_0,dataPointNo_0;        int sampleNo_0,dataPointNo_0;
2561        int numSamples_0 = arg_0_Z.getNumSamples();        int numSamples_0 = arg_0_Z.getNumSamples();
2562        int numDataPointsPerSample_0 = arg_0_Z.getNumDataPointsPerSample();        int numDataPointsPerSample_0 = arg_0_Z.getNumDataPointsPerSample();
2563          res.requireWrite();
2564        #pragma omp parallel for private(sampleNo_0,dataPointNo_0) schedule(static)        #pragma omp parallel for private(sampleNo_0,dataPointNo_0) schedule(static)
2565        for (sampleNo_0 = 0; sampleNo_0 < numSamples_0; sampleNo_0++) {        for (sampleNo_0 = 0; sampleNo_0 < numSamples_0; sampleNo_0++) {
2566          int offset_0 = tmp_0->getPointOffset(sampleNo_0,0); // They're all the same, so just use #0          int offset_0 = tmp_0->getPointOffset(sampleNo_0,0); // They're all the same, so just use #0
2567          double *ptr_0 = &((arg_0_Z.getPointDataView().getData())[offset_0]);          const double *ptr_0 = &(arg_0_Z.getDataAtOffsetRO(offset_0));
2568          for (dataPointNo_0 = 0; dataPointNo_0 < numDataPointsPerSample_0; dataPointNo_0++) {          for (dataPointNo_0 = 0; dataPointNo_0 < numDataPointsPerSample_0; dataPointNo_0++) {
2569            int offset_1 = tmp_1->getPointOffset(sampleNo_0,dataPointNo_0);            int offset_1 = tmp_1->getPointOffset(sampleNo_0,dataPointNo_0);
2570            int offset_2 = tmp_2->getPointOffset(sampleNo_0,dataPointNo_0);            int offset_2 = tmp_2->getPointOffset(sampleNo_0,dataPointNo_0);
2571            double *ptr_1 = &((arg_1_Z.getPointDataView().getData())[offset_1]);            const double *ptr_1 = &(arg_1_Z.getDataAtOffsetRO(offset_1));
2572            double *ptr_2 = &((res.getPointDataView().getData())[offset_2]);            double *ptr_2 = &(res.getDataAtOffsetRW(offset_2));
2573            tensor_binary_operation(size1, ptr_0[0], ptr_1, ptr_2, operation);            tensor_binary_operation(size1, ptr_0[0], ptr_1, ptr_2, operation);
2574          }          }
2575        }        }
2576    
2577      }      }
2578      else if (arg_0_Z.isExpanded()   && arg_1_Z.isConstant()) {      else if (arg_0_Z.isExpanded()   && arg_1_Z.isConstant()) {
   
2579        res = Data(0.0, shape1, arg_1_Z.getFunctionSpace(),true); // DataExpanded output        res = Data(0.0, shape1, arg_1_Z.getFunctionSpace(),true); // DataExpanded output
2580        DataExpanded* tmp_0=dynamic_cast<DataExpanded*>(arg_0_Z.borrowData());        DataExpanded* tmp_0=dynamic_cast<DataExpanded*>(arg_0_Z.borrowData());
2581        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 2585  C_TensorBinaryOperation(Data const &arg_
2585        int numSamples_0 = arg_0_Z.getNumSamples();        int numSamples_0 = arg_0_Z.getNumSamples();
2586        int numDataPointsPerSample_0 = arg_0_Z.getNumDataPointsPerSample();        int numDataPointsPerSample_0 = arg_0_Z.getNumDataPointsPerSample();
2587        int offset_1 = tmp_1->getPointOffset(0,0);        int offset_1 = tmp_1->getPointOffset(0,0);
2588          res.requireWrite();
2589        #pragma omp parallel for private(sampleNo_0,dataPointNo_0) schedule(static)        #pragma omp parallel for private(sampleNo_0,dataPointNo_0) schedule(static)
2590        for (sampleNo_0 = 0; sampleNo_0 < numSamples_0; sampleNo_0++) {        for (sampleNo_0 = 0; sampleNo_0 < numSamples_0; sampleNo_0++) {
2591          for (dataPointNo_0 = 0; dataPointNo_0 < numDataPointsPerSample_0; dataPointNo_0++) {          for (dataPointNo_0 = 0; dataPointNo_0 < numDataPointsPerSample_0; dataPointNo_0++) {
2592            int offset_0 = tmp_0->getPointOffset(sampleNo_0,dataPointNo_0);            int offset_0 = tmp_0->getPointOffset(sampleNo_0,dataPointNo_0);
2593            int offset_2 = tmp_2->getPointOffset(sampleNo_0,dataPointNo_0);            int offset_2 = tmp_2->getPointOffset(sampleNo_0,dataPointNo_0);
2594            double *ptr_0 = &((arg_0_Z.getPointDataView().getData())[offset_0]);            const double *ptr_0 = &(arg_0_Z.getDataAtOffsetRO(offset_0));
2595            double *ptr_1 = &((arg_1_Z.getPointDataView().getData())[offset_1]);            const double *ptr_1 = &(arg_1_Z.getDataAtOffsetRO(offset_1));
2596            double *ptr_2 = &((res.getPointDataView().getData())[offset_2]);            double *ptr_2 = &(res.getDataAtOffsetRW(offset_2));
2597            tensor_binary_operation(size1, ptr_0[0], ptr_1, ptr_2, operation);            tensor_binary_operation(size1, ptr_0[0], ptr_1, ptr_2, operation);
2598          }          }
2599        }        }
# Line 2209  C_TensorBinaryOperation(Data const &arg_ Line 2611  C_TensorBinaryOperation(Data const &arg_
2611        int sampleNo_0,dataPointNo_0;        int sampleNo_0,dataPointNo_0;
2612        int numSamples_0 = arg_0_Z.getNumSamples();        int numSamples_0 = arg_0_Z.getNumSamples();
2613        int numDataPointsPerSample_0 = arg_0_Z.getNumDataPointsPerSample();        int numDataPointsPerSample_0 = arg_0_Z.getNumDataPointsPerSample();
2614          res.requireWrite();
2615        #pragma omp parallel for private(sampleNo_0,dataPointNo_0) schedule(static)        #pragma omp parallel for private(sampleNo_0,dataPointNo_0) schedule(static)
2616        for (sampleNo_0 = 0; sampleNo_0 < numSamples_0; sampleNo_0++) {        for (sampleNo_0 = 0; sampleNo_0 < numSamples_0; sampleNo_0++) {
2617          int offset_1 = tmp_1->getPointOffset(sampleNo_0,0);          int offset_1 = tmp_1->getPointOffset(sampleNo_0,0);
2618          double *ptr_1 = &((arg_1_Z.getPointDataView().getData())[offset_1]);          const double *ptr_1 = &(arg_1_Z.getDataAtOffsetRO(offset_1));
2619          for (dataPointNo_0 = 0; dataPointNo_0 < numDataPointsPerSample_0; dataPointNo_0++) {          for (dataPointNo_0 = 0; dataPointNo_0 < numDataPointsPerSample_0; dataPointNo_0++) {
2620            int offset_0 = tmp_0->getPointOffset(sampleNo_0,dataPointNo_0);            int offset_0 = tmp_0->getPointOffset(sampleNo_0,dataPointNo_0);
2621            int offset_2 = tmp_2->getPointOffset(sampleNo_0,dataPointNo_0);            int offset_2 = tmp_2->getPointOffset(sampleNo_0,dataPointNo_0);
2622            double *ptr_0 = &((arg_0_Z.getPointDataView().getData())[offset_0]);            const double *ptr_0 = &(arg_0_Z.getDataAtOffsetRO(offset_0));
2623            double *ptr_2 = &((res.getPointDataView().getData())[offset_2]);            double *ptr_2 = &(res.getDataAtOffsetRW(offset_2));
2624            tensor_binary_operation(size1, ptr_0[0], ptr_1, ptr_2, operation);            tensor_binary_operation(size1, ptr_0[0], ptr_1, ptr_2, operation);
2625          }          }
2626        }        }
# Line 2234  C_TensorBinaryOperation(Data const &arg_ Line 2637  C_TensorBinaryOperation(Data const &arg_
2637        int sampleNo_0,dataPointNo_0;        int sampleNo_0,dataPointNo_0;
2638        int numSamples_0 = arg_0_Z.getNumSamples();        int numSamples_0 = arg_0_Z.getNumSamples();
2639        int numDataPointsPerSample_0 = arg_0_Z.getNumDataPointsPerSample();        int numDataPointsPerSample_0 = arg_0_Z.getNumDataPointsPerSample();
2640          res.requireWrite();
2641        #pragma omp parallel for private(sampleNo_0,dataPointNo_0) schedule(static)        #pragma omp parallel for private(sampleNo_0,dataPointNo_0) schedule(static)
2642        for (sampleNo_0 = 0; sampleNo_0 < numSamples_0; sampleNo_0++) {        for (sampleNo_0 = 0; sampleNo_0 < numSamples_0; sampleNo_0++) {
2643          for (dataPointNo_0 = 0; dataPointNo_0 < numDataPointsPerSample_0; dataPointNo_0++) {          for (dataPointNo_0 = 0; dataPointNo_0 < numDataPointsPerSample_0; dataPointNo_0++) {
2644            int offset_0 = tmp_0->getPointOffset(sampleNo_0,dataPointNo_0);            int offset_0 = tmp_0->getPointOffset(sampleNo_0,dataPointNo_0);
2645            int offset_1 = tmp_1->getPointOffset(sampleNo_0,dataPointNo_0);            int offset_1 = tmp_1->getPointOffset(sampleNo_0,dataPointNo_0);
2646            int offset_2 = tmp_2->getPointOffset(sampleNo_0,dataPointNo_0);            int offset_2 = tmp_2->getPointOffset(sampleNo_0,dataPointNo_0);
2647            double *ptr_0 = &((arg_0_Z.getPointDataView().getData())[offset_0]);            const double *ptr_0 = &(arg_0_Z.getDataAtOffsetRO(offset_0));
2648            double *ptr_1 = &((arg_1_Z.getPointDataView().getData())[offset_1]);            const double *ptr_1 = &(arg_1_Z.getDataAtOffsetRO(offset_1));
2649            double *ptr_2 = &((res.getPointDataView().getData())[offset_2]);            double *ptr_2 = &(res.getDataAtOffsetRW(offset_2));
2650            tensor_binary_operation(size1, ptr_0[0], ptr_1, ptr_2, operation);            tensor_binary_operation(size1, ptr_0[0], ptr_1, ptr_2, operation);
2651          }          }
2652        }        }
# Line 2253  C_TensorBinaryOperation(Data const &arg_ Line 2657  C_TensorBinaryOperation(Data const &arg_
2657      }      }
2658    
2659    } else if (0 == rank1) {    } else if (0 == rank1) {
   
2660      if (arg_0_Z.isConstant()   && arg_1_Z.isConstant()) {      if (arg_0_Z.isConstant()   && arg_1_Z.isConstant()) {
2661        res = Data(0.0, shape0, arg_1_Z.getFunctionSpace());      // DataConstant output        res = Data(0.0, shape0, arg_1_Z.getFunctionSpace());      // DataConstant output
2662        double *ptr_0 = &((arg_0_Z.getPointDataView().getData())[0]);        const double *ptr_0 = &(arg_0_Z.getDataAtOffsetRO(0));
2663        double *ptr_1 = &((arg_1_Z.getPointDataView().getData())[0]);        const double *ptr_1 = &(arg_1_Z.getDataAtOffsetRO(0));
2664        double *ptr_2 = &((res.getPointDataView().getData())[0]);        double *ptr_2 = &(res.getDataAtOffsetRW(0));
2665        tensor_binary_operation(size0, ptr_0, ptr_1[0], ptr_2, operation);        tensor_binary_operation(size0, ptr_0, ptr_1[0], ptr_2, operation);
2666      }      }
2667      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 2679  C_TensorBinaryOperation(Data const &arg_
2679    
2680        // Prepare offset into DataConstant        // Prepare offset into DataConstant
2681        int offset_0 = tmp_0->getPointOffset(0,0);        int offset_0 = tmp_0->getPointOffset(0,0);
2682        double *ptr_0 = &((arg_0_Z.getPointDataView().getData())[offset_0]);        const double *ptr_0 = &(arg_0_Z.getDataAtOffsetRO(offset_0));
2683        // Get the views  
2684        DataArrayView view_1 = tmp_1->getDefaultValue();        //Get the pointers to the actual data
2685        DataArrayView view_2 = tmp_2->getDefaultValue();        const double *ptr_1 = &(tmp_1->getDefaultValueRO(0));
2686        // Get the pointers to the actual data        double *ptr_2 = &(tmp_2->getDefaultValueRW(0));
2687        double *ptr_1 = &((view_1.getData())[0]);  
       double *ptr_2 = &((view_2.getData())[0]);  
2688        // Compute a result for the default        // Compute a result for the default
2689        tensor_binary_operation(size0, ptr_0, ptr_1[0], ptr_2, operation);        tensor_binary_operation(size0, ptr_0, ptr_1[0], ptr_2, operation);
2690        // Compute a result for each tag        // Compute a result for each tag
2691        const DataTagged::DataMapType& lookup_1=tmp_1->getTagLookup();        const DataTagged::DataMapType& lookup_1=tmp_1->getTagLookup();
2692        DataTagged::DataMapType::const_iterator i; // i->first is a tag, i->second is an offset into memory        DataTagged::DataMapType::const_iterator i; // i->first is a tag, i->second is an offset into memory
2693        for (i=lookup_1.begin();i!=lookup_1.end();i++) {        for (i=lookup_1.begin();i!=lookup_1.end();i++) {
2694          tmp_2->addTaggedValue(i->first,tmp_2->getDefaultValue());          tmp_2->addTag(i->first);
2695          DataArrayView view_1 = tmp_1->getDataPointByTag(i->first);          const double *ptr_1 = &(tmp_1->getDataByTagRO(i->first,0));
2696          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);  
2697          tensor_binary_operation(size0, ptr_0, ptr_1[0], ptr_2, operation);          tensor_binary_operation(size0, ptr_0, ptr_1[0], ptr_2, operation);
2698        }        }
   
2699      }      }
2700      else if (arg_0_Z.isConstant()   && arg_1_Z.isExpanded()) {      else if (arg_0_Z.isConstant()   && arg_1_Z.isExpanded()) {
2701    
# Line 2309  C_TensorBinaryOperation(Data const &arg_ Line 2708  C_TensorBinaryOperation(Data const &arg_
2708        int numSamples_1 = arg_1_Z.getNumSamples();        int numSamples_1 = arg_1_Z.getNumSamples();
2709        int numDataPointsPerSample_1 = arg_1_Z.getNumDataPointsPerSample();        int numDataPointsPerSample_1 = arg_1_Z.getNumDataPointsPerSample();
2710        int offset_0 = tmp_0->getPointOffset(0,0);        int offset_0 = tmp_0->getPointOffset(0,0);
2711          res.requireWrite();
2712        #pragma omp parallel for private(sampleNo_1,dataPointNo_1) schedule(static)        #pragma omp parallel for private(sampleNo_1,dataPointNo_1) schedule(static)
2713        for (sampleNo_1 = 0; sampleNo_1 < numSamples_1; sampleNo_1++) {        for (sampleNo_1 = 0; sampleNo_1 < numSamples_1; sampleNo_1++) {
2714          for (dataPointNo_1 = 0; dataPointNo_1 < numDataPointsPerSample_1; dataPointNo_1++) {          for (dataPointNo_1 = 0; dataPointNo_1 < numDataPointsPerSample_1; dataPointNo_1++) {
2715            int offset_1 = tmp_1->getPointOffset(sampleNo_1,dataPointNo_1);            int offset_1 = tmp_1->getPointOffset(sampleNo_1,dataPointNo_1);
2716            int offset_2 = tmp_2->getPointOffset(sampleNo_1,dataPointNo_1);            int offset_2 = tmp_2->getPointOffset(sampleNo_1,dataPointNo_1);
2717            double *ptr_0 = &((arg_0_Z.getPointDataView().getData())[offset_0]);            const double *ptr_0 = &(arg_0_Z.getDataAtOffsetRO(offset_0));
2718            double *ptr_1 = &((arg_1_Z.getPointDataView().getData())[offset_1]);            const double *ptr_1 = &(arg_1_Z.getDataAtOffsetRO(offset_1));
2719            double *ptr_2 = &((res.getPointDataView().getData())[offset_2]);            double *ptr_2 = &(res.getDataAtOffsetRW(offset_2));
2720            tensor_binary_operation(size0, ptr_0, ptr_1[0], ptr_2, operation);            tensor_binary_operation(size0, ptr_0, ptr_1[0], ptr_2, operation);
2721          }          }
2722        }        }
# Line 2337  C_TensorBinaryOperation(Data const &arg_ Line 2737  C_TensorBinaryOperation(Data const &arg_
2737    
2738        // Prepare offset into DataConstant        // Prepare offset into DataConstant
2739        int offset_1 = tmp_1->getPointOffset(0,0);        int offset_1 = tmp_1->getPointOffset(0,0);
2740        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();  
2741        // Get the pointers to the actual data        // Get the pointers to the actual data
2742        double *ptr_0 = &((view_0.getData())[0]);        const double *ptr_0 = &(tmp_0->getDefaultValueRO(0));
2743        double *ptr_2 = &((view_2.getData())[0]);        double *ptr_2 = &(tmp_2->getDefaultValueRW(0));
2744        // Compute a result for the default        // Compute a result for the default
2745        tensor_binary_operation(size0, ptr_0, ptr_1[0], ptr_2, operation);        tensor_binary_operation(size0, ptr_0, ptr_1[0], ptr_2, operation);
2746        // Compute a result for each tag        // Compute a result for each tag
2747        const DataTagged::DataMapType& lookup_0=tmp_0->getTagLookup();        const DataTagged::DataMapType& lookup_0=tmp_0->getTagLookup();
2748        DataTagged::DataMapType::const_iterator i; // i->first is a tag, i->second is an offset into memory        DataTagged::DataMapType::const_iterator i; // i->first is a tag, i->second is an offset into memory
2749        for (i=lookup_0.begin();i!=lookup_0.end();i++) {        for (i=lookup_0.begin();i!=lookup_0.end();i++) {
2750          tmp_2->addTaggedValue(i->first,tmp_2->getDefaultValue());          tmp_2->addTag(i->first);
2751          DataArrayView view_0 = tmp_0->getDataPointByTag(i->first);          const double *ptr_0 = &(tmp_0->getDataByTagRO(i->first,0));
2752          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);  
2753          tensor_binary_operation(size0, ptr_0, ptr_1[0], ptr_2, operation);          tensor_binary_operation(size0, ptr_0, ptr_1[0], ptr_2, operation);
2754        }        }
2755    
# Line 2372  C_TensorBinaryOperation(Data const &arg_ Line 2767  C_TensorBinaryOperation(Data const &arg_
2767        res.tag();        // DataTagged output        res.tag();        // DataTagged output
2768        DataTagged* tmp_2=dynamic_cast<DataTagged*>(res.borrowData());        DataTagged* tmp_2=dynamic_cast<DataTagged*>(res.borrowData());
2769    
       // Get the views  
       DataArrayView view_0 = tmp_0->getDefaultValue();  
       DataArrayView view_1 = tmp_1->getDefaultValue();  
       DataArrayView view_2 = tmp_2->getDefaultValue();  
2770        // Get the pointers to the actual data        // Get the pointers to the actual data
2771        double *ptr_0 = &((view_0.getData())[0]);        const double *ptr_0 = &(tmp_0->getDefaultValueRO(0));
2772        double *ptr_1 = &((view_1.getData())[0]);        const double *ptr_1 = &(tmp_1->getDefaultValueRO(0));
2773        double *ptr_2 = &((view_2.getData())[0]);        double *ptr_2 = &(tmp_2->getDefaultValueRW(0));
2774    
2775        // Compute a result for the default        // Compute a result for the default
2776        tensor_binary_operation(size0, ptr_0, ptr_1[0], ptr_2, operation);        tensor_binary_operation(size0, ptr_0, ptr_1[0], ptr_2, operation);
2777        // Merge the tags        // Merge the tags
# Line 2387  C_TensorBinaryOperation(Data const &arg_ Line 2779  C_TensorBinaryOperation(Data const &arg_
2779        const DataTagged::DataMapType& lookup_0=tmp_0->getTagLookup();        const DataTagged::DataMapType& lookup_0=tmp_0->getTagLookup();
2780        const DataTagged::DataMapType& lookup_1=tmp_1->getTagLookup();        const DataTagged::DataMapType& lookup_1=tmp_1->getTagLookup();
2781        for (i=lookup_0.begin();i!=lookup_0.end();i++) {        for (i=lookup_0.begin();i!=lookup_0.end();i++) {
2782          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
2783        }        }
2784        for (i=lookup_1.begin();i!=lookup_1.end();i++) {        for (i=lookup_1.begin();i!=lookup_1.end();i++) {
2785          tmp_2->addTaggedValue(i->first,tmp_2->getDefaultValue());          tmp_2->addTag(i->first);
2786        }        }
2787        // Compute a result for each tag        // Compute a result for each tag
2788        const DataTagged::DataMapType& lookup_2=tmp_2->getTagLookup();        const DataTagged::DataMapType& lookup_2=tmp_2->getTagLookup();
2789        for (i=lookup_2.begin();i!=lookup_2.end();i++) {        for (i=lookup_2.begin();i!=lookup_2.end();i++) {
2790          DataArrayView view_0 = tmp_0->getDataPointByTag(i->first);          const double *ptr_0 = &(tmp_0->getDataByTagRO(i->first,0));
2791          DataArrayView view_1 = tmp_1->getDataPointByTag(i->first);          const double *ptr_1 = &(tmp_1->getDataByTagRO(i->first,0));
2792          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);  
2793          tensor_binary_operation(size0, ptr_0, ptr_1[0], ptr_2, operation);          tensor_binary_operation(size0, ptr_0, ptr_1[0], ptr_2, operation);
2794        }        }
2795    
# Line 2416  C_TensorBinaryOperation(Data const &arg_ Line 2805  C_TensorBinaryOperation(Data const &arg_
2805        int sampleNo_0,dataPointNo_0;        int sampleNo_0,dataPointNo_0;
2806        int numSamples_0 = arg_0_Z.getNumSamples();        int numSamples_0 = arg_0_Z.getNumSamples();
2807        int numDataPointsPerSample_0 = arg_0_Z.getNumDataPointsPerSample();        int numDataPointsPerSample_0 = arg_0_Z.getNumDataPointsPerSample();
2808          res.requireWrite();
2809        #pragma omp parallel for private(sampleNo_0,dataPointNo_0) schedule(static)        #pragma omp parallel for private(sampleNo_0,dataPointNo_0) schedule(static)
2810        for (sampleNo_0 = 0; sampleNo_0 < numSamples_0; sampleNo_0++) {        for (sampleNo_0 = 0; sampleNo_0 < numSamples_0; sampleNo_0++) {
2811          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
2812          double *ptr_0 = &((arg_0_Z.getPointDataView().getData())[offset_0]);          const double *ptr_0 = &(arg_0_Z.getDataAtOffsetRO(offset_0));
2813          for (dataPointNo_0 = 0; dataPointNo_0 < numDataPointsPerSample_0; dataPointNo_0++) {          for (dataPointNo_0 = 0; dataPointNo_0 < numDataPointsPerSample_0; dataPointNo_0++) {
2814            int offset_1 = tmp_1->getPointOffset(sampleNo_0,dataPointNo_0);            int offset_1 = tmp_1->getPointOffset(sampleNo_0,dataPointNo_0);
2815            int offset_2 = tmp_2->getPointOffset(sampleNo_0,dataPointNo_0);            int offset_2 = tmp_2->getPointOffset(sampleNo_0,dataPointNo_0);
2816            double *ptr_1 = &((arg_1_Z.getPointDataView().getData())[offset_1]);            const double *ptr_1 = &(arg_1_Z.getDataAtOffsetRO(offset_1));
2817            double *ptr_2 = &((res.getPointDataView().getData())[offset_2]);            double *ptr_2 = &(res.getDataAtOffsetRW(offset_2));
2818            tensor_binary_operation(size0, ptr_0, ptr_1[0], ptr_2, operation);            tensor_binary_operation(size0, ptr_0, ptr_1[0], ptr_2, operation);
2819          }          }
2820        }        }
2821    
2822      }      }
2823      else if (arg_0_Z.isExpanded()   && arg_1_Z.isConstant()) {      else if (arg_0_Z.isExpanded()   && arg_1_Z.isConstant()) {
   
2824        res = Data(0.0, shape0, arg_1_Z.getFunctionSpace(),true); // DataExpanded output        res = Data(0.0, shape0, arg_1_Z.getFunctionSpace(),true); // DataExpanded output
2825        DataExpanded* tmp_0=dynamic_cast<DataExpanded*>(arg_0_Z.borrowData());        DataExpanded* tmp_0=dynamic_cast<DataExpanded*>(arg_0_Z.borrowData());
2826        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 2830  C_TensorBinaryOperation(Data const &arg_
2830        int numSamples_0 = arg_0_Z.getNumSamples();        int numSamples_0 = arg_0_Z.getNumSamples();
2831        int numDataPointsPerSample_0 = arg_0_Z.getNumDataPointsPerSample();        int numDataPointsPerSample_0 = arg_0_Z.getNumDataPointsPerSample();
2832        int offset_1 = tmp_1->getPointOffset(0,0);        int offset_1 = tmp_1->getPointOffset(0,0);
2833          res.requireWrite();
2834        #pragma omp parallel for private(sampleNo_0,dataPointNo_0) schedule(static)        #pragma omp parallel for private(sampleNo_0,dataPointNo_0) schedule(static)
2835        for (sampleNo_0 = 0; sampleNo_0 < numSamples_0; sampleNo_0++) {        for (sampleNo_0 = 0; sampleNo_0 < numSamples_0; sampleNo_0++) {
2836          for (dataPointNo_0 = 0; dataPointNo_0 < numDataPointsPerSample_0; dataPointNo_0++) {          for (dataPointNo_0 = 0; dataPointNo_0 < numDataPointsPerSample_0; dataPointNo_0++) {
2837            int offset_0 = tmp_0->getPointOffset(sampleNo_0,dataPointNo_0);            int offset_0 = tmp_0->getPointOffset(sampleNo_0,dataPointNo_0);
2838            int offset_2 = tmp_2->getPointOffset(sampleNo_0,dataPointNo_0);            int offset_2 = tmp_2->getPointOffset(sampleNo_0,dataPointNo_0);
2839            double *ptr_0 = &((arg_0_Z.getPointDataView().getData())[offset_0]);            const double *ptr_0 = &(arg_0_Z.getDataAtOffsetRO(offset_0));
2840            double *ptr_1 = &((arg_1_Z.getPointDataView().getData())[offset_1]);            const double *ptr_1 = &(arg_1_Z.getDataAtOffsetRO(offset_1));
2841            double *ptr_2 = &((res.getPointDataView().getData())[offset_2]);            double *ptr_2 = &(res.getDataAtOffsetRW(offset_2));
2842            tensor_binary_operation(size0, ptr_0, ptr_1[0], ptr_2, operation);            tensor_binary_operation(size0, ptr_0, ptr_1[0], ptr_2, operation);
2843          }          }
2844        }        }
# Line 2466  C_TensorBinaryOperation(Data const &arg_ Line 2856  C_TensorBinaryOperation(Data const &arg_
2856        int sampleNo_0,dataPointNo_0;        int sampleNo_0,dataPointNo_0;
2857        int numSamples_0 = arg_0_Z.getNumSamples();        int numSamples_0 = arg_0_Z.getNumSamples();
2858        int numDataPointsPerSample_0 = arg_0_Z.getNumDataPointsPerSample();        int numDataPointsPerSample_0 = arg_0_Z.getNumDataPointsPerSample();
2859          res.requireWrite();
2860        #pragma omp parallel for private(sampleNo_0,dataPointNo_0) schedule(static)        #pragma omp parallel for private(sampleNo_0,dataPointNo_0) schedule(static)
2861        for (sampleNo_0 = 0; sampleNo_0 < numSamples_0; sampleNo_0++) {        for (sampleNo_0 = 0; sampleNo_0 < numSamples_0; sampleNo_0++) {
2862          int offset_1 = tmp_1->getPointOffset(sampleNo_0,0);          int offset_1 = tmp_1->getPointOffset(sampleNo_0,0);
2863          double *ptr_1 = &((arg_1_Z.getPointDataView().getData())[offset_1]);          const double *ptr_1 = &(arg_1_Z.getDataAtOffsetRO(offset_1));
2864          for (dataPointNo_0 = 0; dataPointNo_0 < numDataPointsPerSample_0; dataPointNo_0++) {          for (dataPointNo_0 = 0; dataPointNo_0 < numDataPointsPerSample_0; dataPointNo_0++) {
2865            int offset_0 = tmp_0->getPointOffset(sampleNo_0,dataPointNo_0);            int offset_0 = tmp_0->getPointOffset(sampleNo_0,dataPointNo_0);
2866            int offset_2 = tmp_2->getPointOffset(sampleNo_0,dataPointNo_0);            int offset_2 = tmp_2->getPointOffset(sampleNo_0,dataPointNo_0);
2867            double *ptr_0 = &((arg_0_Z.getPointDataView().getData())[offset_0]);            const double *ptr_0 = &(arg_0_Z.getDataAtOffsetRO(offset_0));
2868            double *ptr_2 = &((res.getPointDataView().getData())[offset_2]);            double *ptr_2 = &(res.getDataAtOffsetRW(offset_2));
2869            tensor_binary_operation(size0, ptr_0, ptr_1[0], ptr_2, operation);            tensor_binary_operation(size0, ptr_0, ptr_1[0], ptr_2, operation);
2870          }          }
2871        }        }
# Line 2491  C_TensorBinaryOperation(Data const &arg_ Line 2882  C_TensorBinaryOperation(Data const &arg_
2882        int sampleNo_0,dataPointNo_0;        int sampleNo_0,dataPointNo_0;
2883        int numSamples_0 = arg_0_Z.getNumSamples();        int numSamples_0 = arg_0_Z.getNumSamples();
2884        int numDataPointsPerSample_0 = arg_0_Z.getNumDataPointsPerSample();        int numDataPointsPerSample_0 = arg_0_Z.getNumDataPointsPerSample();
2885          res.requireWrite();
2886        #pragma omp parallel for private(sampleNo_0,dataPointNo_0) schedule(static)        #pragma omp parallel for private(sampleNo_0,dataPointNo_0) schedule(static)
2887        for (sampleNo_0 = 0; sampleNo_0 < numSamples_0; sampleNo_0++) {        for (sampleNo_0 = 0; sampleNo_0 < numSamples_0; sampleNo_0++) {
2888          for (dataPointNo_0 = 0; dataPointNo_0 < numDataPointsPerSample_0; dataPointNo_0++) {          for (dataPointNo_0 = 0; dataPointNo_0 < numDataPointsPerSample_0; dataPointNo_0++) {
2889            int offset_0 = tmp_0->getPointOffset(sampleNo_0,dataPointNo_0);            int offset_0 = tmp_0->getPointOffset(sampleNo_0,dataPointNo_0);
2890            int offset_1 = tmp_1->getPointOffset(sampleNo_0,dataPointNo_0);            int offset_1 = tmp_1->getPointOffset(sampleNo_0,dataPointNo_0);
2891            int offset_2 = tmp_2->getPointOffset(sampleNo_0,dataPointNo_0);            int offset_2 = tmp_2->getPointOffset(sampleNo_0,dataPointNo_0);
2892            double *ptr_0 = &((arg_0_Z.getPointDataView().getData())[offset_0]);            const double *ptr_0 = &(arg_0_Z.getDataAtOffsetRO(offset_0));
2893            double *ptr_1 = &((arg_1_Z.getPointDataView().getData())[offset_1]);            const double *ptr_1 = &(arg_1_Z.getDataAtOffsetRO(offset_1));
2894            double *ptr_2 = &((res.getPointDataView().getData())[offset_2]);            double *ptr_2 = &(res.getDataAtOffsetRW(offset_2));
2895            tensor_binary_operation(size0, ptr_0, ptr_1[0], ptr_2, operation);            tensor_binary_operation(size0, ptr_0, ptr_1[0], ptr_2, operation);
2896          }          }
2897        }        }
# Line 2521  Data Line 2913  Data
2913  C_TensorUnaryOperation(Data const &arg_0,  C_TensorUnaryOperation(Data const &arg_0,
2914                         UnaryFunction operation)                         UnaryFunction operation)
2915  {  {
2916      if (arg_0.isEmpty())  // do this before we attempt to interpolate
2917      {
2918         throw DataException("Error - Operations not permitted on instances of DataEmpty.");
2919      }
2920      if (arg_0.isLazy())
2921      {
2922         throw DataException("Error - Operations not permitted on lazy data.");
2923      }
2924    // Interpolate if necessary and find an appropriate function space    // Interpolate if necessary and find an appropriate function space
2925    Data arg_0_Z = Data(arg_0);    Data arg_0_Z = Data(arg_0);
2926    
2927    // Get rank and shape of inputs    // Get rank and shape of inputs
2928    int rank0 = arg_0_Z.getDataPointRank();    const DataTypes::ShapeType& shape0 = arg_0_Z.getDataPointShape();
   DataArrayView::ShapeType shape0 = arg_0_Z.getDataPointShape();  
2929    int size0 = arg_0_Z.getDataPointSize();    int size0 = arg_0_Z.getDataPointSize();
2930    
2931    // Declare output Data object    // Declare output Data object
# Line 2534  C_TensorUnaryOperation(Data const &arg_0 Line 2933  C_TensorUnaryOperation(Data const &arg_0
2933    
2934    if (arg_0_Z.isConstant()) {    if (arg_0_Z.isConstant()) {
2935      res = Data(0.0, shape0, arg_0_Z.getFunctionSpace());      // DataConstant output      res = Data(0.0, shape0, arg_0_Z.getFunctionSpace());      // DataConstant output
2936      double *ptr_0 = &((arg_0_Z.getPointDataView().getData())[0]);      const double *ptr_0 = &(arg_0_Z.getDataAtOffsetRO(0));
2937      double *ptr_2 = &((res.getPointDataView().getData())[0]);      double *ptr_2 = &(res.getDataAtOffsetRW(0));
2938      tensor_unary_operation(size0, ptr_0, ptr_2, operation);      tensor_unary_operation(size0, ptr_0, ptr_2, operation);
2939    }    }
2940    else if (arg_0_Z.isTagged()) {    else if (arg_0_Z.isTagged()) {
# Line 2548  C_TensorUnaryOperation(Data const &arg_0 Line 2947  C_TensorUnaryOperation(Data const &arg_0
2947      res.tag();      res.tag();
2948      DataTagged* tmp_2=dynamic_cast<DataTagged*>(res.borrowData());      DataTagged* tmp_2=dynamic_cast<DataTagged*>(res.borrowData());
2949    
     // Get the views  
     DataArrayView view_0 = tmp_0->getDefaultValue();  
     DataArrayView view_2 = tmp_2->getDefaultValue();  
2950      // Get the pointers to the actual data      // Get the pointers to the actual data
2951      double *ptr_0 = &((view_0.getData())[0]);      const double *ptr_0 = &(tmp_0->getDefaultValueRO(0));
2952      double *ptr_2 = &((view_2.getData())[0]);      double *ptr_2 = &(tmp_2->getDefaultValueRW(0));
2953      // Compute a result for the default      // Compute a result for the default
2954      tensor_unary_operation(size0, ptr_0, ptr_2, operation);      tensor_unary_operation(size0, ptr_0, ptr_2, operation);
2955      // Compute a result for each tag      // Compute a result for each tag
2956      const DataTagged::DataMapType& lookup_0=tmp_0->getTagLookup();      const DataTagged::DataMapType& lookup_0=tmp_0->getTagLookup();
2957      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
2958      for (i=lookup_0.begin();i!=lookup_0.end();i++) {      for (i=lookup_0.begin();i!=lookup_0.end();i++) {
2959        tmp_2->addTaggedValue(i->first,tmp_2->getDefaultValue());        tmp_2->addTag(i->first);
2960        DataArrayView view_0 = tmp_0->getDataPointByTag(i->first);        const double *ptr_0 = &(tmp_0->getDataByTagRO(i->first,0));
2961        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);  
2962        tensor_unary_operation(size0, ptr_0, ptr_2, operation);        tensor_unary_operation(size0, ptr_0, ptr_2, operation);
2963      }      }
2964    
# Line 2583  C_TensorUnaryOperation(Data const &arg_0 Line 2977  C_TensorUnaryOperation(Data const &arg_0
2977        for (dataPointNo_0 = 0; dataPointNo_0 < numDataPointsPerSample_0; dataPointNo_0++) {        for (dataPointNo_0 = 0; dataPointNo_0 < numDataPointsPerSample_0; dataPointNo_0++) {
2978          int offset_0 = tmp_0->getPointOffset(sampleNo_0,dataPointNo_0);          int offset_0 = tmp_0->getPointOffset(sampleNo_0,dataPointNo_0);
2979          int offset_2 = tmp_2->getPointOffset(sampleNo_0,dataPointNo_0);          int offset_2 = tmp_2->getPointOffset(sampleNo_0,dataPointNo_0);
2980          double *ptr_0 = &((arg_0_Z.getPointDataView().getData())[offset_0]);          const double *ptr_0 = &(arg_0_Z.getDataAtOffsetRO(offset_0));
2981          double *ptr_2 = &((res.getPointDataView().getData())[offset_2]);          double *ptr_2 = &(res.getDataAtOffsetRW(offset_2));
2982          tensor_unary_operation(size0, ptr_0, ptr_2, operation);          tensor_unary_operation(size0, ptr_0, ptr_2, operation);
2983        }        }
2984      }      }
   
2985    }    }
2986    else {    else {
2987      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.2628

  ViewVC Help
Powered by ViewVC 1.1.26