/[escript]/branches/schroedinger/escript/src/Data.h
ViewVC logotype

Annotation of /branches/schroedinger/escript/src/Data.h

Parent Directory Parent Directory | Revision Log Revision Log


Revision 782 - (hide annotations)
Tue Jul 18 00:47:47 2006 UTC (13 years, 4 months ago) by bcumming
Original Path: trunk/escript/src/Data.h
File MIME type: text/plain
File size: 42574 byte(s)
Large number of changes to Finley for meshing in MPI.

- optimisation and neatening up of rectcanglular mesh generation code
- first and second order 1D, 2D and 3D rectangular meshes are now
  available in finley and escript using MPI.
- reduced meshes now generated in MPI, and interpolation to and from 
  reduced data types now supported.  

1 jgs 94 // $Id$
2 jgs 121 /*
3 elspeth 615 ************************************************************
4     * Copyright 2006 by ACcESS MNRF *
5     * *
6     * http://www.access.edu.au *
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     ************************************************************
12 jgs 121 */
13 jgs 94
14 jgs 121 /** \file Data.h */
15 jgs 94
16     #ifndef DATA_H
17     #define DATA_H
18 woo409 757 #include "system_dep.h"
19 jgs 94
20 jgs 474 #include "DataAbstract.h"
21     #include "DataAlgorithm.h"
22     #include "FunctionSpace.h"
23     #include "BinaryOp.h"
24     #include "UnaryOp.h"
25     #include "DataException.h"
26 jgs 94
27     extern "C" {
28 jgs 474 #include "DataC.h"
29 woo409 757 #include "paso/Paso.h"
30 jgs 94 }
31    
32     #include <string>
33     #include <algorithm>
34    
35     #include <boost/shared_ptr.hpp>
36     #include <boost/python/object.hpp>
37     #include <boost/python/tuple.hpp>
38     #include <boost/python/numeric.hpp>
39    
40 jgs 121 namespace escript {
41    
42     //
43 jgs 149 // Forward declaration for various implementations of Data.
44 jgs 121 class DataConstant;
45     class DataTagged;
46     class DataExpanded;
47    
48 jgs 94 /**
49     \brief
50 jgs 121 Data creates the appropriate Data object for the given construction
51     arguments.
52 jgs 94
53     Description:
54     Data is essentially a factory class which creates the appropriate Data
55     object for the given construction arguments. It retains control over
56     the object created for the lifetime of the object.
57     The type of Data object referred to may change during the lifetime of
58     the Data object.
59     */
60     class Data {
61    
62     public:
63    
64 jgs 110 // These typedefs allow function names to be cast to pointers
65     // to functions of the appropriate type when calling unaryOp etc.
66 jgs 94 typedef double (*UnaryDFunPtr)(double);
67     typedef double (*BinaryDFunPtr)(double,double);
68    
69     /**
70 jgs 102 Constructors.
71     */
72    
73     /**
74 jgs 94 \brief
75     Default constructor.
76     Creates a DataEmpty object.
77     */
78 woo409 757 ESCRIPT_DLL_API
79 jgs 94 Data();
80    
81     /**
82     \brief
83     Copy constructor.
84     WARNING: Only performs a shallow copy.
85     */
86 woo409 757 ESCRIPT_DLL_API
87 jgs 94 Data(const Data& inData);
88    
89     /**
90     \brief
91     Constructor from another Data object. If "what" is different from the
92 jgs 102 function space of inData the inData are tried to be interpolated to what,
93 jgs 94 otherwise a shallow copy of inData is returned.
94     */
95 woo409 757 ESCRIPT_DLL_API
96 jgs 94 Data(const Data& inData,
97     const FunctionSpace& what);
98    
99     /**
100     \brief
101     Constructor which copies data from a DataArrayView.
102    
103     \param value - Input - Data value for a single point.
104     \param what - Input - A description of what this data represents.
105     \param expanded - Input - Flag, if true fill the entire container with
106     the value. Otherwise a more efficient storage
107     mechanism will be used.
108     */
109 woo409 757 ESCRIPT_DLL_API
110 jgs 94 Data(const DataArrayView& value,
111     const FunctionSpace& what=FunctionSpace(),
112     bool expanded=false);
113    
114     /**
115     \brief
116     Constructor which creates a Data from a DataArrayView shape.
117    
118     \param value - Input - Single value applied to all Data.
119     \param dataPointShape - Input - The shape of each data point.
120     \param what - Input - A description of what this data represents.
121     \param expanded - Input - Flag, if true fill the entire container with
122     the given value. Otherwise a more efficient storage
123     mechanism will be used.
124     */
125 woo409 757 ESCRIPT_DLL_API
126 jgs 94 Data(double value,
127     const DataArrayView::ShapeType& dataPointShape=DataArrayView::ShapeType(),
128     const FunctionSpace& what=FunctionSpace(),
129     bool expanded=false);
130    
131     /**
132     \brief
133     Constructor which performs a deep copy of a region from another Data object.
134    
135     \param inData - Input - Input Data object.
136     \param region - Input - Region to copy.
137     */
138 woo409 757 ESCRIPT_DLL_API
139 jgs 94 Data(const Data& inData,
140     const DataArrayView::RegionType& region);
141    
142     /**
143     \brief
144     Constructor which will create Tagged data if expanded is false.
145     No attempt is made to ensure the tag keys match the tag keys
146     within the function space.
147    
148     \param tagKeys - Input - List of tag values.
149     \param values - Input - List of values, one for each tag.
150     \param defaultValue - Input - A default value, used if tag doesn't exist.
151     \param what - Input - A description of what this data represents.
152     \param expanded - Input - Flag, if true fill the entire container with
153     the appropriate values.
154 jgs 544 ==>*
155 jgs 94 */
156 woo409 757 ESCRIPT_DLL_API
157 jgs 94 Data(const DataTagged::TagListType& tagKeys,
158     const DataTagged::ValueListType& values,
159     const DataArrayView& defaultValue,
160     const FunctionSpace& what=FunctionSpace(),
161     bool expanded=false);
162    
163     /**
164     \brief
165     Constructor which copies data from a python numarray.
166    
167     \param value - Input - Data value for a single point.
168     \param what - Input - A description of what this data represents.
169     \param expanded - Input - Flag, if true fill the entire container with
170     the value. Otherwise a more efficient storage
171     mechanism will be used.
172     */
173 woo409 757 ESCRIPT_DLL_API
174 jgs 94 Data(const boost::python::numeric::array& value,
175     const FunctionSpace& what=FunctionSpace(),
176     bool expanded=false);
177    
178     /**
179     \brief
180     Constructor which copies data from any object that can be converted into
181     a python numarray.
182    
183     \param value - Input - Input data.
184     \param what - Input - A description of what this data represents.
185     \param expanded - Input - Flag, if true fill the entire container with
186     the value. Otherwise a more efficient storage
187     mechanism will be used.
188     */
189 woo409 757 ESCRIPT_DLL_API
190 jgs 94 Data(const boost::python::object& value,
191     const FunctionSpace& what=FunctionSpace(),
192     bool expanded=false);
193    
194     /**
195     \brief
196     Constructor which creates a DataConstant.
197     Copies data from any object that can be converted
198     into a numarray. All other parameters are copied from other.
199    
200     \param value - Input - Input data.
201     \param other - Input - contains all other parameters.
202     */
203 woo409 757 ESCRIPT_DLL_API
204 jgs 94 Data(const boost::python::object& value,
205     const Data& other);
206    
207     /**
208     \brief
209     Constructor which creates a DataConstant of "shape" with constant value.
210     */
211 woo409 757 ESCRIPT_DLL_API
212 jgs 94 Data(double value,
213     const boost::python::tuple& shape=boost::python::make_tuple(),
214     const FunctionSpace& what=FunctionSpace(),
215     bool expanded=false);
216 jgs 151 /**
217     \brief
218     Destructor
219     */
220 woo409 757 ESCRIPT_DLL_API
221 jgs 151 ~Data();
222 jgs 94
223     /**
224     \brief
225 jgs 102 Perform a deep copy.
226 jgs 94 */
227 woo409 757 ESCRIPT_DLL_API
228 jgs 94 void
229 jgs 102 copy(const Data& other);
230 jgs 94
231     /**
232 jgs 102 Member access methods.
233 jgs 94 */
234    
235     /**
236     \brief
237 jgs 121 Return the values of all data-points as a single python numarray object.
238 jgs 117 */
239 woo409 757 ESCRIPT_DLL_API
240 jgs 121 const boost::python::numeric::array
241 jgs 117 convertToNumArray();
242    
243     /**
244     \brief
245 jgs 121 Return the values of all data-points for the given sample as a single python numarray object.
246     */
247 woo409 757 ESCRIPT_DLL_API
248 jgs 121 const boost::python::numeric::array
249     convertToNumArrayFromSampleNo(int sampleNo);
250    
251     /**
252     \brief
253     Return the value of the specified data-point as a single python numarray object.
254     */
255 woo409 757 ESCRIPT_DLL_API
256 jgs 121 const boost::python::numeric::array
257     convertToNumArrayFromDPNo(int sampleNo,
258     int dataPointNo);
259    
260     /**
261     \brief
262 jgs 149 Fills the expanded Data object from values of a python numarray object.
263     */
264 woo409 757 ESCRIPT_DLL_API
265 jgs 149 void
266     fillFromNumArray(const boost::python::numeric::array);
267    
268     /**
269     \brief
270     Return the tag number associated with the given data-point.
271    
272     The data-point number here corresponds to the data-point number in the
273     numarray returned by convertToNumArray.
274     */
275 woo409 757 ESCRIPT_DLL_API
276 jgs 149 int
277     getTagNumber(int dpno);
278    
279     /**
280     \brief
281 jgs 94 Return the C wrapper for the Data object.
282     */
283 woo409 757 ESCRIPT_DLL_API
284 jgs 102 escriptDataC
285     getDataC();
286 jgs 94
287     /**
288     \brief
289     Return the C wrapper for the Data object - const version.
290     */
291 woo409 757 ESCRIPT_DLL_API
292 jgs 102 escriptDataC
293     getDataC() const;
294 jgs 94
295     /**
296     \brief
297     Write the data as a string.
298     */
299 woo409 757 ESCRIPT_DLL_API
300 jgs 102 inline
301     std::string
302     toString() const
303     {
304     return m_data->toString();
305     }
306 jgs 94
307     /**
308     \brief
309     Return the DataArrayView of the point data. This essentially contains
310     the shape information for each data point although it also may be used
311     to manipulate the point data.
312     */
313 woo409 757 ESCRIPT_DLL_API
314 jgs 94 inline
315     const DataArrayView&
316     getPointDataView() const
317     {
318     return m_data->getPointDataView();
319     }
320    
321     /**
322     \brief
323 jgs 102 Whatever the current Data type make this into a DataExpanded.
324 jgs 94 */
325 woo409 757 ESCRIPT_DLL_API
326 jgs 94 void
327     expand();
328    
329     /**
330     \brief
331 jgs 102 If possible convert this Data to DataTagged. This will only allow
332 jgs 94 Constant data to be converted to tagged. An attempt to convert
333     Expanded data to tagged will throw an exception.
334 jgs 544 ==>*
335 jgs 94 */
336 woo409 757 ESCRIPT_DLL_API
337 jgs 94 void
338     tag();
339    
340     /**
341     \brief
342     Return true if this Data is expanded.
343     */
344 woo409 757 ESCRIPT_DLL_API
345 jgs 94 bool
346     isExpanded() const;
347    
348     /**
349     \brief
350     Return true if this Data is tagged.
351     */
352 woo409 757 ESCRIPT_DLL_API
353 jgs 94 bool
354     isTagged() const;
355    
356     /**
357     \brief
358 jgs 102 Return true if this Data is constant.
359 jgs 94 */
360 woo409 757 ESCRIPT_DLL_API
361 jgs 94 bool
362 jgs 102 isConstant() const;
363 jgs 94
364     /**
365     \brief
366 jgs 102 Return true if this Data is empty.
367 jgs 94 */
368 woo409 757 ESCRIPT_DLL_API
369 jgs 94 bool
370 jgs 102 isEmpty() const;
371 jgs 94
372     /**
373     \brief
374     Return the function space.
375     */
376 woo409 757 ESCRIPT_DLL_API
377 jgs 94 inline
378     const FunctionSpace&
379     getFunctionSpace() const
380     {
381     return m_data->getFunctionSpace();
382     }
383    
384     /**
385     \brief
386     Return a copy of the function space.
387     */
388 woo409 757 ESCRIPT_DLL_API
389 jgs 94 const FunctionSpace
390     getCopyOfFunctionSpace() const;
391    
392     /**
393     \brief
394     Return the domain.
395     */
396 woo409 757 ESCRIPT_DLL_API
397 jgs 94 inline
398     const AbstractDomain&
399     getDomain() const
400     {
401     return getFunctionSpace().getDomain();
402     }
403    
404     /**
405     \brief
406     Return a copy of the domain.
407     */
408 woo409 757 ESCRIPT_DLL_API
409 jgs 94 const AbstractDomain
410     getCopyOfDomain() const;
411    
412     /**
413     \brief
414     Return the rank of the point data.
415     */
416 woo409 757 ESCRIPT_DLL_API
417 jgs 94 inline
418     int
419     getDataPointRank() const
420     {
421     return m_data->getPointDataView().getRank();
422     }
423    
424     /**
425     \brief
426 jgs 102 Return the number of samples.
427 jgs 94 */
428 woo409 757 ESCRIPT_DLL_API
429 jgs 94 inline
430     int
431 jgs 102 getNumSamples() const
432 jgs 94 {
433 jgs 102 return m_data->getNumSamples();
434 jgs 94 }
435    
436     /**
437     \brief
438 jgs 102 Return the number of data points per sample.
439 jgs 94 */
440 woo409 757 ESCRIPT_DLL_API
441 jgs 102 inline
442 jgs 94 int
443 jgs 102 getNumDataPointsPerSample() const
444     {
445     return m_data->getNumDPPSample();
446     }
447 jgs 94
448     /**
449     \brief
450     Return the sample data for the given sample no. This is not the
451     preferred interface but is provided for use by C code.
452     \param sampleNo - Input - the given sample no.
453     */
454 woo409 757 ESCRIPT_DLL_API
455 jgs 102 inline
456 jgs 94 DataAbstract::ValueType::value_type*
457 jgs 102 getSampleData(DataAbstract::ValueType::size_type sampleNo)
458     {
459     return m_data->getSampleData(sampleNo);
460     }
461 jgs 94
462     /**
463     \brief
464     Return the sample data for the given tag. If an attempt is made to
465     access data that isn't tagged an exception will be thrown.
466     \param tag - Input - the tag key.
467     */
468 woo409 757 ESCRIPT_DLL_API
469 jgs 102 inline
470 jgs 94 DataAbstract::ValueType::value_type*
471 jgs 102 getSampleDataByTag(int tag)
472     {
473     return m_data->getSampleDataByTag(tag);
474     }
475 jgs 94
476     /**
477     \brief
478 jgs 110 Assign the given value to the data-points referenced by the given
479     reference number.
480    
481     The value supplied is a python numarray object. The data from this numarray
482     is unpacked into a DataArray, and this is used to set the corresponding
483     data-points in the underlying Data object.
484    
485     If the underlying Data object cannot be accessed via reference numbers, an
486     exception will be thrown.
487    
488     \param ref - Input - reference number.
489     \param value - Input - value to assign to data-points associated with
490     the given reference number.
491     */
492 woo409 757 ESCRIPT_DLL_API
493 jgs 110 void
494     setRefValue(int ref,
495     const boost::python::numeric::array& value);
496    
497     /**
498     \brief
499     Return the values associated with the data-points referenced by the given
500     reference number.
501    
502     The value supplied is a python numarray object. The data from the corresponding
503     data-points in this Data object are packed into the given numarray object.
504    
505     If the underlying Data object cannot be accessed via reference numbers, an
506     exception will be thrown.
507    
508     \param ref - Input - reference number.
509     \param value - Output - object to receive values from data-points
510     associated with the given reference number.
511     */
512 woo409 757 ESCRIPT_DLL_API
513 jgs 110 void
514     getRefValue(int ref,
515     boost::python::numeric::array& value);
516    
517     /**
518     \brief
519 jgs 94 Return a view into the data for the data point specified.
520     NOTE: Construction of the DataArrayView is a relatively expensive
521     operation.
522     \param sampleNo - Input -
523     \param dataPointNo - Input -
524     */
525 woo409 757 ESCRIPT_DLL_API
526 jgs 94 inline
527     DataArrayView
528     getDataPoint(int sampleNo,
529     int dataPointNo)
530     {
531     return m_data->getDataPoint(sampleNo,dataPointNo);
532     }
533    
534     /**
535     \brief
536     Return a reference to the data point shape.
537     */
538 woo409 757 ESCRIPT_DLL_API
539 jgs 94 const DataArrayView::ShapeType&
540     getDataPointShape() const;
541    
542     /**
543     \brief
544 jgs 102 Return the data point shape as a tuple of integers.
545 jgs 94 */
546 woo409 757 ESCRIPT_DLL_API
547 jgs 121 const boost::python::tuple
548 jgs 94 getShapeTuple() const;
549    
550     /**
551     \brief
552     Return the size of the data point. It is the product of the
553     data point shape dimensions.
554     */
555 woo409 757 ESCRIPT_DLL_API
556 jgs 94 int
557     getDataPointSize() const;
558    
559     /**
560     \brief
561 jgs 102 Return the number of doubles stored for this Data.
562 jgs 94 */
563 woo409 757 ESCRIPT_DLL_API
564 jgs 94 DataArrayView::ValueType::size_type
565     getLength() const;
566    
567     /**
568     \brief
569 jgs 102 Assign the given value to the tag. Implicitly converts this
570     object to type DataTagged. Throws an exception if this object
571     cannot be converted to a DataTagged object.
572     \param tagKey - Input - Integer key.
573     \param value - Input - Value to associate with given key.
574 jgs 544 ==>*
575 jgs 94 */
576 woo409 757 ESCRIPT_DLL_API
577 jgs 102 void
578     setTaggedValue(int tagKey,
579     const boost::python::object& value);
580    
581     /**
582     \brief
583     Assign the given value to the tag. Implicitly converts this
584     object to type DataTagged. Throws an exception if this object
585     cannot be converted to a DataTagged object.
586     \param tagKey - Input - Integer key.
587     \param value - Input - Value to associate with given key.
588 jgs 544 ==>*
589 jgs 102 */
590 woo409 757 ESCRIPT_DLL_API
591 jgs 102 void
592 jgs 121 setTaggedValueFromCPP(int tagKey,
593     const DataArrayView& value);
594 jgs 102
595     /**
596     \brief
597     Copy other Data object into this Data object where mask is positive.
598     */
599 woo409 757 ESCRIPT_DLL_API
600 jgs 102 void
601     copyWithMask(const Data& other,
602     const Data& mask);
603    
604     /**
605     Data object operation methods and operators.
606     */
607    
608     /**
609     \brief
610     Interpolates this onto the given functionspace and returns
611     the result as a Data object.
612 jgs 123 *
613 jgs 102 */
614 woo409 757 ESCRIPT_DLL_API
615 jgs 94 Data
616     interpolate(const FunctionSpace& functionspace) const;
617    
618     /**
619     \brief
620     Calculates the gradient of the data at the data points of functionspace.
621     If functionspace is not present the function space of Function(getDomain()) is used.
622 jgs 123 *
623 jgs 94 */
624 woo409 757 ESCRIPT_DLL_API
625 jgs 94 Data
626     gradOn(const FunctionSpace& functionspace) const;
627    
628 woo409 757 ESCRIPT_DLL_API
629 jgs 94 Data
630     grad() const;
631    
632     /**
633     \brief
634     Calculate the integral over the function space domain.
635 jgs 123 *
636 jgs 94 */
637 woo409 757 ESCRIPT_DLL_API
638 jgs 94 boost::python::numeric::array
639     integrate() const;
640    
641     /**
642     \brief
643     Return a Data with a 1 for +ive values and a 0 for 0 or -ive values.
644 jgs 123 *
645 jgs 94 */
646 woo409 757 ESCRIPT_DLL_API
647 jgs 94 Data
648 gross 698 wherePositive() const;
649 jgs 94
650     /**
651     \brief
652 jgs 102 Return a Data with a 1 for -ive values and a 0 for +ive or 0 values.
653 jgs 123 *
654 jgs 94 */
655 woo409 757 ESCRIPT_DLL_API
656 jgs 94 Data
657 gross 698 whereNegative() const;
658 jgs 102
659     /**
660     \brief
661     Return a Data with a 1 for +ive or 0 values and a 0 for -ive values.
662 jgs 123 *
663 jgs 102 */
664 woo409 757 ESCRIPT_DLL_API
665 jgs 102 Data
666 gross 698 whereNonNegative() const;
667 jgs 94
668     /**
669     \brief
670 jgs 102 Return a Data with a 1 for -ive or 0 values and a 0 for +ive values.
671 jgs 123 *
672 jgs 94 */
673 woo409 757 ESCRIPT_DLL_API
674 jgs 94 Data
675 gross 698 whereNonPositive() const;
676 jgs 94
677     /**
678     \brief
679 jgs 102 Return a Data with a 1 for 0 values and a 0 for +ive or -ive values.
680 jgs 123 *
681 jgs 94 */
682 woo409 757 ESCRIPT_DLL_API
683 jgs 94 Data
684 jgs 571 whereZero(double tol=0.0) const;
685 jgs 94
686     /**
687     \brief
688 jgs 102 Return a Data with a 0 for 0 values and a 1 for +ive or -ive values.
689 jgs 123 *
690 jgs 94 */
691 woo409 757 ESCRIPT_DLL_API
692 jgs 94 Data
693 jgs 571 whereNonZero(double tol=0.0) const;
694 jgs 102
695     /**
696     \brief
697     Return the maximum absolute value of this Data object.
698 jgs 123 *
699 jgs 94 */
700 woo409 757 ESCRIPT_DLL_API
701 jgs 102 double
702     Lsup() const;
703 jgs 94
704     /**
705     \brief
706 jgs 117 Return the minimum absolute value of this Data object.
707 jgs 123 *
708 jgs 117 */
709 woo409 757 ESCRIPT_DLL_API
710 jgs 117 double
711     Linf() const;
712    
713     /**
714     \brief
715 jgs 102 Return the maximum value of this Data object.
716 jgs 123 *
717 jgs 94 */
718 woo409 757 ESCRIPT_DLL_API
719 jgs 102 double
720     sup() const;
721 jgs 94
722     /**
723     \brief
724 jgs 102 Return the minimum value of this Data object.
725 jgs 123 *
726 jgs 94 */
727 woo409 757 ESCRIPT_DLL_API
728 jgs 94 double
729 jgs 102 inf() const;
730 jgs 94
731     /**
732     \brief
733 jgs 102 Return the absolute value of each data point of this Data object.
734 jgs 123 *
735 jgs 94 */
736 woo409 757 ESCRIPT_DLL_API
737 jgs 102 Data
738     abs() const;
739 jgs 94
740     /**
741     \brief
742 jgs 102 Return the maximum value of each data point of this Data object.
743 jgs 123 *
744 jgs 94 */
745 woo409 757 ESCRIPT_DLL_API
746 jgs 102 Data
747     maxval() const;
748 jgs 94
749     /**
750     \brief
751 jgs 102 Return the minimum value of each data point of this Data object.
752 jgs 123 *
753 jgs 94 */
754 woo409 757 ESCRIPT_DLL_API
755 jgs 94 Data
756 jgs 102 minval() const;
757 jgs 94
758     /**
759     \brief
760 jgs 121 Return the (sample number, data-point number) of the data point with
761     the minimum value in this Data object.
762     */
763 woo409 757 ESCRIPT_DLL_API
764 jgs 121 const boost::python::tuple
765     mindp() const;
766    
767 woo409 757 ESCRIPT_DLL_API
768 jgs 148 void
769 bcumming 782 #ifndef PASO_MPI
770 jgs 148 calc_mindp(int& SampleNo,
771     int& DataPointNo) const;
772 bcumming 782 #else
773     calc_mindp(int& ProcNo,
774     int& SampleNo,
775     int& DataPointNo) const;
776     #endif
777 jgs 121 /**
778     \brief
779 jgs 102 Return the sign of each data point of this Data object.
780     -1 for negative values, zero for zero values, 1 for positive values.
781 jgs 123 *
782 jgs 94 */
783 woo409 757 ESCRIPT_DLL_API
784 jgs 102 Data
785     sign() const;
786 jgs 94
787     /**
788 jgs 123 \brief
789 ksteube 775 Return the symmetric part of a matrix which is half the matrix plus its transpose.
790     *
791     */
792     ESCRIPT_DLL_API
793     Data
794     symmetric() const;
795    
796     /**
797     \brief
798     Return the nonsymmetric part of a matrix which is half the matrix minus its transpose.
799     *
800     */
801     ESCRIPT_DLL_API
802     Data
803     nonsymmetric() const;
804    
805     /**
806     \brief
807     Return the trace of a matrix
808     *
809     */
810     ESCRIPT_DLL_API
811     Data
812     matrixtrace(int axis_offset) const;
813    
814     /**
815     \brief
816     Transpose each data point of this Data object around the given axis.
817     *
818     */
819     ESCRIPT_DLL_API
820     Data
821     transpose(int axis_offset) const;
822    
823     /**
824     \brief
825 gross 576 Return the eigenvalues of the symmetric part at each data point of this Data object in increasing values.
826     Currently this function is restricted to rank 2, square shape, and dimension 3.
827     *
828     */
829 woo409 757 ESCRIPT_DLL_API
830 gross 576 Data
831     eigenvalues() const;
832    
833     /**
834     \brief
835     Return the eigenvalues and corresponding eigenvcetors of the symmetric part at each data point of this Data object.
836     the eigenvalues are ordered in increasing size where eigenvalues with relative difference less than
837     tol are treated as equal. The eigenvectors are orthogonal, normalized and the sclaed such that the
838     first non-zero entry is positive.
839     Currently this function is restricted to rank 2, square shape, and dimension 3
840     *
841     */
842 woo409 757 ESCRIPT_DLL_API
843 gross 576 const boost::python::tuple
844     eigenvalues_and_eigenvectors(const double tol=1.e-12) const;
845    
846     /**
847     \brief
848 jgs 123 Calculate the trace of each data point of this Data object.
849     *
850 jgs 102 */
851 woo409 757 ESCRIPT_DLL_API
852 jgs 102 Data
853     trace() const;
854    
855     /**
856 jgs 123 \brief
857     Return the sin of each data point of this Data object.
858     *
859 jgs 102 */
860 woo409 757 ESCRIPT_DLL_API
861 jgs 102 Data
862 jgs 123 sin() const;
863    
864     /**
865     \brief
866     Return the cos of each data point of this Data object.
867     *
868     */
869 woo409 757 ESCRIPT_DLL_API
870 jgs 123 Data
871     cos() const;
872    
873     /**
874     \brief
875     Return the tan of each data point of this Data object.
876     *
877     */
878 woo409 757 ESCRIPT_DLL_API
879 jgs 123 Data
880     tan() const;
881    
882     /**
883     \brief
884 jgs 150 Return the asin of each data point of this Data object.
885     *
886     */
887 woo409 757 ESCRIPT_DLL_API
888 jgs 150 Data
889     asin() const;
890    
891     /**
892     \brief
893     Return the acos of each data point of this Data object.
894     *
895     */
896 woo409 757 ESCRIPT_DLL_API
897 jgs 150 Data
898     acos() const;
899    
900     /**
901     \brief
902     Return the atan of each data point of this Data object.
903     *
904     */
905 woo409 757 ESCRIPT_DLL_API
906 jgs 150 Data
907     atan() const;
908    
909     /**
910     \brief
911     Return the sinh of each data point of this Data object.
912     *
913     */
914 woo409 757 ESCRIPT_DLL_API
915 jgs 150 Data
916     sinh() const;
917    
918     /**
919     \brief
920     Return the cosh of each data point of this Data object.
921     *
922     */
923 woo409 757 ESCRIPT_DLL_API
924 jgs 150 Data
925     cosh() const;
926    
927     /**
928     \brief
929     Return the tanh of each data point of this Data object.
930     *
931     */
932 woo409 757 ESCRIPT_DLL_API
933 jgs 150 Data
934     tanh() const;
935    
936     /**
937     \brief
938     Return the asinh of each data point of this Data object.
939     *
940     */
941 woo409 757 ESCRIPT_DLL_API
942 jgs 150 Data
943     asinh() const;
944    
945     /**
946     \brief
947     Return the acosh of each data point of this Data object.
948     *
949     */
950 woo409 757 ESCRIPT_DLL_API
951 jgs 150 Data
952     acosh() const;
953    
954     /**
955     \brief
956     Return the atanh of each data point of this Data object.
957     *
958     */
959 woo409 757 ESCRIPT_DLL_API
960 jgs 150 Data
961     atanh() const;
962    
963     /**
964     \brief
965 jgs 123 Return the log to base 10 of each data point of this Data object.
966     *
967     */
968 woo409 757 ESCRIPT_DLL_API
969 jgs 123 Data
970 gross 286 log10() const;
971 jgs 123
972     /**
973     \brief
974     Return the natural log of each data point of this Data object.
975     *
976     */
977 woo409 757 ESCRIPT_DLL_API
978 jgs 123 Data
979 gross 286 log() const;
980 jgs 123
981     /**
982     \brief
983     Return the exponential function of each data point of this Data object.
984     *
985     */
986 woo409 757 ESCRIPT_DLL_API
987 jgs 123 Data
988 jgs 102 exp() const;
989    
990     /**
991 jgs 123 \brief
992     Return the square root of each data point of this Data object.
993     *
994 jgs 102 */
995 woo409 757 ESCRIPT_DLL_API
996 jgs 102 Data
997     sqrt() const;
998    
999     /**
1000 jgs 123 \brief
1001     Return the negation of each data point of this Data object.
1002     *
1003 jgs 121 */
1004 woo409 757 ESCRIPT_DLL_API
1005 jgs 121 Data
1006     neg() const;
1007    
1008     /**
1009 jgs 123 \brief
1010     Return the identity of each data point of this Data object.
1011     Simply returns this object unmodified.
1012     *
1013 jgs 121 */
1014 woo409 757 ESCRIPT_DLL_API
1015 jgs 121 Data
1016     pos() const;
1017    
1018     /**
1019 jgs 94 \brief
1020 jgs 102 Return the given power of each data point of this Data object.
1021 jgs 121
1022     \param right Input - the power to raise the object to.
1023 jgs 123 *
1024 jgs 102 */
1025 woo409 757 ESCRIPT_DLL_API
1026 jgs 102 Data
1027     powD(const Data& right) const;
1028    
1029 jgs 121 /**
1030 jgs 123 \brief
1031     Return the given power of each data point of this boost python object.
1032    
1033     \param right Input - the power to raise the object to.
1034     *
1035 jgs 121 */
1036 woo409 757 ESCRIPT_DLL_API
1037 jgs 102 Data
1038     powO(const boost::python::object& right) const;
1039    
1040     /**
1041 jgs 123 \brief
1042 gross 699 Return the given power of each data point of this boost python object.
1043    
1044     \param left Input - the bases
1045     *
1046     */
1047    
1048 woo409 757 ESCRIPT_DLL_API
1049 gross 699 Data
1050     rpowO(const boost::python::object& left) const;
1051    
1052     /**
1053     \brief
1054 jgs 123 writes the object to a file in the DX file format
1055 jgs 104 */
1056 woo409 757 ESCRIPT_DLL_API
1057 jgs 104 void
1058     saveDX(std::string fileName) const;
1059    
1060     /**
1061 jgs 123 \brief
1062     writes the object to a file in the VTK file format
1063 jgs 110 */
1064 woo409 757 ESCRIPT_DLL_API
1065 jgs 110 void
1066     saveVTK(std::string fileName) const;
1067    
1068     /**
1069 jgs 102 \brief
1070     Overloaded operator +=
1071     \param right - Input - The right hand side.
1072 jgs 123 *
1073 jgs 102 */
1074 woo409 757 ESCRIPT_DLL_API
1075 jgs 102 Data& operator+=(const Data& right);
1076 woo409 757 ESCRIPT_DLL_API
1077 jgs 102 Data& operator+=(const boost::python::object& right);
1078    
1079     /**
1080     \brief
1081     Overloaded operator -=
1082     \param right - Input - The right hand side.
1083 jgs 123 *
1084 jgs 102 */
1085 woo409 757 ESCRIPT_DLL_API
1086 jgs 102 Data& operator-=(const Data& right);
1087 woo409 757 ESCRIPT_DLL_API
1088 jgs 102 Data& operator-=(const boost::python::object& right);
1089    
1090     /**
1091     \brief
1092     Overloaded operator *=
1093     \param right - Input - The right hand side.
1094 jgs 123 *
1095 jgs 102 */
1096 woo409 757 ESCRIPT_DLL_API
1097 jgs 102 Data& operator*=(const Data& right);
1098 woo409 757 ESCRIPT_DLL_API
1099 jgs 102 Data& operator*=(const boost::python::object& right);
1100    
1101     /**
1102     \brief
1103     Overloaded operator /=
1104     \param right - Input - The right hand side.
1105 jgs 123 *
1106 jgs 102 */
1107 woo409 757 ESCRIPT_DLL_API
1108 jgs 102 Data& operator/=(const Data& right);
1109 woo409 757 ESCRIPT_DLL_API
1110 jgs 102 Data& operator/=(const boost::python::object& right);
1111    
1112     /**
1113     \brief
1114 jgs 94 Returns true if this can be interpolated to functionspace.
1115     */
1116 woo409 757 ESCRIPT_DLL_API
1117 jgs 94 bool
1118     probeInterpolation(const FunctionSpace& functionspace) const;
1119    
1120     /**
1121 jgs 102 Data object slicing methods.
1122     */
1123    
1124     /**
1125 jgs 94 \brief
1126 jgs 102 Returns a slice from this Data object.
1127    
1128     /description
1129     Implements the [] get operator in python.
1130     Calls getSlice.
1131    
1132     \param key - Input - python slice tuple specifying
1133     slice to return.
1134 jgs 94 */
1135 woo409 757 ESCRIPT_DLL_API
1136 jgs 102 Data
1137     getItem(const boost::python::object& key) const;
1138    
1139     /**
1140     \brief
1141     Copies slice from value into this Data object.
1142    
1143     Implements the [] set operator in python.
1144     Calls setSlice.
1145    
1146     \param key - Input - python slice tuple specifying
1147     slice to copy from value.
1148     \param value - Input - Data object to copy from.
1149     */
1150 woo409 757 ESCRIPT_DLL_API
1151 jgs 94 void
1152 jgs 102 setItemD(const boost::python::object& key,
1153     const Data& value);
1154 jgs 94
1155 woo409 757 ESCRIPT_DLL_API
1156 jgs 102 void
1157     setItemO(const boost::python::object& key,
1158     const boost::python::object& value);
1159    
1160     // These following public methods should be treated as private.
1161    
1162 jgs 94 /**
1163     \brief
1164 jgs 102 Perform the given unary operation on every element of every data point in
1165     this Data object.
1166 jgs 94 */
1167 jgs 102 template <class UnaryFunction>
1168 woo409 757 ESCRIPT_DLL_API
1169 jgs 102 inline
1170 jgs 94 void
1171 jgs 102 unaryOp(UnaryFunction operation);
1172    
1173     /**
1174     \brief
1175     Return a Data object containing the specified slice of
1176     this Data object.
1177     \param region - Input - Region to copy.
1178 jgs 123 *
1179 jgs 94 */
1180 woo409 757 ESCRIPT_DLL_API
1181 jgs 102 Data
1182     getSlice(const DataArrayView::RegionType& region) const;
1183 jgs 94
1184 jgs 102 /**
1185     \brief
1186     Copy the specified slice from the given value into this
1187     Data object.
1188     \param value - Input - Data to copy from.
1189     \param region - Input - Region to copy.
1190 jgs 123 *
1191 jgs 102 */
1192 woo409 757 ESCRIPT_DLL_API
1193 jgs 102 void
1194     setSlice(const Data& value,
1195     const DataArrayView::RegionType& region);
1196    
1197 jgs 119 /**
1198     \brief
1199     Archive the current Data object to the given file.
1200     \param fileName - Input - file to archive to.
1201     */
1202 woo409 757 ESCRIPT_DLL_API
1203 jgs 119 void
1204     archiveData(const std::string fileName);
1205    
1206     /**
1207     \brief
1208     Extract the Data object archived in the given file, overwriting
1209     the current Data object.
1210     Note - the current object must be of type DataEmpty.
1211     \param fileName - Input - file to extract from.
1212 jgs 121 \param fspace - Input - a suitable FunctionSpace descibing the data.
1213 jgs 119 */
1214 woo409 757 ESCRIPT_DLL_API
1215 jgs 119 void
1216     extractData(const std::string fileName,
1217     const FunctionSpace& fspace);
1218    
1219 bcumming 751
1220 bcumming 782 #ifdef PASO_MPI
1221 bcumming 751 /**
1222     \brief
1223     print the data values to stdout. Used for debugging
1224     */
1225 bcumming 782 ESCRIPT_DLL_API
1226     void
1227     print(void);
1228 bcumming 751
1229 bcumming 782 /**
1230     \brief
1231     return the MPI rank number of the local data
1232     MPI_COMM_WORLD is assumed and the result of MPI_Comm_size()
1233     is returned
1234     */
1235     ESCRIPT_DLL_API
1236     int
1237     get_MPIRank(void) const;
1238    
1239     /**
1240     \brief
1241     return the MPI rank number of the local data
1242     MPI_COMM_WORLD is assumed and the result of MPI_Comm_rank()
1243     is returned
1244     */
1245     ESCRIPT_DLL_API
1246     int
1247     get_MPISize(void) const;
1248    
1249     /**
1250     \brief
1251     return the MPI rank number of the local data
1252     MPI_COMM_WORLD is assumed and returned.
1253     */
1254     ESCRIPT_DLL_API
1255     MPI_Comm
1256     get_MPIComm(void) const;
1257     #endif
1258    
1259 jgs 102 protected:
1260    
1261 jgs 94 private:
1262    
1263     /**
1264     \brief
1265 jgs 102 Check *this and the right operand are compatible. Throws
1266     an exception if they aren't.
1267     \param right - Input - The right hand side.
1268     */
1269     inline
1270     void
1271     operandCheck(const Data& right) const
1272     {
1273     return m_data->operandCheck(*(right.m_data.get()));
1274     }
1275    
1276     /**
1277     \brief
1278     Perform the specified reduction algorithm on every element of every data point in
1279 jgs 113 this Data object according to the given function and return the single value result.
1280 jgs 102 */
1281 jgs 147 template <class BinaryFunction>
1282 jgs 102 inline
1283     double
1284 jgs 147 algorithm(BinaryFunction operation,
1285     double initial_value) const;
1286 jgs 102
1287 jgs 113 /**
1288     \brief
1289     Reduce each data-point in this Data object using the given operation. Return a Data
1290     object with the same number of data-points, but with each data-point containing only
1291     one value - the result of the reduction operation on the corresponding data-point in
1292     this Data object
1293     */
1294 jgs 147 template <class BinaryFunction>
1295 jgs 106 inline
1296     Data
1297 jgs 147 dp_algorithm(BinaryFunction operation,
1298     double initial_value) const;
1299 jgs 106
1300 jgs 102 /**
1301     \brief
1302     Perform the given binary operation on all of the data's elements.
1303     The underlying type of the right hand side (right) determines the final
1304     type of *this after the operation. For example if the right hand side
1305     is expanded *this will be expanded if necessary.
1306     RHS is a Data object.
1307     */
1308     template <class BinaryFunction>
1309     inline
1310     void
1311     binaryOp(const Data& right,
1312     BinaryFunction operation);
1313    
1314     /**
1315     \brief
1316     Perform the given binary operation on all of the data's elements.
1317     RHS is a boost::python object.
1318     */
1319     template <class BinaryFunction>
1320     inline
1321     void
1322     binaryOp(const boost::python::object& right,
1323     BinaryFunction operation);
1324    
1325     /**
1326     \brief
1327     Convert the data type of the RHS to match this.
1328     \param right - Input - data type to match.
1329     */
1330     void
1331     typeMatchLeft(Data& right) const;
1332    
1333     /**
1334     \brief
1335     Convert the data type of this to match the RHS.
1336     \param right - Input - data type to match.
1337     */
1338     void
1339     typeMatchRight(const Data& right);
1340    
1341     /**
1342     \brief
1343 jgs 94 Construct a Data object of the appropriate type.
1344     */
1345     template <class IValueType>
1346     void
1347     initialise(const IValueType& value,
1348     const FunctionSpace& what,
1349     bool expanded);
1350    
1351     /**
1352     \brief
1353     Reshape the data point if the data point is currently rank 0.
1354     Will throw an exception if the data points are not rank 0.
1355     The original data point value is used for all values of the new
1356     data point.
1357     */
1358     void
1359     reshapeDataPoint(const DataArrayView::ShapeType& shape);
1360    
1361     //
1362 jgs 102 // pointer to the actual data object
1363 jgs 94 boost::shared_ptr<DataAbstract> m_data;
1364    
1365 jgs 123 //
1366     // pointer to the internal profiling data
1367     struct profDataEntry *profData;
1368    
1369 jgs 94 };
1370    
1371     template <class IValueType>
1372     void
1373     Data::initialise(const IValueType& value,
1374     const FunctionSpace& what,
1375     bool expanded)
1376     {
1377     //
1378     // Construct a Data object of the appropriate type.
1379     // Construct the object first as there seems to be a bug which causes
1380     // undefined behaviour if an exception is thrown during construction
1381     // within the shared_ptr constructor.
1382     if (expanded) {
1383     DataAbstract* temp=new DataExpanded(value,what);
1384 jgs 102 boost::shared_ptr<DataAbstract> temp_data(temp);
1385     m_data=temp_data;
1386 jgs 94 } else {
1387     DataAbstract* temp=new DataConstant(value,what);
1388 jgs 102 boost::shared_ptr<DataAbstract> temp_data(temp);
1389     m_data=temp_data;
1390 jgs 94 }
1391     }
1392    
1393 jgs 102 /**
1394     Binary Data object operators.
1395     */
1396 jgs 94
1397     /**
1398     \brief
1399     Operator+
1400     Takes two Data objects.
1401     */
1402 woo409 757 ESCRIPT_DLL_API Data operator+(const Data& left, const Data& right);
1403 jgs 94
1404     /**
1405     \brief
1406     Operator-
1407     Takes two Data objects.
1408     */
1409 woo409 757 ESCRIPT_DLL_API Data operator-(const Data& left, const Data& right);
1410 jgs 94
1411     /**
1412     \brief
1413     Operator*
1414     Takes two Data objects.
1415     */
1416 woo409 757 ESCRIPT_DLL_API Data operator*(const Data& left, const Data& right);
1417 jgs 94
1418     /**
1419     \brief
1420     Operator/
1421     Takes two Data objects.
1422     */
1423 woo409 757 ESCRIPT_DLL_API Data operator/(const Data& left, const Data& right);
1424 jgs 94
1425     /**
1426     \brief
1427     Operator+
1428     Takes LHS Data object and RHS python::object.
1429     python::object must be convertable to Data type.
1430     */
1431 woo409 757 ESCRIPT_DLL_API Data operator+(const Data& left, const boost::python::object& right);
1432 jgs 94
1433     /**
1434     \brief
1435     Operator-
1436     Takes LHS Data object and RHS python::object.
1437     python::object must be convertable to Data type.
1438     */
1439 woo409 757 ESCRIPT_DLL_API Data operator-(const Data& left, const boost::python::object& right);
1440 jgs 94
1441     /**
1442     \brief
1443     Operator*
1444     Takes LHS Data object and RHS python::object.
1445     python::object must be convertable to Data type.
1446     */
1447 woo409 757 ESCRIPT_DLL_API Data operator*(const Data& left, const boost::python::object& right);
1448 jgs 94
1449     /**
1450     \brief
1451     Operator/
1452     Takes LHS Data object and RHS python::object.
1453     python::object must be convertable to Data type.
1454     */
1455 woo409 757 ESCRIPT_DLL_API Data operator/(const Data& left, const boost::python::object& right);
1456 jgs 94
1457     /**
1458     \brief
1459     Operator+
1460     Takes LHS python::object and RHS Data object.
1461     python::object must be convertable to Data type.
1462     */
1463 woo409 757 ESCRIPT_DLL_API Data operator+(const boost::python::object& left, const Data& right);
1464 jgs 94
1465     /**
1466     \brief
1467     Operator-
1468     Takes LHS python::object and RHS Data object.
1469     python::object must be convertable to Data type.
1470     */
1471 woo409 757 ESCRIPT_DLL_API Data operator-(const boost::python::object& left, const Data& right);
1472 jgs 94
1473     /**
1474     \brief
1475     Operator*
1476     Takes LHS python::object and RHS Data object.
1477     python::object must be convertable to Data type.
1478     */
1479 woo409 757 ESCRIPT_DLL_API Data operator*(const boost::python::object& left, const Data& right);
1480 jgs 94
1481     /**
1482     \brief
1483     Operator/
1484     Takes LHS python::object and RHS Data object.
1485     python::object must be convertable to Data type.
1486     */
1487 woo409 757 ESCRIPT_DLL_API Data operator/(const boost::python::object& left, const Data& right);
1488 jgs 94
1489     /**
1490     \brief
1491     Output operator
1492     */
1493 woo409 757 ESCRIPT_DLL_API std::ostream& operator<<(std::ostream& o, const Data& data);
1494 jgs 94
1495     /**
1496     \brief
1497     Return true if operands are equivalent, else return false.
1498     NB: this operator does very little at this point, and isn't to
1499     be relied on. Requires further implementation.
1500     */
1501 woo409 757 //ESCRIPT_DLL_API bool operator==(const Data& left, const Data& right);
1502 jgs 94
1503 jgs 102 /**
1504     \brief
1505     Perform the given binary operation with this and right as operands.
1506     Right is a Data object.
1507     */
1508 jgs 94 template <class BinaryFunction>
1509     inline
1510     void
1511     Data::binaryOp(const Data& right,
1512     BinaryFunction operation)
1513     {
1514     //
1515     // if this has a rank of zero promote it to the rank of the RHS
1516     if (getPointDataView().getRank()==0 && right.getPointDataView().getRank()!=0) {
1517     reshapeDataPoint(right.getPointDataView().getShape());
1518     }
1519     //
1520     // initially make the temporary a shallow copy
1521     Data tempRight(right);
1522     if (getFunctionSpace()!=right.getFunctionSpace()) {
1523     if (right.probeInterpolation(getFunctionSpace())) {
1524     //
1525     // an interpolation is required so create a new Data
1526     tempRight=Data(right,this->getFunctionSpace());
1527     } else if (probeInterpolation(right.getFunctionSpace())) {
1528     //
1529     // interpolate onto the RHS function space
1530     Data tempLeft(*this,right.getFunctionSpace());
1531     m_data=tempLeft.m_data;
1532     }
1533     }
1534     operandCheck(tempRight);
1535     //
1536     // ensure this has the right type for the RHS
1537 jgs 102 typeMatchRight(tempRight);
1538 jgs 94 //
1539     // Need to cast to the concrete types so that the correct binaryOp
1540     // is called.
1541     if (isExpanded()) {
1542     //
1543     // Expanded data will be done in parallel, the right hand side can be
1544     // of any data type
1545     DataExpanded* leftC=dynamic_cast<DataExpanded*>(m_data.get());
1546     EsysAssert((leftC!=0), "Programming error - casting to DataExpanded.");
1547     escript::binaryOp(*leftC,*(tempRight.m_data.get()),operation);
1548     } else if (isTagged()) {
1549     //
1550     // Tagged data is operated on serially, the right hand side can be
1551     // either DataConstant or DataTagged
1552     DataTagged* leftC=dynamic_cast<DataTagged*>(m_data.get());
1553     EsysAssert((leftC!=0), "Programming error - casting to DataTagged.");
1554     if (right.isTagged()) {
1555     DataTagged* rightC=dynamic_cast<DataTagged*>(tempRight.m_data.get());
1556     EsysAssert((rightC!=0), "Programming error - casting to DataTagged.");
1557     escript::binaryOp(*leftC,*rightC,operation);
1558     } else {
1559     DataConstant* rightC=dynamic_cast<DataConstant*>(tempRight.m_data.get());
1560     EsysAssert((rightC!=0), "Programming error - casting to DataConstant.");
1561     escript::binaryOp(*leftC,*rightC,operation);
1562     }
1563 jgs 102 } else if (isConstant()) {
1564 jgs 94 DataConstant* leftC=dynamic_cast<DataConstant*>(m_data.get());
1565     DataConstant* rightC=dynamic_cast<DataConstant*>(tempRight.m_data.get());
1566     EsysAssert((leftC!=0 && rightC!=0), "Programming error - casting to DataConstant.");
1567     escript::binaryOp(*leftC,*rightC,operation);
1568     }
1569     }
1570    
1571 jgs 102 /**
1572     \brief
1573     Perform the given binary operation with this and right as operands.
1574     Right is a boost::python object.
1575     */
1576 jgs 94 template <class BinaryFunction>
1577     inline
1578     void
1579     Data::binaryOp(const boost::python::object& right,
1580     BinaryFunction operation)
1581     {
1582     DataArray temp(right);
1583     //
1584     // if this has a rank of zero promote it to the rank of the RHS.
1585     if (getPointDataView().getRank()==0 && temp.getView().getRank()!=0) {
1586     reshapeDataPoint(temp.getView().getShape());
1587     }
1588     //
1589     // Always allow scalar values for the RHS but check other shapes
1590     if (temp.getView().getRank()!=0) {
1591     if (!getPointDataView().checkShape(temp.getView().getShape())) {
1592     throw DataException(getPointDataView().createShapeErrorMessage(
1593     "Error - RHS shape doesn't match LHS shape.",temp.getView().getShape()));
1594     }
1595     }
1596     if (isExpanded()) {
1597     DataExpanded* leftC=dynamic_cast<DataExpanded*>(m_data.get());
1598     EsysAssert((leftC!=0),"Programming error - casting to DataExpanded.");
1599     escript::binaryOp(*leftC,temp.getView(),operation);
1600     } else if (isTagged()) {
1601     DataTagged* leftC=dynamic_cast<DataTagged*>(m_data.get());
1602     EsysAssert((leftC!=0), "Programming error - casting to DataTagged.");
1603     escript::binaryOp(*leftC,temp.getView(),operation);
1604 jgs 102 } else if (isConstant()) {
1605 jgs 94 DataConstant* leftC=dynamic_cast<DataConstant*>(m_data.get());
1606     EsysAssert((leftC!=0),"Programming error - casting to DataConstant.");
1607     escript::binaryOp(*leftC,temp.getView(),operation);
1608     }
1609     }
1610    
1611 jgs 102 /**
1612     \brief
1613     Perform the given unary operation on other and return the result.
1614     Given operation is performed on each element of each data point, thus
1615     argument object is a rank n Data object, and returned object is a rank n
1616     Data object.
1617     Calls Data::unaryOp.
1618     */
1619 jgs 94 template <class UnaryFunction>
1620     inline
1621 jgs 102 Data
1622     unaryOp(const Data& other,
1623     UnaryFunction operation)
1624     {
1625     Data result;
1626     result.copy(other);
1627     result.unaryOp(operation);
1628     return result;
1629     }
1630    
1631     /**
1632     \brief
1633     Perform the given unary operation on this.
1634     Given operation is performed on each element of each data point.
1635     Calls escript::unaryOp.
1636     */
1637     template <class UnaryFunction>
1638     inline
1639 jgs 94 void
1640     Data::unaryOp(UnaryFunction operation)
1641     {
1642     if (isExpanded()) {
1643     DataExpanded* leftC=dynamic_cast<DataExpanded*>(m_data.get());
1644     EsysAssert((leftC!=0), "Programming error - casting to DataExpanded.");
1645     escript::unaryOp(*leftC,operation);
1646     } else if (isTagged()) {
1647     DataTagged* leftC=dynamic_cast<DataTagged*>(m_data.get());
1648     EsysAssert((leftC!=0), "Programming error - casting to DataTagged.");
1649     escript::unaryOp(*leftC,operation);
1650 jgs 102 } else if (isConstant()) {
1651 jgs 94 DataConstant* leftC=dynamic_cast<DataConstant*>(m_data.get());
1652     EsysAssert((leftC!=0), "Programming error - casting to DataConstant.");
1653     escript::unaryOp(*leftC,operation);
1654     }
1655     }
1656    
1657 jgs 102 /**
1658     \brief
1659     Perform the given Data object reduction algorithm on this and return the result.
1660     Given operation combines each element of each data point, thus argument
1661     object (*this) is a rank n Data object, and returned object is a scalar.
1662     Calls escript::algorithm.
1663     */
1664 jgs 147 template <class BinaryFunction>
1665 jgs 94 inline
1666     double
1667 jgs 147 Data::algorithm(BinaryFunction operation, double initial_value) const
1668 jgs 94 {
1669     if (isExpanded()) {
1670     DataExpanded* leftC=dynamic_cast<DataExpanded*>(m_data.get());
1671     EsysAssert((leftC!=0), "Programming error - casting to DataExpanded.");
1672 jgs 147 return escript::algorithm(*leftC,operation,initial_value);
1673 jgs 102 } else if (isTagged()) {
1674 jgs 94 DataTagged* leftC=dynamic_cast<DataTagged*>(m_data.get());
1675     EsysAssert((leftC!=0), "Programming error - casting to DataTagged.");
1676 jgs 147 return escript::algorithm(*leftC,operation,initial_value);
1677 jgs 102 } else if (isConstant()) {
1678 jgs 94 DataConstant* leftC=dynamic_cast<DataConstant*>(m_data.get());
1679     EsysAssert((leftC!=0), "Programming error - casting to DataConstant.");
1680 jgs 147 return escript::algorithm(*leftC,operation,initial_value);
1681 jgs 94 }
1682 jgs 102 return 0;
1683 jgs 94 }
1684    
1685 jgs 102 /**
1686     \brief
1687     Perform the given data point reduction algorithm on data and return the result.
1688     Given operation combines each element within each data point into a scalar,
1689     thus argument object is a rank n Data object, and returned object is a
1690     rank 0 Data object.
1691     Calls escript::dp_algorithm.
1692     */
1693 jgs 147 template <class BinaryFunction>
1694 jgs 94 inline
1695     Data
1696 jgs 147 Data::dp_algorithm(BinaryFunction operation, double initial_value) const
1697 jgs 94 {
1698 jgs 106 if (isExpanded()) {
1699 jgs 559 Data result(0,DataArrayView::ShapeType(),getFunctionSpace(),isExpanded());
1700 jgs 106 DataExpanded* dataE=dynamic_cast<DataExpanded*>(m_data.get());
1701 jgs 102 DataExpanded* resultE=dynamic_cast<DataExpanded*>(result.m_data.get());
1702     EsysAssert((dataE!=0), "Programming error - casting data to DataExpanded.");
1703     EsysAssert((resultE!=0), "Programming error - casting result to DataExpanded.");
1704 jgs 147 escript::dp_algorithm(*dataE,*resultE,operation,initial_value);
1705 jgs 559 return result;
1706 jgs 106 } else if (isTagged()) {
1707     DataTagged* dataT=dynamic_cast<DataTagged*>(m_data.get());
1708 jgs 562 DataArrayView::ShapeType viewShape;
1709     DataArrayView::ValueType viewData(1);
1710     viewData[0]=0;
1711     DataArrayView defaultValue(viewData,viewShape);
1712 jgs 559 DataTagged::TagListType keys;
1713 jgs 562 DataTagged::ValueListType values;
1714 jgs 559 DataTagged::DataMapType::const_iterator i;
1715     for (i=dataT->getTagLookup().begin();i!=dataT->getTagLookup().end();i++) {
1716     keys.push_back(i->first);
1717 jgs 562 values.push_back(defaultValue);
1718 jgs 559 }
1719     Data result(keys,values,defaultValue,getFunctionSpace());
1720 jgs 102 DataTagged* resultT=dynamic_cast<DataTagged*>(result.m_data.get());
1721     EsysAssert((dataT!=0), "Programming error - casting data to DataTagged.");
1722     EsysAssert((resultT!=0), "Programming error - casting result to DataTagged.");
1723 jgs 147 escript::dp_algorithm(*dataT,*resultT,operation,initial_value);
1724 jgs 559 return result;
1725 jgs 106 } else if (isConstant()) {
1726 jgs 559 Data result(0,DataArrayView::ShapeType(),getFunctionSpace(),isExpanded());
1727 jgs 106 DataConstant* dataC=dynamic_cast<DataConstant*>(m_data.get());
1728 jgs 102 DataConstant* resultC=dynamic_cast<DataConstant*>(result.m_data.get());
1729     EsysAssert((dataC!=0), "Programming error - casting data to DataConstant.");
1730     EsysAssert((resultC!=0), "Programming error - casting result to DataConstant.");
1731 jgs 147 escript::dp_algorithm(*dataC,*resultC,operation,initial_value);
1732 jgs 559 return result;
1733 jgs 102 }
1734 jgs 559 Data falseRetVal; // to keep compiler quiet
1735     return falseRetVal;
1736 jgs 94 }
1737    
1738     }
1739     #endif

Properties

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

  ViewVC Help
Powered by ViewVC 1.1.26