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

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

Parent Directory Parent Directory | Revision Log Revision Log


Revision 3992 - (show annotations)
Wed Sep 26 02:50:48 2012 UTC (7 years ago) by caltinay
File MIME type: text/plain
File size: 96215 byte(s)
Data.cpp: Removed some namespace ambiguities which fixes doxygen warnings, some
cleanup and reindent.
*: doxygen fixes.

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

Properties

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

  ViewVC Help
Powered by ViewVC 1.1.26