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

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

Parent Directory Parent Directory | Revision Log Revision Log


Revision 1358 - (hide annotations)
Wed Dec 5 03:41:06 2007 UTC (11 years, 7 months ago) by gross
File MIME type: text/plain
File size: 84097 byte(s)
setTaggedValue works now for expanded data.
1 jgs 94
2 ksteube 1312 /* $Id$ */
3    
4     /*******************************************************
5     *
6     * Copyright 2003-2007 by ACceSS MNRF
7     * Copyright 2007 by University of Queensland
8     *
9     * http://esscc.uq.edu.au
10     * Primary Business: Queensland, Australia
11     * Licensed under the Open Software License version 3.0
12     * http://www.opensource.org/licenses/osl-3.0.php
13     *
14     *******************************************************/
15    
16 jgs 121 /** \file Data.h */
17 jgs 94
18     #ifndef DATA_H
19     #define DATA_H
20 woo409 757 #include "system_dep.h"
21 jgs 94
22 jgs 474 #include "DataAbstract.h"
23     #include "DataAlgorithm.h"
24     #include "FunctionSpace.h"
25     #include "BinaryOp.h"
26     #include "UnaryOp.h"
27     #include "DataException.h"
28 jgs 94
29     extern "C" {
30 jgs 474 #include "DataC.h"
31 ksteube 1312 /* #include "paso/Paso.h" doesn't belong in this file...causes trouble for BruceFactory.cpp */
32 jgs 94 }
33    
34 ksteube 1312 #include "esysmpi.h"
35 jgs 94 #include <string>
36     #include <algorithm>
37    
38     #include <boost/shared_ptr.hpp>
39     #include <boost/python/object.hpp>
40     #include <boost/python/tuple.hpp>
41     #include <boost/python/numeric.hpp>
42    
43 jgs 121 namespace escript {
44    
45     //
46 jgs 149 // Forward declaration for various implementations of Data.
47 jgs 121 class DataConstant;
48     class DataTagged;
49     class DataExpanded;
50    
51 jgs 94 /**
52     \brief
53 matt 1327 Data creates the appropriate Data object for the given construction
54     arguments.
55 jgs 94
56     Description:
57     Data is essentially a factory class which creates the appropriate Data
58     object for the given construction arguments. It retains control over
59     the object created for the lifetime of the object.
60     The type of Data object referred to may change during the lifetime of
61     the Data object.
62     */
63     class Data {
64    
65     public:
66    
67 jgs 110 // These typedefs allow function names to be cast to pointers
68     // to functions of the appropriate type when calling unaryOp etc.
69 jgs 94 typedef double (*UnaryDFunPtr)(double);
70     typedef double (*BinaryDFunPtr)(double,double);
71    
72 gross 854
73 jgs 94 /**
74 jgs 102 Constructors.
75     */
76    
77     /**
78 jgs 94 \brief
79     Default constructor.
80     Creates a DataEmpty object.
81     */
82 woo409 757 ESCRIPT_DLL_API
83 jgs 94 Data();
84    
85     /**
86     \brief
87     Copy constructor.
88     WARNING: Only performs a shallow copy.
89     */
90 woo409 757 ESCRIPT_DLL_API
91 jgs 94 Data(const Data& inData);
92    
93     /**
94     \brief
95     Constructor from another Data object. If "what" is different from the
96 jgs 102 function space of inData the inData are tried to be interpolated to what,
97 jgs 94 otherwise a shallow copy of inData is returned.
98     */
99 woo409 757 ESCRIPT_DLL_API
100 jgs 94 Data(const Data& inData,
101     const FunctionSpace& what);
102    
103     /**
104     \brief
105     Constructor which copies data from a DataArrayView.
106    
107     \param value - Input - Data value for a single point.
108     \param what - Input - A description of what this data represents.
109     \param expanded - Input - Flag, if true fill the entire container with
110     the value. Otherwise a more efficient storage
111     mechanism will be used.
112     */
113 woo409 757 ESCRIPT_DLL_API
114 jgs 94 Data(const DataArrayView& value,
115     const FunctionSpace& what=FunctionSpace(),
116     bool expanded=false);
117    
118     /**
119     \brief
120     Constructor which creates a Data from a DataArrayView shape.
121    
122     \param value - Input - Single value applied to all Data.
123     \param dataPointShape - Input - The shape of each data point.
124     \param what - Input - A description of what this data represents.
125     \param expanded - Input - Flag, if true fill the entire container with
126     the given value. Otherwise a more efficient storage
127     mechanism will be used.
128     */
129 woo409 757 ESCRIPT_DLL_API
130 jgs 94 Data(double value,
131     const DataArrayView::ShapeType& dataPointShape=DataArrayView::ShapeType(),
132     const FunctionSpace& what=FunctionSpace(),
133     bool expanded=false);
134    
135     /**
136     \brief
137     Constructor which performs a deep copy of a region from another Data object.
138    
139     \param inData - Input - Input Data object.
140     \param region - Input - Region to copy.
141     */
142 woo409 757 ESCRIPT_DLL_API
143 jgs 94 Data(const Data& inData,
144     const DataArrayView::RegionType& region);
145    
146     /**
147     \brief
148     Constructor which will create Tagged data if expanded is false.
149     No attempt is made to ensure the tag keys match the tag keys
150     within the function space.
151    
152     \param tagKeys - Input - List of tag values.
153     \param values - Input - List of values, one for each tag.
154     \param defaultValue - Input - A default value, used if tag doesn't exist.
155     \param what - Input - A description of what this data represents.
156     \param expanded - Input - Flag, if true fill the entire container with
157     the appropriate values.
158 jgs 544 ==>*
159 jgs 94 */
160 woo409 757 ESCRIPT_DLL_API
161 jgs 94 Data(const DataTagged::TagListType& tagKeys,
162     const DataTagged::ValueListType& values,
163     const DataArrayView& defaultValue,
164     const FunctionSpace& what=FunctionSpace(),
165     bool expanded=false);
166    
167     /**
168     \brief
169     Constructor which copies data from a python numarray.
170    
171     \param value - Input - Data value for a single point.
172     \param what - Input - A description of what this data represents.
173     \param expanded - Input - Flag, if true fill the entire container with
174     the value. Otherwise a more efficient storage
175     mechanism will be used.
176     */
177 woo409 757 ESCRIPT_DLL_API
178 jgs 94 Data(const boost::python::numeric::array& value,
179     const FunctionSpace& what=FunctionSpace(),
180     bool expanded=false);
181    
182     /**
183     \brief
184     Constructor which copies data from any object that can be converted into
185     a python numarray.
186    
187     \param value - Input - Input data.
188     \param what - Input - A description of what this data represents.
189     \param expanded - Input - Flag, if true fill the entire container with
190     the value. Otherwise a more efficient storage
191     mechanism will be used.
192     */
193 woo409 757 ESCRIPT_DLL_API
194 jgs 94 Data(const boost::python::object& value,
195     const FunctionSpace& what=FunctionSpace(),
196     bool expanded=false);
197    
198     /**
199     \brief
200     Constructor which creates a DataConstant.
201     Copies data from any object that can be converted
202     into a numarray. All other parameters are copied from other.
203    
204     \param value - Input - Input data.
205     \param other - Input - contains all other parameters.
206     */
207 woo409 757 ESCRIPT_DLL_API
208 jgs 94 Data(const boost::python::object& value,
209     const Data& other);
210    
211     /**
212     \brief
213     Constructor which creates a DataConstant of "shape" with constant value.
214     */
215 woo409 757 ESCRIPT_DLL_API
216 matt 1327 Data(double value,
217     const boost::python::tuple& shape=boost::python::make_tuple(),
218 jgs 94 const FunctionSpace& what=FunctionSpace(),
219     bool expanded=false);
220 jgs 151 /**
221     \brief
222     Destructor
223     */
224 woo409 757 ESCRIPT_DLL_API
225 jgs 151 ~Data();
226 jgs 94
227     /**
228     \brief
229 jgs 102 Perform a deep copy.
230 jgs 94 */
231 woo409 757 ESCRIPT_DLL_API
232 jgs 94 void
233 jgs 102 copy(const Data& other);
234 jgs 94
235     /**
236 jgs 102 Member access methods.
237 jgs 94 */
238    
239     /**
240     \brief
241 matt 1327 switches on update protection
242 gross 783
243     */
244     ESCRIPT_DLL_API
245     void
246     setProtection();
247    
248     /**
249     \brief
250     Returns trueif the data object is protected against update
251    
252     */
253     ESCRIPT_DLL_API
254     bool
255     isProtected() const;
256 gross 1034
257 gross 783 /**
258     \brief
259 matt 1327 Return the values of a data point on this process
260 jgs 117 */
261 woo409 757 ESCRIPT_DLL_API
262 jgs 121 const boost::python::numeric::array
263 gross 1034 getValueOfDataPoint(int dataPointNo);
264 jgs 117
265     /**
266     \brief
267 gross 1034 sets the values of a data-point from a python object on this process
268 jgs 121 */
269 woo409 757 ESCRIPT_DLL_API
270 gross 921 void
271 gross 1034 setValueOfDataPointToPyObject(int dataPointNo, const boost::python::object& py_object);
272 jgs 121
273     /**
274     \brief
275 gross 1034 sets the values of a data-point from a numarray object on this process
276 jgs 121 */
277 woo409 757 ESCRIPT_DLL_API
278 jgs 149 void
279 gross 1034 setValueOfDataPointToArray(int dataPointNo, const boost::python::numeric::array&);
280 jgs 149
281     /**
282     \brief
283 gross 922 sets the values of a data-point on this process
284     */
285     ESCRIPT_DLL_API
286     void
287     setValueOfDataPoint(int dataPointNo, const double);
288    
289     /**
290     \brief
291 gross 921 Return the value of the specified data-point across all processors
292     */
293     ESCRIPT_DLL_API
294     const boost::python::numeric::array
295     getValueOfGlobalDataPoint(int procNo, int dataPointNo);
296    
297     /**
298     \brief
299 jgs 149 Return the tag number associated with the given data-point.
300    
301     The data-point number here corresponds to the data-point number in the
302     numarray returned by convertToNumArray.
303     */
304 woo409 757 ESCRIPT_DLL_API
305 jgs 149 int
306     getTagNumber(int dpno);
307    
308     /**
309     \brief
310 jgs 94 Return the C wrapper for the Data object.
311     */
312 woo409 757 ESCRIPT_DLL_API
313 jgs 102 escriptDataC
314     getDataC();
315 jgs 94
316     /**
317     \brief
318     Return the C wrapper for the Data object - const version.
319     */
320 woo409 757 ESCRIPT_DLL_API
321 jgs 102 escriptDataC
322     getDataC() const;
323 jgs 94
324     /**
325     \brief
326     Write the data as a string.
327     */
328 woo409 757 ESCRIPT_DLL_API
329 jgs 102 inline
330     std::string
331     toString() const
332     {
333     return m_data->toString();
334     }
335 jgs 94
336     /**
337     \brief
338     Return the DataArrayView of the point data. This essentially contains
339     the shape information for each data point although it also may be used
340     to manipulate the point data.
341     */
342 woo409 757 ESCRIPT_DLL_API
343 jgs 94 inline
344     const DataArrayView&
345     getPointDataView() const
346     {
347     return m_data->getPointDataView();
348     }
349    
350     /**
351     \brief
352 jgs 102 Whatever the current Data type make this into a DataExpanded.
353 jgs 94 */
354 woo409 757 ESCRIPT_DLL_API
355 jgs 94 void
356     expand();
357    
358     /**
359     \brief
360 jgs 102 If possible convert this Data to DataTagged. This will only allow
361 jgs 94 Constant data to be converted to tagged. An attempt to convert
362     Expanded data to tagged will throw an exception.
363 jgs 544 ==>*
364 jgs 94 */
365 woo409 757 ESCRIPT_DLL_API
366 jgs 94 void
367     tag();
368    
369     /**
370     \brief
371     Return true if this Data is expanded.
372     */
373 woo409 757 ESCRIPT_DLL_API
374 jgs 94 bool
375     isExpanded() const;
376    
377     /**
378     \brief
379     Return true if this Data is tagged.
380     */
381 woo409 757 ESCRIPT_DLL_API
382 jgs 94 bool
383     isTagged() const;
384    
385     /**
386     \brief
387 jgs 102 Return true if this Data is constant.
388 jgs 94 */
389 woo409 757 ESCRIPT_DLL_API
390 jgs 94 bool
391 jgs 102 isConstant() const;
392 jgs 94
393     /**
394     \brief
395 jgs 102 Return true if this Data is empty.
396 jgs 94 */
397 woo409 757 ESCRIPT_DLL_API
398 jgs 94 bool
399 jgs 102 isEmpty() const;
400 jgs 94
401     /**
402     \brief
403     Return the function space.
404     */
405 woo409 757 ESCRIPT_DLL_API
406 jgs 94 inline
407     const FunctionSpace&
408     getFunctionSpace() const
409     {
410     return m_data->getFunctionSpace();
411     }
412    
413     /**
414     \brief
415     Return a copy of the function space.
416     */
417 woo409 757 ESCRIPT_DLL_API
418 jgs 94 const FunctionSpace
419     getCopyOfFunctionSpace() const;
420    
421     /**
422     \brief
423     Return the domain.
424     */
425 woo409 757 ESCRIPT_DLL_API
426 jgs 94 inline
427     const AbstractDomain&
428     getDomain() const
429     {
430     return getFunctionSpace().getDomain();
431     }
432    
433     /**
434     \brief
435     Return a copy of the domain.
436     */
437 woo409 757 ESCRIPT_DLL_API
438 jgs 94 const AbstractDomain
439     getCopyOfDomain() const;
440    
441     /**
442     \brief
443     Return the rank of the point data.
444     */
445 woo409 757 ESCRIPT_DLL_API
446 jgs 94 inline
447     int
448     getDataPointRank() const
449     {
450     return m_data->getPointDataView().getRank();
451     }
452    
453     /**
454     \brief
455 gross 921 Return the number of data points
456     */
457     ESCRIPT_DLL_API
458     inline
459     int
460     getNumDataPoints() const
461     {
462     return getNumSamples() * getNumDataPointsPerSample();
463     }
464     /**
465     \brief
466 jgs 102 Return the number of samples.
467 jgs 94 */
468 woo409 757 ESCRIPT_DLL_API
469 jgs 94 inline
470     int
471 jgs 102 getNumSamples() const
472 jgs 94 {
473 jgs 102 return m_data->getNumSamples();
474 jgs 94 }
475    
476     /**
477     \brief
478 jgs 102 Return the number of data points per sample.
479 jgs 94 */
480 woo409 757 ESCRIPT_DLL_API
481 jgs 102 inline
482 jgs 94 int
483 jgs 102 getNumDataPointsPerSample() const
484     {
485     return m_data->getNumDPPSample();
486     }
487 gross 950 /**
488     \brief
489 matt 1327 dumps the object into a netCDF file
490 gross 950 */
491     ESCRIPT_DLL_API
492     void
493 gross 1141 dump(const std::string fileName) const;
494 jgs 94 /**
495     \brief
496     Return the sample data for the given sample no. This is not the
497     preferred interface but is provided for use by C code.
498     \param sampleNo - Input - the given sample no.
499     */
500 woo409 757 ESCRIPT_DLL_API
501 jgs 102 inline
502 jgs 94 DataAbstract::ValueType::value_type*
503 jgs 102 getSampleData(DataAbstract::ValueType::size_type sampleNo)
504     {
505     return m_data->getSampleData(sampleNo);
506     }
507 jgs 94
508     /**
509     \brief
510     Return the sample data for the given tag. If an attempt is made to
511     access data that isn't tagged an exception will be thrown.
512     \param tag - Input - the tag key.
513     */
514 woo409 757 ESCRIPT_DLL_API
515 jgs 102 inline
516 jgs 94 DataAbstract::ValueType::value_type*
517 jgs 102 getSampleDataByTag(int tag)
518     {
519     return m_data->getSampleDataByTag(tag);
520     }
521 jgs 94
522     /**
523     \brief
524     Return a view into the data for the data point specified.
525     NOTE: Construction of the DataArrayView is a relatively expensive
526     operation.
527     \param sampleNo - Input -
528     \param dataPointNo - Input -
529     */
530 woo409 757 ESCRIPT_DLL_API
531 jgs 94 inline
532     DataArrayView
533     getDataPoint(int sampleNo,
534     int dataPointNo)
535     {
536 matt 1332 return m_data->getDataPoint(sampleNo,dataPointNo);
537 jgs 94 }
538    
539     /**
540     \brief
541     Return a reference to the data point shape.
542     */
543 woo409 757 ESCRIPT_DLL_API
544 jgs 94 const DataArrayView::ShapeType&
545     getDataPointShape() const;
546    
547     /**
548     \brief
549 jgs 102 Return the data point shape as a tuple of integers.
550 jgs 94 */
551 woo409 757 ESCRIPT_DLL_API
552 jgs 121 const boost::python::tuple
553 jgs 94 getShapeTuple() const;
554    
555     /**
556     \brief
557     Return the size of the data point. It is the product of the
558     data point shape dimensions.
559     */
560 woo409 757 ESCRIPT_DLL_API
561 jgs 94 int
562     getDataPointSize() const;
563    
564     /**
565     \brief
566 jgs 102 Return the number of doubles stored for this Data.
567 jgs 94 */
568 woo409 757 ESCRIPT_DLL_API
569 jgs 94 DataArrayView::ValueType::size_type
570     getLength() const;
571    
572 gross 1044
573    
574 jgs 94 /**
575     \brief
576 gross 1044 Assign the given value to the tag assocciated with name. Implicitly converts this
577     object to type DataTagged. Throws an exception if this object
578     cannot be converted to a DataTagged object or name cannot be mapped onto a tag key.
579     \param tagKey - Input - Integer key.
580     \param value - Input - Value to associate with given key.
581     ==>*
582     */
583     ESCRIPT_DLL_API
584     void
585     setTaggedValueByName(std::string name,
586     const boost::python::object& value);
587    
588     /**
589     \brief
590 jgs 102 Assign the given value to the tag. Implicitly converts this
591 gross 1358 object to type DataTagged if it is constant.
592    
593 jgs 102 \param tagKey - Input - Integer key.
594     \param value - Input - Value to associate with given key.
595 jgs 544 ==>*
596 jgs 94 */
597 woo409 757 ESCRIPT_DLL_API
598 jgs 102 void
599     setTaggedValue(int tagKey,
600     const boost::python::object& value);
601    
602     /**
603     \brief
604     Assign the given value to the tag. Implicitly converts this
605 gross 1358 object to type DataTagged if it is constant.
606    
607 jgs 102 \param tagKey - Input - Integer key.
608     \param value - Input - Value to associate with given key.
609 jgs 544 ==>*
610 jgs 102 */
611 woo409 757 ESCRIPT_DLL_API
612 jgs 102 void
613 jgs 121 setTaggedValueFromCPP(int tagKey,
614     const DataArrayView& value);
615 jgs 102
616     /**
617     \brief
618     Copy other Data object into this Data object where mask is positive.
619     */
620 woo409 757 ESCRIPT_DLL_API
621 jgs 102 void
622     copyWithMask(const Data& other,
623     const Data& mask);
624    
625     /**
626     Data object operation methods and operators.
627     */
628    
629     /**
630     \brief
631 gross 1118 set all values to zero
632     *
633     */
634     ESCRIPT_DLL_API
635     void
636     setToZero();
637    
638     /**
639     \brief
640 jgs 102 Interpolates this onto the given functionspace and returns
641     the result as a Data object.
642 jgs 123 *
643 jgs 102 */
644 woo409 757 ESCRIPT_DLL_API
645 jgs 94 Data
646     interpolate(const FunctionSpace& functionspace) const;
647     /**
648     \brief
649     Calculates the gradient of the data at the data points of functionspace.
650     If functionspace is not present the function space of Function(getDomain()) is used.
651 jgs 123 *
652 jgs 94 */
653 woo409 757 ESCRIPT_DLL_API
654 jgs 94 Data
655     gradOn(const FunctionSpace& functionspace) const;
656    
657 woo409 757 ESCRIPT_DLL_API
658 jgs 94 Data
659     grad() const;
660    
661     /**
662     \brief
663     Calculate the integral over the function space domain.
664 jgs 123 *
665 jgs 94 */
666 woo409 757 ESCRIPT_DLL_API
667 jgs 94 boost::python::numeric::array
668     integrate() const;
669    
670     /**
671     \brief
672 gross 854 Returns 1./ Data object
673     *
674     */
675     ESCRIPT_DLL_API
676     Data
677     oneOver() const;
678     /**
679     \brief
680 jgs 94 Return a Data with a 1 for +ive values and a 0 for 0 or -ive values.
681 jgs 123 *
682 jgs 94 */
683 woo409 757 ESCRIPT_DLL_API
684 jgs 94 Data
685 gross 698 wherePositive() const;
686 jgs 94
687     /**
688     \brief
689 jgs 102 Return a Data with a 1 for -ive values and a 0 for +ive or 0 values.
690 jgs 123 *
691 jgs 94 */
692 woo409 757 ESCRIPT_DLL_API
693 jgs 94 Data
694 gross 698 whereNegative() const;
695 jgs 102
696     /**
697     \brief
698     Return a Data with a 1 for +ive or 0 values and a 0 for -ive values.
699 jgs 123 *
700 jgs 102 */
701 woo409 757 ESCRIPT_DLL_API
702 jgs 102 Data
703 gross 698 whereNonNegative() const;
704 jgs 94
705     /**
706     \brief
707 jgs 102 Return a Data with a 1 for -ive or 0 values and a 0 for +ive values.
708 jgs 123 *
709 jgs 94 */
710 woo409 757 ESCRIPT_DLL_API
711 jgs 94 Data
712 gross 698 whereNonPositive() const;
713 jgs 94
714     /**
715     \brief
716 jgs 102 Return a Data with a 1 for 0 values and a 0 for +ive or -ive values.
717 jgs 123 *
718 jgs 94 */
719 woo409 757 ESCRIPT_DLL_API
720 jgs 94 Data
721 jgs 571 whereZero(double tol=0.0) const;
722 jgs 94
723     /**
724     \brief
725 jgs 102 Return a Data with a 0 for 0 values and a 1 for +ive or -ive values.
726 jgs 123 *
727 jgs 94 */
728 woo409 757 ESCRIPT_DLL_API
729 jgs 94 Data
730 jgs 571 whereNonZero(double tol=0.0) const;
731 jgs 102
732     /**
733     \brief
734     Return the maximum absolute value of this Data object.
735 jgs 123 *
736 jgs 94 */
737 woo409 757 ESCRIPT_DLL_API
738 jgs 102 double
739     Lsup() const;
740 jgs 94
741     /**
742     \brief
743 jgs 102 Return the maximum value of this Data object.
744 jgs 123 *
745 jgs 94 */
746 woo409 757 ESCRIPT_DLL_API
747 jgs 102 double
748     sup() const;
749 jgs 94
750     /**
751     \brief
752 jgs 102 Return the minimum value of this Data object.
753 jgs 123 *
754 jgs 94 */
755 woo409 757 ESCRIPT_DLL_API
756 jgs 94 double
757 jgs 102 inf() const;
758 jgs 94
759     /**
760     \brief
761 jgs 102 Return the absolute value of each data point of this Data object.
762 jgs 123 *
763 jgs 94 */
764 woo409 757 ESCRIPT_DLL_API
765 jgs 102 Data
766     abs() const;
767 jgs 94
768     /**
769     \brief
770 jgs 102 Return the maximum value of each data point of this Data object.
771 jgs 123 *
772 jgs 94 */
773 woo409 757 ESCRIPT_DLL_API
774 jgs 102 Data
775     maxval() const;
776 jgs 94
777     /**
778     \brief
779 jgs 102 Return the minimum value of each data point of this Data object.
780 jgs 123 *
781 jgs 94 */
782 woo409 757 ESCRIPT_DLL_API
783 jgs 94 Data
784 jgs 102 minval() const;
785 jgs 94
786     /**
787     \brief
788 jgs 121 Return the (sample number, data-point number) of the data point with
789     the minimum value in this Data object.
790     */
791 woo409 757 ESCRIPT_DLL_API
792 jgs 121 const boost::python::tuple
793 gross 921 minGlobalDataPoint() const;
794 jgs 121
795 woo409 757 ESCRIPT_DLL_API
796 jgs 148 void
797 gross 921 calc_minGlobalDataPoint(int& ProcNo, int& DataPointNo) const;
798 jgs 121 /**
799     \brief
800 jgs 102 Return the sign of each data point of this Data object.
801     -1 for negative values, zero for zero values, 1 for positive values.
802 jgs 123 *
803 jgs 94 */
804 woo409 757 ESCRIPT_DLL_API
805 jgs 102 Data
806     sign() const;
807 jgs 94
808     /**
809 jgs 123 \brief
810 ksteube 775 Return the symmetric part of a matrix which is half the matrix plus its transpose.
811     *
812     */
813     ESCRIPT_DLL_API
814     Data
815     symmetric() const;
816    
817     /**
818     \brief
819     Return the nonsymmetric part of a matrix which is half the matrix minus its transpose.
820     *
821     */
822     ESCRIPT_DLL_API
823     Data
824     nonsymmetric() const;
825    
826     /**
827     \brief
828     Return the trace of a matrix
829     *
830     */
831     ESCRIPT_DLL_API
832     Data
833 gross 800 trace(int axis_offset) const;
834 ksteube 775
835     /**
836     \brief
837     Transpose each data point of this Data object around the given axis.
838     *
839     */
840     ESCRIPT_DLL_API
841     Data
842     transpose(int axis_offset) const;
843    
844     /**
845     \brief
846 gross 576 Return the eigenvalues of the symmetric part at each data point of this Data object in increasing values.
847     Currently this function is restricted to rank 2, square shape, and dimension 3.
848     *
849     */
850 woo409 757 ESCRIPT_DLL_API
851 gross 576 Data
852     eigenvalues() const;
853    
854     /**
855     \brief
856     Return the eigenvalues and corresponding eigenvcetors of the symmetric part at each data point of this Data object.
857 matt 1327 the eigenvalues are ordered in increasing size where eigenvalues with relative difference less than
858     tol are treated as equal. The eigenvectors are orthogonal, normalized and the sclaed such that the
859 gross 576 first non-zero entry is positive.
860     Currently this function is restricted to rank 2, square shape, and dimension 3
861     *
862     */
863 woo409 757 ESCRIPT_DLL_API
864 gross 576 const boost::python::tuple
865     eigenvalues_and_eigenvectors(const double tol=1.e-12) const;
866    
867     /**
868     \brief
869 gross 804 swaps the components axis0 and axis1
870 jgs 123 *
871 jgs 102 */
872 woo409 757 ESCRIPT_DLL_API
873 jgs 102 Data
874 gross 804 swapaxes(const int axis0, const int axis1) const;
875 jgs 102
876     /**
877 jgs 123 \brief
878 ksteube 876 Return the error function erf of each data point of this Data object.
879     *
880     */
881     ESCRIPT_DLL_API
882     Data
883     erf() const;
884    
885     /**
886     \brief
887 jgs 123 Return the sin of each data point of this Data object.
888     *
889 jgs 102 */
890 woo409 757 ESCRIPT_DLL_API
891 jgs 102 Data
892 jgs 123 sin() const;
893    
894     /**
895     \brief
896     Return the cos of each data point of this Data object.
897     *
898     */
899 woo409 757 ESCRIPT_DLL_API
900 jgs 123 Data
901     cos() const;
902    
903     /**
904     \brief
905     Return the tan of each data point of this Data object.
906     *
907     */
908 woo409 757 ESCRIPT_DLL_API
909 jgs 123 Data
910     tan() const;
911    
912     /**
913     \brief
914 jgs 150 Return the asin of each data point of this Data object.
915     *
916     */
917 woo409 757 ESCRIPT_DLL_API
918 jgs 150 Data
919     asin() const;
920    
921     /**
922     \brief
923     Return the acos of each data point of this Data object.
924     *
925     */
926 woo409 757 ESCRIPT_DLL_API
927 jgs 150 Data
928     acos() const;
929    
930     /**
931     \brief
932     Return the atan of each data point of this Data object.
933     *
934     */
935 woo409 757 ESCRIPT_DLL_API
936 jgs 150 Data
937     atan() const;
938    
939     /**
940     \brief
941     Return the sinh of each data point of this Data object.
942     *
943     */
944 woo409 757 ESCRIPT_DLL_API
945 jgs 150 Data
946     sinh() const;
947    
948     /**
949     \brief
950     Return the cosh of each data point of this Data object.
951     *
952     */
953 woo409 757 ESCRIPT_DLL_API
954 jgs 150 Data
955     cosh() const;
956    
957     /**
958     \brief
959     Return the tanh of each data point of this Data object.
960     *
961     */
962 woo409 757 ESCRIPT_DLL_API
963 jgs 150 Data
964     tanh() const;
965    
966     /**
967     \brief
968     Return the asinh of each data point of this Data object.
969     *
970     */
971 woo409 757 ESCRIPT_DLL_API
972 jgs 150 Data
973     asinh() const;
974    
975     /**
976     \brief
977     Return the acosh of each data point of this Data object.
978     *
979     */
980 woo409 757 ESCRIPT_DLL_API
981 jgs 150 Data
982     acosh() const;
983    
984     /**
985     \brief
986     Return the atanh of each data point of this Data object.
987     *
988     */
989 woo409 757 ESCRIPT_DLL_API
990 jgs 150 Data
991     atanh() const;
992    
993     /**
994     \brief
995 jgs 123 Return the log to base 10 of each data point of this Data object.
996     *
997     */
998 woo409 757 ESCRIPT_DLL_API
999 jgs 123 Data
1000 gross 286 log10() const;
1001 jgs 123
1002     /**
1003     \brief
1004     Return the natural log of each data point of this Data object.
1005     *
1006     */
1007 woo409 757 ESCRIPT_DLL_API
1008 jgs 123 Data
1009 gross 286 log() const;
1010 jgs 123
1011     /**
1012     \brief
1013     Return the exponential function of each data point of this Data object.
1014     *
1015     */
1016 woo409 757 ESCRIPT_DLL_API
1017 jgs 123 Data
1018 jgs 102 exp() const;
1019    
1020     /**
1021 jgs 123 \brief
1022     Return the square root of each data point of this Data object.
1023     *
1024 jgs 102 */
1025 woo409 757 ESCRIPT_DLL_API
1026 jgs 102 Data
1027     sqrt() const;
1028    
1029     /**
1030 jgs 123 \brief
1031     Return the negation of each data point of this Data object.
1032     *
1033 jgs 121 */
1034 woo409 757 ESCRIPT_DLL_API
1035 jgs 121 Data
1036     neg() const;
1037    
1038     /**
1039 jgs 123 \brief
1040     Return the identity of each data point of this Data object.
1041     Simply returns this object unmodified.
1042     *
1043 jgs 121 */
1044 woo409 757 ESCRIPT_DLL_API
1045 jgs 121 Data
1046     pos() const;
1047    
1048     /**
1049 jgs 94 \brief
1050 jgs 102 Return the given power of each data point of this Data object.
1051 jgs 121
1052     \param right Input - the power to raise the object to.
1053 jgs 123 *
1054 jgs 102 */
1055 woo409 757 ESCRIPT_DLL_API
1056 jgs 102 Data
1057     powD(const Data& right) const;
1058    
1059 jgs 121 /**
1060 jgs 123 \brief
1061     Return the given power of each data point of this boost python object.
1062 matt 1327
1063 jgs 123 \param right Input - the power to raise the object to.
1064     *
1065 jgs 121 */
1066 woo409 757 ESCRIPT_DLL_API
1067 jgs 102 Data
1068     powO(const boost::python::object& right) const;
1069    
1070     /**
1071 jgs 123 \brief
1072 gross 699 Return the given power of each data point of this boost python object.
1073 matt 1327
1074 gross 699 \param left Input - the bases
1075     *
1076     */
1077    
1078 woo409 757 ESCRIPT_DLL_API
1079 gross 699 Data
1080     rpowO(const boost::python::object& left) const;
1081    
1082     /**
1083     \brief
1084 jgs 123 writes the object to a file in the DX file format
1085 jgs 104 */
1086 woo409 757 ESCRIPT_DLL_API
1087 jgs 104 void
1088     saveDX(std::string fileName) const;
1089    
1090     /**
1091 jgs 123 \brief
1092     writes the object to a file in the VTK file format
1093 jgs 110 */
1094 woo409 757 ESCRIPT_DLL_API
1095 jgs 110 void
1096     saveVTK(std::string fileName) const;
1097    
1098     /**
1099 jgs 102 \brief
1100     Overloaded operator +=
1101     \param right - Input - The right hand side.
1102 jgs 123 *
1103 jgs 102 */
1104 woo409 757 ESCRIPT_DLL_API
1105 jgs 102 Data& operator+=(const Data& right);
1106 woo409 757 ESCRIPT_DLL_API
1107 jgs 102 Data& operator+=(const boost::python::object& right);
1108    
1109 ksteube 1312 ESCRIPT_DLL_API
1110     Data& operator=(const Data& other);
1111    
1112 jgs 102 /**
1113     \brief
1114     Overloaded operator -=
1115     \param right - Input - The right hand side.
1116 jgs 123 *
1117 jgs 102 */
1118 woo409 757 ESCRIPT_DLL_API
1119 jgs 102 Data& operator-=(const Data& right);
1120 woo409 757 ESCRIPT_DLL_API
1121 jgs 102 Data& operator-=(const boost::python::object& right);
1122    
1123     /**
1124     \brief
1125     Overloaded operator *=
1126     \param right - Input - The right hand side.
1127 jgs 123 *
1128 jgs 102 */
1129 woo409 757 ESCRIPT_DLL_API
1130 jgs 102 Data& operator*=(const Data& right);
1131 woo409 757 ESCRIPT_DLL_API
1132 jgs 102 Data& operator*=(const boost::python::object& right);
1133    
1134     /**
1135     \brief
1136     Overloaded operator /=
1137     \param right - Input - The right hand side.
1138 jgs 123 *
1139 jgs 102 */
1140 woo409 757 ESCRIPT_DLL_API
1141 jgs 102 Data& operator/=(const Data& right);
1142 woo409 757 ESCRIPT_DLL_API
1143 jgs 102 Data& operator/=(const boost::python::object& right);
1144    
1145     /**
1146     \brief
1147 jgs 94 Returns true if this can be interpolated to functionspace.
1148     */
1149 woo409 757 ESCRIPT_DLL_API
1150 jgs 94 bool
1151     probeInterpolation(const FunctionSpace& functionspace) const;
1152    
1153     /**
1154 jgs 102 Data object slicing methods.
1155     */
1156    
1157     /**
1158 jgs 94 \brief
1159 jgs 102 Returns a slice from this Data object.
1160    
1161     /description
1162     Implements the [] get operator in python.
1163     Calls getSlice.
1164    
1165     \param key - Input - python slice tuple specifying
1166     slice to return.
1167 jgs 94 */
1168 woo409 757 ESCRIPT_DLL_API
1169 jgs 102 Data
1170     getItem(const boost::python::object& key) const;
1171    
1172     /**
1173     \brief
1174     Copies slice from value into this Data object.
1175    
1176     Implements the [] set operator in python.
1177     Calls setSlice.
1178    
1179     \param key - Input - python slice tuple specifying
1180     slice to copy from value.
1181     \param value - Input - Data object to copy from.
1182     */
1183 woo409 757 ESCRIPT_DLL_API
1184 jgs 94 void
1185 jgs 102 setItemD(const boost::python::object& key,
1186     const Data& value);
1187 jgs 94
1188 woo409 757 ESCRIPT_DLL_API
1189 jgs 102 void
1190     setItemO(const boost::python::object& key,
1191     const boost::python::object& value);
1192    
1193     // These following public methods should be treated as private.
1194    
1195 jgs 94 /**
1196     \brief
1197 jgs 102 Perform the given unary operation on every element of every data point in
1198     this Data object.
1199 jgs 94 */
1200 jgs 102 template <class UnaryFunction>
1201 woo409 757 ESCRIPT_DLL_API
1202 jgs 102 inline
1203 jgs 94 void
1204 matt 1334 unaryOp2(UnaryFunction operation);
1205 jgs 102
1206     /**
1207     \brief
1208     Return a Data object containing the specified slice of
1209     this Data object.
1210     \param region - Input - Region to copy.
1211 jgs 123 *
1212 jgs 94 */
1213 woo409 757 ESCRIPT_DLL_API
1214 jgs 102 Data
1215     getSlice(const DataArrayView::RegionType& region) const;
1216 jgs 94
1217 jgs 102 /**
1218     \brief
1219     Copy the specified slice from the given value into this
1220     Data object.
1221     \param value - Input - Data to copy from.
1222     \param region - Input - Region to copy.
1223 jgs 123 *
1224 jgs 102 */
1225 woo409 757 ESCRIPT_DLL_API
1226 jgs 102 void
1227     setSlice(const Data& value,
1228     const DataArrayView::RegionType& region);
1229    
1230 jgs 119 /**
1231     \brief
1232     Archive the current Data object to the given file.
1233     \param fileName - Input - file to archive to.
1234     */
1235 woo409 757 ESCRIPT_DLL_API
1236 jgs 119 void
1237     archiveData(const std::string fileName);
1238    
1239     /**
1240     \brief
1241     Extract the Data object archived in the given file, overwriting
1242     the current Data object.
1243     Note - the current object must be of type DataEmpty.
1244     \param fileName - Input - file to extract from.
1245 jgs 121 \param fspace - Input - a suitable FunctionSpace descibing the data.
1246 jgs 119 */
1247 woo409 757 ESCRIPT_DLL_API
1248 jgs 119 void
1249     extractData(const std::string fileName,
1250     const FunctionSpace& fspace);
1251    
1252 bcumming 751
1253     /**
1254     \brief
1255 matt 1327 print the data values to stdout. Used for debugging
1256 bcumming 751 */
1257 bcumming 782 ESCRIPT_DLL_API
1258 matt 1327 void
1259 matt 1332 print(void);
1260 bcumming 751
1261 bcumming 782 /**
1262     \brief
1263     return the MPI rank number of the local data
1264 matt 1332 MPI_COMM_WORLD is assumed and the result of MPI_Comm_size()
1265     is returned
1266 bcumming 782 */
1267     ESCRIPT_DLL_API
1268 matt 1332 int
1269     get_MPIRank(void) const;
1270 bcumming 782
1271     /**
1272     \brief
1273     return the MPI rank number of the local data
1274 matt 1332 MPI_COMM_WORLD is assumed and the result of MPI_Comm_rank()
1275     is returned
1276 bcumming 782 */
1277     ESCRIPT_DLL_API
1278 matt 1332 int
1279     get_MPISize(void) const;
1280 bcumming 782
1281     /**
1282     \brief
1283     return the MPI rank number of the local data
1284 matt 1332 MPI_COMM_WORLD is assumed and returned.
1285 bcumming 782 */
1286     ESCRIPT_DLL_API
1287 matt 1332 MPI_Comm
1288     get_MPIComm(void) const;
1289 ksteube 813
1290     /**
1291     \brief
1292     return the object produced by the factory, which is a DataConstant or DataExpanded
1293     */
1294     ESCRIPT_DLL_API
1295 matt 1332 DataAbstract*
1296     borrowData(void) const;
1297 ksteube 813
1298 jgs 102 protected:
1299    
1300 jgs 94 private:
1301    
1302     /**
1303     \brief
1304 jgs 102 Check *this and the right operand are compatible. Throws
1305     an exception if they aren't.
1306     \param right - Input - The right hand side.
1307     */
1308     inline
1309     void
1310     operandCheck(const Data& right) const
1311     {
1312     return m_data->operandCheck(*(right.m_data.get()));
1313     }
1314    
1315     /**
1316     \brief
1317     Perform the specified reduction algorithm on every element of every data point in
1318 jgs 113 this Data object according to the given function and return the single value result.
1319 jgs 102 */
1320 jgs 147 template <class BinaryFunction>
1321 jgs 102 inline
1322     double
1323 jgs 147 algorithm(BinaryFunction operation,
1324     double initial_value) const;
1325 jgs 102
1326 jgs 113 /**
1327     \brief
1328     Reduce each data-point in this Data object using the given operation. Return a Data
1329     object with the same number of data-points, but with each data-point containing only
1330     one value - the result of the reduction operation on the corresponding data-point in
1331     this Data object
1332     */
1333 jgs 147 template <class BinaryFunction>
1334 jgs 106 inline
1335     Data
1336 jgs 147 dp_algorithm(BinaryFunction operation,
1337     double initial_value) const;
1338 jgs 106
1339 jgs 102 /**
1340     \brief
1341     Perform the given binary operation on all of the data's elements.
1342     The underlying type of the right hand side (right) determines the final
1343     type of *this after the operation. For example if the right hand side
1344     is expanded *this will be expanded if necessary.
1345     RHS is a Data object.
1346     */
1347     template <class BinaryFunction>
1348     inline
1349     void
1350     binaryOp(const Data& right,
1351     BinaryFunction operation);
1352    
1353     /**
1354     \brief
1355     Convert the data type of the RHS to match this.
1356     \param right - Input - data type to match.
1357     */
1358     void
1359     typeMatchLeft(Data& right) const;
1360    
1361     /**
1362     \brief
1363     Convert the data type of this to match the RHS.
1364     \param right - Input - data type to match.
1365     */
1366     void
1367     typeMatchRight(const Data& right);
1368    
1369     /**
1370     \brief
1371 jgs 94 Construct a Data object of the appropriate type.
1372     */
1373     template <class IValueType>
1374     void
1375     initialise(const IValueType& value,
1376     const FunctionSpace& what,
1377     bool expanded);
1378    
1379     //
1380 gross 783 // flag to protect the data object against any update
1381     bool m_protected;
1382    
1383     //
1384 jgs 102 // pointer to the actual data object
1385 jgs 94 boost::shared_ptr<DataAbstract> m_data;
1386    
1387     };
1388    
1389     template <class IValueType>
1390     void
1391     Data::initialise(const IValueType& value,
1392     const FunctionSpace& what,
1393     bool expanded)
1394     {
1395     //
1396     // Construct a Data object of the appropriate type.
1397     // Construct the object first as there seems to be a bug which causes
1398     // undefined behaviour if an exception is thrown during construction
1399     // within the shared_ptr constructor.
1400     if (expanded) {
1401     DataAbstract* temp=new DataExpanded(value,what);
1402 jgs 102 boost::shared_ptr<DataAbstract> temp_data(temp);
1403     m_data=temp_data;
1404 jgs 94 } else {
1405     DataAbstract* temp=new DataConstant(value,what);
1406 jgs 102 boost::shared_ptr<DataAbstract> temp_data(temp);
1407     m_data=temp_data;
1408 jgs 94 }
1409     }
1410    
1411 jgs 102 /**
1412     Binary Data object operators.
1413     */
1414 matt 1327 inline double rpow(double x,double y)
1415 gross 854 {
1416     return pow(y,x);
1417 gross 1028 }
1418 jgs 94
1419     /**
1420     \brief
1421     Operator+
1422     Takes two Data objects.
1423     */
1424 woo409 757 ESCRIPT_DLL_API Data operator+(const Data& left, const Data& right);
1425 jgs 94
1426     /**
1427     \brief
1428     Operator-
1429     Takes two Data objects.
1430     */
1431 woo409 757 ESCRIPT_DLL_API Data operator-(const Data& left, const Data& right);
1432 jgs 94
1433     /**
1434     \brief
1435     Operator*
1436     Takes two Data objects.
1437     */
1438 woo409 757 ESCRIPT_DLL_API Data operator*(const Data& left, const Data& right);
1439 jgs 94
1440     /**
1441     \brief
1442     Operator/
1443     Takes two Data objects.
1444     */
1445 woo409 757 ESCRIPT_DLL_API Data operator/(const Data& left, const Data& right);
1446 jgs 94
1447     /**
1448     \brief
1449     Operator+
1450     Takes LHS Data object and RHS python::object.
1451     python::object must be convertable to Data type.
1452     */
1453 woo409 757 ESCRIPT_DLL_API Data operator+(const Data& left, const boost::python::object& right);
1454 jgs 94
1455     /**
1456     \brief
1457     Operator-
1458     Takes LHS Data object and RHS python::object.
1459     python::object must be convertable to Data type.
1460     */
1461 woo409 757 ESCRIPT_DLL_API Data operator-(const Data& left, const boost::python::object& right);
1462 jgs 94
1463     /**
1464     \brief
1465     Operator*
1466     Takes LHS Data object and RHS python::object.
1467     python::object must be convertable to Data type.
1468     */
1469 woo409 757 ESCRIPT_DLL_API Data operator*(const Data& left, const boost::python::object& right);
1470 jgs 94
1471     /**
1472     \brief
1473     Operator/
1474     Takes LHS Data object and RHS python::object.
1475     python::object must be convertable to Data type.
1476     */
1477 woo409 757 ESCRIPT_DLL_API Data operator/(const Data& left, const boost::python::object& right);
1478 jgs 94
1479     /**
1480     \brief
1481     Operator+
1482     Takes LHS python::object and RHS Data object.
1483     python::object must be convertable to Data type.
1484     */
1485 woo409 757 ESCRIPT_DLL_API Data operator+(const boost::python::object& left, const Data& right);
1486 jgs 94
1487     /**
1488     \brief
1489     Operator-
1490     Takes LHS python::object and RHS Data object.
1491     python::object must be convertable to Data type.
1492     */
1493 woo409 757 ESCRIPT_DLL_API Data operator-(const boost::python::object& left, const Data& right);
1494 jgs 94
1495     /**
1496     \brief
1497     Operator*
1498     Takes LHS python::object and RHS Data object.
1499     python::object must be convertable to Data type.
1500     */
1501 woo409 757 ESCRIPT_DLL_API Data operator*(const boost::python::object& left, const Data& right);
1502 jgs 94
1503     /**
1504     \brief
1505     Operator/
1506     Takes LHS python::object and RHS Data object.
1507     python::object must be convertable to Data type.
1508     */
1509 woo409 757 ESCRIPT_DLL_API Data operator/(const boost::python::object& left, const Data& right);
1510 jgs 94
1511 ksteube 1312
1512    
1513 jgs 94 /**
1514     \brief
1515     Output operator
1516     */
1517 woo409 757 ESCRIPT_DLL_API std::ostream& operator<<(std::ostream& o, const Data& data);
1518 jgs 94
1519     /**
1520     \brief
1521 ksteube 813 Compute a tensor product of two Data objects
1522     \param arg0 - Input - Data object
1523     \param arg1 - Input - Data object
1524     \param axis_offset - Input - axis offset
1525     \param transpose - Input - 0: transpose neither, 1: transpose arg0, 2: transpose arg1
1526     */
1527     ESCRIPT_DLL_API
1528     Data
1529     C_GeneralTensorProduct(Data& arg0,
1530     Data& arg1,
1531     int axis_offset=0,
1532     int transpose=0);
1533    
1534 matt 1327
1535 ksteube 813 /**
1536     \brief
1537 matt 1327 Compute a tensor operation with two Data objects
1538     \param arg0 - Input - Data object
1539     \param arg1 - Input - Data object
1540     \param operation - Input - Binary op functor
1541     */
1542     template <typename BinaryFunction>
1543     ESCRIPT_DLL_API
1544     Data
1545     C_TensorBinaryOperation(Data const &arg0,
1546 matt 1332 Data const &arg1,
1547     BinaryFunction operation);
1548 matt 1327
1549     /**
1550     \brief
1551 jgs 94 Return true if operands are equivalent, else return false.
1552 matt 1327 NB: this operator does very little at this point, and isn't to
1553 jgs 94 be relied on. Requires further implementation.
1554     */
1555 ksteube 1312 // ESCRIPT_DLL_API bool operator==(const Data& left, const Data& right);
1556 jgs 94
1557 jgs 102 /**
1558     \brief
1559     Perform the given binary operation with this and right as operands.
1560     Right is a Data object.
1561     */
1562 jgs 94 template <class BinaryFunction>
1563     inline
1564     void
1565     Data::binaryOp(const Data& right,
1566     BinaryFunction operation)
1567     {
1568     //
1569     // if this has a rank of zero promote it to the rank of the RHS
1570     if (getPointDataView().getRank()==0 && right.getPointDataView().getRank()!=0) {
1571 gross 854 throw DataException("Error - attempt to update rank zero object with object with rank bigger than zero.");
1572 jgs 94 }
1573     //
1574     // initially make the temporary a shallow copy
1575     Data tempRight(right);
1576     if (getFunctionSpace()!=right.getFunctionSpace()) {
1577     if (right.probeInterpolation(getFunctionSpace())) {
1578     //
1579 matt 1327 // an interpolation is required so create a new Data
1580 jgs 94 tempRight=Data(right,this->getFunctionSpace());
1581     } else if (probeInterpolation(right.getFunctionSpace())) {
1582     //
1583     // interpolate onto the RHS function space
1584     Data tempLeft(*this,right.getFunctionSpace());
1585     m_data=tempLeft.m_data;
1586     }
1587     }
1588     operandCheck(tempRight);
1589     //
1590     // ensure this has the right type for the RHS
1591 jgs 102 typeMatchRight(tempRight);
1592 jgs 94 //
1593     // Need to cast to the concrete types so that the correct binaryOp
1594     // is called.
1595     if (isExpanded()) {
1596     //
1597     // Expanded data will be done in parallel, the right hand side can be
1598     // of any data type
1599     DataExpanded* leftC=dynamic_cast<DataExpanded*>(m_data.get());
1600     EsysAssert((leftC!=0), "Programming error - casting to DataExpanded.");
1601     escript::binaryOp(*leftC,*(tempRight.m_data.get()),operation);
1602     } else if (isTagged()) {
1603     //
1604     // Tagged data is operated on serially, the right hand side can be
1605     // either DataConstant or DataTagged
1606     DataTagged* leftC=dynamic_cast<DataTagged*>(m_data.get());
1607     EsysAssert((leftC!=0), "Programming error - casting to DataTagged.");
1608     if (right.isTagged()) {
1609     DataTagged* rightC=dynamic_cast<DataTagged*>(tempRight.m_data.get());
1610     EsysAssert((rightC!=0), "Programming error - casting to DataTagged.");
1611     escript::binaryOp(*leftC,*rightC,operation);
1612     } else {
1613     DataConstant* rightC=dynamic_cast<DataConstant*>(tempRight.m_data.get());
1614     EsysAssert((rightC!=0), "Programming error - casting to DataConstant.");
1615     escript::binaryOp(*leftC,*rightC,operation);
1616     }
1617 jgs 102 } else if (isConstant()) {
1618 jgs 94 DataConstant* leftC=dynamic_cast<DataConstant*>(m_data.get());
1619     DataConstant* rightC=dynamic_cast<DataConstant*>(tempRight.m_data.get());
1620     EsysAssert((leftC!=0 && rightC!=0), "Programming error - casting to DataConstant.");
1621     escript::binaryOp(*leftC,*rightC,operation);
1622     }
1623     }
1624    
1625 jgs 102 /**
1626     \brief
1627     Perform the given Data object reduction algorithm on this and return the result.
1628     Given operation combines each element of each data point, thus argument
1629     object (*this) is a rank n Data object, and returned object is a scalar.
1630     Calls escript::algorithm.
1631     */
1632 jgs 147 template <class BinaryFunction>
1633 jgs 94 inline
1634     double
1635 jgs 147 Data::algorithm(BinaryFunction operation, double initial_value) const
1636 jgs 94 {
1637     if (isExpanded()) {
1638     DataExpanded* leftC=dynamic_cast<DataExpanded*>(m_data.get());
1639     EsysAssert((leftC!=0), "Programming error - casting to DataExpanded.");
1640 jgs 147 return escript::algorithm(*leftC,operation,initial_value);
1641 jgs 102 } else if (isTagged()) {
1642 jgs 94 DataTagged* leftC=dynamic_cast<DataTagged*>(m_data.get());
1643     EsysAssert((leftC!=0), "Programming error - casting to DataTagged.");
1644 jgs 147 return escript::algorithm(*leftC,operation,initial_value);
1645 jgs 102 } else if (isConstant()) {
1646 jgs 94 DataConstant* leftC=dynamic_cast<DataConstant*>(m_data.get());
1647     EsysAssert((leftC!=0), "Programming error - casting to DataConstant.");
1648 jgs 147 return escript::algorithm(*leftC,operation,initial_value);
1649 jgs 94 }
1650 jgs 102 return 0;
1651 jgs 94 }
1652    
1653 jgs 102 /**
1654     \brief
1655     Perform the given data point reduction algorithm on data and return the result.
1656     Given operation combines each element within each data point into a scalar,
1657 matt 1327 thus argument object is a rank n Data object, and returned object is a
1658 jgs 102 rank 0 Data object.
1659     Calls escript::dp_algorithm.
1660     */
1661 jgs 147 template <class BinaryFunction>
1662 jgs 94 inline
1663     Data
1664 jgs 147 Data::dp_algorithm(BinaryFunction operation, double initial_value) const
1665 jgs 94 {
1666 jgs 106 if (isExpanded()) {
1667 jgs 559 Data result(0,DataArrayView::ShapeType(),getFunctionSpace(),isExpanded());
1668 jgs 106 DataExpanded* dataE=dynamic_cast<DataExpanded*>(m_data.get());
1669 jgs 102 DataExpanded* resultE=dynamic_cast<DataExpanded*>(result.m_data.get());
1670     EsysAssert((dataE!=0), "Programming error - casting data to DataExpanded.");
1671     EsysAssert((resultE!=0), "Programming error - casting result to DataExpanded.");
1672 jgs 147 escript::dp_algorithm(*dataE,*resultE,operation,initial_value);
1673 jgs 559 return result;
1674 jgs 106 } else if (isTagged()) {
1675     DataTagged* dataT=dynamic_cast<DataTagged*>(m_data.get());
1676 jgs 562 DataArrayView::ShapeType viewShape;
1677     DataArrayView::ValueType viewData(1);
1678     viewData[0]=0;
1679     DataArrayView defaultValue(viewData,viewShape);
1680 jgs 559 DataTagged::TagListType keys;
1681 jgs 562 DataTagged::ValueListType values;
1682 jgs 559 DataTagged::DataMapType::const_iterator i;
1683     for (i=dataT->getTagLookup().begin();i!=dataT->getTagLookup().end();i++) {
1684     keys.push_back(i->first);
1685 jgs 562 values.push_back(defaultValue);
1686 jgs 559 }
1687     Data result(keys,values,defaultValue,getFunctionSpace());
1688 jgs 102 DataTagged* resultT=dynamic_cast<DataTagged*>(result.m_data.get());
1689     EsysAssert((dataT!=0), "Programming error - casting data to DataTagged.");
1690     EsysAssert((resultT!=0), "Programming error - casting result to DataTagged.");
1691 jgs 147 escript::dp_algorithm(*dataT,*resultT,operation,initial_value);
1692 jgs 559 return result;
1693 jgs 106 } else if (isConstant()) {
1694 jgs 559 Data result(0,DataArrayView::ShapeType(),getFunctionSpace(),isExpanded());
1695 jgs 106 DataConstant* dataC=dynamic_cast<DataConstant*>(m_data.get());
1696 jgs 102 DataConstant* resultC=dynamic_cast<DataConstant*>(result.m_data.get());
1697     EsysAssert((dataC!=0), "Programming error - casting data to DataConstant.");
1698     EsysAssert((resultC!=0), "Programming error - casting result to DataConstant.");
1699 jgs 147 escript::dp_algorithm(*dataC,*resultC,operation,initial_value);
1700 jgs 559 return result;
1701 jgs 102 }
1702 jgs 559 Data falseRetVal; // to keep compiler quiet
1703     return falseRetVal;
1704 jgs 94 }
1705    
1706 matt 1327 template <typename BinaryFunction>
1707     Data
1708     C_TensorBinaryOperation(Data const &arg_0,
1709 matt 1332 Data const &arg_1,
1710     BinaryFunction operation)
1711 matt 1327 {
1712     // Interpolate if necessary and find an appropriate function space
1713     Data arg_0_Z, arg_1_Z;
1714     if (arg_0.getFunctionSpace()!=arg_1.getFunctionSpace()) {
1715     if (arg_0.probeInterpolation(arg_1.getFunctionSpace())) {
1716     arg_0_Z = arg_0.interpolate(arg_1.getFunctionSpace());
1717     arg_1_Z = Data(arg_1);
1718     }
1719     else if (arg_1.probeInterpolation(arg_0.getFunctionSpace())) {
1720     arg_1_Z=arg_1.interpolate(arg_0.getFunctionSpace());
1721     arg_0_Z =Data(arg_0);
1722     }
1723     else {
1724     throw DataException("Error - C_TensorBinaryOperation: arguments have incompatible function spaces.");
1725     }
1726     } else {
1727     arg_0_Z = Data(arg_0);
1728     arg_1_Z = Data(arg_1);
1729     }
1730     // Get rank and shape of inputs
1731     int rank0 = arg_0_Z.getDataPointRank();
1732     int rank1 = arg_1_Z.getDataPointRank();
1733     DataArrayView::ShapeType shape0 = arg_0_Z.getDataPointShape();
1734     DataArrayView::ShapeType shape1 = arg_1_Z.getDataPointShape();
1735     int size0 = arg_0_Z.getDataPointSize();
1736     int size1 = arg_1_Z.getDataPointSize();
1737    
1738     // Declare output Data object
1739     Data res;
1740    
1741     if (shape0 == shape1) {
1742    
1743     if (arg_0_Z.isConstant() && arg_1_Z.isConstant()) {
1744 matt 1332 res = Data(0.0, shape0, arg_1_Z.getFunctionSpace()); // DataConstant output
1745 matt 1327 double *ptr_0 = &((arg_0_Z.getPointDataView().getData())[0]);
1746     double *ptr_1 = &((arg_1_Z.getPointDataView().getData())[0]);
1747     double *ptr_2 = &((res.getPointDataView().getData())[0]);
1748     tensor_binary_operation(size0, ptr_0, ptr_1, ptr_2, operation);
1749     }
1750     else if (arg_0_Z.isConstant() && arg_1_Z.isTagged()) {
1751    
1752     // Prepare the DataConstant input
1753     DataConstant* tmp_0=dynamic_cast<DataConstant*>(arg_0_Z.borrowData());
1754    
1755     // Borrow DataTagged input from Data object
1756     DataTagged* tmp_1=dynamic_cast<DataTagged*>(arg_1_Z.borrowData());
1757    
1758     // Prepare a DataTagged output 2
1759 matt 1332 res = Data(0.0, shape0, arg_1_Z.getFunctionSpace()); // DataTagged output
1760 matt 1327 res.tag();
1761     DataTagged* tmp_2=dynamic_cast<DataTagged*>(res.borrowData());
1762    
1763     // Prepare offset into DataConstant
1764     int offset_0 = tmp_0->getPointOffset(0,0);
1765     double *ptr_0 = &((arg_0_Z.getPointDataView().getData())[offset_0]);
1766     // Get the views
1767     DataArrayView view_1 = tmp_1->getDefaultValue();
1768     DataArrayView view_2 = tmp_2->getDefaultValue();
1769     // Get the pointers to the actual data
1770     double *ptr_1 = &((view_1.getData())[0]);
1771     double *ptr_2 = &((view_2.getData())[0]);
1772     // Compute a result for the default
1773     tensor_binary_operation(size0, ptr_0, ptr_1, ptr_2, operation);
1774     // Compute a result for each tag
1775     const DataTagged::DataMapType& lookup_1=tmp_1->getTagLookup();
1776     DataTagged::DataMapType::const_iterator i; // i->first is a tag, i->second is an offset into memory
1777     for (i=lookup_1.begin();i!=lookup_1.end();i++) {
1778 matt 1332 tmp_2->addTaggedValue(i->first,tmp_2->getDefaultValue());
1779     DataArrayView view_1 = tmp_1->getDataPointByTag(i->first);
1780     DataArrayView view_2 = tmp_2->getDataPointByTag(i->first);
1781     double *ptr_1 = &view_1.getData(0);
1782     double *ptr_2 = &view_2.getData(0);
1783     tensor_binary_operation(size0, ptr_0, ptr_1, ptr_2, operation);
1784 matt 1327 }
1785    
1786     }
1787     else if (arg_0_Z.isConstant() && arg_1_Z.isExpanded()) {
1788    
1789     res = Data(0.0, shape0, arg_1_Z.getFunctionSpace(),true); // DataExpanded output
1790     DataConstant* tmp_0=dynamic_cast<DataConstant*>(arg_0_Z.borrowData());
1791     DataExpanded* tmp_1=dynamic_cast<DataExpanded*>(arg_1_Z.borrowData());
1792     DataExpanded* tmp_2=dynamic_cast<DataExpanded*>(res.borrowData());
1793    
1794     int sampleNo_1,dataPointNo_1;
1795     int numSamples_1 = arg_1_Z.getNumSamples();
1796     int numDataPointsPerSample_1 = arg_1_Z.getNumDataPointsPerSample();
1797     int offset_0 = tmp_0->getPointOffset(0,0);
1798     #pragma omp parallel for private(sampleNo_1,dataPointNo_1) schedule(static)
1799     for (sampleNo_1 = 0; sampleNo_1 < numSamples_1; sampleNo_1++) {
1800 matt 1332 for (dataPointNo_1 = 0; dataPointNo_1 < numDataPointsPerSample_1; dataPointNo_1++) {
1801     int offset_1 = tmp_1->getPointOffset(sampleNo_1,dataPointNo_1);
1802     int offset_2 = tmp_2->getPointOffset(sampleNo_1,dataPointNo_1);
1803     double *ptr_0 = &((arg_0_Z.getPointDataView().getData())[offset_0]);
1804     double *ptr_1 = &((arg_1_Z.getPointDataView().getData())[offset_1]);
1805     double *ptr_2 = &((res.getPointDataView().getData())[offset_2]);
1806     tensor_binary_operation(size0, ptr_0, ptr_1, ptr_2, operation);
1807     }
1808 matt 1327 }
1809    
1810     }
1811     else if (arg_0_Z.isTagged() && arg_1_Z.isConstant()) {
1812    
1813     // Borrow DataTagged input from Data object
1814     DataTagged* tmp_0=dynamic_cast<DataTagged*>(arg_0_Z.borrowData());
1815    
1816     // Prepare the DataConstant input
1817     DataConstant* tmp_1=dynamic_cast<DataConstant*>(arg_1_Z.borrowData());
1818    
1819     // Prepare a DataTagged output 2
1820 matt 1332 res = Data(0.0, shape0, arg_0_Z.getFunctionSpace()); // DataTagged output
1821 matt 1327 res.tag();
1822     DataTagged* tmp_2=dynamic_cast<DataTagged*>(res.borrowData());
1823    
1824     // Prepare offset into DataConstant
1825     int offset_1 = tmp_1->getPointOffset(0,0);
1826     double *ptr_1 = &((arg_1_Z.getPointDataView().getData())[offset_1]);
1827     // Get the views
1828     DataArrayView view_0 = tmp_0->getDefaultValue();
1829     DataArrayView view_2 = tmp_2->getDefaultValue();
1830     // Get the pointers to the actual data
1831     double *ptr_0 = &((view_0.getData())[0]);
1832     double *ptr_2 = &((view_2.getData())[0]);
1833     // Compute a result for the default
1834     tensor_binary_operation(size0, ptr_0, ptr_1, ptr_2, operation);
1835     // Compute a result for each tag
1836     const DataTagged::DataMapType& lookup_0=tmp_0->getTagLookup();
1837     DataTagged::DataMapType::const_iterator i; // i->first is a tag, i->second is an offset into memory
1838     for (i=lookup_0.begin();i!=lookup_0.end();i++) {
1839 matt 1332 tmp_2->addTaggedValue(i->first,tmp_2->getDefaultValue());
1840     DataArrayView view_0 = tmp_0->getDataPointByTag(i->first);
1841     DataArrayView view_2 = tmp_2->getDataPointByTag(i->first);
1842     double *ptr_0 = &view_0.getData(0);
1843     double *ptr_2 = &view_2.getData(0);
1844     tensor_binary_operation(size0, ptr_0, ptr_1, ptr_2, operation);
1845 matt 1327 }
1846    
1847     }
1848     else if (arg_0_Z.isTagged() && arg_1_Z.isTagged()) {
1849    
1850     // Borrow DataTagged input from Data object
1851     DataTagged* tmp_0=dynamic_cast<DataTagged*>(arg_0_Z.borrowData());
1852    
1853     // Borrow DataTagged input from Data object
1854     DataTagged* tmp_1=dynamic_cast<DataTagged*>(arg_1_Z.borrowData());
1855    
1856     // Prepare a DataTagged output 2
1857     res = Data(0.0, shape0, arg_1_Z.getFunctionSpace());
1858 matt 1332 res.tag(); // DataTagged output
1859 matt 1327 DataTagged* tmp_2=dynamic_cast<DataTagged*>(res.borrowData());
1860    
1861     // Get the views
1862     DataArrayView view_0 = tmp_0->getDefaultValue();
1863     DataArrayView view_1 = tmp_1->getDefaultValue();
1864     DataArrayView view_2 = tmp_2->getDefaultValue();
1865     // Get the pointers to the actual data
1866     double *ptr_0 = &((view_0.getData())[0]);
1867     double *ptr_1 = &((view_1.getData())[0]);
1868     double *ptr_2 = &((view_2.getData())[0]);
1869     // Compute a result for the default
1870     tensor_binary_operation(size0, ptr_0, ptr_1, ptr_2, operation);
1871     // Merge the tags
1872     DataTagged::DataMapType::const_iterator i; // i->first is a tag, i->second is an offset into memory
1873     const DataTagged::DataMapType& lookup_0=tmp_0->getTagLookup();
1874     const DataTagged::DataMapType& lookup_1=tmp_1->getTagLookup();
1875     for (i=lookup_0.begin();i!=lookup_0.end();i++) {
1876 matt 1332 tmp_2->addTaggedValue(i->first,tmp_2->getDefaultValue()); // use tmp_2 to get correct shape
1877 matt 1327 }
1878     for (i=lookup_1.begin();i!=lookup_1.end();i++) {
1879 matt 1332 tmp_2->addTaggedValue(i->first,tmp_2->getDefaultValue());
1880 matt 1327 }
1881     // Compute a result for each tag
1882     const DataTagged::DataMapType& lookup_2=tmp_2->getTagLookup();
1883     for (i=lookup_2.begin();i!=lookup_2.end();i++) {
1884 matt 1332 DataArrayView view_0 = tmp_0->getDataPointByTag(i->first);
1885     DataArrayView view_1 = tmp_1->getDataPointByTag(i->first);
1886     DataArrayView view_2 = tmp_2->getDataPointByTag(i->first);
1887     double *ptr_0 = &view_0.getData(0);
1888     double *ptr_1 = &view_1.getData(0);
1889     double *ptr_2 = &view_2.getData(0);
1890     tensor_binary_operation(size0, ptr_0, ptr_1, ptr_2, operation);
1891 matt 1327 }
1892    
1893     }
1894     else if (arg_0_Z.isTagged() && arg_1_Z.isExpanded()) {
1895    
1896     // After finding a common function space above the two inputs have the same numSamples and num DPPS
1897     res = Data(0.0, shape0, arg_1_Z.getFunctionSpace(),true); // DataExpanded output
1898     DataTagged* tmp_0=dynamic_cast<DataTagged*>(arg_0_Z.borrowData());
1899     DataExpanded* tmp_1=dynamic_cast<DataExpanded*>(arg_1_Z.borrowData());
1900     DataExpanded* tmp_2=dynamic_cast<DataExpanded*>(res.borrowData());
1901    
1902     int sampleNo_0,dataPointNo_0;
1903     int numSamples_0 = arg_0_Z.getNumSamples();
1904     int numDataPointsPerSample_0 = arg_0_Z.getNumDataPointsPerSample();
1905     #pragma omp parallel for private(sampleNo_0,dataPointNo_0) schedule(static)
1906     for (sampleNo_0 = 0; sampleNo_0 < numSamples_0; sampleNo_0++) {
1907 matt 1332 int offset_0 = tmp_0->getPointOffset(sampleNo_0,0); // They're all the same, so just use #0
1908     double *ptr_0 = &((arg_0_Z.getPointDataView().getData())[offset_0]);
1909     for (dataPointNo_0 = 0; dataPointNo_0 < numDataPointsPerSample_0; dataPointNo_0++) {
1910     int offset_1 = tmp_1->getPointOffset(sampleNo_0,dataPointNo_0);
1911     int offset_2 = tmp_2->getPointOffset(sampleNo_0,dataPointNo_0);
1912     double *ptr_1 = &((arg_1_Z.getPointDataView().getData())[offset_1]);
1913     double *ptr_2 = &((res.getPointDataView().getData())[offset_2]);
1914     tensor_binary_operation(size0, ptr_0, ptr_1, ptr_2, operation);
1915     }
1916 matt 1327 }
1917    
1918     }
1919     else if (arg_0_Z.isExpanded() && arg_1_Z.isConstant()) {
1920    
1921     res = Data(0.0, shape0, arg_1_Z.getFunctionSpace(),true); // DataExpanded output
1922     DataExpanded* tmp_0=dynamic_cast<DataExpanded*>(arg_0_Z.borrowData());
1923     DataConstant* tmp_1=dynamic_cast<DataConstant*>(arg_1_Z.borrowData());
1924     DataExpanded* tmp_2=dynamic_cast<DataExpanded*>(res.borrowData());
1925    
1926     int sampleNo_0,dataPointNo_0;
1927     int numSamples_0 = arg_0_Z.getNumSamples();
1928     int numDataPointsPerSample_0 = arg_0_Z.getNumDataPointsPerSample();
1929     int offset_1 = tmp_1->getPointOffset(0,0);
1930     #pragma omp parallel for private(sampleNo_0,dataPointNo_0) schedule(static)
1931     for (sampleNo_0 = 0; sampleNo_0 < numSamples_0; sampleNo_0++) {
1932 matt 1332 for (dataPointNo_0 = 0; dataPointNo_0 < numDataPointsPerSample_0; dataPointNo_0++) {
1933     int offset_0 = tmp_0->getPointOffset(sampleNo_0,dataPointNo_0);
1934     int offset_2 = tmp_2->getPointOffset(sampleNo_0,dataPointNo_0);
1935     double *ptr_0 = &((arg_0_Z.getPointDataView().getData())[offset_0]);
1936     double *ptr_1 = &((arg_1_Z.getPointDataView().getData())[offset_1]);
1937     double *ptr_2 = &((res.getPointDataView().getData())[offset_2]);
1938     tensor_binary_operation(size0, ptr_0, ptr_1, ptr_2, operation);
1939     }
1940 matt 1327 }
1941    
1942     }
1943     else if (arg_0_Z.isExpanded() && arg_1_Z.isTagged()) {
1944    
1945     // After finding a common function space above the two inputs have the same numSamples and num DPPS
1946     res = Data(0.0, shape0, arg_1_Z.getFunctionSpace(),true); // DataExpanded output
1947     DataExpanded* tmp_0=dynamic_cast<DataExpanded*>(arg_0_Z.borrowData());
1948     DataTagged* tmp_1=dynamic_cast<DataTagged*>(arg_1_Z.borrowData());
1949     DataExpanded* tmp_2=dynamic_cast<DataExpanded*>(res.borrowData());
1950    
1951     int sampleNo_0,dataPointNo_0;
1952     int numSamples_0 = arg_0_Z.getNumSamples();
1953     int numDataPointsPerSample_0 = arg_0_Z.getNumDataPointsPerSample();
1954     #pragma omp parallel for private(sampleNo_0,dataPointNo_0) schedule(static)
1955     for (sampleNo_0 = 0; sampleNo_0 < numSamples_0; sampleNo_0++) {
1956 matt 1332 int offset_1 = tmp_1->getPointOffset(sampleNo_0,0);
1957     double *ptr_1 = &((arg_1_Z.getPointDataView().getData())[offset_1]);
1958     for (dataPointNo_0 = 0; dataPointNo_0 < numDataPointsPerSample_0; dataPointNo_0++) {
1959     int offset_0 = tmp_0->getPointOffset(sampleNo_0,dataPointNo_0);
1960     int offset_2 = tmp_2->getPointOffset(sampleNo_0,dataPointNo_0);
1961     double *ptr_0 = &((arg_0_Z.getPointDataView().getData())[offset_0]);
1962     double *ptr_2 = &((res.getPointDataView().getData())[offset_2]);
1963     tensor_binary_operation(size0, ptr_0, ptr_1, ptr_2, operation);
1964     }
1965 matt 1327 }
1966    
1967     }
1968     else if (arg_0_Z.isExpanded() && arg_1_Z.isExpanded()) {
1969    
1970     // After finding a common function space above the two inputs have the same numSamples and num DPPS
1971     res = Data(0.0, shape0, arg_1_Z.getFunctionSpace(),true); // DataExpanded output
1972     DataExpanded* tmp_0=dynamic_cast<DataExpanded*>(arg_0_Z.borrowData());
1973     DataExpanded* tmp_1=dynamic_cast<DataExpanded*>(arg_1_Z.borrowData());
1974     DataExpanded* tmp_2=dynamic_cast<DataExpanded*>(res.borrowData());
1975    
1976     int sampleNo_0,dataPointNo_0;
1977     int numSamples_0 = arg_0_Z.getNumSamples();
1978     int numDataPointsPerSample_0 = arg_0_Z.getNumDataPointsPerSample();
1979     #pragma omp parallel for private(sampleNo_0,dataPointNo_0) schedule(static)
1980     for (sampleNo_0 = 0; sampleNo_0 < numSamples_0; sampleNo_0++) {
1981 matt 1332 for (dataPointNo_0 = 0; dataPointNo_0 < numDataPointsPerSample_0; dataPointNo_0++) {
1982     int offset_0 = tmp_0->getPointOffset(sampleNo_0,dataPointNo_0);
1983     int offset_1 = tmp_1->getPointOffset(sampleNo_0,dataPointNo_0);
1984     int offset_2 = tmp_2->getPointOffset(sampleNo_0,dataPointNo_0);
1985     double *ptr_0 = &((arg_0_Z.getPointDataView().getData())[offset_0]);
1986     double *ptr_1 = &((arg_1_Z.getPointDataView().getData())[offset_1]);
1987     double *ptr_2 = &((res.getPointDataView().getData())[offset_2]);
1988     tensor_binary_operation(size0, ptr_0, ptr_1, ptr_2, operation);
1989     }
1990 matt 1327 }
1991    
1992     }
1993     else {
1994     throw DataException("Error - C_TensorBinaryOperation: unknown combination of inputs");
1995     }
1996    
1997     } else if (0 == rank0) {
1998    
1999     if (arg_0_Z.isConstant() && arg_1_Z.isConstant()) {
2000 matt 1332 res = Data(0.0, shape1, arg_1_Z.getFunctionSpace()); // DataConstant output
2001 matt 1327 double *ptr_0 = &((arg_0_Z.getPointDataView().getData())[0]);
2002     double *ptr_1 = &((arg_1_Z.getPointDataView().getData())[0]);
2003     double *ptr_2 = &((res.getPointDataView().getData())[0]);
2004     tensor_binary_operation(size1, ptr_0[0], ptr_1, ptr_2, operation);
2005     }
2006     else if (arg_0_Z.isConstant() && arg_1_Z.isTagged()) {
2007    
2008     // Prepare the DataConstant input
2009     DataConstant* tmp_0=dynamic_cast<DataConstant*>(arg_0_Z.borrowData());
2010    
2011     // Borrow DataTagged input from Data object
2012     DataTagged* tmp_1=dynamic_cast<DataTagged*>(arg_1_Z.borrowData());
2013    
2014     // Prepare a DataTagged output 2
2015 matt 1332 res = Data(0.0, shape1, arg_1_Z.getFunctionSpace()); // DataTagged output
2016 matt 1327 res.tag();
2017     DataTagged* tmp_2=dynamic_cast<DataTagged*>(res.borrowData());
2018    
2019     // Prepare offset into DataConstant
2020     int offset_0 = tmp_0->getPointOffset(0,0);
2021     double *ptr_0 = &((arg_0_Z.getPointDataView().getData())[offset_0]);
2022     // Get the views
2023     DataArrayView view_1 = tmp_1->getDefaultValue();
2024     DataArrayView view_2 = tmp_2->getDefaultValue();
2025     // Get the pointers to the actual data
2026     double *ptr_1 = &((view_1.getData())[0]);
2027     double *ptr_2 = &((view_2.getData())[0]);
2028     // Compute a result for the default
2029     tensor_binary_operation(size1, ptr_0[0], ptr_1, ptr_2, operation);
2030     // Compute a result for each tag
2031     const DataTagged::DataMapType& lookup_1=tmp_1->getTagLookup();
2032     DataTagged::DataMapType::const_iterator i; // i->first is a tag, i->second is an offset into memory
2033     for (i=lookup_1.begin();i!=lookup_1.end();i++) {
2034 matt 1332 tmp_2->addTaggedValue(i->first,tmp_2->getDefaultValue());
2035     DataArrayView view_1 = tmp_1->getDataPointByTag(i->first);
2036     DataArrayView view_2 = tmp_2->getDataPointByTag(i->first);
2037     double *ptr_1 = &view_1.getData(0);
2038     double *ptr_2 = &view_2.getData(0);
2039     tensor_binary_operation(size1, ptr_0[0], ptr_1, ptr_2, operation);
2040 matt 1327 }
2041    
2042     }
2043     else if (arg_0_Z.isConstant() && arg_1_Z.isExpanded()) {
2044    
2045     res = Data(0.0, shape1, arg_1_Z.getFunctionSpace(),true); // DataExpanded output
2046     DataConstant* tmp_0=dynamic_cast<DataConstant*>(arg_0_Z.borrowData());
2047     DataExpanded* tmp_1=dynamic_cast<DataExpanded*>(arg_1_Z.borrowData());
2048     DataExpanded* tmp_2=dynamic_cast<DataExpanded*>(res.borrowData());
2049    
2050     int sampleNo_1,dataPointNo_1;
2051     int numSamples_1 = arg_1_Z.getNumSamples();
2052     int numDataPointsPerSample_1 = arg_1_Z.getNumDataPointsPerSample();
2053     int offset_0 = tmp_0->getPointOffset(0,0);
2054     #pragma omp parallel for private(sampleNo_1,dataPointNo_1) schedule(static)
2055     for (sampleNo_1 = 0; sampleNo_1 < numSamples_1; sampleNo_1++) {
2056 matt 1332 for (dataPointNo_1 = 0; dataPointNo_1 < numDataPointsPerSample_1; dataPointNo_1++) {
2057     int offset_1 = tmp_1->getPointOffset(sampleNo_1,dataPointNo_1);
2058     int offset_2 = tmp_2->getPointOffset(sampleNo_1,dataPointNo_1);
2059     double *ptr_0 = &((arg_0_Z.getPointDataView().getData())[offset_0]);
2060     double *ptr_1 = &((arg_1_Z.getPointDataView().getData())[offset_1]);
2061     double *ptr_2 = &((res.getPointDataView().getData())[offset_2]);
2062     tensor_binary_operation(size1, ptr_0[0], ptr_1, ptr_2, operation);
2063 matt 1327
2064 matt 1332 }
2065 matt 1327 }
2066    
2067     }
2068     else if (arg_0_Z.isTagged() && arg_1_Z.isConstant()) {
2069    
2070     // Borrow DataTagged input from Data object
2071     DataTagged* tmp_0=dynamic_cast<DataTagged*>(arg_0_Z.borrowData());
2072    
2073     // Prepare the DataConstant input
2074     DataConstant* tmp_1=dynamic_cast<DataConstant*>(arg_1_Z.borrowData());
2075    
2076     // Prepare a DataTagged output 2
2077 matt 1332 res = Data(0.0, shape1, arg_0_Z.getFunctionSpace()); // DataTagged output
2078 matt 1327 res.tag();
2079     DataTagged* tmp_2=dynamic_cast<DataTagged*>(res.borrowData());
2080    
2081     // Prepare offset into DataConstant
2082     int offset_1 = tmp_1->getPointOffset(0,0);
2083     double *ptr_1 = &((arg_1_Z.getPointDataView().getData())[offset_1]);
2084     // Get the views
2085     DataArrayView view_0 = tmp_0->getDefaultValue();
2086     DataArrayView view_2 = tmp_2->getDefaultValue();
2087     // Get the pointers to the actual data
2088     double *ptr_0 = &((view_0.getData())[0]);
2089     double *ptr_2 = &((view_2.getData())[0]);
2090     // Compute a result for the default
2091     tensor_binary_operation(size1, ptr_0[0], ptr_1, ptr_2, operation);
2092     // Compute a result for each tag
2093     const DataTagged::DataMapType& lookup_0=tmp_0->getTagLookup();
2094     DataTagged::DataMapType::const_iterator i; // i->first is a tag, i->second is an offset into memory
2095     for (i=lookup_0.begin();i!=lookup_0.end();i++) {
2096 matt 1332 tmp_2->addTaggedValue(i->first,tmp_2->getDefaultValue());
2097     DataArrayView view_0 = tmp_0->getDataPointByTag(i->first);
2098     DataArrayView view_2 = tmp_2->getDataPointByTag(i->first);
2099     double *ptr_0 = &view_0.getData(0);
2100     double *ptr_2 = &view_2.getData(0);
2101     tensor_binary_operation(size1, ptr_0[0], ptr_1, ptr_2, operation);
2102 matt 1327 }
2103    
2104     }
2105     else if (arg_0_Z.isTagged() && arg_1_Z.isTagged()) {
2106    
2107     // Borrow DataTagged input from Data object
2108     DataTagged* tmp_0=dynamic_cast<DataTagged*>(arg_0_Z.borrowData());
2109    
2110     // Borrow DataTagged input from Data object
2111     DataTagged* tmp_1=dynamic_cast<DataTagged*>(arg_1_Z.borrowData());
2112    
2113     // Prepare a DataTagged output 2
2114     res = Data(0.0, shape1, arg_1_Z.getFunctionSpace());
2115 matt 1332 res.tag(); // DataTagged output
2116 matt 1327 DataTagged* tmp_2=dynamic_cast<DataTagged*>(res.borrowData());
2117    
2118     // Get the views
2119     DataArrayView view_0 = tmp_0->getDefaultValue();
2120     DataArrayView view_1 = tmp_1->getDefaultValue();
2121     DataArrayView view_2 = tmp_2->getDefaultValue();
2122     // Get the pointers to the actual data
2123     double *ptr_0 = &((view_0.getData())[0]);
2124     double *ptr_1 = &((view_1.getData())[0]);
2125     double *ptr_2 = &((view_2.getData())[0]);
2126     // Compute a result for the default
2127     tensor_binary_operation(size1, ptr_0[0], ptr_1, ptr_2, operation);
2128     // Merge the tags
2129     DataTagged::DataMapType::const_iterator i; // i->first is a tag, i->second is an offset into memory
2130     const DataTagged::DataMapType& lookup_0=tmp_0->getTagLookup();
2131     const DataTagged::DataMapType& lookup_1=tmp_1->getTagLookup();
2132     for (i=lookup_0.begin();i!=lookup_0.end();i++) {
2133 matt 1332 tmp_2->addTaggedValue(i->first,tmp_2->getDefaultValue()); // use tmp_2 to get correct shape
2134 matt 1327 }
2135     for (i=lookup_1.begin();i!=lookup_1.end();i++) {
2136 matt 1332 tmp_2->addTaggedValue(i->first,tmp_2->getDefaultValue());
2137 matt 1327 }
2138     // Compute a result for each tag
2139     const DataTagged::DataMapType& lookup_2=tmp_2->getTagLookup();
2140     for (i=lookup_2.begin();i!=lookup_2.end();i++) {
2141 matt 1332 DataArrayView view_0 = tmp_0->getDataPointByTag(i->first);
2142     DataArrayView view_1 = tmp_1->getDataPointByTag(i->first);
2143     DataArrayView view_2 = tmp_2->getDataPointByTag(i->first);
2144     double *ptr_0 = &view_0.getData(0);
2145     double *ptr_1 = &view_1.getData(0);
2146     double *ptr_2 = &view_2.getData(0);
2147     tensor_binary_operation(size1, ptr_0[0], ptr_1, ptr_2, operation);
2148 matt 1327 }
2149    
2150     }
2151     else if (arg_0_Z.isTagged() && arg_1_Z.isExpanded()) {
2152    
2153     // After finding a common function space above the two inputs have the same numSamples and num DPPS
2154     res = Data(0.0, shape1, arg_1_Z.getFunctionSpace(),true); // DataExpanded output
2155     DataTagged* tmp_0=dynamic_cast<DataTagged*>(arg_0_Z.borrowData());
2156     DataExpanded* tmp_1=dynamic_cast<DataExpanded*>(arg_1_Z.borrowData());
2157     DataExpanded* tmp_2=dynamic_cast<DataExpanded*>(res.borrowData());
2158    
2159     int sampleNo_0,dataPointNo_0;
2160     int numSamples_0 = arg_0_Z.getNumSamples();
2161     int numDataPointsPerSample_0 = arg_0_Z.getNumDataPointsPerSample();
2162     #pragma omp parallel for private(sampleNo_0,dataPointNo_0) schedule(static)
2163     for (sampleNo_0 = 0; sampleNo_0 < numSamples_0; sampleNo_0++) {
2164 matt 1332 int offset_0 = tmp_0->getPointOffset(sampleNo_0,0); // They're all the same, so just use #0
2165     double *ptr_0 = &((arg_0_Z.getPointDataView().getData())[offset_0]);
2166     for (dataPointNo_0 = 0; dataPointNo_0 < numDataPointsPerSample_0; dataPointNo_0++) {
2167     int offset_1 = tmp_1->getPointOffset(sampleNo_0,dataPointNo_0);
2168     int offset_2 = tmp_2->getPointOffset(sampleNo_0,dataPointNo_0);
2169     double *ptr_1 = &((arg_1_Z.getPointDataView().getData())[offset_1]);
2170     double *ptr_2 = &((res.getPointDataView().getData())[offset_2]);
2171     tensor_binary_operation(size1, ptr_0[0], ptr_1, ptr_2, operation);
2172     }
2173 matt 1327 }
2174    
2175     }
2176     else if (arg_0_Z.isExpanded() && arg_1_Z.isConstant()) {
2177    
2178     res = Data(0.0, shape1, arg_1_Z.getFunctionSpace(),true); // DataExpanded output
2179     DataExpanded* tmp_0=dynamic_cast<DataExpanded*>(arg_0_Z.borrowData());
2180     DataConstant* tmp_1=dynamic_cast<DataConstant*>(arg_1_Z.borrowData());
2181     DataExpanded* tmp_2=dynamic_cast<DataExpanded*>(res.borrowData());
2182    
2183     int sampleNo_0,dataPointNo_0;
2184     int numSamples_0 = arg_0_Z.getNumSamples();
2185     int numDataPointsPerSample_0 = arg_0_Z.getNumDataPointsPerSample();
2186     int offset_1 = tmp_1->getPointOffset(0,0);
2187     #pragma omp parallel for private(sampleNo_0,dataPointNo_0) schedule(static)
2188     for (sampleNo_0 = 0; sampleNo_0 < numSamples_0; sampleNo_0++) {
2189 matt 1332 for (dataPointNo_0 = 0; dataPointNo_0 < numDataPointsPerSample_0; dataPointNo_0++) {
2190     int offset_0 = tmp_0->getPointOffset(sampleNo_0,dataPointNo_0);
2191     int offset_2 = tmp_2->getPointOffset(sampleNo_0,dataPointNo_0);
2192     double *ptr_0 = &((arg_0_Z.getPointDataView().getData())[offset_0]);
2193     double *ptr_1 = &((arg_1_Z.getPointDataView().getData())[offset_1]);
2194     double *ptr_2 = &((res.getPointDataView().getData())[offset_2]);
2195     tensor_binary_operation(size1, ptr_0[0], ptr_1, ptr_2, operation);
2196     }
2197 matt 1327 }
2198    
2199    
2200     }
2201     else if (arg_0_Z.isExpanded() && arg_1_Z.isTagged()) {
2202    
2203     // After finding a common function space above the two inputs have the same numSamples and num DPPS
2204     res = Data(0.0, shape1, arg_1_Z.getFunctionSpace(),true); // DataExpanded output
2205     DataExpanded* tmp_0=dynamic_cast<DataExpanded*>(arg_0_Z.borrowData());
2206     DataTagged* tmp_1=dynamic_cast<DataTagged*>(arg_1_Z.borrowData());
2207     DataExpanded* tmp_2=dynamic_cast<DataExpanded*>(res.borrowData());
2208    
2209     int sampleNo_0,dataPointNo_0;
2210     int numSamples_0 = arg_0_Z.getNumSamples();
2211     int numDataPointsPerSample_0 = arg_0_Z.getNumDataPointsPerSample();
2212     #pragma omp parallel for private(sampleNo_0,dataPointNo_0) schedule(static)
2213     for (sampleNo_0 = 0; sampleNo_0 < numSamples_0; sampleNo_0++) {
2214 matt 1332 int offset_1 = tmp_1->getPointOffset(sampleNo_0,0);
2215     double *ptr_1 = &((arg_1_Z.getPointDataView().getData())[offset_1]);
2216     for (dataPointNo_0 = 0; dataPointNo_0 < numDataPointsPerSample_0; dataPointNo_0++) {
2217     int offset_0 = tmp_0->getPointOffset(sampleNo_0,dataPointNo_0);
2218     int offset_2 = tmp_2->getPointOffset(sampleNo_0,dataPointNo_0);
2219     double *ptr_0 = &((arg_0_Z.getPointDataView().getData())[offset_0]);
2220     double *ptr_2 = &((res.getPointDataView().getData())[offset_2]);
2221     tensor_binary_operation(size1, ptr_0[0], ptr_1, ptr_2, operation);
2222     }
2223 matt 1327 }
2224    
2225     }
2226     else if (arg_0_Z.isExpanded() && arg_1_Z.isExpanded()) {
2227    
2228     // After finding a common function space above the two inputs have the same numSamples and num DPPS
2229     res = Data(0.0, shape1, arg_1_Z.getFunctionSpace(),true); // DataExpanded output
2230     DataExpanded* tmp_0=dynamic_cast<DataExpanded*>(arg_0_Z.borrowData());
2231     DataExpanded* tmp_1=dynamic_cast<DataExpanded*>(arg_1_Z.borrowData());
2232     DataExpanded* tmp_2=dynamic_cast<DataExpanded*>(res.borrowData());
2233    
2234     int sampleNo_0,dataPointNo_0;
2235     int numSamples_0 = arg_0_Z.getNumSamples();
2236     int numDataPointsPerSample_0 = arg_0_Z.getNumDataPointsPerSample();
2237     #pragma omp parallel for private(sampleNo_0,dataPointNo_0) schedule(static)
2238     for (sampleNo_0 = 0; sampleNo_0 < numSamples_0; sampleNo_0++) {
2239 matt 1332 for (dataPointNo_0 = 0; dataPointNo_0 < numDataPointsPerSample_0; dataPointNo_0++) {
2240     int offset_0 = tmp_0->getPointOffset(sampleNo_0,dataPointNo_0);
2241     int offset_1 = tmp_1->getPointOffset(sampleNo_0,dataPointNo_0);
2242     int offset_2 = tmp_2->getPointOffset(sampleNo_0,dataPointNo_0);
2243     double *ptr_0 = &((arg_0_Z.getPointDataView().getData())[offset_0]);
2244     double *ptr_1 = &((arg_1_Z.getPointDataView().getData())[offset_1]);
2245     double *ptr_2 = &((res.getPointDataView().getData())[offset_2]);
2246     tensor_binary_operation(size1, ptr_0[0], ptr_1, ptr_2, operation);
2247     }
2248 matt 1327 }
2249    
2250     }
2251     else {
2252     throw DataException("Error - C_TensorBinaryOperation: unknown combination of inputs");
2253     }
2254    
2255     } else if (0 == rank1) {
2256    
2257     if (arg_0_Z.isConstant() && arg_1_Z.isConstant()) {
2258 matt 1332 res = Data(0.0, shape0, arg_1_Z.getFunctionSpace()); // DataConstant output
2259 matt 1327 double *ptr_0 = &((arg_0_Z.getPointDataView().getData())[0]);
2260     double *ptr_1 = &((arg_1_Z.getPointDataView().getData())[0]);
2261     double *ptr_2 = &((res.getPointDataView().getData())[0]);
2262     tensor_binary_operation(size0, ptr_0, ptr_1[0], ptr_2, operation);
2263     }
2264     else if (arg_0_Z.isConstant() && arg_1_Z.isTagged()) {
2265    
2266     // Prepare the DataConstant input
2267     DataConstant* tmp_0=dynamic_cast<DataConstant*>(arg_0_Z.borrowData());
2268    
2269     // Borrow DataTagged input from Data object
2270     DataTagged* tmp_1=dynamic_cast<DataTagged*>(arg_1_Z.borrowData());
2271    
2272     // Prepare a DataTagged output 2
2273 matt 1332 res = Data(0.0, shape0, arg_1_Z.getFunctionSpace()); // DataTagged output
2274 matt 1327 res.tag();
2275     DataTagged* tmp_2=dynamic_cast<DataTagged*>(res.borrowData());
2276    
2277     // Prepare offset into DataConstant
2278     int offset_0 = tmp_0->getPointOffset(0,0);
2279     double *ptr_0 = &((arg_0_Z.getPointDataView().getData())[offset_0]);
2280     // Get the views
2281     DataArrayView view_1 = tmp_1->getDefaultValue();
2282     DataArrayView view_2 = tmp_2->getDefaultValue();
2283     // Get the pointers to the actual data
2284     double *ptr_1 = &((view_1.getData())[0]);
2285     double *ptr_2 = &((view_2.getData())[0]);
2286     // Compute a result for the default
2287     tensor_binary_operation(size0, ptr_0, ptr_1[0], ptr_2, operation);
2288     // Compute a result for each tag
2289     const DataTagged::DataMapType& lookup_1=tmp_1->getTagLookup();
2290     DataTagged::DataMapType::const_iterator i; // i->first is a tag, i->second is an offset into memory
2291     for (i=lookup_1.begin();i!=lookup_1.end();i++) {
2292 matt 1332 tmp_2->addTaggedValue(i->first,tmp_2->getDefaultValue());
2293     DataArrayView view_1 = tmp_1->getDataPointByTag(i->first);
2294     DataArrayView view_2 = tmp_2->getDataPointByTag(i->first);
2295     double *ptr_1 = &view_1.getData(0);
2296     double *ptr_2 = &view_2.getData(0);
2297     tensor_binary_operation(size0, ptr_0, ptr_1[0], ptr_2, operation);
2298 matt 1327 }
2299    
2300     }
2301     else if (arg_0_Z.isConstant() && arg_1_Z.isExpanded()) {
2302    
2303     res = Data(0.0, shape0, arg_1_Z.getFunctionSpace(),true); // DataExpanded output
2304     DataConstant* tmp_0=dynamic_cast<DataConstant*>(arg_0_Z.borrowData());
2305     DataExpanded* tmp_1=dynamic_cast<DataExpanded*>(arg_1_Z.borrowData());
2306     DataExpanded* tmp_2=dynamic_cast<DataExpanded*>(res.borrowData());
2307    
2308     int sampleNo_1,dataPointNo_1;
2309     int numSamples_1 = arg_1_Z.getNumSamples();
2310     int numDataPointsPerSample_1 = arg_1_Z.getNumDataPointsPerSample();
2311     int offset_0 = tmp_0->getPointOffset(0,0);
2312     #pragma omp parallel for private(sampleNo_1,dataPointNo_1) schedule(static)
2313     for (sampleNo_1 = 0; sampleNo_1 < numSamples_1; sampleNo_1++) {
2314 matt 1332 for (dataPointNo_1 = 0; dataPointNo_1 < numDataPointsPerSample_1; dataPointNo_1++) {
2315     int offset_1 = tmp_1->getPointOffset(sampleNo_1,dataPointNo_1);
2316     int offset_2 = tmp_2->getPointOffset(sampleNo_1,dataPointNo_1);
2317     double *ptr_0 = &((arg_0_Z.getPointDataView().getData())[offset_0]);
2318     double *ptr_1 = &((arg_1_Z.getPointDataView().getData())[offset_1]);
2319     double *ptr_2 = &((res.getPointDataView().getData())[offset_2]);
2320     tensor_binary_operation(size0, ptr_0, ptr_1[0], ptr_2, operation);
2321     }
2322 matt 1327 }
2323    
2324     }
2325     else if (arg_0_Z.isTagged() && arg_1_Z.isConstant()) {
2326    
2327     // Borrow DataTagged input from Data object
2328     DataTagged* tmp_0=dynamic_cast<DataTagged*>(arg_0_Z.borrowData());
2329    
2330     // Prepare the DataConstant input
2331     DataConstant* tmp_1=dynamic_cast<DataConstant*>(arg_1_Z.borrowData());
2332    
2333     // Prepare a DataTagged output 2
2334 matt 1332 res = Data(0.0, shape0, arg_0_Z.getFunctionSpace()); // DataTagged output
2335 matt 1327 res.tag();
2336     DataTagged* tmp_2=dynamic_cast<DataTagged*>(res.borrowData());
2337    
2338     // Prepare offset into DataConstant
2339     int offset_1 = tmp_1->getPointOffset(0,0);
2340     double *ptr_1 = &((arg_1_Z.getPointDataView().getData())[offset_1]);
2341     // Get the views
2342     DataArrayView view_0 = tmp_0->getDefaultValue();
2343     DataArrayView view_2 = tmp_2->getDefaultValue();
2344     // Get the pointers to the actual data
2345     double *ptr_0 = &((view_0.getData())[0]);
2346     double *ptr_2 = &((view_2.getData())[0]);
2347     // Compute a result for the default
2348     tensor_binary_operation(size0, ptr_0, ptr_1[0], ptr_2, operation);
2349     // Compute a result for each tag
2350     const DataTagged::DataMapType& lookup_0=tmp_0->getTagLookup();
2351     DataTagged::DataMapType::const_iterator i; // i->first is a tag, i->second is an offset into memory
2352     for (i=lookup_0.begin();i!=lookup_0.end();i++) {
2353 matt 1332 tmp_2->addTaggedValue(i->first,tmp_2->getDefaultValue());
2354     DataArrayView view_0 = tmp_0->getDataPointByTag(i->first);
2355     DataArrayView view_2 = tmp_2->getDataPointByTag(i->first);
2356     double *ptr_0 = &view_0.getData(0);
2357     double *ptr_2 = &view_2.getData(0);
2358     tensor_binary_operation(size0, ptr_0, ptr_1[0], ptr_2, operation);
2359 matt 1327 }
2360    
2361     }
2362     else if (arg_0_Z.isTagged() && arg_1_Z.isTagged()) {
2363    
2364     // Borrow DataTagged input from Data object
2365     DataTagged* tmp_0=dynamic_cast<DataTagged*>(arg_0_Z.borrowData());
2366    
2367     // Borrow DataTagged input from Data object
2368     DataTagged* tmp_1=dynamic_cast<DataTagged*>(arg_1_Z.borrowData());
2369    
2370     // Prepare a DataTagged output 2
2371     res = Data(0.0, shape0, arg_1_Z.getFunctionSpace());
2372 matt 1332 res.tag(); // DataTagged output
2373 matt 1327 DataTagged* tmp_2=dynamic_cast<DataTagged*>(res.borrowData());
2374    
2375     // Get the views
2376     DataArrayView view_0 = tmp_0->getDefaultValue();
2377     DataArrayView view_1 = tmp_1->getDefaultValue();
2378     DataArrayView view_2 = tmp_2->getDefaultValue();
2379     // Get the pointers to the actual data
2380     double *ptr_0 = &((view_0.getData())[0]);
2381     double *ptr_1 = &((view_1.getData())[0]);
2382     double *ptr_2 = &((view_2.getData())[0]);
2383     // Compute a result for the default
2384     tensor_binary_operation(size0, ptr_0, ptr_1[0], ptr_2, operation);
2385     // Merge the tags
2386     DataTagged::DataMapType::const_iterator i; // i->first is a tag, i->second is an offset into memory
2387     const DataTagged::DataMapType& lookup_0=tmp_0->getTagLookup();
2388     const DataTagged::DataMapType& lookup_1=tmp_1->getTagLookup();
2389     for (i=lookup_0.begin();i!=lookup_0.end();i++) {
2390 matt 1332 tmp_2->addTaggedValue(i->first,tmp_2->getDefaultValue()); // use tmp_2 to get correct shape
2391 matt 1327 }
2392     for (i=lookup_1.begin();i!=lookup_1.end();i++) {
2393 matt 1332 tmp_2->addTaggedValue(i->first,tmp_2->getDefaultValue());
2394 matt 1327 }
2395     // Compute a result for each tag
2396     const DataTagged::DataMapType& lookup_2=tmp_2->getTagLookup();
2397     for (i=lookup_2.begin();i!=lookup_2.end();i++) {
2398 matt 1332 DataArrayView view_0 = tmp_0->getDataPointByTag(i->first);
2399     DataArrayView view_1 = tmp_1->getDataPointByTag(i->first);
2400     DataArrayView view_2 = tmp_2->getDataPointByTag(i->first);
2401     double *ptr_0 = &view_0.getData(0);
2402     double *ptr_1 = &view_1.getData(0);
2403     double *ptr_2 = &view_2.getData(0);
2404     tensor_binary_operation(size0, ptr_0, ptr_1[0], ptr_2, operation);
2405 matt 1327 }
2406    
2407     }
2408     else if (arg_0_Z.isTagged() && arg_1_Z.isExpanded()) {
2409    
2410     // After finding a common function space above the two inputs have the same numSamples and num DPPS
2411     res = Data(0.0, shape0, arg_1_Z.getFunctionSpace(),true); // DataExpanded output
2412     DataTagged* tmp_0=dynamic_cast<DataTagged*>(arg_0_Z.borrowData());
2413     DataExpanded* tmp_1=dynamic_cast<DataExpanded*>(arg_1_Z.borrowData());
2414     DataExpanded* tmp_2=dynamic_cast<DataExpanded*>(res.borrowData());
2415    
2416     int sampleNo_0,dataPointNo_0;
2417     int numSamples_0 = arg_0_Z.getNumSamples();
2418     int numDataPointsPerSample_0 = arg_0_Z.getNumDataPointsPerSample();
2419     #pragma omp parallel for private(sampleNo_0,dataPointNo_0) schedule(static)
2420     for (sampleNo_0 = 0; sampleNo_0 < numSamples_0; sampleNo_0++) {
2421 matt 1332 int offset_0 = tmp_0->getPointOffset(sampleNo_0,0); // They're all the same, so just use #0
2422     double *ptr_0 = &((arg_0_Z.getPointDataView().getData())[offset_0]);
2423     for (dataPointNo_0 = 0; dataPointNo_0 < numDataPointsPerSample_0; dataPointNo_0++) {
2424     int offset_1 = tmp_1->getPointOffset(sampleNo_0,dataPointNo_0);
2425     int offset_2 = tmp_2->getPointOffset(sampleNo_0,dataPointNo_0);
2426     double *ptr_1 = &((arg_1_Z.getPointDataView().getData())[offset_1]);
2427     double *ptr_2 = &((res.getPointDataView().getData())[offset_2]);
2428     tensor_binary_operation(size0, ptr_0, ptr_1[0], ptr_2, operation);
2429     }
2430 matt 1327 }
2431    
2432     }
2433     else if (arg_0_Z.isExpanded() && arg_1_Z.isConstant()) {
2434    
2435     res = Data(0.0, shape0, arg_1_Z.getFunctionSpace(),true); // DataExpanded output
2436     DataExpanded* tmp_0=dynamic_cast<DataExpanded*>(arg_0_Z.borrowData());
2437     DataConstant* tmp_1=dynamic_cast<DataConstant*>(arg_1_Z.borrowData());
2438     DataExpanded* tmp_2=dynamic_cast<DataExpanded*>(res.borrowData());
2439    
2440     int sampleNo_0,dataPointNo_0;
2441     int numSamples_0 = arg_0_Z.getNumSamples();
2442     int numDataPointsPerSample_0 = arg_0_Z.getNumDataPointsPerSample();
2443     int offset_1 = tmp_1->getPointOffset(0,0);
2444     #pragma omp parallel for private(sampleNo_0,dataPointNo_0) schedule(static)
2445     for (sampleNo_0 = 0; sampleNo_0 < numSamples_0; sampleNo_0++) {
2446 matt 1332 for (dataPointNo_0 = 0; dataPointNo_0 < numDataPointsPerSample_0; dataPointNo_0++) {
2447     int offset_0 = tmp_0->getPointOffset(sampleNo_0,dataPointNo_0);
2448     int offset_2 = tmp_2->getPointOffset(sampleNo_0,dataPointNo_0);
2449     double *ptr_0 = &((arg_0_Z.getPointDataView().getData())[offset_0]);
2450     double *ptr_1 = &((arg_1_Z.getPointDataView().getData())[offset_1]);
2451     double *ptr_2 = &((res.getPointDataView().getData())[offset_2]);
2452     tensor_binary_operation(size0, ptr_0, ptr_1[0], ptr_2, operation);
2453     }
2454 matt 1327 }
2455    
2456    
2457     }
2458     else if (arg_0_Z.isExpanded() && arg_1_Z.isTagged()) {
2459    
2460     // After finding a common function space above the two inputs have the same numSamples and num DPPS
2461     res = Data(0.0, shape0, arg_1_Z.getFunctionSpace(),true); // DataExpanded output
2462     DataExpanded* tmp_0=dynamic_cast<DataExpanded*>(arg_0_Z.borrowData());
2463     DataTagged* tmp_1=dynamic_cast<DataTagged*>(arg_1_Z.borrowData());
2464     DataExpanded* tmp_2=dynamic_cast<DataExpanded*>(res.borrowData());
2465    
2466     int sampleNo_0,dataPointNo_0;
2467     int numSamples_0 = arg_0_Z.getNumSamples();
2468     int numDataPointsPerSample_0 = arg_0_Z.getNumDataPointsPerSample();
2469     #pragma omp parallel for private(sampleNo_0,dataPointNo_0) schedule(static)
2470     for (sampleNo_0 = 0; sampleNo_0 < numSamples_0; sampleNo_0++) {
2471 matt 1332 int offset_1 = tmp_1->getPointOffset(sampleNo_0,0);
2472     double *ptr_1 = &((arg_1_Z.getPointDataView().getData())[offset_1]);
2473     for (dataPointNo_0 = 0; dataPointNo_0 < numDataPointsPerSample_0; dataPointNo_0++) {
2474     int offset_0 = tmp_0->getPointOffset(sampleNo_0,dataPointNo_0);
2475     int offset_2 = tmp_2->getPointOffset(sampleNo_0,dataPointNo_0);
2476     double *ptr_0 = &((arg_0_Z.getPointDataView().getData())[offset_0]);
2477     double *ptr_2 = &((res.getPointDataView().getData())[offset_2]);
2478     tensor_binary_operation(size0, ptr_0, ptr_1[0], ptr_2, operation);
2479     }
2480 matt 1327 }
2481    
2482     }
2483     else if (arg_0_Z.isExpanded() && arg_1_Z.isExpanded()) {
2484    
2485     // After finding a common function space above the two inputs have the same numSamples and num DPPS
2486     res = Data(0.0, shape0, arg_1_Z.getFunctionSpace(),true); // DataExpanded output
2487     DataExpanded* tmp_0=dynamic_cast<DataExpanded*>(arg_0_Z.borrowData());
2488     DataExpanded* tmp_1=dynamic_cast<DataExpanded*>(arg_1_Z.borrowData());
2489     DataExpanded* tmp_2=dynamic_cast<DataExpanded*>(res.borrowData());
2490    
2491     int sampleNo_0,dataPointNo_0;
2492     int numSamples_0 = arg_0_Z.getNumSamples();
2493     int numDataPointsPerSample_0 = arg_0_Z.getNumDataPointsPerSample();
2494     #pragma omp parallel for private(sampleNo_0,dataPointNo_0) schedule(static)
2495     for (sampleNo_0 = 0; sampleNo_0 < numSamples_0; sampleNo_0++) {
2496 matt 1332 for (dataPointNo_0 = 0; dataPointNo_0 < numDataPointsPerSample_0; dataPointNo_0++) {
2497     int offset_0 = tmp_0->getPointOffset(sampleNo_0,dataPointNo_0);
2498     int offset_1 = tmp_1->getPointOffset(sampleNo_0,dataPointNo_0);
2499     int offset_2 = tmp_2->getPointOffset(sampleNo_0,dataPointNo_0);
2500     double *ptr_0 = &((arg_0_Z.getPointDataView().getData())[offset_0]);
2501     double *ptr_1 = &((arg_1_Z.getPointDataView().getData())[offset_1]);
2502     double *ptr_2 = &((res.getPointDataView().getData())[offset_2]);
2503     tensor_binary_operation(size0, ptr_0, ptr_1[0], ptr_2, operation);
2504     }
2505 matt 1327 }
2506    
2507     }
2508     else {
2509     throw DataException("Error - C_TensorBinaryOperation: unknown combination of inputs");
2510     }
2511    
2512     } else {
2513     throw DataException("Error - C_TensorBinaryOperation: arguments have incompatible shapes");
2514     }
2515    
2516     return res;
2517 jgs 94 }
2518 matt 1327
2519 matt 1334 template <typename UnaryFunction>
2520     Data
2521     C_TensorUnaryOperation(Data const &arg_0,
2522     UnaryFunction operation)
2523     {
2524     // Interpolate if necessary and find an appropriate function space
2525     Data arg_0_Z = Data(arg_0);
2526    
2527     // Get rank and shape of inputs
2528     int rank0 = arg_0_Z.getDataPointRank();
2529     DataArrayView::ShapeType shape0 = arg_0_Z.getDataPointShape();
2530     int size0 = arg_0_Z.getDataPointSize();
2531    
2532     // Declare output Data object
2533     Data res;
2534    
2535     if (arg_0_Z.isConstant()) {
2536     res = Data(0.0, shape0, arg_0_Z.getFunctionSpace()); // DataConstant output
2537     double *ptr_0 = &((arg_0_Z.getPointDataView().getData())[0]);
2538     double *ptr_2 = &((res.getPointDataView().getData())[0]);
2539     tensor_unary_operation(size0, ptr_0, ptr_2, operation);
2540     }
2541     else if (arg_0_Z.isTagged()) {
2542    
2543     // Borrow DataTagged input from Data object
2544     DataTagged* tmp_0=dynamic_cast<DataTagged*>(arg_0_Z.borrowData());
2545    
2546     // Prepare a DataTagged output 2
2547     res = Data(0.0, shape0, arg_0_Z.getFunctionSpace()); // DataTagged output
2548     res.tag();
2549     DataTagged* tmp_2=dynamic_cast<DataTagged*>(res.borrowData());
2550    
2551     // Get the views
2552     DataArrayView view_0 = tmp_0->getDefaultValue();
2553     DataArrayView view_2 = tmp_2->getDefaultValue();
2554     // Get the pointers to the actual data
2555     double *ptr_0 = &((view_0.getData())[0]);
2556     double *ptr_2 = &((view_2.getData())[0]);
2557     // Compute a result for the default
2558     tensor_unary_operation(size0, ptr_0, ptr_2, operation);
2559     // Compute a result for each tag
2560     const DataTagged::DataMapType& lookup_0=tmp_0->getTagLookup();
2561     DataTagged::DataMapType::const_iterator i; // i->first is a tag, i->second is an offset into memory
2562     for (i=lookup_0.begin();i!=lookup_0.end();i++) {
2563     tmp_2->addTaggedValue(i->first,tmp_2->getDefaultValue());
2564     DataArrayView view_0 = tmp_0->getDataPointByTag(i->first);
2565     DataArrayView view_2 = tmp_2->getDataPointByTag(i->first);
2566     double *ptr_0 = &view_0.getData(0);
2567     double *ptr_2 = &view_2.getData(0);
2568     tensor_unary_operation(size0, ptr_0, ptr_2, operation);
2569     }
2570    
2571     }
2572     else if (arg_0_Z.isExpanded()) {
2573    
2574     res = Data(0.0, shape0, arg_0_Z.getFunctionSpace(),true); // DataExpanded output
2575     DataExpanded* tmp_0=dynamic_cast<DataExpanded*>(arg_0_Z.borrowData());
2576     DataExpanded* tmp_2=dynamic_cast<DataExpanded*>(res.borrowData());
2577    
2578     int sampleNo_0,dataPointNo_0;
2579     int numSamples_0 = arg_0_Z.getNumSamples();
2580     int numDataPointsPerSample_0 = arg_0_Z.getNumDataPointsPerSample();
2581     #pragma omp parallel for private(sampleNo_0,dataPointNo_0) schedule(static)
2582     for (sampleNo_0 = 0; sampleNo_0 < numSamples_0; sampleNo_0++) {
2583     for (dataPointNo_0 = 0; dataPointNo_0 < numDataPointsPerSample_0; dataPointNo_0++) {
2584     int offset_0 = tmp_0->getPointOffset(sampleNo_0,dataPointNo_0);
2585     int offset_2 = tmp_2->getPointOffset(sampleNo_0,dataPointNo_0);
2586     double *ptr_0 = &((arg_0_Z.getPointDataView().getData())[offset_0]);
2587     double *ptr_2 = &((res.getPointDataView().getData())[offset_2]);
2588     tensor_unary_operation(size0, ptr_0, ptr_2, operation);
2589     }
2590     }
2591    
2592     }
2593     else {
2594     throw DataException("Error - C_TensorUnaryOperation: unknown combination of inputs");
2595     }
2596    
2597     return res;
2598 matt 1327 }
2599 matt 1334
2600     }
2601 jgs 94 #endif

Properties

Name Value
svn:eol-style native
svn:keywords Author Date Id Revision

  ViewVC Help
Powered by ViewVC 1.1.26