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

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

Parent Directory Parent Directory | Revision Log Revision Log


Revision 1332 - (hide annotations)
Tue Oct 23 03:28:51 2007 UTC (12 years, 1 month ago) by matt
File MIME type: text/plain
File size: 82298 byte(s)
Pow now uses the new binary function interface of C_TensorBinaryOperation.

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     object to type DataTagged. Throws an exception if this object
592     cannot be converted to a DataTagged object.
593     \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     object to type DataTagged. Throws an exception if this object
606     cannot be converted to a DataTagged object.
607     \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 jgs 102 unaryOp(UnaryFunction operation);
1205    
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 unary operation on other and return the result.
1628     Given operation is performed on each element of each data point, thus
1629     argument object is a rank n Data object, and returned object is a rank n
1630     Data object.
1631     Calls Data::unaryOp.
1632     */
1633 jgs 94 template <class UnaryFunction>
1634     inline
1635 jgs 102 Data
1636     unaryOp(const Data& other,
1637     UnaryFunction operation)
1638     {
1639     Data result;
1640     result.copy(other);
1641     result.unaryOp(operation);
1642     return result;
1643     }
1644    
1645     /**
1646     \brief
1647     Perform the given unary operation on this.
1648     Given operation is performed on each element of each data point.
1649     Calls escript::unaryOp.
1650     */
1651     template <class UnaryFunction>
1652     inline
1653 jgs 94 void
1654     Data::unaryOp(UnaryFunction operation)
1655     {
1656     if (isExpanded()) {
1657     DataExpanded* leftC=dynamic_cast<DataExpanded*>(m_data.get());
1658     EsysAssert((leftC!=0), "Programming error - casting to DataExpanded.");
1659     escript::unaryOp(*leftC,operation);
1660     } else if (isTagged()) {
1661     DataTagged* leftC=dynamic_cast<DataTagged*>(m_data.get());
1662     EsysAssert((leftC!=0), "Programming error - casting to DataTagged.");
1663     escript::unaryOp(*leftC,operation);
1664 jgs 102 } else if (isConstant()) {
1665 jgs 94 DataConstant* leftC=dynamic_cast<DataConstant*>(m_data.get());
1666     EsysAssert((leftC!=0), "Programming error - casting to DataConstant.");
1667     escript::unaryOp(*leftC,operation);
1668     }
1669     }
1670    
1671 jgs 102 /**
1672     \brief
1673     Perform the given Data object reduction algorithm on this and return the result.
1674     Given operation combines each element of each data point, thus argument
1675     object (*this) is a rank n Data object, and returned object is a scalar.
1676     Calls escript::algorithm.
1677     */
1678 jgs 147 template <class BinaryFunction>
1679 jgs 94 inline
1680     double
1681 jgs 147 Data::algorithm(BinaryFunction operation, double initial_value) const
1682 jgs 94 {
1683     if (isExpanded()) {
1684     DataExpanded* leftC=dynamic_cast<DataExpanded*>(m_data.get());
1685     EsysAssert((leftC!=0), "Programming error - casting to DataExpanded.");
1686 jgs 147 return escript::algorithm(*leftC,operation,initial_value);
1687 jgs 102 } else if (isTagged()) {
1688 jgs 94 DataTagged* leftC=dynamic_cast<DataTagged*>(m_data.get());
1689     EsysAssert((leftC!=0), "Programming error - casting to DataTagged.");
1690 jgs 147 return escript::algorithm(*leftC,operation,initial_value);
1691 jgs 102 } else if (isConstant()) {
1692 jgs 94 DataConstant* leftC=dynamic_cast<DataConstant*>(m_data.get());
1693     EsysAssert((leftC!=0), "Programming error - casting to DataConstant.");
1694 jgs 147 return escript::algorithm(*leftC,operation,initial_value);
1695 jgs 94 }
1696 jgs 102 return 0;
1697 jgs 94 }
1698    
1699 jgs 102 /**
1700     \brief
1701     Perform the given data point reduction algorithm on data and return the result.
1702     Given operation combines each element within each data point into a scalar,
1703 matt 1327 thus argument object is a rank n Data object, and returned object is a
1704 jgs 102 rank 0 Data object.
1705     Calls escript::dp_algorithm.
1706     */
1707 jgs 147 template <class BinaryFunction>
1708 jgs 94 inline
1709     Data
1710 jgs 147 Data::dp_algorithm(BinaryFunction operation, double initial_value) const
1711 jgs 94 {
1712 jgs 106 if (isExpanded()) {
1713 jgs 559 Data result(0,DataArrayView::ShapeType(),getFunctionSpace(),isExpanded());
1714 jgs 106 DataExpanded* dataE=dynamic_cast<DataExpanded*>(m_data.get());
1715 jgs 102 DataExpanded* resultE=dynamic_cast<DataExpanded*>(result.m_data.get());
1716     EsysAssert((dataE!=0), "Programming error - casting data to DataExpanded.");
1717     EsysAssert((resultE!=0), "Programming error - casting result to DataExpanded.");
1718 jgs 147 escript::dp_algorithm(*dataE,*resultE,operation,initial_value);
1719 jgs 559 return result;
1720 jgs 106 } else if (isTagged()) {
1721     DataTagged* dataT=dynamic_cast<DataTagged*>(m_data.get());
1722 jgs 562 DataArrayView::ShapeType viewShape;
1723     DataArrayView::ValueType viewData(1);
1724     viewData[0]=0;
1725     DataArrayView defaultValue(viewData,viewShape);
1726 jgs 559 DataTagged::TagListType keys;
1727 jgs 562 DataTagged::ValueListType values;
1728 jgs 559 DataTagged::DataMapType::const_iterator i;
1729     for (i=dataT->getTagLookup().begin();i!=dataT->getTagLookup().end();i++) {
1730     keys.push_back(i->first);
1731 jgs 562 values.push_back(defaultValue);
1732 jgs 559 }
1733     Data result(keys,values,defaultValue,getFunctionSpace());
1734 jgs 102 DataTagged* resultT=dynamic_cast<DataTagged*>(result.m_data.get());
1735     EsysAssert((dataT!=0), "Programming error - casting data to DataTagged.");
1736     EsysAssert((resultT!=0), "Programming error - casting result to DataTagged.");
1737 jgs 147 escript::dp_algorithm(*dataT,*resultT,operation,initial_value);
1738 jgs 559 return result;
1739 jgs 106 } else if (isConstant()) {
1740 jgs 559 Data result(0,DataArrayView::ShapeType(),getFunctionSpace(),isExpanded());
1741 jgs 106 DataConstant* dataC=dynamic_cast<DataConstant*>(m_data.get());
1742 jgs 102 DataConstant* resultC=dynamic_cast<DataConstant*>(result.m_data.get());
1743     EsysAssert((dataC!=0), "Programming error - casting data to DataConstant.");
1744     EsysAssert((resultC!=0), "Programming error - casting result to DataConstant.");
1745 jgs 147 escript::dp_algorithm(*dataC,*resultC,operation,initial_value);
1746 jgs 559 return result;
1747 jgs 102 }
1748 jgs 559 Data falseRetVal; // to keep compiler quiet
1749     return falseRetVal;
1750 jgs 94 }
1751    
1752 matt 1327 template <typename BinaryFunction>
1753     Data
1754     C_TensorBinaryOperation(Data const &arg_0,
1755 matt 1332 Data const &arg_1,
1756     BinaryFunction operation)
1757 matt 1327 {
1758     // Interpolate if necessary and find an appropriate function space
1759     Data arg_0_Z, arg_1_Z;
1760     if (arg_0.getFunctionSpace()!=arg_1.getFunctionSpace()) {
1761     if (arg_0.probeInterpolation(arg_1.getFunctionSpace())) {
1762     arg_0_Z = arg_0.interpolate(arg_1.getFunctionSpace());
1763     arg_1_Z = Data(arg_1);
1764     }
1765     else if (arg_1.probeInterpolation(arg_0.getFunctionSpace())) {
1766     arg_1_Z=arg_1.interpolate(arg_0.getFunctionSpace());
1767     arg_0_Z =Data(arg_0);
1768     }
1769     else {
1770     throw DataException("Error - C_TensorBinaryOperation: arguments have incompatible function spaces.");
1771     }
1772     } else {
1773     arg_0_Z = Data(arg_0);
1774     arg_1_Z = Data(arg_1);
1775     }
1776     // Get rank and shape of inputs
1777     int rank0 = arg_0_Z.getDataPointRank();
1778     int rank1 = arg_1_Z.getDataPointRank();
1779     DataArrayView::ShapeType shape0 = arg_0_Z.getDataPointShape();
1780     DataArrayView::ShapeType shape1 = arg_1_Z.getDataPointShape();
1781     int size0 = arg_0_Z.getDataPointSize();
1782     int size1 = arg_1_Z.getDataPointSize();
1783    
1784     // Declare output Data object
1785     Data res;
1786    
1787     if (shape0 == shape1) {
1788    
1789     if (arg_0_Z.isConstant() && arg_1_Z.isConstant()) {
1790 matt 1332 res = Data(0.0, shape0, arg_1_Z.getFunctionSpace()); // DataConstant output
1791 matt 1327 double *ptr_0 = &((arg_0_Z.getPointDataView().getData())[0]);
1792     double *ptr_1 = &((arg_1_Z.getPointDataView().getData())[0]);
1793     double *ptr_2 = &((res.getPointDataView().getData())[0]);
1794     tensor_binary_operation(size0, ptr_0, ptr_1, ptr_2, operation);
1795     }
1796     else if (arg_0_Z.isConstant() && arg_1_Z.isTagged()) {
1797    
1798     // Prepare the DataConstant input
1799     DataConstant* tmp_0=dynamic_cast<DataConstant*>(arg_0_Z.borrowData());
1800    
1801     // Borrow DataTagged input from Data object
1802     DataTagged* tmp_1=dynamic_cast<DataTagged*>(arg_1_Z.borrowData());
1803    
1804     // Prepare a DataTagged output 2
1805 matt 1332 res = Data(0.0, shape0, arg_1_Z.getFunctionSpace()); // DataTagged output
1806 matt 1327 res.tag();
1807     DataTagged* tmp_2=dynamic_cast<DataTagged*>(res.borrowData());
1808    
1809     // Prepare offset into DataConstant
1810     int offset_0 = tmp_0->getPointOffset(0,0);
1811     double *ptr_0 = &((arg_0_Z.getPointDataView().getData())[offset_0]);
1812     // Get the views
1813     DataArrayView view_1 = tmp_1->getDefaultValue();
1814     DataArrayView view_2 = tmp_2->getDefaultValue();
1815     // Get the pointers to the actual data
1816     double *ptr_1 = &((view_1.getData())[0]);
1817     double *ptr_2 = &((view_2.getData())[0]);
1818     // Compute a result for the default
1819     tensor_binary_operation(size0, ptr_0, ptr_1, ptr_2, operation);
1820     // Compute a result for each tag
1821     const DataTagged::DataMapType& lookup_1=tmp_1->getTagLookup();
1822     DataTagged::DataMapType::const_iterator i; // i->first is a tag, i->second is an offset into memory
1823     for (i=lookup_1.begin();i!=lookup_1.end();i++) {
1824 matt 1332 tmp_2->addTaggedValue(i->first,tmp_2->getDefaultValue());
1825     DataArrayView view_1 = tmp_1->getDataPointByTag(i->first);
1826     DataArrayView view_2 = tmp_2->getDataPointByTag(i->first);
1827     double *ptr_1 = &view_1.getData(0);
1828     double *ptr_2 = &view_2.getData(0);
1829     tensor_binary_operation(size0, ptr_0, ptr_1, ptr_2, operation);
1830 matt 1327 }
1831    
1832     }
1833     else if (arg_0_Z.isConstant() && arg_1_Z.isExpanded()) {
1834    
1835     res = Data(0.0, shape0, arg_1_Z.getFunctionSpace(),true); // DataExpanded output
1836     DataConstant* tmp_0=dynamic_cast<DataConstant*>(arg_0_Z.borrowData());
1837     DataExpanded* tmp_1=dynamic_cast<DataExpanded*>(arg_1_Z.borrowData());
1838     DataExpanded* tmp_2=dynamic_cast<DataExpanded*>(res.borrowData());
1839    
1840     int sampleNo_1,dataPointNo_1;
1841     int numSamples_1 = arg_1_Z.getNumSamples();
1842     int numDataPointsPerSample_1 = arg_1_Z.getNumDataPointsPerSample();
1843     int offset_0 = tmp_0->getPointOffset(0,0);
1844     #pragma omp parallel for private(sampleNo_1,dataPointNo_1) schedule(static)
1845     for (sampleNo_1 = 0; sampleNo_1 < numSamples_1; sampleNo_1++) {
1846 matt 1332 for (dataPointNo_1 = 0; dataPointNo_1 < numDataPointsPerSample_1; dataPointNo_1++) {
1847     int offset_1 = tmp_1->getPointOffset(sampleNo_1,dataPointNo_1);
1848     int offset_2 = tmp_2->getPointOffset(sampleNo_1,dataPointNo_1);
1849     double *ptr_0 = &((arg_0_Z.getPointDataView().getData())[offset_0]);
1850     double *ptr_1 = &((arg_1_Z.getPointDataView().getData())[offset_1]);
1851     double *ptr_2 = &((res.getPointDataView().getData())[offset_2]);
1852     tensor_binary_operation(size0, ptr_0, ptr_1, ptr_2, operation);
1853     }
1854 matt 1327 }
1855    
1856     }
1857     else if (arg_0_Z.isTagged() && arg_1_Z.isConstant()) {
1858    
1859     // Borrow DataTagged input from Data object
1860     DataTagged* tmp_0=dynamic_cast<DataTagged*>(arg_0_Z.borrowData());
1861    
1862     // Prepare the DataConstant input
1863     DataConstant* tmp_1=dynamic_cast<DataConstant*>(arg_1_Z.borrowData());
1864    
1865     // Prepare a DataTagged output 2
1866 matt 1332 res = Data(0.0, shape0, arg_0_Z.getFunctionSpace()); // DataTagged output
1867 matt 1327 res.tag();
1868     DataTagged* tmp_2=dynamic_cast<DataTagged*>(res.borrowData());
1869    
1870     // Prepare offset into DataConstant
1871     int offset_1 = tmp_1->getPointOffset(0,0);
1872     double *ptr_1 = &((arg_1_Z.getPointDataView().getData())[offset_1]);
1873     // Get the views
1874     DataArrayView view_0 = tmp_0->getDefaultValue();
1875     DataArrayView view_2 = tmp_2->getDefaultValue();
1876     // Get the pointers to the actual data
1877     double *ptr_0 = &((view_0.getData())[0]);
1878     double *ptr_2 = &((view_2.getData())[0]);
1879     // Compute a result for the default
1880     tensor_binary_operation(size0, ptr_0, ptr_1, ptr_2, operation);
1881     // Compute a result for each tag
1882     const DataTagged::DataMapType& lookup_0=tmp_0->getTagLookup();
1883     DataTagged::DataMapType::const_iterator i; // i->first is a tag, i->second is an offset into memory
1884     for (i=lookup_0.begin();i!=lookup_0.end();i++) {
1885 matt 1332 tmp_2->addTaggedValue(i->first,tmp_2->getDefaultValue());
1886     DataArrayView view_0 = tmp_0->getDataPointByTag(i->first);
1887     DataArrayView view_2 = tmp_2->getDataPointByTag(i->first);
1888     double *ptr_0 = &view_0.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.isTagged()) {
1895    
1896     // Borrow DataTagged input from Data object
1897     DataTagged* tmp_0=dynamic_cast<DataTagged*>(arg_0_Z.borrowData());
1898    
1899     // Borrow DataTagged input from Data object
1900     DataTagged* tmp_1=dynamic_cast<DataTagged*>(arg_1_Z.borrowData());
1901    
1902     // Prepare a DataTagged output 2
1903     res = Data(0.0, shape0, arg_1_Z.getFunctionSpace());
1904 matt 1332 res.tag(); // DataTagged output
1905 matt 1327 DataTagged* tmp_2=dynamic_cast<DataTagged*>(res.borrowData());
1906    
1907     // Get the views
1908     DataArrayView view_0 = tmp_0->getDefaultValue();
1909     DataArrayView view_1 = tmp_1->getDefaultValue();
1910     DataArrayView view_2 = tmp_2->getDefaultValue();
1911     // Get the pointers to the actual data
1912     double *ptr_0 = &((view_0.getData())[0]);
1913     double *ptr_1 = &((view_1.getData())[0]);
1914     double *ptr_2 = &((view_2.getData())[0]);
1915     // Compute a result for the default
1916     tensor_binary_operation(size0, ptr_0, ptr_1, ptr_2, operation);
1917     // Merge the tags
1918     DataTagged::DataMapType::const_iterator i; // i->first is a tag, i->second is an offset into memory
1919     const DataTagged::DataMapType& lookup_0=tmp_0->getTagLookup();
1920     const DataTagged::DataMapType& lookup_1=tmp_1->getTagLookup();
1921     for (i=lookup_0.begin();i!=lookup_0.end();i++) {
1922 matt 1332 tmp_2->addTaggedValue(i->first,tmp_2->getDefaultValue()); // use tmp_2 to get correct shape
1923 matt 1327 }
1924     for (i=lookup_1.begin();i!=lookup_1.end();i++) {
1925 matt 1332 tmp_2->addTaggedValue(i->first,tmp_2->getDefaultValue());
1926 matt 1327 }
1927     // Compute a result for each tag
1928     const DataTagged::DataMapType& lookup_2=tmp_2->getTagLookup();
1929     for (i=lookup_2.begin();i!=lookup_2.end();i++) {
1930 matt 1332 DataArrayView view_0 = tmp_0->getDataPointByTag(i->first);
1931     DataArrayView view_1 = tmp_1->getDataPointByTag(i->first);
1932     DataArrayView view_2 = tmp_2->getDataPointByTag(i->first);
1933     double *ptr_0 = &view_0.getData(0);
1934     double *ptr_1 = &view_1.getData(0);
1935     double *ptr_2 = &view_2.getData(0);
1936     tensor_binary_operation(size0, ptr_0, ptr_1, ptr_2, operation);
1937 matt 1327 }
1938    
1939     }
1940     else if (arg_0_Z.isTagged() && arg_1_Z.isExpanded()) {
1941    
1942     // After finding a common function space above the two inputs have the same numSamples and num DPPS
1943     res = Data(0.0, shape0, arg_1_Z.getFunctionSpace(),true); // DataExpanded output
1944     DataTagged* tmp_0=dynamic_cast<DataTagged*>(arg_0_Z.borrowData());
1945     DataExpanded* tmp_1=dynamic_cast<DataExpanded*>(arg_1_Z.borrowData());
1946     DataExpanded* tmp_2=dynamic_cast<DataExpanded*>(res.borrowData());
1947    
1948     int sampleNo_0,dataPointNo_0;
1949     int numSamples_0 = arg_0_Z.getNumSamples();
1950     int numDataPointsPerSample_0 = arg_0_Z.getNumDataPointsPerSample();
1951     #pragma omp parallel for private(sampleNo_0,dataPointNo_0) schedule(static)
1952     for (sampleNo_0 = 0; sampleNo_0 < numSamples_0; sampleNo_0++) {
1953 matt 1332 int offset_0 = tmp_0->getPointOffset(sampleNo_0,0); // They're all the same, so just use #0
1954     double *ptr_0 = &((arg_0_Z.getPointDataView().getData())[offset_0]);
1955     for (dataPointNo_0 = 0; dataPointNo_0 < numDataPointsPerSample_0; dataPointNo_0++) {
1956     int offset_1 = tmp_1->getPointOffset(sampleNo_0,dataPointNo_0);
1957     int offset_2 = tmp_2->getPointOffset(sampleNo_0,dataPointNo_0);
1958     double *ptr_1 = &((arg_1_Z.getPointDataView().getData())[offset_1]);
1959     double *ptr_2 = &((res.getPointDataView().getData())[offset_2]);
1960     tensor_binary_operation(size0, ptr_0, ptr_1, ptr_2, operation);
1961     }
1962 matt 1327 }
1963    
1964     }
1965     else if (arg_0_Z.isExpanded() && arg_1_Z.isConstant()) {
1966    
1967     res = Data(0.0, shape0, arg_1_Z.getFunctionSpace(),true); // DataExpanded output
1968     DataExpanded* tmp_0=dynamic_cast<DataExpanded*>(arg_0_Z.borrowData());
1969     DataConstant* tmp_1=dynamic_cast<DataConstant*>(arg_1_Z.borrowData());
1970     DataExpanded* tmp_2=dynamic_cast<DataExpanded*>(res.borrowData());
1971    
1972     int sampleNo_0,dataPointNo_0;
1973     int numSamples_0 = arg_0_Z.getNumSamples();
1974     int numDataPointsPerSample_0 = arg_0_Z.getNumDataPointsPerSample();
1975     int offset_1 = tmp_1->getPointOffset(0,0);
1976     #pragma omp parallel for private(sampleNo_0,dataPointNo_0) schedule(static)
1977     for (sampleNo_0 = 0; sampleNo_0 < numSamples_0; sampleNo_0++) {
1978 matt 1332 for (dataPointNo_0 = 0; dataPointNo_0 < numDataPointsPerSample_0; dataPointNo_0++) {
1979     int offset_0 = tmp_0->getPointOffset(sampleNo_0,dataPointNo_0);
1980     int offset_2 = tmp_2->getPointOffset(sampleNo_0,dataPointNo_0);
1981     double *ptr_0 = &((arg_0_Z.getPointDataView().getData())[offset_0]);
1982     double *ptr_1 = &((arg_1_Z.getPointDataView().getData())[offset_1]);
1983     double *ptr_2 = &((res.getPointDataView().getData())[offset_2]);
1984     tensor_binary_operation(size0, ptr_0, ptr_1, ptr_2, operation);
1985     }
1986 matt 1327 }
1987    
1988     }
1989     else if (arg_0_Z.isExpanded() && arg_1_Z.isTagged()) {
1990    
1991     // After finding a common function space above the two inputs have the same numSamples and num DPPS
1992     res = Data(0.0, shape0, arg_1_Z.getFunctionSpace(),true); // DataExpanded output
1993     DataExpanded* tmp_0=dynamic_cast<DataExpanded*>(arg_0_Z.borrowData());
1994     DataTagged* tmp_1=dynamic_cast<DataTagged*>(arg_1_Z.borrowData());
1995     DataExpanded* tmp_2=dynamic_cast<DataExpanded*>(res.borrowData());
1996    
1997     int sampleNo_0,dataPointNo_0;
1998     int numSamples_0 = arg_0_Z.getNumSamples();
1999     int numDataPointsPerSample_0 = arg_0_Z.getNumDataPointsPerSample();
2000     #pragma omp parallel for private(sampleNo_0,dataPointNo_0) schedule(static)
2001     for (sampleNo_0 = 0; sampleNo_0 < numSamples_0; sampleNo_0++) {
2002 matt 1332 int offset_1 = tmp_1->getPointOffset(sampleNo_0,0);
2003     double *ptr_1 = &((arg_1_Z.getPointDataView().getData())[offset_1]);
2004     for (dataPointNo_0 = 0; dataPointNo_0 < numDataPointsPerSample_0; dataPointNo_0++) {
2005     int offset_0 = tmp_0->getPointOffset(sampleNo_0,dataPointNo_0);
2006     int offset_2 = tmp_2->getPointOffset(sampleNo_0,dataPointNo_0);
2007     double *ptr_0 = &((arg_0_Z.getPointDataView().getData())[offset_0]);
2008     double *ptr_2 = &((res.getPointDataView().getData())[offset_2]);
2009     tensor_binary_operation(size0, ptr_0, ptr_1, ptr_2, operation);
2010     }
2011 matt 1327 }
2012    
2013     }
2014     else if (arg_0_Z.isExpanded() && arg_1_Z.isExpanded()) {
2015    
2016     // After finding a common function space above the two inputs have the same numSamples and num DPPS
2017     res = Data(0.0, shape0, arg_1_Z.getFunctionSpace(),true); // DataExpanded output
2018     DataExpanded* tmp_0=dynamic_cast<DataExpanded*>(arg_0_Z.borrowData());
2019     DataExpanded* tmp_1=dynamic_cast<DataExpanded*>(arg_1_Z.borrowData());
2020     DataExpanded* tmp_2=dynamic_cast<DataExpanded*>(res.borrowData());
2021    
2022     int sampleNo_0,dataPointNo_0;
2023     int numSamples_0 = arg_0_Z.getNumSamples();
2024     int numDataPointsPerSample_0 = arg_0_Z.getNumDataPointsPerSample();
2025     #pragma omp parallel for private(sampleNo_0,dataPointNo_0) schedule(static)
2026     for (sampleNo_0 = 0; sampleNo_0 < numSamples_0; sampleNo_0++) {
2027 matt 1332 for (dataPointNo_0 = 0; dataPointNo_0 < numDataPointsPerSample_0; dataPointNo_0++) {
2028     int offset_0 = tmp_0->getPointOffset(sampleNo_0,dataPointNo_0);
2029     int offset_1 = tmp_1->getPointOffset(sampleNo_0,dataPointNo_0);
2030     int offset_2 = tmp_2->getPointOffset(sampleNo_0,dataPointNo_0);
2031     double *ptr_0 = &((arg_0_Z.getPointDataView().getData())[offset_0]);
2032     double *ptr_1 = &((arg_1_Z.getPointDataView().getData())[offset_1]);
2033     double *ptr_2 = &((res.getPointDataView().getData())[offset_2]);
2034     tensor_binary_operation(size0, ptr_0, ptr_1, ptr_2, operation);
2035     }
2036 matt 1327 }
2037    
2038     }
2039     else {
2040     throw DataException("Error - C_TensorBinaryOperation: unknown combination of inputs");
2041     }
2042    
2043     } else if (0 == rank0) {
2044    
2045     if (arg_0_Z.isConstant() && arg_1_Z.isConstant()) {
2046 matt 1332 res = Data(0.0, shape1, arg_1_Z.getFunctionSpace()); // DataConstant output
2047 matt 1327 double *ptr_0 = &((arg_0_Z.getPointDataView().getData())[0]);
2048     double *ptr_1 = &((arg_1_Z.getPointDataView().getData())[0]);
2049     double *ptr_2 = &((res.getPointDataView().getData())[0]);
2050     tensor_binary_operation(size1, ptr_0[0], ptr_1, ptr_2, operation);
2051     }
2052     else if (arg_0_Z.isConstant() && arg_1_Z.isTagged()) {
2053    
2054     // Prepare the DataConstant input
2055     DataConstant* tmp_0=dynamic_cast<DataConstant*>(arg_0_Z.borrowData());
2056    
2057     // Borrow DataTagged input from Data object
2058     DataTagged* tmp_1=dynamic_cast<DataTagged*>(arg_1_Z.borrowData());
2059    
2060     // Prepare a DataTagged output 2
2061 matt 1332 res = Data(0.0, shape1, arg_1_Z.getFunctionSpace()); // DataTagged output
2062 matt 1327 res.tag();
2063     DataTagged* tmp_2=dynamic_cast<DataTagged*>(res.borrowData());
2064    
2065     // Prepare offset into DataConstant
2066     int offset_0 = tmp_0->getPointOffset(0,0);
2067     double *ptr_0 = &((arg_0_Z.getPointDataView().getData())[offset_0]);
2068     // Get the views
2069     DataArrayView view_1 = tmp_1->getDefaultValue();
2070     DataArrayView view_2 = tmp_2->getDefaultValue();
2071     // Get the pointers to the actual data
2072     double *ptr_1 = &((view_1.getData())[0]);
2073     double *ptr_2 = &((view_2.getData())[0]);
2074     // Compute a result for the default
2075     tensor_binary_operation(size1, ptr_0[0], ptr_1, ptr_2, operation);
2076     // Compute a result for each tag
2077     const DataTagged::DataMapType& lookup_1=tmp_1->getTagLookup();
2078     DataTagged::DataMapType::const_iterator i; // i->first is a tag, i->second is an offset into memory
2079     for (i=lookup_1.begin();i!=lookup_1.end();i++) {
2080 matt 1332 tmp_2->addTaggedValue(i->first,tmp_2->getDefaultValue());
2081     DataArrayView view_1 = tmp_1->getDataPointByTag(i->first);
2082     DataArrayView view_2 = tmp_2->getDataPointByTag(i->first);
2083     double *ptr_1 = &view_1.getData(0);
2084     double *ptr_2 = &view_2.getData(0);
2085     tensor_binary_operation(size1, ptr_0[0], ptr_1, ptr_2, operation);
2086 matt 1327 }
2087    
2088     }
2089     else if (arg_0_Z.isConstant() && arg_1_Z.isExpanded()) {
2090    
2091     res = Data(0.0, shape1, arg_1_Z.getFunctionSpace(),true); // DataExpanded output
2092     DataConstant* tmp_0=dynamic_cast<DataConstant*>(arg_0_Z.borrowData());
2093     DataExpanded* tmp_1=dynamic_cast<DataExpanded*>(arg_1_Z.borrowData());
2094     DataExpanded* tmp_2=dynamic_cast<DataExpanded*>(res.borrowData());
2095    
2096     int sampleNo_1,dataPointNo_1;
2097     int numSamples_1 = arg_1_Z.getNumSamples();
2098     int numDataPointsPerSample_1 = arg_1_Z.getNumDataPointsPerSample();
2099     int offset_0 = tmp_0->getPointOffset(0,0);
2100     #pragma omp parallel for private(sampleNo_1,dataPointNo_1) schedule(static)
2101     for (sampleNo_1 = 0; sampleNo_1 < numSamples_1; sampleNo_1++) {
2102 matt 1332 for (dataPointNo_1 = 0; dataPointNo_1 < numDataPointsPerSample_1; dataPointNo_1++) {
2103     int offset_1 = tmp_1->getPointOffset(sampleNo_1,dataPointNo_1);
2104     int offset_2 = tmp_2->getPointOffset(sampleNo_1,dataPointNo_1);
2105     double *ptr_0 = &((arg_0_Z.getPointDataView().getData())[offset_0]);
2106     double *ptr_1 = &((arg_1_Z.getPointDataView().getData())[offset_1]);
2107     double *ptr_2 = &((res.getPointDataView().getData())[offset_2]);
2108     tensor_binary_operation(size1, ptr_0[0], ptr_1, ptr_2, operation);
2109 matt 1327
2110 matt 1332 }
2111 matt 1327 }
2112    
2113     }
2114     else if (arg_0_Z.isTagged() && arg_1_Z.isConstant()) {
2115    
2116     // Borrow DataTagged input from Data object
2117     DataTagged* tmp_0=dynamic_cast<DataTagged*>(arg_0_Z.borrowData());
2118    
2119     // Prepare the DataConstant input
2120     DataConstant* tmp_1=dynamic_cast<DataConstant*>(arg_1_Z.borrowData());
2121    
2122     // Prepare a DataTagged output 2
2123 matt 1332 res = Data(0.0, shape1, arg_0_Z.getFunctionSpace()); // DataTagged output
2124 matt 1327 res.tag();
2125     DataTagged* tmp_2=dynamic_cast<DataTagged*>(res.borrowData());
2126    
2127     // Prepare offset into DataConstant
2128     int offset_1 = tmp_1->getPointOffset(0,0);
2129     double *ptr_1 = &((arg_1_Z.getPointDataView().getData())[offset_1]);
2130     // Get the views
2131     DataArrayView view_0 = tmp_0->getDefaultValue();
2132     DataArrayView view_2 = tmp_2->getDefaultValue();
2133     // Get the pointers to the actual data
2134     double *ptr_0 = &((view_0.getData())[0]);
2135     double *ptr_2 = &((view_2.getData())[0]);
2136     // Compute a result for the default
2137     tensor_binary_operation(size1, ptr_0[0], ptr_1, ptr_2, operation);
2138     // Compute a result for each tag
2139     const DataTagged::DataMapType& lookup_0=tmp_0->getTagLookup();
2140     DataTagged::DataMapType::const_iterator i; // i->first is a tag, i->second is an offset into memory
2141     for (i=lookup_0.begin();i!=lookup_0.end();i++) {
2142 matt 1332 tmp_2->addTaggedValue(i->first,tmp_2->getDefaultValue());
2143     DataArrayView view_0 = tmp_0->getDataPointByTag(i->first);
2144     DataArrayView view_2 = tmp_2->getDataPointByTag(i->first);
2145     double *ptr_0 = &view_0.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.isTagged()) {
2152    
2153     // Borrow DataTagged input from Data object
2154     DataTagged* tmp_0=dynamic_cast<DataTagged*>(arg_0_Z.borrowData());
2155    
2156     // Borrow DataTagged input from Data object
2157     DataTagged* tmp_1=dynamic_cast<DataTagged*>(arg_1_Z.borrowData());
2158    
2159     // Prepare a DataTagged output 2
2160     res = Data(0.0, shape1, arg_1_Z.getFunctionSpace());
2161 matt 1332 res.tag(); // DataTagged output
2162 matt 1327 DataTagged* tmp_2=dynamic_cast<DataTagged*>(res.borrowData());
2163    
2164     // Get the views
2165     DataArrayView view_0 = tmp_0->getDefaultValue();
2166     DataArrayView view_1 = tmp_1->getDefaultValue();
2167     DataArrayView view_2 = tmp_2->getDefaultValue();
2168     // Get the pointers to the actual data
2169     double *ptr_0 = &((view_0.getData())[0]);
2170     double *ptr_1 = &((view_1.getData())[0]);
2171     double *ptr_2 = &((view_2.getData())[0]);
2172     // Compute a result for the default
2173     tensor_binary_operation(size1, ptr_0[0], ptr_1, ptr_2, operation);
2174     // Merge the tags
2175     DataTagged::DataMapType::const_iterator i; // i->first is a tag, i->second is an offset into memory
2176     const DataTagged::DataMapType& lookup_0=tmp_0->getTagLookup();
2177     const DataTagged::DataMapType& lookup_1=tmp_1->getTagLookup();
2178     for (i=lookup_0.begin();i!=lookup_0.end();i++) {
2179 matt 1332 tmp_2->addTaggedValue(i->first,tmp_2->getDefaultValue()); // use tmp_2 to get correct shape
2180 matt 1327 }
2181     for (i=lookup_1.begin();i!=lookup_1.end();i++) {
2182 matt 1332 tmp_2->addTaggedValue(i->first,tmp_2->getDefaultValue());
2183 matt 1327 }
2184     // Compute a result for each tag
2185     const DataTagged::DataMapType& lookup_2=tmp_2->getTagLookup();
2186     for (i=lookup_2.begin();i!=lookup_2.end();i++) {
2187 matt 1332 DataArrayView view_0 = tmp_0->getDataPointByTag(i->first);
2188     DataArrayView view_1 = tmp_1->getDataPointByTag(i->first);
2189     DataArrayView view_2 = tmp_2->getDataPointByTag(i->first);
2190     double *ptr_0 = &view_0.getData(0);
2191     double *ptr_1 = &view_1.getData(0);
2192     double *ptr_2 = &view_2.getData(0);
2193     tensor_binary_operation(size1, ptr_0[0], ptr_1, ptr_2, operation);
2194 matt 1327 }
2195    
2196     }
2197     else if (arg_0_Z.isTagged() && arg_1_Z.isExpanded()) {
2198    
2199     // After finding a common function space above the two inputs have the same numSamples and num DPPS
2200     res = Data(0.0, shape1, arg_1_Z.getFunctionSpace(),true); // DataExpanded output
2201     DataTagged* tmp_0=dynamic_cast<DataTagged*>(arg_0_Z.borrowData());
2202     DataExpanded* tmp_1=dynamic_cast<DataExpanded*>(arg_1_Z.borrowData());
2203     DataExpanded* tmp_2=dynamic_cast<DataExpanded*>(res.borrowData());
2204    
2205     int sampleNo_0,dataPointNo_0;
2206     int numSamples_0 = arg_0_Z.getNumSamples();
2207     int numDataPointsPerSample_0 = arg_0_Z.getNumDataPointsPerSample();
2208     #pragma omp parallel for private(sampleNo_0,dataPointNo_0) schedule(static)
2209     for (sampleNo_0 = 0; sampleNo_0 < numSamples_0; sampleNo_0++) {
2210 matt 1332 int offset_0 = tmp_0->getPointOffset(sampleNo_0,0); // They're all the same, so just use #0
2211     double *ptr_0 = &((arg_0_Z.getPointDataView().getData())[offset_0]);
2212     for (dataPointNo_0 = 0; dataPointNo_0 < numDataPointsPerSample_0; dataPointNo_0++) {
2213     int offset_1 = tmp_1->getPointOffset(sampleNo_0,dataPointNo_0);
2214     int offset_2 = tmp_2->getPointOffset(sampleNo_0,dataPointNo_0);
2215     double *ptr_1 = &((arg_1_Z.getPointDataView().getData())[offset_1]);
2216     double *ptr_2 = &((res.getPointDataView().getData())[offset_2]);
2217     tensor_binary_operation(size1, ptr_0[0], ptr_1, ptr_2, operation);
2218     }
2219 matt 1327 }
2220    
2221     }
2222     else if (arg_0_Z.isExpanded() && arg_1_Z.isConstant()) {
2223    
2224     res = Data(0.0, shape1, arg_1_Z.getFunctionSpace(),true); // DataExpanded output
2225     DataExpanded* tmp_0=dynamic_cast<DataExpanded*>(arg_0_Z.borrowData());
2226     DataConstant* tmp_1=dynamic_cast<DataConstant*>(arg_1_Z.borrowData());
2227     DataExpanded* tmp_2=dynamic_cast<DataExpanded*>(res.borrowData());
2228    
2229     int sampleNo_0,dataPointNo_0;
2230     int numSamples_0 = arg_0_Z.getNumSamples();
2231     int numDataPointsPerSample_0 = arg_0_Z.getNumDataPointsPerSample();
2232     int offset_1 = tmp_1->getPointOffset(0,0);
2233     #pragma omp parallel for private(sampleNo_0,dataPointNo_0) schedule(static)
2234     for (sampleNo_0 = 0; sampleNo_0 < numSamples_0; sampleNo_0++) {
2235 matt 1332 for (dataPointNo_0 = 0; dataPointNo_0 < numDataPointsPerSample_0; dataPointNo_0++) {
2236     int offset_0 = tmp_0->getPointOffset(sampleNo_0,dataPointNo_0);
2237     int offset_2 = tmp_2->getPointOffset(sampleNo_0,dataPointNo_0);
2238     double *ptr_0 = &((arg_0_Z.getPointDataView().getData())[offset_0]);
2239     double *ptr_1 = &((arg_1_Z.getPointDataView().getData())[offset_1]);
2240     double *ptr_2 = &((res.getPointDataView().getData())[offset_2]);
2241     tensor_binary_operation(size1, ptr_0[0], ptr_1, ptr_2, operation);
2242     }
2243 matt 1327 }
2244    
2245    
2246     }
2247     else if (arg_0_Z.isExpanded() && arg_1_Z.isTagged()) {
2248    
2249     // After finding a common function space above the two inputs have the same numSamples and num DPPS
2250     res = Data(0.0, shape1, arg_1_Z.getFunctionSpace(),true); // DataExpanded output
2251     DataExpanded* tmp_0=dynamic_cast<DataExpanded*>(arg_0_Z.borrowData());
2252     DataTagged* tmp_1=dynamic_cast<DataTagged*>(arg_1_Z.borrowData());
2253     DataExpanded* tmp_2=dynamic_cast<DataExpanded*>(res.borrowData());
2254    
2255     int sampleNo_0,dataPointNo_0;
2256     int numSamples_0 = arg_0_Z.getNumSamples();
2257     int numDataPointsPerSample_0 = arg_0_Z.getNumDataPointsPerSample();
2258     #pragma omp parallel for private(sampleNo_0,dataPointNo_0) schedule(static)
2259     for (sampleNo_0 = 0; sampleNo_0 < numSamples_0; sampleNo_0++) {
2260 matt 1332 int offset_1 = tmp_1->getPointOffset(sampleNo_0,0);
2261     double *ptr_1 = &((arg_1_Z.getPointDataView().getData())[offset_1]);
2262     for (dataPointNo_0 = 0; dataPointNo_0 < numDataPointsPerSample_0; dataPointNo_0++) {
2263     int offset_0 = tmp_0->getPointOffset(sampleNo_0,dataPointNo_0);
2264     int offset_2 = tmp_2->getPointOffset(sampleNo_0,dataPointNo_0);
2265     double *ptr_0 = &((arg_0_Z.getPointDataView().getData())[offset_0]);
2266     double *ptr_2 = &((res.getPointDataView().getData())[offset_2]);
2267     tensor_binary_operation(size1, ptr_0[0], ptr_1, ptr_2, operation);
2268     }
2269 matt 1327 }
2270    
2271     }
2272     else if (arg_0_Z.isExpanded() && arg_1_Z.isExpanded()) {
2273    
2274     // After finding a common function space above the two inputs have the same numSamples and num DPPS
2275     res = Data(0.0, shape1, arg_1_Z.getFunctionSpace(),true); // DataExpanded output
2276     DataExpanded* tmp_0=dynamic_cast<DataExpanded*>(arg_0_Z.borrowData());
2277     DataExpanded* tmp_1=dynamic_cast<DataExpanded*>(arg_1_Z.borrowData());
2278     DataExpanded* tmp_2=dynamic_cast<DataExpanded*>(res.borrowData());
2279    
2280     int sampleNo_0,dataPointNo_0;
2281     int numSamples_0 = arg_0_Z.getNumSamples();
2282     int numDataPointsPerSample_0 = arg_0_Z.getNumDataPointsPerSample();
2283     #pragma omp parallel for private(sampleNo_0,dataPointNo_0) schedule(static)
2284     for (sampleNo_0 = 0; sampleNo_0 < numSamples_0; sampleNo_0++) {
2285 matt 1332 for (dataPointNo_0 = 0; dataPointNo_0 < numDataPointsPerSample_0; dataPointNo_0++) {
2286     int offset_0 = tmp_0->getPointOffset(sampleNo_0,dataPointNo_0);
2287     int offset_1 = tmp_1->getPointOffset(sampleNo_0,dataPointNo_0);
2288     int offset_2 = tmp_2->getPointOffset(sampleNo_0,dataPointNo_0);
2289     double *ptr_0 = &((arg_0_Z.getPointDataView().getData())[offset_0]);
2290     double *ptr_1 = &((arg_1_Z.getPointDataView().getData())[offset_1]);
2291     double *ptr_2 = &((res.getPointDataView().getData())[offset_2]);
2292     tensor_binary_operation(size1, ptr_0[0], ptr_1, ptr_2, operation);
2293     }
2294 matt 1327 }
2295    
2296     }
2297     else {
2298     throw DataException("Error - C_TensorBinaryOperation: unknown combination of inputs");
2299     }
2300    
2301     } else if (0 == rank1) {
2302    
2303     if (arg_0_Z.isConstant() && arg_1_Z.isConstant()) {
2304 matt 1332 res = Data(0.0, shape0, arg_1_Z.getFunctionSpace()); // DataConstant output
2305 matt 1327 double *ptr_0 = &((arg_0_Z.getPointDataView().getData())[0]);
2306     double *ptr_1 = &((arg_1_Z.getPointDataView().getData())[0]);
2307     double *ptr_2 = &((res.getPointDataView().getData())[0]);
2308     tensor_binary_operation(size0, ptr_0, ptr_1[0], ptr_2, operation);
2309     }
2310     else if (arg_0_Z.isConstant() && arg_1_Z.isTagged()) {
2311    
2312     // Prepare the DataConstant input
2313     DataConstant* tmp_0=dynamic_cast<DataConstant*>(arg_0_Z.borrowData());
2314    
2315     // Borrow DataTagged input from Data object
2316     DataTagged* tmp_1=dynamic_cast<DataTagged*>(arg_1_Z.borrowData());
2317    
2318     // Prepare a DataTagged output 2
2319 matt 1332 res = Data(0.0, shape0, arg_1_Z.getFunctionSpace()); // DataTagged output
2320 matt 1327 res.tag();
2321     DataTagged* tmp_2=dynamic_cast<DataTagged*>(res.borrowData());
2322    
2323     // Prepare offset into DataConstant
2324     int offset_0 = tmp_0->getPointOffset(0,0);
2325     double *ptr_0 = &((arg_0_Z.getPointDataView().getData())[offset_0]);
2326     // Get the views
2327     DataArrayView view_1 = tmp_1->getDefaultValue();
2328     DataArrayView view_2 = tmp_2->getDefaultValue();
2329     // Get the pointers to the actual data
2330     double *ptr_1 = &((view_1.getData())[0]);
2331     double *ptr_2 = &((view_2.getData())[0]);
2332     // Compute a result for the default
2333     tensor_binary_operation(size0, ptr_0, ptr_1[0], ptr_2, operation);
2334     // Compute a result for each tag
2335     const DataTagged::DataMapType& lookup_1=tmp_1->getTagLookup();
2336     DataTagged::DataMapType::const_iterator i; // i->first is a tag, i->second is an offset into memory
2337     for (i=lookup_1.begin();i!=lookup_1.end();i++) {
2338 matt 1332 tmp_2->addTaggedValue(i->first,tmp_2->getDefaultValue());
2339     DataArrayView view_1 = tmp_1->getDataPointByTag(i->first);
2340     DataArrayView view_2 = tmp_2->getDataPointByTag(i->first);
2341     double *ptr_1 = &view_1.getData(0);
2342     double *ptr_2 = &view_2.getData(0);
2343     tensor_binary_operation(size0, ptr_0, ptr_1[0], ptr_2, operation);
2344 matt 1327 }
2345    
2346     }
2347     else if (arg_0_Z.isConstant() && arg_1_Z.isExpanded()) {
2348    
2349     res = Data(0.0, shape0, arg_1_Z.getFunctionSpace(),true); // DataExpanded output
2350     DataConstant* tmp_0=dynamic_cast<DataConstant*>(arg_0_Z.borrowData());
2351     DataExpanded* tmp_1=dynamic_cast<DataExpanded*>(arg_1_Z.borrowData());
2352     DataExpanded* tmp_2=dynamic_cast<DataExpanded*>(res.borrowData());
2353    
2354     int sampleNo_1,dataPointNo_1;
2355     int numSamples_1 = arg_1_Z.getNumSamples();
2356     int numDataPointsPerSample_1 = arg_1_Z.getNumDataPointsPerSample();
2357     int offset_0 = tmp_0->getPointOffset(0,0);
2358     #pragma omp parallel for private(sampleNo_1,dataPointNo_1) schedule(static)
2359     for (sampleNo_1 = 0; sampleNo_1 < numSamples_1; sampleNo_1++) {
2360 matt 1332 for (dataPointNo_1 = 0; dataPointNo_1 < numDataPointsPerSample_1; dataPointNo_1++) {
2361     int offset_1 = tmp_1->getPointOffset(sampleNo_1,dataPointNo_1);
2362     int offset_2 = tmp_2->getPointOffset(sampleNo_1,dataPointNo_1);
2363     double *ptr_0 = &((arg_0_Z.getPointDataView().getData())[offset_0]);
2364     double *ptr_1 = &((arg_1_Z.getPointDataView().getData())[offset_1]);
2365     double *ptr_2 = &((res.getPointDataView().getData())[offset_2]);
2366     tensor_binary_operation(size0, ptr_0, ptr_1[0], ptr_2, operation);
2367     }
2368 matt 1327 }
2369    
2370     }
2371     else if (arg_0_Z.isTagged() && arg_1_Z.isConstant()) {
2372    
2373     // Borrow DataTagged input from Data object
2374     DataTagged* tmp_0=dynamic_cast<DataTagged*>(arg_0_Z.borrowData());
2375    
2376     // Prepare the DataConstant input
2377     DataConstant* tmp_1=dynamic_cast<DataConstant*>(arg_1_Z.borrowData());
2378    
2379     // Prepare a DataTagged output 2
2380 matt 1332 res = Data(0.0, shape0, arg_0_Z.getFunctionSpace()); // DataTagged output
2381 matt 1327 res.tag();
2382     DataTagged* tmp_2=dynamic_cast<DataTagged*>(res.borrowData());
2383    
2384     // Prepare offset into DataConstant
2385     int offset_1 = tmp_1->getPointOffset(0,0);
2386     double *ptr_1 = &((arg_1_Z.getPointDataView().getData())[offset_1]);
2387     // Get the views
2388     DataArrayView view_0 = tmp_0->getDefaultValue();
2389     DataArrayView view_2 = tmp_2->getDefaultValue();
2390     // Get the pointers to the actual data
2391     double *ptr_0 = &((view_0.getData())[0]);
2392     double *ptr_2 = &((view_2.getData())[0]);
2393     // Compute a result for the default
2394     tensor_binary_operation(size0, ptr_0, ptr_1[0], ptr_2, operation);
2395     // Compute a result for each tag
2396     const DataTagged::DataMapType& lookup_0=tmp_0->getTagLookup();
2397     DataTagged::DataMapType::const_iterator i; // i->first is a tag, i->second is an offset into memory
2398     for (i=lookup_0.begin();i!=lookup_0.end();i++) {
2399 matt 1332 tmp_2->addTaggedValue(i->first,tmp_2->getDefaultValue());
2400     DataArrayView view_0 = tmp_0->getDataPointByTag(i->first);
2401     DataArrayView view_2 = tmp_2->getDataPointByTag(i->first);
2402     double *ptr_0 = &view_0.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.isTagged()) {
2409    
2410     // Borrow DataTagged input from Data object
2411     DataTagged* tmp_0=dynamic_cast<DataTagged*>(arg_0_Z.borrowData());
2412    
2413     // Borrow DataTagged input from Data object
2414     DataTagged* tmp_1=dynamic_cast<DataTagged*>(arg_1_Z.borrowData());
2415    
2416     // Prepare a DataTagged output 2
2417     res = Data(0.0, shape0, arg_1_Z.getFunctionSpace());
2418 matt 1332 res.tag(); // DataTagged output
2419 matt 1327 DataTagged* tmp_2=dynamic_cast<DataTagged*>(res.borrowData());
2420    
2421     // Get the views
2422     DataArrayView view_0 = tmp_0->getDefaultValue();
2423     DataArrayView view_1 = tmp_1->getDefaultValue();
2424     DataArrayView view_2 = tmp_2->getDefaultValue();
2425     // Get the pointers to the actual data
2426     double *ptr_0 = &((view_0.getData())[0]);
2427     double *ptr_1 = &((view_1.getData())[0]);
2428     double *ptr_2 = &((view_2.getData())[0]);
2429     // Compute a result for the default
2430     tensor_binary_operation(size0, ptr_0, ptr_1[0], ptr_2, operation);
2431     // Merge the tags
2432     DataTagged::DataMapType::const_iterator i; // i->first is a tag, i->second is an offset into memory
2433     const DataTagged::DataMapType& lookup_0=tmp_0->getTagLookup();
2434     const DataTagged::DataMapType& lookup_1=tmp_1->getTagLookup();
2435     for (i=lookup_0.begin();i!=lookup_0.end();i++) {
2436 matt 1332 tmp_2->addTaggedValue(i->first,tmp_2->getDefaultValue()); // use tmp_2 to get correct shape
2437 matt 1327 }
2438     for (i=lookup_1.begin();i!=lookup_1.end();i++) {
2439 matt 1332 tmp_2->addTaggedValue(i->first,tmp_2->getDefaultValue());
2440 matt 1327 }
2441     // Compute a result for each tag
2442     const DataTagged::DataMapType& lookup_2=tmp_2->getTagLookup();
2443     for (i=lookup_2.begin();i!=lookup_2.end();i++) {
2444 matt 1332 DataArrayView view_0 = tmp_0->getDataPointByTag(i->first);
2445     DataArrayView view_1 = tmp_1->getDataPointByTag(i->first);
2446     DataArrayView view_2 = tmp_2->getDataPointByTag(i->first);
2447     double *ptr_0 = &view_0.getData(0);
2448     double *ptr_1 = &view_1.getData(0);
2449     double *ptr_2 = &view_2.getData(0);
2450     tensor_binary_operation(size0, ptr_0, ptr_1[0], ptr_2, operation);
2451 matt 1327 }
2452    
2453     }
2454     else if (arg_0_Z.isTagged() && arg_1_Z.isExpanded()) {
2455    
2456     // After finding a common function space above the two inputs have the same numSamples and num DPPS
2457     res = Data(0.0, shape0, arg_1_Z.getFunctionSpace(),true); // DataExpanded output
2458     DataTagged* tmp_0=dynamic_cast<DataTagged*>(arg_0_Z.borrowData());
2459     DataExpanded* tmp_1=dynamic_cast<DataExpanded*>(arg_1_Z.borrowData());
2460     DataExpanded* tmp_2=dynamic_cast<DataExpanded*>(res.borrowData());
2461    
2462     int sampleNo_0,dataPointNo_0;
2463     int numSamples_0 = arg_0_Z.getNumSamples();
2464     int numDataPointsPerSample_0 = arg_0_Z.getNumDataPointsPerSample();
2465     #pragma omp parallel for private(sampleNo_0,dataPointNo_0) schedule(static)
2466     for (sampleNo_0 = 0; sampleNo_0 < numSamples_0; sampleNo_0++) {
2467 matt 1332 int offset_0 = tmp_0->getPointOffset(sampleNo_0,0); // They're all the same, so just use #0
2468     double *ptr_0 = &((arg_0_Z.getPointDataView().getData())[offset_0]);
2469     for (dataPointNo_0 = 0; dataPointNo_0 < numDataPointsPerSample_0; dataPointNo_0++) {
2470     int offset_1 = tmp_1->getPointOffset(sampleNo_0,dataPointNo_0);
2471     int offset_2 = tmp_2->getPointOffset(sampleNo_0,dataPointNo_0);
2472     double *ptr_1 = &((arg_1_Z.getPointDataView().getData())[offset_1]);
2473     double *ptr_2 = &((res.getPointDataView().getData())[offset_2]);
2474     tensor_binary_operation(size0, ptr_0, ptr_1[0], ptr_2, operation);
2475     }
2476 matt 1327 }
2477    
2478     }
2479     else if (arg_0_Z.isExpanded() && arg_1_Z.isConstant()) {
2480    
2481     res = Data(0.0, shape0, arg_1_Z.getFunctionSpace(),true); // DataExpanded output
2482     DataExpanded* tmp_0=dynamic_cast<DataExpanded*>(arg_0_Z.borrowData());
2483     DataConstant* tmp_1=dynamic_cast<DataConstant*>(arg_1_Z.borrowData());
2484     DataExpanded* tmp_2=dynamic_cast<DataExpanded*>(res.borrowData());
2485    
2486     int sampleNo_0,dataPointNo_0;
2487     int numSamples_0 = arg_0_Z.getNumSamples();
2488     int numDataPointsPerSample_0 = arg_0_Z.getNumDataPointsPerSample();
2489     int offset_1 = tmp_1->getPointOffset(0,0);
2490     #pragma omp parallel for private(sampleNo_0,dataPointNo_0) schedule(static)
2491     for (sampleNo_0 = 0; sampleNo_0 < numSamples_0; sampleNo_0++) {
2492 matt 1332 for (dataPointNo_0 = 0; dataPointNo_0 < numDataPointsPerSample_0; dataPointNo_0++) {
2493     int offset_0 = tmp_0->getPointOffset(sampleNo_0,dataPointNo_0);
2494     int offset_2 = tmp_2->getPointOffset(sampleNo_0,dataPointNo_0);
2495     double *ptr_0 = &((arg_0_Z.getPointDataView().getData())[offset_0]);
2496     double *ptr_1 = &((arg_1_Z.getPointDataView().getData())[offset_1]);
2497     double *ptr_2 = &((res.getPointDataView().getData())[offset_2]);
2498     tensor_binary_operation(size0, ptr_0, ptr_1[0], ptr_2, operation);
2499     }
2500 matt 1327 }
2501    
2502    
2503     }
2504     else if (arg_0_Z.isExpanded() && arg_1_Z.isTagged()) {
2505    
2506     // After finding a common function space above the two inputs have the same numSamples and num DPPS
2507     res = Data(0.0, shape0, arg_1_Z.getFunctionSpace(),true); // DataExpanded output
2508     DataExpanded* tmp_0=dynamic_cast<DataExpanded*>(arg_0_Z.borrowData());
2509     DataTagged* tmp_1=dynamic_cast<DataTagged*>(arg_1_Z.borrowData());
2510     DataExpanded* tmp_2=dynamic_cast<DataExpanded*>(res.borrowData());
2511    
2512     int sampleNo_0,dataPointNo_0;
2513     int numSamples_0 = arg_0_Z.getNumSamples();
2514     int numDataPointsPerSample_0 = arg_0_Z.getNumDataPointsPerSample();
2515     #pragma omp parallel for private(sampleNo_0,dataPointNo_0) schedule(static)
2516     for (sampleNo_0 = 0; sampleNo_0 < numSamples_0; sampleNo_0++) {
2517 matt 1332 int offset_1 = tmp_1->getPointOffset(sampleNo_0,0);
2518     double *ptr_1 = &((arg_1_Z.getPointDataView().getData())[offset_1]);
2519     for (dataPointNo_0 = 0; dataPointNo_0 < numDataPointsPerSample_0; dataPointNo_0++) {
2520     int offset_0 = tmp_0->getPointOffset(sampleNo_0,dataPointNo_0);
2521     int offset_2 = tmp_2->getPointOffset(sampleNo_0,dataPointNo_0);
2522     double *ptr_0 = &((arg_0_Z.getPointDataView().getData())[offset_0]);
2523     double *ptr_2 = &((res.getPointDataView().getData())[offset_2]);
2524     tensor_binary_operation(size0, ptr_0, ptr_1[0], ptr_2, operation);
2525     }
2526 matt 1327 }
2527    
2528     }
2529     else if (arg_0_Z.isExpanded() && arg_1_Z.isExpanded()) {
2530    
2531     // After finding a common function space above the two inputs have the same numSamples and num DPPS
2532     res = Data(0.0, shape0, arg_1_Z.getFunctionSpace(),true); // DataExpanded output
2533     DataExpanded* tmp_0=dynamic_cast<DataExpanded*>(arg_0_Z.borrowData());
2534     DataExpanded* tmp_1=dynamic_cast<DataExpanded*>(arg_1_Z.borrowData());
2535     DataExpanded* tmp_2=dynamic_cast<DataExpanded*>(res.borrowData());
2536    
2537     int sampleNo_0,dataPointNo_0;
2538     int numSamples_0 = arg_0_Z.getNumSamples();
2539     int numDataPointsPerSample_0 = arg_0_Z.getNumDataPointsPerSample();
2540     #pragma omp parallel for private(sampleNo_0,dataPointNo_0) schedule(static)
2541     for (sampleNo_0 = 0; sampleNo_0 < numSamples_0; sampleNo_0++) {
2542 matt 1332 for (dataPointNo_0 = 0; dataPointNo_0 < numDataPointsPerSample_0; dataPointNo_0++) {
2543     int offset_0 = tmp_0->getPointOffset(sampleNo_0,dataPointNo_0);
2544     int offset_1 = tmp_1->getPointOffset(sampleNo_0,dataPointNo_0);
2545     int offset_2 = tmp_2->getPointOffset(sampleNo_0,dataPointNo_0);
2546     double *ptr_0 = &((arg_0_Z.getPointDataView().getData())[offset_0]);
2547     double *ptr_1 = &((arg_1_Z.getPointDataView().getData())[offset_1]);
2548     double *ptr_2 = &((res.getPointDataView().getData())[offset_2]);
2549     tensor_binary_operation(size0, ptr_0, ptr_1[0], ptr_2, operation);
2550     }
2551 matt 1327 }
2552    
2553     }
2554     else {
2555     throw DataException("Error - C_TensorBinaryOperation: unknown combination of inputs");
2556     }
2557    
2558     } else {
2559     throw DataException("Error - C_TensorBinaryOperation: arguments have incompatible shapes");
2560     }
2561    
2562     return res;
2563 jgs 94 }
2564 matt 1327
2565     }
2566 jgs 94 #endif

Properties

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

  ViewVC Help
Powered by ViewVC 1.1.26