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

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

Parent Directory Parent Directory | Revision Log Revision Log


Revision 775 - (show annotations)
Mon Jul 10 04:00:08 2006 UTC (13 years, 3 months ago) by ksteube
File size: 59356 byte(s)
Modified the following python methods in escript/py_src/util.py to
call faster C++ methods:
	escript_trace
	escript_transpose
	escript_symmetric
	escript_nonsymmetric

1 // $Id$
2
3 /*
4 ************************************************************
5 * Copyright 2006 by ACcESS MNRF *
6 * *
7 * http://www.access.edu.au *
8 * Primary Business: Queensland, Australia *
9 * Licensed under the Open Software License version 3.0 *
10 * http://www.opensource.org/licenses/osl-3.0.php *
11 * *
12 ************************************************************
13 */
14 #include "Data.h"
15
16 #include "DataExpanded.h"
17 #include "DataConstant.h"
18 #include "DataTagged.h"
19 #include "DataEmpty.h"
20 #include "DataArray.h"
21 #include "DataArrayView.h"
22 #include "DataProf.h"
23 #include "FunctionSpaceFactory.h"
24 #include "AbstractContinuousDomain.h"
25 #include "UnaryFuncs.h"
26
27 #include <fstream>
28 #include <algorithm>
29 #include <vector>
30 #include <functional>
31
32 #include <boost/python/dict.hpp>
33 #include <boost/python/extract.hpp>
34 #include <boost/python/long.hpp>
35
36 using namespace std;
37 using namespace boost::python;
38 using namespace boost;
39 using namespace escript;
40
41 #if defined DOPROF
42 //
43 // global table of profiling data for all Data objects
44 DataProf dataProfTable;
45 #endif
46
47 Data::Data()
48 {
49 //
50 // Default data is type DataEmpty
51 DataAbstract* temp=new DataEmpty();
52 shared_ptr<DataAbstract> temp_data(temp);
53 m_data=temp_data;
54 #if defined DOPROF
55 // create entry in global profiling table for this object
56 profData = dataProfTable.newData();
57 #endif
58 }
59
60 Data::Data(double value,
61 const tuple& shape,
62 const FunctionSpace& what,
63 bool expanded)
64 {
65 DataArrayView::ShapeType dataPointShape;
66 for (int i = 0; i < shape.attr("__len__")(); ++i) {
67 dataPointShape.push_back(extract<const int>(shape[i]));
68 }
69 DataArray temp(dataPointShape,value);
70 initialise(temp.getView(),what,expanded);
71 #if defined DOPROF
72 // create entry in global profiling table for this object
73 profData = dataProfTable.newData();
74 #endif
75 }
76
77 Data::Data(double value,
78 const DataArrayView::ShapeType& dataPointShape,
79 const FunctionSpace& what,
80 bool expanded)
81 {
82 DataArray temp(dataPointShape,value);
83 pair<int,int> dataShape=what.getDataShape();
84 initialise(temp.getView(),what,expanded);
85 #if defined DOPROF
86 // create entry in global profiling table for this object
87 profData = dataProfTable.newData();
88 #endif
89 }
90
91 Data::Data(const Data& inData)
92 {
93 m_data=inData.m_data;
94 #if defined DOPROF
95 // create entry in global profiling table for this object
96 profData = dataProfTable.newData();
97 #endif
98 }
99
100 Data::Data(const Data& inData,
101 const DataArrayView::RegionType& region)
102 {
103 //
104 // Create Data which is a slice of another Data
105 DataAbstract* tmp = inData.m_data->getSlice(region);
106 shared_ptr<DataAbstract> temp_data(tmp);
107 m_data=temp_data;
108 #if defined DOPROF
109 // create entry in global profiling table for this object
110 profData = dataProfTable.newData();
111 #endif
112 }
113
114 Data::Data(const Data& inData,
115 const FunctionSpace& functionspace)
116 {
117 #if defined DOPROF
118 // create entry in global profiling table for this object
119 profData = dataProfTable.newData();
120 #endif
121 if (inData.getFunctionSpace()==functionspace) {
122 m_data=inData.m_data;
123 } else {
124 #if defined DOPROF
125 profData->interpolate++;
126 #endif
127 Data tmp(0,inData.getPointDataView().getShape(),functionspace,true);
128 // Note: Must use a reference or pointer to a derived object
129 // in order to get polymorphic behaviour. Shouldn't really
130 // be able to create an instance of AbstractDomain but that was done
131 // as a boost:python work around which may no longer be required.
132 const AbstractDomain& inDataDomain=inData.getDomain();
133 if (inDataDomain==functionspace.getDomain()) {
134 inDataDomain.interpolateOnDomain(tmp,inData);
135 } else {
136 inDataDomain.interpolateACross(tmp,inData);
137 }
138 m_data=tmp.m_data;
139 }
140 }
141
142 Data::Data(const DataTagged::TagListType& tagKeys,
143 const DataTagged::ValueListType & values,
144 const DataArrayView& defaultValue,
145 const FunctionSpace& what,
146 bool expanded)
147 {
148 DataAbstract* temp=new DataTagged(tagKeys,values,defaultValue,what);
149 shared_ptr<DataAbstract> temp_data(temp);
150 m_data=temp_data;
151 if (expanded) {
152 expand();
153 }
154 #if defined DOPROF
155 // create entry in global profiling table for this object
156 profData = dataProfTable.newData();
157 #endif
158 }
159
160 Data::Data(const numeric::array& value,
161 const FunctionSpace& what,
162 bool expanded)
163 {
164 initialise(value,what,expanded);
165 #if defined DOPROF
166 // create entry in global profiling table for this object
167 profData = dataProfTable.newData();
168 #endif
169 }
170
171 Data::Data(const DataArrayView& value,
172 const FunctionSpace& what,
173 bool expanded)
174 {
175 initialise(value,what,expanded);
176 #if defined DOPROF
177 // create entry in global profiling table for this object
178 profData = dataProfTable.newData();
179 #endif
180 }
181
182 Data::Data(const object& value,
183 const FunctionSpace& what,
184 bool expanded)
185 {
186 numeric::array asNumArray(value);
187 initialise(asNumArray,what,expanded);
188 #if defined DOPROF
189 // create entry in global profiling table for this object
190 profData = dataProfTable.newData();
191 #endif
192 }
193
194 Data::Data(const object& value,
195 const Data& other)
196 {
197 //
198 // Create DataConstant using the given value and all other parameters
199 // copied from other. If value is a rank 0 object this Data
200 // will assume the point data shape of other.
201 DataArray temp(value);
202 if (temp.getView().getRank()==0) {
203 //
204 // Create a DataArray with the scalar value for all elements
205 DataArray temp2(other.getPointDataView().getShape(),temp.getView()());
206 initialise(temp2.getView(),other.getFunctionSpace(),false);
207 } else {
208 //
209 // Create a DataConstant with the same sample shape as other
210 initialise(temp.getView(),other.getFunctionSpace(),false);
211 }
212 #if defined DOPROF
213 // create entry in global profiling table for this object
214 profData = dataProfTable.newData();
215 #endif
216 }
217
218 Data::~Data()
219 {
220
221 }
222
223 escriptDataC
224 Data::getDataC()
225 {
226 escriptDataC temp;
227 temp.m_dataPtr=(void*)this;
228 return temp;
229 }
230
231 escriptDataC
232 Data::getDataC() const
233 {
234 escriptDataC temp;
235 temp.m_dataPtr=(void*)this;
236 return temp;
237 }
238
239 const boost::python::tuple
240 Data::getShapeTuple() const
241 {
242 const DataArrayView::ShapeType& shape=getDataPointShape();
243 switch(getDataPointRank()) {
244 case 0:
245 return make_tuple();
246 case 1:
247 return make_tuple(long_(shape[0]));
248 case 2:
249 return make_tuple(long_(shape[0]),long_(shape[1]));
250 case 3:
251 return make_tuple(long_(shape[0]),long_(shape[1]),long_(shape[2]));
252 case 4:
253 return make_tuple(long_(shape[0]),long_(shape[1]),long_(shape[2]),long_(shape[3]));
254 default:
255 throw DataException("Error - illegal Data rank.");
256 }
257 }
258
259 void
260 Data::copy(const Data& other)
261 {
262 //
263 // Perform a deep copy
264 {
265 DataExpanded* temp=dynamic_cast<DataExpanded*>(other.m_data.get());
266 if (temp!=0) {
267 //
268 // Construct a DataExpanded copy
269 DataAbstract* newData=new DataExpanded(*temp);
270 shared_ptr<DataAbstract> temp_data(newData);
271 m_data=temp_data;
272 return;
273 }
274 }
275 {
276 DataTagged* temp=dynamic_cast<DataTagged*>(other.m_data.get());
277 if (temp!=0) {
278 //
279 // Construct a DataTagged copy
280 DataAbstract* newData=new DataTagged(*temp);
281 shared_ptr<DataAbstract> temp_data(newData);
282 m_data=temp_data;
283 return;
284 }
285 }
286 {
287 DataConstant* temp=dynamic_cast<DataConstant*>(other.m_data.get());
288 if (temp!=0) {
289 //
290 // Construct a DataConstant copy
291 DataAbstract* newData=new DataConstant(*temp);
292 shared_ptr<DataAbstract> temp_data(newData);
293 m_data=temp_data;
294 return;
295 }
296 }
297 {
298 DataEmpty* temp=dynamic_cast<DataEmpty*>(other.m_data.get());
299 if (temp!=0) {
300 //
301 // Construct a DataEmpty copy
302 DataAbstract* newData=new DataEmpty();
303 shared_ptr<DataAbstract> temp_data(newData);
304 m_data=temp_data;
305 return;
306 }
307 }
308 throw DataException("Error - Copy not implemented for this Data type.");
309 }
310
311 void
312 Data::copyWithMask(const Data& other,
313 const Data& mask)
314 {
315 Data mask1;
316 Data mask2;
317
318 mask1 = mask.wherePositive();
319 mask2.copy(mask1);
320
321 mask1 *= other;
322 mask2 *= *this;
323 mask2 = *this - mask2;
324
325 *this = mask1 + mask2;
326 }
327
328 bool
329 Data::isExpanded() const
330 {
331 DataExpanded* temp=dynamic_cast<DataExpanded*>(m_data.get());
332 return (temp!=0);
333 }
334
335 bool
336 Data::isTagged() const
337 {
338 DataTagged* temp=dynamic_cast<DataTagged*>(m_data.get());
339 return (temp!=0);
340 }
341
342 /* TODO */
343 /* global reduction -- the local data being empty does not imply that it is empty on other processers*/
344 bool
345 Data::isEmpty() const
346 {
347 DataEmpty* temp=dynamic_cast<DataEmpty*>(m_data.get());
348 return (temp!=0);
349 }
350
351 bool
352 Data::isConstant() const
353 {
354 DataConstant* temp=dynamic_cast<DataConstant*>(m_data.get());
355 return (temp!=0);
356 }
357
358 void
359 Data::expand()
360 {
361 if (isConstant()) {
362 DataConstant* tempDataConst=dynamic_cast<DataConstant*>(m_data.get());
363 DataAbstract* temp=new DataExpanded(*tempDataConst);
364 shared_ptr<DataAbstract> temp_data(temp);
365 m_data=temp_data;
366 } else if (isTagged()) {
367 DataTagged* tempDataTag=dynamic_cast<DataTagged*>(m_data.get());
368 DataAbstract* temp=new DataExpanded(*tempDataTag);
369 shared_ptr<DataAbstract> temp_data(temp);
370 m_data=temp_data;
371 } else if (isExpanded()) {
372 //
373 // do nothing
374 } else if (isEmpty()) {
375 throw DataException("Error - Expansion of DataEmpty not possible.");
376 } else {
377 throw DataException("Error - Expansion not implemented for this Data type.");
378 }
379 }
380
381 void
382 Data::tag()
383 {
384 if (isConstant()) {
385 DataConstant* tempDataConst=dynamic_cast<DataConstant*>(m_data.get());
386 DataAbstract* temp=new DataTagged(*tempDataConst);
387 shared_ptr<DataAbstract> temp_data(temp);
388 m_data=temp_data;
389 } else if (isTagged()) {
390 // do nothing
391 } else if (isExpanded()) {
392 throw DataException("Error - Creating tag data from DataExpanded not possible.");
393 } else if (isEmpty()) {
394 throw DataException("Error - Creating tag data from DataEmpty not possible.");
395 } else {
396 throw DataException("Error - Tagging not implemented for this Data type.");
397 }
398 }
399
400 void
401 Data::reshapeDataPoint(const DataArrayView::ShapeType& shape)
402 {
403 m_data->reshapeDataPoint(shape);
404 }
405
406 Data
407 Data::wherePositive() const
408 {
409 #if defined DOPROF
410 profData->where++;
411 #endif
412 return escript::unaryOp(*this,bind2nd(greater<double>(),0.0));
413 }
414
415 Data
416 Data::whereNegative() const
417 {
418 #if defined DOPROF
419 profData->where++;
420 #endif
421 return escript::unaryOp(*this,bind2nd(less<double>(),0.0));
422 }
423
424 Data
425 Data::whereNonNegative() const
426 {
427 #if defined DOPROF
428 profData->where++;
429 #endif
430 return escript::unaryOp(*this,bind2nd(greater_equal<double>(),0.0));
431 }
432
433 Data
434 Data::whereNonPositive() const
435 {
436 #if defined DOPROF
437 profData->where++;
438 #endif
439 return escript::unaryOp(*this,bind2nd(less_equal<double>(),0.0));
440 }
441
442 Data
443 Data::whereZero(double tol) const
444 {
445 #if defined DOPROF
446 profData->where++;
447 #endif
448 Data dataAbs=abs();
449 return escript::unaryOp(dataAbs,bind2nd(less_equal<double>(),tol));
450 }
451
452 Data
453 Data::whereNonZero(double tol) const
454 {
455 #if defined DOPROF
456 profData->where++;
457 #endif
458 Data dataAbs=abs();
459 return escript::unaryOp(dataAbs,bind2nd(greater<double>(),tol));
460 }
461
462 Data
463 Data::interpolate(const FunctionSpace& functionspace) const
464 {
465 #if defined DOPROF
466 profData->interpolate++;
467 #endif
468 return Data(*this,functionspace);
469 }
470
471 bool
472 Data::probeInterpolation(const FunctionSpace& functionspace) const
473 {
474 if (getFunctionSpace()==functionspace) {
475 return true;
476 } else {
477 const AbstractDomain& domain=getDomain();
478 if (domain==functionspace.getDomain()) {
479 return domain.probeInterpolationOnDomain(getFunctionSpace().getTypeCode(),functionspace.getTypeCode());
480 } else {
481 return domain.probeInterpolationACross(getFunctionSpace().getTypeCode(),functionspace.getDomain(),functionspace.getTypeCode());
482 }
483 }
484 }
485
486 Data
487 Data::gradOn(const FunctionSpace& functionspace) const
488 {
489 #if defined DOPROF
490 profData->grad++;
491 #endif
492 if (functionspace.getDomain()!=getDomain())
493 throw DataException("Error - gradient cannot be calculated on different domains.");
494 DataArrayView::ShapeType grad_shape=getPointDataView().getShape();
495 grad_shape.push_back(functionspace.getDim());
496 Data out(0.0,grad_shape,functionspace,true);
497 getDomain().setToGradient(out,*this);
498 return out;
499 }
500
501 Data
502 Data::grad() const
503 {
504 return gradOn(escript::function(getDomain()));
505 }
506
507 int
508 Data::getDataPointSize() const
509 {
510 return getPointDataView().noValues();
511 }
512
513 DataArrayView::ValueType::size_type
514 Data::getLength() const
515 {
516 return m_data->getLength();
517 }
518
519 const DataArrayView::ShapeType&
520 Data::getDataPointShape() const
521 {
522 return getPointDataView().getShape();
523 }
524
525 void
526 Data::fillFromNumArray(const boost::python::numeric::array num_array)
527 {
528 //
529 // check rank
530 if (num_array.getrank()<getDataPointRank())
531 throw DataException("Rank of numarray does not match Data object rank");
532
533 //
534 // check shape of num_array
535 for (int i=0; i<getDataPointRank(); i++) {
536 if (extract<int>(num_array.getshape()[i+1])!=getDataPointShape()[i])
537 throw DataException("Shape of numarray does not match Data object rank");
538 }
539
540 //
541 // make sure data is expanded:
542 if (!isExpanded()) {
543 expand();
544 }
545
546 //
547 // and copy over
548 m_data->copyAll(num_array);
549 }
550
551 const
552 boost::python::numeric::array
553 Data::convertToNumArray()
554 {
555 //
556 // determine the total number of data points
557 int numSamples = getNumSamples();
558 int numDataPointsPerSample = getNumDataPointsPerSample();
559 int numDataPoints = numSamples * numDataPointsPerSample;
560
561 //
562 // determine the rank and shape of each data point
563 int dataPointRank = getDataPointRank();
564 DataArrayView::ShapeType dataPointShape = getDataPointShape();
565
566 //
567 // create the numeric array to be returned
568 boost::python::numeric::array numArray(0.0);
569
570 //
571 // the rank of the returned numeric array will be the rank of
572 // the data points, plus one. Where the rank of the array is n,
573 // the last n-1 dimensions will be equal to the shape of the
574 // data points, whilst the first dimension will be equal to the
575 // total number of data points. Thus the array will consist of
576 // a serial vector of the data points.
577 int arrayRank = dataPointRank + 1;
578 DataArrayView::ShapeType arrayShape;
579 arrayShape.push_back(numDataPoints);
580 for (int d=0; d<dataPointRank; d++) {
581 arrayShape.push_back(dataPointShape[d]);
582 }
583
584 //
585 // resize the numeric array to the shape just calculated
586 if (arrayRank==1) {
587 numArray.resize(arrayShape[0]);
588 }
589 if (arrayRank==2) {
590 numArray.resize(arrayShape[0],arrayShape[1]);
591 }
592 if (arrayRank==3) {
593 numArray.resize(arrayShape[0],arrayShape[1],arrayShape[2]);
594 }
595 if (arrayRank==4) {
596 numArray.resize(arrayShape[0],arrayShape[1],arrayShape[2],arrayShape[3]);
597 }
598 if (arrayRank==5) {
599 numArray.resize(arrayShape[0],arrayShape[1],arrayShape[2],arrayShape[3],arrayShape[4]);
600 }
601
602 //
603 // loop through each data point in turn, loading the values for that data point
604 // into the numeric array.
605 int dataPoint = 0;
606 for (int sampleNo = 0; sampleNo < numSamples; sampleNo++) {
607 for (int dataPointNo = 0; dataPointNo < numDataPointsPerSample; dataPointNo++) {
608 DataArrayView dataPointView = getDataPoint(sampleNo, dataPointNo);
609 if (dataPointRank==0) {
610 numArray[dataPoint]=dataPointView();
611 }
612 if (dataPointRank==1) {
613 for (int i=0; i<dataPointShape[0]; i++) {
614 numArray[dataPoint][i]=dataPointView(i);
615 }
616 }
617 if (dataPointRank==2) {
618 for (int i=0; i<dataPointShape[0]; i++) {
619 for (int j=0; j<dataPointShape[1]; j++) {
620 numArray[dataPoint][i][j] = dataPointView(i,j);
621 }
622 }
623 }
624 if (dataPointRank==3) {
625 for (int i=0; i<dataPointShape[0]; i++) {
626 for (int j=0; j<dataPointShape[1]; j++) {
627 for (int k=0; k<dataPointShape[2]; k++) {
628 numArray[dataPoint][i][j][k]=dataPointView(i,j,k);
629 }
630 }
631 }
632 }
633 if (dataPointRank==4) {
634 for (int i=0; i<dataPointShape[0]; i++) {
635 for (int j=0; j<dataPointShape[1]; j++) {
636 for (int k=0; k<dataPointShape[2]; k++) {
637 for (int l=0; l<dataPointShape[3]; l++) {
638 numArray[dataPoint][i][j][k][l]=dataPointView(i,j,k,l);
639 }
640 }
641 }
642 }
643 }
644 dataPoint++;
645 }
646 }
647
648 //
649 // return the loaded array
650 return numArray;
651 }
652
653 const
654 boost::python::numeric::array
655 Data::convertToNumArrayFromSampleNo(int sampleNo)
656 {
657 //
658 // Check a valid sample number has been supplied
659 if (sampleNo >= getNumSamples()) {
660 throw DataException("Error - Data::convertToNumArray: invalid sampleNo.");
661 }
662
663 //
664 // determine the number of data points per sample
665 int numDataPointsPerSample = getNumDataPointsPerSample();
666
667 //
668 // determine the rank and shape of each data point
669 int dataPointRank = getDataPointRank();
670 DataArrayView::ShapeType dataPointShape = getDataPointShape();
671
672 //
673 // create the numeric array to be returned
674 boost::python::numeric::array numArray(0.0);
675
676 //
677 // the rank of the returned numeric array will be the rank of
678 // the data points, plus one. Where the rank of the array is n,
679 // the last n-1 dimensions will be equal to the shape of the
680 // data points, whilst the first dimension will be equal to the
681 // total number of data points. Thus the array will consist of
682 // a serial vector of the data points.
683 int arrayRank = dataPointRank + 1;
684 DataArrayView::ShapeType arrayShape;
685 arrayShape.push_back(numDataPointsPerSample);
686 for (int d=0; d<dataPointRank; d++) {
687 arrayShape.push_back(dataPointShape[d]);
688 }
689
690 //
691 // resize the numeric array to the shape just calculated
692 if (arrayRank==1) {
693 numArray.resize(arrayShape[0]);
694 }
695 if (arrayRank==2) {
696 numArray.resize(arrayShape[0],arrayShape[1]);
697 }
698 if (arrayRank==3) {
699 numArray.resize(arrayShape[0],arrayShape[1],arrayShape[2]);
700 }
701 if (arrayRank==4) {
702 numArray.resize(arrayShape[0],arrayShape[1],arrayShape[2],arrayShape[3]);
703 }
704 if (arrayRank==5) {
705 numArray.resize(arrayShape[0],arrayShape[1],arrayShape[2],arrayShape[3],arrayShape[4]);
706 }
707
708 //
709 // loop through each data point in turn, loading the values for that data point
710 // into the numeric array.
711 for (int dataPoint = 0; dataPoint < numDataPointsPerSample; dataPoint++) {
712 DataArrayView dataPointView = getDataPoint(sampleNo, dataPoint);
713 if (dataPointRank==0) {
714 numArray[dataPoint]=dataPointView();
715 }
716 if (dataPointRank==1) {
717 for (int i=0; i<dataPointShape[0]; i++) {
718 numArray[dataPoint][i]=dataPointView(i);
719 }
720 }
721 if (dataPointRank==2) {
722 for (int i=0; i<dataPointShape[0]; i++) {
723 for (int j=0; j<dataPointShape[1]; j++) {
724 numArray[dataPoint][i][j] = dataPointView(i,j);
725 }
726 }
727 }
728 if (dataPointRank==3) {
729 for (int i=0; i<dataPointShape[0]; i++) {
730 for (int j=0; j<dataPointShape[1]; j++) {
731 for (int k=0; k<dataPointShape[2]; k++) {
732 numArray[dataPoint][i][j][k]=dataPointView(i,j,k);
733 }
734 }
735 }
736 }
737 if (dataPointRank==4) {
738 for (int i=0; i<dataPointShape[0]; i++) {
739 for (int j=0; j<dataPointShape[1]; j++) {
740 for (int k=0; k<dataPointShape[2]; k++) {
741 for (int l=0; l<dataPointShape[3]; l++) {
742 numArray[dataPoint][i][j][k][l]=dataPointView(i,j,k,l);
743 }
744 }
745 }
746 }
747 }
748 }
749
750 //
751 // return the loaded array
752 return numArray;
753 }
754
755 const
756 boost::python::numeric::array
757 Data::convertToNumArrayFromDPNo(int sampleNo,
758 int dataPointNo)
759 {
760 //
761 // Check a valid sample number has been supplied
762 if (sampleNo >= getNumSamples()) {
763 throw DataException("Error - Data::convertToNumArray: invalid sampleNo.");
764 }
765
766 //
767 // Check a valid data point number has been supplied
768 if (dataPointNo >= getNumDataPointsPerSample()) {
769 throw DataException("Error - Data::convertToNumArray: invalid dataPointNo.");
770 }
771
772 //
773 // determine the rank and shape of each data point
774 int dataPointRank = getDataPointRank();
775 DataArrayView::ShapeType dataPointShape = getDataPointShape();
776
777 //
778 // create the numeric array to be returned
779 boost::python::numeric::array numArray(0.0);
780
781 //
782 // the shape of the returned numeric array will be the same
783 // as that of the data point
784 int arrayRank = dataPointRank;
785 DataArrayView::ShapeType arrayShape = dataPointShape;
786
787 //
788 // resize the numeric array to the shape just calculated
789 if (arrayRank==0) {
790 numArray.resize(1);
791 }
792 if (arrayRank==1) {
793 numArray.resize(arrayShape[0]);
794 }
795 if (arrayRank==2) {
796 numArray.resize(arrayShape[0],arrayShape[1]);
797 }
798 if (arrayRank==3) {
799 numArray.resize(arrayShape[0],arrayShape[1],arrayShape[2]);
800 }
801 if (arrayRank==4) {
802 numArray.resize(arrayShape[0],arrayShape[1],arrayShape[2],arrayShape[3]);
803 }
804
805 //
806 // load the values for the data point into the numeric array.
807 DataArrayView dataPointView = getDataPoint(sampleNo, dataPointNo);
808 if (dataPointRank==0) {
809 numArray[0]=dataPointView();
810 }
811 if (dataPointRank==1) {
812 for (int i=0; i<dataPointShape[0]; i++) {
813 numArray[i]=dataPointView(i);
814 }
815 }
816 if (dataPointRank==2) {
817 for (int i=0; i<dataPointShape[0]; i++) {
818 for (int j=0; j<dataPointShape[1]; j++) {
819 numArray[i][j] = dataPointView(i,j);
820 }
821 }
822 }
823 if (dataPointRank==3) {
824 for (int i=0; i<dataPointShape[0]; i++) {
825 for (int j=0; j<dataPointShape[1]; j++) {
826 for (int k=0; k<dataPointShape[2]; k++) {
827 numArray[i][j][k]=dataPointView(i,j,k);
828 }
829 }
830 }
831 }
832 if (dataPointRank==4) {
833 for (int i=0; i<dataPointShape[0]; i++) {
834 for (int j=0; j<dataPointShape[1]; j++) {
835 for (int k=0; k<dataPointShape[2]; k++) {
836 for (int l=0; l<dataPointShape[3]; l++) {
837 numArray[i][j][k][l]=dataPointView(i,j,k,l);
838 }
839 }
840 }
841 }
842 }
843
844 //
845 // return the loaded array
846 return numArray;
847 }
848
849 boost::python::numeric::array
850 Data::integrate() const
851 {
852 int index;
853 int rank = getDataPointRank();
854 DataArrayView::ShapeType shape = getDataPointShape();
855
856 #if defined DOPROF
857 profData->integrate++;
858 #endif
859
860 //
861 // calculate the integral values
862 vector<double> integrals(getDataPointSize());
863 AbstractContinuousDomain::asAbstractContinuousDomain(getDomain()).setToIntegrals(integrals,*this);
864
865 //
866 // create the numeric array to be returned
867 // and load the array with the integral values
868 boost::python::numeric::array bp_array(1.0);
869 if (rank==0) {
870 bp_array.resize(1);
871 index = 0;
872 bp_array[0] = integrals[index];
873 }
874 if (rank==1) {
875 bp_array.resize(shape[0]);
876 for (int i=0; i<shape[0]; i++) {
877 index = i;
878 bp_array[i] = integrals[index];
879 }
880 }
881 if (rank==2) {
882 bp_array.resize(shape[0],shape[1]);
883 for (int i=0; i<shape[0]; i++) {
884 for (int j=0; j<shape[1]; j++) {
885 index = i + shape[0] * j;
886 bp_array[make_tuple(i,j)] = integrals[index];
887 }
888 }
889 }
890 if (rank==3) {
891 bp_array.resize(shape[0],shape[1],shape[2]);
892 for (int i=0; i<shape[0]; i++) {
893 for (int j=0; j<shape[1]; j++) {
894 for (int k=0; k<shape[2]; k++) {
895 index = i + shape[0] * ( j + shape[1] * k );
896 bp_array[make_tuple(i,j,k)] = integrals[index];
897 }
898 }
899 }
900 }
901 if (rank==4) {
902 bp_array.resize(shape[0],shape[1],shape[2],shape[3]);
903 for (int i=0; i<shape[0]; i++) {
904 for (int j=0; j<shape[1]; j++) {
905 for (int k=0; k<shape[2]; k++) {
906 for (int l=0; l<shape[3]; l++) {
907 index = i + shape[0] * ( j + shape[1] * ( k + shape[2] * l ) );
908 bp_array[make_tuple(i,j,k,l)] = integrals[index];
909 }
910 }
911 }
912 }
913 }
914
915 //
916 // return the loaded array
917 return bp_array;
918 }
919
920 Data
921 Data::sin() const
922 {
923 #if defined DOPROF
924 profData->unary++;
925 #endif
926 return escript::unaryOp(*this,(Data::UnaryDFunPtr)::sin);
927 }
928
929 Data
930 Data::cos() const
931 {
932 #if defined DOPROF
933 profData->unary++;
934 #endif
935 return escript::unaryOp(*this,(Data::UnaryDFunPtr)::cos);
936 }
937
938 Data
939 Data::tan() const
940 {
941 #if defined DOPROF
942 profData->unary++;
943 #endif
944 return escript::unaryOp(*this,(Data::UnaryDFunPtr)::tan);
945 }
946
947 Data
948 Data::asin() const
949 {
950 #if defined DOPROF
951 profData->unary++;
952 #endif
953 return escript::unaryOp(*this,(Data::UnaryDFunPtr)::asin);
954 }
955
956 Data
957 Data::acos() const
958 {
959 #if defined DOPROF
960 profData->unary++;
961 #endif
962 return escript::unaryOp(*this,(Data::UnaryDFunPtr)::acos);
963 }
964
965 Data
966 Data::atan() const
967 {
968 #if defined DOPROF
969 profData->unary++;
970 #endif
971 return escript::unaryOp(*this,(Data::UnaryDFunPtr)::atan);
972 }
973
974 Data
975 Data::sinh() const
976 {
977 #if defined DOPROF
978 profData->unary++;
979 #endif
980 return escript::unaryOp(*this,(Data::UnaryDFunPtr)::sinh);
981 }
982
983 Data
984 Data::cosh() const
985 {
986 #if defined DOPROF
987 profData->unary++;
988 #endif
989 return escript::unaryOp(*this,(Data::UnaryDFunPtr)::cosh);
990 }
991
992 Data
993 Data::tanh() const
994 {
995 #if defined DOPROF
996 profData->unary++;
997 #endif
998 return escript::unaryOp(*this,(Data::UnaryDFunPtr)::tanh);
999 }
1000
1001 Data
1002 Data::asinh() const
1003 {
1004 #if defined DOPROF
1005 profData->unary++;
1006 #endif
1007 return escript::unaryOp(*this,(Data::UnaryDFunPtr)::asinh);
1008 }
1009
1010 Data
1011 Data::acosh() const
1012 {
1013 #if defined DOPROF
1014 profData->unary++;
1015 #endif
1016 return escript::unaryOp(*this,(Data::UnaryDFunPtr)::acosh);
1017 }
1018
1019 Data
1020 Data::atanh() const
1021 {
1022 #if defined DOPROF
1023 profData->unary++;
1024 #endif
1025 return escript::unaryOp(*this,(Data::UnaryDFunPtr)::atanh);
1026 }
1027
1028 Data
1029 Data::log10() const
1030 {
1031 #if defined DOPROF
1032 profData->unary++;
1033 #endif
1034 return escript::unaryOp(*this,(Data::UnaryDFunPtr)::log10);
1035 }
1036
1037 Data
1038 Data::log() const
1039 {
1040 #if defined DOPROF
1041 profData->unary++;
1042 #endif
1043 return escript::unaryOp(*this,(Data::UnaryDFunPtr)::log);
1044 }
1045
1046 Data
1047 Data::sign() const
1048 {
1049 #if defined DOPROF
1050 profData->unary++;
1051 #endif
1052 return escript::unaryOp(*this,escript::fsign);
1053 }
1054
1055 Data
1056 Data::abs() const
1057 {
1058 #if defined DOPROF
1059 profData->unary++;
1060 #endif
1061 return escript::unaryOp(*this,(Data::UnaryDFunPtr)::fabs);
1062 }
1063
1064 Data
1065 Data::neg() const
1066 {
1067 #if defined DOPROF
1068 profData->unary++;
1069 #endif
1070 return escript::unaryOp(*this,negate<double>());
1071 }
1072
1073 Data
1074 Data::pos() const
1075 {
1076 #if defined DOPROF
1077 profData->unary++;
1078 #endif
1079 Data result;
1080 // perform a deep copy
1081 result.copy(*this);
1082 return result;
1083 }
1084
1085 Data
1086 Data::exp() const
1087 {
1088 #if defined DOPROF
1089 profData->unary++;
1090 #endif
1091 return escript::unaryOp(*this,(Data::UnaryDFunPtr)::exp);
1092 }
1093
1094 Data
1095 Data::sqrt() const
1096 {
1097 #if defined DOPROF
1098 profData->unary++;
1099 #endif
1100 return escript::unaryOp(*this,(Data::UnaryDFunPtr)::sqrt);
1101 }
1102
1103 double
1104 Data::Lsup() const
1105 {
1106 double localValue, globalValue;
1107 #if defined DOPROF
1108 profData->reduction1++;
1109 #endif
1110 //
1111 // set the initial absolute maximum value to zero
1112
1113 AbsMax abs_max_func;
1114 localValue = algorithm(abs_max_func,0);
1115 #ifdef PASO_MPI
1116 MPI_Allreduce( &localValue, &globalValue, 1, MPI_DOUBLE, MPI_MAX, MPI_COMM_WORLD );
1117 return globalValue;
1118 #else
1119 return localValue;
1120 #endif
1121 }
1122
1123 double
1124 Data::Linf() const
1125 {
1126 double localValue, globalValue;
1127 #if defined DOPROF
1128 profData->reduction1++;
1129 #endif
1130 //
1131 // set the initial absolute minimum value to max double
1132 AbsMin abs_min_func;
1133 localValue = algorithm(abs_min_func,numeric_limits<double>::max());
1134
1135 #ifdef PASO_MPI
1136 MPI_Allreduce( &localValue, &globalValue, 1, MPI_DOUBLE, MPI_MIN, MPI_COMM_WORLD );
1137 return globalValue;
1138 #else
1139 return localValue;
1140 #endif
1141 }
1142
1143 double
1144 Data::sup() const
1145 {
1146 double localValue, globalValue;
1147 #if defined DOPROF
1148 profData->reduction1++;
1149 #endif
1150 //
1151 // set the initial maximum value to min possible double
1152 FMax fmax_func;
1153 localValue = algorithm(fmax_func,numeric_limits<double>::max()*-1);
1154 #ifdef PASO_MPI
1155 MPI_Allreduce( &localValue, &globalValue, 1, MPI_DOUBLE, MPI_MAX, MPI_COMM_WORLD );
1156 return globalValue;
1157 #else
1158 return localValue;
1159 #endif
1160 }
1161
1162 double
1163 Data::inf() const
1164 {
1165 double localValue, globalValue;
1166 #if defined DOPROF
1167 profData->reduction1++;
1168 #endif
1169 //
1170 // set the initial minimum value to max possible double
1171 FMin fmin_func;
1172 localValue = algorithm(fmin_func,numeric_limits<double>::max());
1173 #ifdef PASO_MPI
1174 MPI_Allreduce( &localValue, &globalValue, 1, MPI_DOUBLE, MPI_MIN, MPI_COMM_WORLD );
1175 return globalValue;
1176 #else
1177 return localValue;
1178 #endif
1179 }
1180
1181 /* TODO */
1182 /* global reduction */
1183 Data
1184 Data::maxval() const
1185 {
1186 #if defined DOPROF
1187 profData->reduction2++;
1188 #endif
1189 //
1190 // set the initial maximum value to min possible double
1191 FMax fmax_func;
1192 return dp_algorithm(fmax_func,numeric_limits<double>::max()*-1);
1193 }
1194
1195 Data
1196 Data::minval() const
1197 {
1198 #if defined DOPROF
1199 profData->reduction2++;
1200 #endif
1201 //
1202 // set the initial minimum value to max possible double
1203 FMin fmin_func;
1204 return dp_algorithm(fmin_func,numeric_limits<double>::max());
1205 }
1206
1207 Data
1208 Data::trace() const
1209 {
1210 #if defined DOPROF
1211 profData->reduction2++;
1212 #endif
1213 Trace trace_func;
1214 return dp_algorithm(trace_func,0);
1215 }
1216
1217 Data
1218 Data::symmetric() const
1219 {
1220 #if defined DOPROF
1221 profData->unary++;
1222 #endif
1223 // check input
1224 DataArrayView::ShapeType s=getDataPointShape();
1225 if (getDataPointRank()==2) {
1226 if(s[0] != s[1])
1227 throw DataException("Error - Data::symmetric can only be calculated for rank 2 object with equal first and second dimension.");
1228 }
1229 else if (getDataPointRank()==4) {
1230 if(!(s[0] == s[2] && s[1] == s[3]))
1231 throw DataException("Error - Data::symmetric can only be calculated for rank 4 object with dim0==dim2 and dim1==dim3.");
1232 }
1233 else {
1234 throw DataException("Error - Data::symmetric can only be calculated for rank 2 or 4 object.");
1235 }
1236 Data ev(0.,getDataPointShape(),getFunctionSpace());
1237 ev.typeMatchRight(*this);
1238 m_data->symmetric(ev.m_data.get());
1239 return ev;
1240 }
1241
1242 Data
1243 Data::nonsymmetric() const
1244 {
1245 #if defined DOPROF
1246 profData->unary++;
1247 #endif
1248 // check input
1249 DataArrayView::ShapeType s=getDataPointShape();
1250 if (getDataPointRank()==2) {
1251 if(s[0] != s[1])
1252 throw DataException("Error - Data::nonsymmetric can only be calculated for rank 2 object with equal first and second dimension.");
1253 DataArrayView::ShapeType ev_shape;
1254 ev_shape.push_back(s[0]);
1255 ev_shape.push_back(s[1]);
1256 Data ev(0.,ev_shape,getFunctionSpace());
1257 ev.typeMatchRight(*this);
1258 m_data->nonsymmetric(ev.m_data.get());
1259 return ev;
1260 }
1261 else if (getDataPointRank()==4) {
1262 if(!(s[0] == s[2] && s[1] == s[3]))
1263 throw DataException("Error - Data::nonsymmetric can only be calculated for rank 4 object with dim0==dim2 and dim1==dim3.");
1264 DataArrayView::ShapeType ev_shape;
1265 ev_shape.push_back(s[0]);
1266 ev_shape.push_back(s[1]);
1267 ev_shape.push_back(s[2]);
1268 ev_shape.push_back(s[3]);
1269 Data ev(0.,ev_shape,getFunctionSpace());
1270 ev.typeMatchRight(*this);
1271 m_data->nonsymmetric(ev.m_data.get());
1272 return ev;
1273 }
1274 else {
1275 throw DataException("Error - Data::nonsymmetric can only be calculated for rank 2 or 4 object.");
1276 }
1277 }
1278
1279 Data
1280 Data::matrixtrace(int axis_offset) const
1281 {
1282 #if defined DOPROF
1283 profData->unary++;
1284 #endif
1285 DataArrayView::ShapeType s=getDataPointShape();
1286 if (getDataPointRank()==2) {
1287 DataArrayView::ShapeType ev_shape;
1288 Data ev(0.,ev_shape,getFunctionSpace());
1289 ev.typeMatchRight(*this);
1290 m_data->matrixtrace(ev.m_data.get(), axis_offset);
1291 return ev;
1292 }
1293 if (getDataPointRank()==3) {
1294 DataArrayView::ShapeType ev_shape;
1295 if (axis_offset==0) {
1296 int s2=s[2];
1297 ev_shape.push_back(s2);
1298 }
1299 else if (axis_offset==1) {
1300 int s0=s[0];
1301 ev_shape.push_back(s0);
1302 }
1303 Data ev(0.,ev_shape,getFunctionSpace());
1304 ev.typeMatchRight(*this);
1305 m_data->matrixtrace(ev.m_data.get(), axis_offset);
1306 return ev;
1307 }
1308 if (getDataPointRank()==4) {
1309 DataArrayView::ShapeType ev_shape;
1310 if (axis_offset==0) {
1311 ev_shape.push_back(s[2]);
1312 ev_shape.push_back(s[3]);
1313 }
1314 else if (axis_offset==1) {
1315 ev_shape.push_back(s[0]);
1316 ev_shape.push_back(s[3]);
1317 }
1318 else if (axis_offset==2) {
1319 ev_shape.push_back(s[0]);
1320 ev_shape.push_back(s[1]);
1321 }
1322 Data ev(0.,ev_shape,getFunctionSpace());
1323 ev.typeMatchRight(*this);
1324 m_data->matrixtrace(ev.m_data.get(), axis_offset);
1325 return ev;
1326 }
1327 else {
1328 throw DataException("Error - Data::matrixtrace can only be calculated for rank 2, 3 or 4 object.");
1329 }
1330 }
1331
1332 Data
1333 Data::transpose(int axis_offset) const
1334 {
1335 #if defined DOPROF
1336 profData->reduction2++;
1337 #endif
1338 DataArrayView::ShapeType s=getDataPointShape();
1339 DataArrayView::ShapeType ev_shape;
1340 // Here's the equivalent of python s_out=s[axis_offset:]+s[:axis_offset]
1341 // which goes thru all shape vector elements starting with axis_offset (at index=rank wrap around to 0)
1342 int rank=getDataPointRank();
1343 if (axis_offset<0 || axis_offset>rank) {
1344 throw DataException("Error - Data::transpose must have 0 <= axis_offset <= rank=" + rank);
1345 }
1346 for (int i=0; i<rank; i++) {
1347 int index = (axis_offset+i)%rank;
1348 ev_shape.push_back(s[index]); // Append to new shape
1349 }
1350 Data ev(0.,ev_shape,getFunctionSpace());
1351 ev.typeMatchRight(*this);
1352 m_data->transpose(ev.m_data.get(), axis_offset);
1353 return ev;
1354 }
1355
1356 Data
1357 Data::eigenvalues() const
1358 {
1359 #if defined DOPROF
1360 profData->unary++;
1361 #endif
1362 // check input
1363 DataArrayView::ShapeType s=getDataPointShape();
1364 if (getDataPointRank()!=2)
1365 throw DataException("Error - Data::eigenvalues can only be calculated for rank 2 object.");
1366 if(s[0] != s[1])
1367 throw DataException("Error - Data::eigenvalues can only be calculated for object with equal first and second dimension.");
1368 // create return
1369 DataArrayView::ShapeType ev_shape(1,s[0]);
1370 Data ev(0.,ev_shape,getFunctionSpace());
1371 ev.typeMatchRight(*this);
1372 m_data->eigenvalues(ev.m_data.get());
1373 return ev;
1374 }
1375
1376 const boost::python::tuple
1377 Data::eigenvalues_and_eigenvectors(const double tol) const
1378 {
1379 #if defined DOPROF
1380 profData->unary++;
1381 #endif
1382 DataArrayView::ShapeType s=getDataPointShape();
1383 if (getDataPointRank()!=2)
1384 throw DataException("Error - Data::eigenvalues and eigenvectors can only be calculated for rank 2 object.");
1385 if(s[0] != s[1])
1386 throw DataException("Error - Data::eigenvalues and eigenvectors can only be calculated for object with equal first and second dimension.");
1387 // create return
1388 DataArrayView::ShapeType ev_shape(1,s[0]);
1389 Data ev(0.,ev_shape,getFunctionSpace());
1390 ev.typeMatchRight(*this);
1391 DataArrayView::ShapeType V_shape(2,s[0]);
1392 Data V(0.,V_shape,getFunctionSpace());
1393 V.typeMatchRight(*this);
1394 m_data->eigenvalues_and_eigenvectors(ev.m_data.get(),V.m_data.get(),tol);
1395 return make_tuple(boost::python::object(ev),boost::python::object(V));
1396 }
1397
1398 const boost::python::tuple
1399 Data::mindp() const
1400 {
1401 // NB: calc_mindp had to be split off from mindp as boost::make_tuple causes an
1402 // abort (for unknown reasons) if there are openmp directives with it in the
1403 // surrounding function
1404
1405 int SampleNo;
1406 int DataPointNo;
1407
1408 calc_mindp(SampleNo,DataPointNo);
1409
1410 return make_tuple(SampleNo,DataPointNo);
1411 }
1412
1413 void
1414 Data::calc_mindp(int& SampleNo,
1415 int& DataPointNo) const
1416 {
1417 int i,j;
1418 int lowi=0,lowj=0;
1419 double min=numeric_limits<double>::max();
1420
1421 Data temp=minval();
1422
1423 int numSamples=temp.getNumSamples();
1424 int numDPPSample=temp.getNumDataPointsPerSample();
1425
1426 double next,local_min;
1427 int local_lowi,local_lowj;
1428
1429 #pragma omp parallel private(next,local_min,local_lowi,local_lowj)
1430 {
1431 local_min=min;
1432 #pragma omp for private(i,j) schedule(static)
1433 for (i=0; i<numSamples; i++) {
1434 for (j=0; j<numDPPSample; j++) {
1435 next=temp.getDataPoint(i,j)();
1436 if (next<local_min) {
1437 local_min=next;
1438 local_lowi=i;
1439 local_lowj=j;
1440 }
1441 }
1442 }
1443 #pragma omp critical
1444 if (local_min<min) {
1445 min=local_min;
1446 lowi=local_lowi;
1447 lowj=local_lowj;
1448 }
1449 }
1450
1451 SampleNo = lowi;
1452 DataPointNo = lowj;
1453 }
1454
1455 void
1456 Data::saveDX(std::string fileName) const
1457 {
1458 boost::python::dict args;
1459 args["data"]=boost::python::object(this);
1460 getDomain().saveDX(fileName,args);
1461 return;
1462 }
1463
1464 void
1465 Data::saveVTK(std::string fileName) const
1466 {
1467 boost::python::dict args;
1468 args["data"]=boost::python::object(this);
1469 getDomain().saveVTK(fileName,args);
1470 return;
1471 }
1472
1473 Data&
1474 Data::operator+=(const Data& right)
1475 {
1476 #if defined DOPROF
1477 profData->binary++;
1478 #endif
1479 binaryOp(right,plus<double>());
1480 return (*this);
1481 }
1482
1483 Data&
1484 Data::operator+=(const boost::python::object& right)
1485 {
1486 #if defined DOPROF
1487 profData->binary++;
1488 #endif
1489 binaryOp(right,plus<double>());
1490 return (*this);
1491 }
1492
1493 Data&
1494 Data::operator-=(const Data& right)
1495 {
1496 #if defined DOPROF
1497 profData->binary++;
1498 #endif
1499 binaryOp(right,minus<double>());
1500 return (*this);
1501 }
1502
1503 Data&
1504 Data::operator-=(const boost::python::object& right)
1505 {
1506 #if defined DOPROF
1507 profData->binary++;
1508 #endif
1509 binaryOp(right,minus<double>());
1510 return (*this);
1511 }
1512
1513 Data&
1514 Data::operator*=(const Data& right)
1515 {
1516 #if defined DOPROF
1517 profData->binary++;
1518 #endif
1519 binaryOp(right,multiplies<double>());
1520 return (*this);
1521 }
1522
1523 Data&
1524 Data::operator*=(const boost::python::object& right)
1525 {
1526 #if defined DOPROF
1527 profData->binary++;
1528 #endif
1529 binaryOp(right,multiplies<double>());
1530 return (*this);
1531 }
1532
1533 Data&
1534 Data::operator/=(const Data& right)
1535 {
1536 #if defined DOPROF
1537 profData->binary++;
1538 #endif
1539 binaryOp(right,divides<double>());
1540 return (*this);
1541 }
1542
1543 Data&
1544 Data::operator/=(const boost::python::object& right)
1545 {
1546 #if defined DOPROF
1547 profData->binary++;
1548 #endif
1549 binaryOp(right,divides<double>());
1550 return (*this);
1551 }
1552
1553 Data
1554 Data::rpowO(const boost::python::object& left) const
1555 {
1556 #if defined DOPROF
1557 profData->binary++;
1558 #endif
1559 Data left_d(left,*this);
1560 return left_d.powD(*this);
1561 }
1562
1563 Data
1564 Data::powO(const boost::python::object& right) const
1565 {
1566 #if defined DOPROF
1567 profData->binary++;
1568 #endif
1569 Data result;
1570 result.copy(*this);
1571 result.binaryOp(right,(Data::BinaryDFunPtr)::pow);
1572 return result;
1573 }
1574
1575 Data
1576 Data::powD(const Data& right) const
1577 {
1578 #if defined DOPROF
1579 profData->binary++;
1580 #endif
1581 Data result;
1582 result.copy(*this);
1583 result.binaryOp(right,(Data::BinaryDFunPtr)::pow);
1584 return result;
1585 }
1586
1587 void
1588 Data::print()
1589 {
1590 int i,j;
1591
1592 printf( "Data is %dX%d\n", getNumSamples(), getNumDataPointsPerSample() );
1593 for( i=0; i<getNumSamples(); i++ )
1594 {
1595 printf( "[%6d]", i );
1596 for( j=0; j<getNumDataPointsPerSample(); j++ )
1597 printf( "\t%10.7g", (getSampleData(i))[j] );
1598 printf( "\n" );
1599 }
1600 }
1601
1602 //
1603 // NOTE: It is essential to specify the namespace this operator belongs to
1604 Data
1605 escript::operator+(const Data& left, const Data& right)
1606 {
1607 Data result;
1608 //
1609 // perform a deep copy
1610 result.copy(left);
1611 result+=right;
1612 return result;
1613 }
1614
1615 //
1616 // NOTE: It is essential to specify the namespace this operator belongs to
1617 Data
1618 escript::operator-(const Data& left, const Data& right)
1619 {
1620 Data result;
1621 //
1622 // perform a deep copy
1623 result.copy(left);
1624 result-=right;
1625 return result;
1626 }
1627
1628 //
1629 // NOTE: It is essential to specify the namespace this operator belongs to
1630 Data
1631 escript::operator*(const Data& left, const Data& right)
1632 {
1633 Data result;
1634 //
1635 // perform a deep copy
1636 result.copy(left);
1637 result*=right;
1638 return result;
1639 }
1640
1641 //
1642 // NOTE: It is essential to specify the namespace this operator belongs to
1643 Data
1644 escript::operator/(const Data& left, const Data& right)
1645 {
1646 Data result;
1647 //
1648 // perform a deep copy
1649 result.copy(left);
1650 result/=right;
1651 return result;
1652 }
1653
1654 //
1655 // NOTE: It is essential to specify the namespace this operator belongs to
1656 Data
1657 escript::operator+(const Data& left, const boost::python::object& right)
1658 {
1659 //
1660 // Convert to DataArray format if possible
1661 DataArray temp(right);
1662 Data result;
1663 //
1664 // perform a deep copy
1665 result.copy(left);
1666 result+=right;
1667 return result;
1668 }
1669
1670 //
1671 // NOTE: It is essential to specify the namespace this operator belongs to
1672 Data
1673 escript::operator-(const Data& left, const boost::python::object& right)
1674 {
1675 //
1676 // Convert to DataArray format if possible
1677 DataArray temp(right);
1678 Data result;
1679 //
1680 // perform a deep copy
1681 result.copy(left);
1682 result-=right;
1683 return result;
1684 }
1685
1686 //
1687 // NOTE: It is essential to specify the namespace this operator belongs to
1688 Data
1689 escript::operator*(const Data& left, const boost::python::object& right)
1690 {
1691 //
1692 // Convert to DataArray format if possible
1693 DataArray temp(right);
1694 Data result;
1695 //
1696 // perform a deep copy
1697 result.copy(left);
1698 result*=right;
1699 return result;
1700 }
1701
1702 //
1703 // NOTE: It is essential to specify the namespace this operator belongs to
1704 Data
1705 escript::operator/(const Data& left, const boost::python::object& right)
1706 {
1707 //
1708 // Convert to DataArray format if possible
1709 DataArray temp(right);
1710 Data result;
1711 //
1712 // perform a deep copy
1713 result.copy(left);
1714 result/=right;
1715 return result;
1716 }
1717
1718 //
1719 // NOTE: It is essential to specify the namespace this operator belongs to
1720 Data
1721 escript::operator+(const boost::python::object& left, const Data& right)
1722 {
1723 //
1724 // Construct the result using the given value and the other parameters
1725 // from right
1726 Data result(left,right);
1727 result+=right;
1728 return result;
1729 }
1730
1731 //
1732 // NOTE: It is essential to specify the namespace this operator belongs to
1733 Data
1734 escript::operator-(const boost::python::object& left, const Data& right)
1735 {
1736 //
1737 // Construct the result using the given value and the other parameters
1738 // from right
1739 Data result(left,right);
1740 result-=right;
1741 return result;
1742 }
1743
1744 //
1745 // NOTE: It is essential to specify the namespace this operator belongs to
1746 Data
1747 escript::operator*(const boost::python::object& left, const Data& right)
1748 {
1749 //
1750 // Construct the result using the given value and the other parameters
1751 // from right
1752 Data result(left,right);
1753 result*=right;
1754 return result;
1755 }
1756
1757 //
1758 // NOTE: It is essential to specify the namespace this operator belongs to
1759 Data
1760 escript::operator/(const boost::python::object& left, const Data& right)
1761 {
1762 //
1763 // Construct the result using the given value and the other parameters
1764 // from right
1765 Data result(left,right);
1766 result/=right;
1767 return result;
1768 }
1769
1770 //
1771 //bool escript::operator==(const Data& left, const Data& right)
1772 //{
1773 // /*
1774 // NB: this operator does very little at this point, and isn't to
1775 // be relied on. Requires further implementation.
1776 // */
1777 //
1778 // bool ret;
1779 //
1780 // if (left.isEmpty()) {
1781 // if(!right.isEmpty()) {
1782 // ret = false;
1783 // } else {
1784 // ret = true;
1785 // }
1786 // }
1787 //
1788 // if (left.isConstant()) {
1789 // if(!right.isConstant()) {
1790 // ret = false;
1791 // } else {
1792 // ret = true;
1793 // }
1794 // }
1795 //
1796 // if (left.isTagged()) {
1797 // if(!right.isTagged()) {
1798 // ret = false;
1799 // } else {
1800 // ret = true;
1801 // }
1802 // }
1803 //
1804 // if (left.isExpanded()) {
1805 // if(!right.isExpanded()) {
1806 // ret = false;
1807 // } else {
1808 // ret = true;
1809 // }
1810 // }
1811 //
1812 // return ret;
1813 //}
1814
1815 /* TODO */
1816 /* global reduction */
1817 Data
1818 Data::getItem(const boost::python::object& key) const
1819 {
1820 const DataArrayView& view=getPointDataView();
1821
1822 DataArrayView::RegionType slice_region=view.getSliceRegion(key);
1823
1824 if (slice_region.size()!=view.getRank()) {
1825 throw DataException("Error - slice size does not match Data rank.");
1826 }
1827
1828 return getSlice(slice_region);
1829 }
1830
1831 /* TODO */
1832 /* global reduction */
1833 Data
1834 Data::getSlice(const DataArrayView::RegionType& region) const
1835 {
1836 #if defined DOPROF
1837 profData->slicing++;
1838 #endif
1839 return Data(*this,region);
1840 }
1841
1842 /* TODO */
1843 /* global reduction */
1844 void
1845 Data::setItemO(const boost::python::object& key,
1846 const boost::python::object& value)
1847 {
1848 Data tempData(value,getFunctionSpace());
1849 setItemD(key,tempData);
1850 }
1851
1852 /* TODO */
1853 /* global reduction */
1854 void
1855 Data::setItemD(const boost::python::object& key,
1856 const Data& value)
1857 {
1858 const DataArrayView& view=getPointDataView();
1859
1860 DataArrayView::RegionType slice_region=view.getSliceRegion(key);
1861 if (slice_region.size()!=view.getRank()) {
1862 throw DataException("Error - slice size does not match Data rank.");
1863 }
1864 if (getFunctionSpace()!=value.getFunctionSpace()) {
1865 setSlice(Data(value,getFunctionSpace()),slice_region);
1866 } else {
1867 setSlice(value,slice_region);
1868 }
1869 }
1870
1871 /* TODO */
1872 /* global reduction */
1873 void
1874 Data::setSlice(const Data& value,
1875 const DataArrayView::RegionType& region)
1876 {
1877 #if defined DOPROF
1878 profData->slicing++;
1879 #endif
1880 Data tempValue(value);
1881 typeMatchLeft(tempValue);
1882 typeMatchRight(tempValue);
1883 m_data->setSlice(tempValue.m_data.get(),region);
1884 }
1885
1886 void
1887 Data::typeMatchLeft(Data& right) const
1888 {
1889 if (isExpanded()){
1890 right.expand();
1891 } else if (isTagged()) {
1892 if (right.isConstant()) {
1893 right.tag();
1894 }
1895 }
1896 }
1897
1898 void
1899 Data::typeMatchRight(const Data& right)
1900 {
1901 if (isTagged()) {
1902 if (right.isExpanded()) {
1903 expand();
1904 }
1905 } else if (isConstant()) {
1906 if (right.isExpanded()) {
1907 expand();
1908 } else if (right.isTagged()) {
1909 tag();
1910 }
1911 }
1912 }
1913
1914 /* TODO */
1915 /* global reduction */
1916 void
1917 Data::setTaggedValue(int tagKey,
1918 const boost::python::object& value)
1919 {
1920 //
1921 // Ensure underlying data object is of type DataTagged
1922 tag();
1923
1924 if (!isTagged()) {
1925 throw DataException("Error - DataTagged conversion failed!!");
1926 }
1927
1928 //
1929 // Construct DataArray from boost::python::object input value
1930 DataArray valueDataArray(value);
1931
1932 //
1933 // Call DataAbstract::setTaggedValue
1934 m_data->setTaggedValue(tagKey,valueDataArray.getView());
1935 }
1936
1937 /* TODO */
1938 /* global reduction */
1939 void
1940 Data::setTaggedValueFromCPP(int tagKey,
1941 const DataArrayView& value)
1942 {
1943 //
1944 // Ensure underlying data object is of type DataTagged
1945 tag();
1946
1947 if (!isTagged()) {
1948 throw DataException("Error - DataTagged conversion failed!!");
1949 }
1950
1951 //
1952 // Call DataAbstract::setTaggedValue
1953 m_data->setTaggedValue(tagKey,value);
1954 }
1955
1956 /* TODO */
1957 /* global reduction */
1958 int
1959 Data::getTagNumber(int dpno)
1960 {
1961 return m_data->getTagNumber(dpno);
1962 }
1963
1964 /* TODO */
1965 /* global reduction */
1966 void
1967 Data::setRefValue(int ref,
1968 const boost::python::numeric::array& value)
1969 {
1970 //
1971 // Construct DataArray from boost::python::object input value
1972 DataArray valueDataArray(value);
1973
1974 //
1975 // Call DataAbstract::setRefValue
1976 m_data->setRefValue(ref,valueDataArray);
1977 }
1978
1979 /* TODO */
1980 /* global reduction */
1981 void
1982 Data::getRefValue(int ref,
1983 boost::python::numeric::array& value)
1984 {
1985 //
1986 // Construct DataArray for boost::python::object return value
1987 DataArray valueDataArray(value);
1988
1989 //
1990 // Load DataArray with values from data-points specified by ref
1991 m_data->getRefValue(ref,valueDataArray);
1992
1993 //
1994 // Load values from valueDataArray into return numarray
1995
1996 // extract the shape of the numarray
1997 int rank = value.getrank();
1998 DataArrayView::ShapeType shape;
1999 for (int i=0; i < rank; i++) {
2000 shape.push_back(extract<int>(value.getshape()[i]));
2001 }
2002
2003 // and load the numarray with the data from the DataArray
2004 DataArrayView valueView = valueDataArray.getView();
2005
2006 if (rank==0) {
2007 boost::python::numeric::array temp_numArray(valueView());
2008 value = temp_numArray;
2009 }
2010 if (rank==1) {
2011 for (int i=0; i < shape[0]; i++) {
2012 value[i] = valueView(i);
2013 }
2014 }
2015 if (rank==2) {
2016 for (int i=0; i < shape[0]; i++) {
2017 for (int j=0; j < shape[1]; j++) {
2018 value[i][j] = valueView(i,j);
2019 }
2020 }
2021 }
2022 if (rank==3) {
2023 for (int i=0; i < shape[0]; i++) {
2024 for (int j=0; j < shape[1]; j++) {
2025 for (int k=0; k < shape[2]; k++) {
2026 value[i][j][k] = valueView(i,j,k);
2027 }
2028 }
2029 }
2030 }
2031 if (rank==4) {
2032 for (int i=0; i < shape[0]; i++) {
2033 for (int j=0; j < shape[1]; j++) {
2034 for (int k=0; k < shape[2]; k++) {
2035 for (int l=0; l < shape[3]; l++) {
2036 value[i][j][k][l] = valueView(i,j,k,l);
2037 }
2038 }
2039 }
2040 }
2041 }
2042
2043 }
2044
2045 void
2046 Data::archiveData(const std::string fileName)
2047 {
2048 cout << "Archiving Data object to: " << fileName << endl;
2049
2050 //
2051 // Determine type of this Data object
2052 int dataType = -1;
2053
2054 if (isEmpty()) {
2055 dataType = 0;
2056 cout << "\tdataType: DataEmpty" << endl;
2057 }
2058 if (isConstant()) {
2059 dataType = 1;
2060 cout << "\tdataType: DataConstant" << endl;
2061 }
2062 if (isTagged()) {
2063 dataType = 2;
2064 cout << "\tdataType: DataTagged" << endl;
2065 }
2066 if (isExpanded()) {
2067 dataType = 3;
2068 cout << "\tdataType: DataExpanded" << endl;
2069 }
2070
2071 if (dataType == -1) {
2072 throw DataException("archiveData Error: undefined dataType");
2073 }
2074
2075 //
2076 // Collect data items common to all Data types
2077 int noSamples = getNumSamples();
2078 int noDPPSample = getNumDataPointsPerSample();
2079 int functionSpaceType = getFunctionSpace().getTypeCode();
2080 int dataPointRank = getDataPointRank();
2081 int dataPointSize = getDataPointSize();
2082 int dataLength = getLength();
2083 DataArrayView::ShapeType dataPointShape = getDataPointShape();
2084 vector<int> referenceNumbers(noSamples);
2085 for (int sampleNo=0; sampleNo<noSamples; sampleNo++) {
2086 referenceNumbers[sampleNo] = getFunctionSpace().getReferenceNoFromSampleNo(sampleNo);
2087 }
2088 vector<int> tagNumbers(noSamples);
2089 if (isTagged()) {
2090 for (int sampleNo=0; sampleNo<noSamples; sampleNo++) {
2091 tagNumbers[sampleNo] = getFunctionSpace().getTagFromSampleNo(sampleNo);
2092 }
2093 }
2094
2095 cout << "\tnoSamples: " << noSamples << " noDPPSample: " << noDPPSample << endl;
2096 cout << "\tfunctionSpaceType: " << functionSpaceType << endl;
2097 cout << "\trank: " << dataPointRank << " size: " << dataPointSize << " length: " << dataLength << endl;
2098
2099 //
2100 // Flatten Shape to an array of integers suitable for writing to file
2101 int flatShape[4] = {0,0,0,0};
2102 cout << "\tshape: < ";
2103 for (int dim=0; dim<dataPointRank; dim++) {
2104 flatShape[dim] = dataPointShape[dim];
2105 cout << dataPointShape[dim] << " ";
2106 }
2107 cout << ">" << endl;
2108
2109 //
2110 // Open archive file
2111 ofstream archiveFile;
2112 archiveFile.open(fileName.data(), ios::out);
2113
2114 if (!archiveFile.good()) {
2115 throw DataException("archiveData Error: problem opening archive file");
2116 }
2117
2118 //
2119 // Write common data items to archive file
2120 archiveFile.write(reinterpret_cast<char *>(&dataType),sizeof(int));
2121 archiveFile.write(reinterpret_cast<char *>(&noSamples),sizeof(int));
2122 archiveFile.write(reinterpret_cast<char *>(&noDPPSample),sizeof(int));
2123 archiveFile.write(reinterpret_cast<char *>(&functionSpaceType),sizeof(int));
2124 archiveFile.write(reinterpret_cast<char *>(&dataPointRank),sizeof(int));
2125 archiveFile.write(reinterpret_cast<char *>(&dataPointSize),sizeof(int));
2126 archiveFile.write(reinterpret_cast<char *>(&dataLength),sizeof(int));
2127 for (int dim = 0; dim < 4; dim++) {
2128 archiveFile.write(reinterpret_cast<char *>(&flatShape[dim]),sizeof(int));
2129 }
2130 for (int sampleNo=0; sampleNo<noSamples; sampleNo++) {
2131 archiveFile.write(reinterpret_cast<char *>(&referenceNumbers[sampleNo]),sizeof(int));
2132 }
2133 if (isTagged()) {
2134 for (int sampleNo=0; sampleNo<noSamples; sampleNo++) {
2135 archiveFile.write(reinterpret_cast<char *>(&tagNumbers[sampleNo]),sizeof(int));
2136 }
2137 }
2138
2139 if (!archiveFile.good()) {
2140 throw DataException("archiveData Error: problem writing to archive file");
2141 }
2142
2143 //
2144 // Archive underlying data values for each Data type
2145 int noValues;
2146 switch (dataType) {
2147 case 0:
2148 // DataEmpty
2149 noValues = 0;
2150 archiveFile.write(reinterpret_cast<char *>(&noValues),sizeof(int));
2151 cout << "\tnoValues: " << noValues << endl;
2152 break;
2153 case 1:
2154 // DataConstant
2155 noValues = m_data->getLength();
2156 archiveFile.write(reinterpret_cast<char *>(&noValues),sizeof(int));
2157 cout << "\tnoValues: " << noValues << endl;
2158 if (m_data->archiveData(archiveFile,noValues)) {
2159 throw DataException("archiveData Error: problem writing data to archive file");
2160 }
2161 break;
2162 case 2:
2163 // DataTagged
2164 noValues = m_data->getLength();
2165 archiveFile.write(reinterpret_cast<char *>(&noValues),sizeof(int));
2166 cout << "\tnoValues: " << noValues << endl;
2167 if (m_data->archiveData(archiveFile,noValues)) {
2168 throw DataException("archiveData Error: problem writing data to archive file");
2169 }
2170 break;
2171 case 3:
2172 // DataExpanded
2173 noValues = m_data->getLength();
2174 archiveFile.write(reinterpret_cast<char *>(&noValues),sizeof(int));
2175 cout << "\tnoValues: " << noValues << endl;
2176 if (m_data->archiveData(archiveFile,noValues)) {
2177 throw DataException("archiveData Error: problem writing data to archive file");
2178 }
2179 break;
2180 }
2181
2182 if (!archiveFile.good()) {
2183 throw DataException("archiveData Error: problem writing data to archive file");
2184 }
2185
2186 //
2187 // Close archive file
2188 archiveFile.close();
2189
2190 if (!archiveFile.good()) {
2191 throw DataException("archiveData Error: problem closing archive file");
2192 }
2193
2194 }
2195
2196 void
2197 Data::extractData(const std::string fileName,
2198 const FunctionSpace& fspace)
2199 {
2200 //
2201 // Can only extract Data to an object which is initially DataEmpty
2202 if (!isEmpty()) {
2203 throw DataException("extractData Error: can only extract to DataEmpty object");
2204 }
2205
2206 cout << "Extracting Data object from: " << fileName << endl;
2207
2208 int dataType;
2209 int noSamples;
2210 int noDPPSample;
2211 int functionSpaceType;
2212 int dataPointRank;
2213 int dataPointSize;
2214 int dataLength;
2215 DataArrayView::ShapeType dataPointShape;
2216 int flatShape[4];
2217
2218 //
2219 // Open the archive file
2220 ifstream archiveFile;
2221 archiveFile.open(fileName.data(), ios::in);
2222
2223 if (!archiveFile.good()) {
2224 throw DataException("extractData Error: problem opening archive file");
2225 }
2226
2227 //
2228 // Read common data items from archive file
2229 archiveFile.read(reinterpret_cast<char *>(&dataType),sizeof(int));
2230 archiveFile.read(reinterpret_cast<char *>(&noSamples),sizeof(int));
2231 archiveFile.read(reinterpret_cast<char *>(&noDPPSample),sizeof(int));
2232 archiveFile.read(reinterpret_cast<char *>(&functionSpaceType),sizeof(int));
2233 archiveFile.read(reinterpret_cast<char *>(&dataPointRank),sizeof(int));
2234 archiveFile.read(reinterpret_cast<char *>(&dataPointSize),sizeof(int));
2235 archiveFile.read(reinterpret_cast<char *>(&dataLength),sizeof(int));
2236 for (int dim = 0; dim < 4; dim++) {
2237 archiveFile.read(reinterpret_cast<char *>(&flatShape[dim]),sizeof(int));
2238 if (flatShape[dim]>0) {
2239 dataPointShape.push_back(flatShape[dim]);
2240 }
2241 }
2242 vector<int> referenceNumbers(noSamples);
2243 for (int sampleNo=0; sampleNo<noSamples; sampleNo++) {
2244 archiveFile.read(reinterpret_cast<char *>(&referenceNumbers[sampleNo]),sizeof(int));
2245 }
2246 vector<int> tagNumbers(noSamples);
2247 if (dataType==2) {
2248 for (int sampleNo=0; sampleNo<noSamples; sampleNo++) {
2249 archiveFile.read(reinterpret_cast<char *>(&tagNumbers[sampleNo]),sizeof(int));
2250 }
2251 }
2252
2253 if (!archiveFile.good()) {
2254 throw DataException("extractData Error: problem reading from archive file");
2255 }
2256
2257 //
2258 // Verify the values just read from the archive file
2259 switch (dataType) {
2260 case 0:
2261 cout << "\tdataType: DataEmpty" << endl;
2262 break;
2263 case 1:
2264 cout << "\tdataType: DataConstant" << endl;
2265 break;
2266 case 2:
2267 cout << "\tdataType: DataTagged" << endl;
2268 break;
2269 case 3:
2270 cout << "\tdataType: DataExpanded" << endl;
2271 break;
2272 default:
2273 throw DataException("extractData Error: undefined dataType read from archive file");
2274 break;
2275 }
2276
2277 cout << "\tnoSamples: " << noSamples << " noDPPSample: " << noDPPSample << endl;
2278 cout << "\tfunctionSpaceType: " << functionSpaceType << endl;
2279 cout << "\trank: " << dataPointRank << " size: " << dataPointSize << " length: " << dataLength << endl;
2280 cout << "\tshape: < ";
2281 for (int dim = 0; dim < dataPointRank; dim++) {
2282 cout << dataPointShape[dim] << " ";
2283 }
2284 cout << ">" << endl;
2285
2286 //
2287 // Verify that supplied FunctionSpace object is compatible with this Data object.
2288 if ( (fspace.getTypeCode()!=functionSpaceType) ||
2289 (fspace.getNumSamples()!=noSamples) ||
2290 (fspace.getNumDPPSample()!=noDPPSample)
2291 ) {
2292 throw DataException("extractData Error: incompatible FunctionSpace");
2293 }
2294 for (int sampleNo=0; sampleNo<noSamples; sampleNo++) {
2295 if (referenceNumbers[sampleNo] != fspace.getReferenceNoFromSampleNo(sampleNo)) {
2296 throw DataException("extractData Error: incompatible FunctionSpace");
2297 }
2298 }
2299 if (dataType==2) {
2300 for (int sampleNo=0; sampleNo<noSamples; sampleNo++) {
2301 if (tagNumbers[sampleNo] != fspace.getTagFromSampleNo(sampleNo)) {
2302 throw DataException("extractData Error: incompatible FunctionSpace");
2303 }
2304 }
2305 }
2306
2307 //
2308 // Construct a DataVector to hold underlying data values
2309 DataVector dataVec(dataLength);
2310
2311 //
2312 // Load this DataVector with the appropriate values
2313 int noValues;
2314 archiveFile.read(reinterpret_cast<char *>(&noValues),sizeof(int));
2315 cout << "\tnoValues: " << noValues << endl;
2316 switch (dataType) {
2317 case 0:
2318 // DataEmpty
2319 if (noValues != 0) {
2320 throw DataException("extractData Error: problem reading data from archive file");
2321 }
2322 break;
2323 case 1:
2324 // DataConstant
2325 if (dataVec.extractData(archiveFile,noValues)) {
2326 throw DataException("extractData Error: problem reading data from archive file");
2327 }
2328 break;
2329 case 2:
2330 // DataTagged
2331 if (dataVec.extractData(archiveFile,noValues)) {
2332 throw DataException("extractData Error: problem reading data from archive file");
2333 }
2334 break;
2335 case 3:
2336 // DataExpanded
2337 if (dataVec.extractData(archiveFile,noValues)) {
2338 throw DataException("extractData Error: problem reading data from archive file");
2339 }
2340 break;
2341 }
2342
2343 if (!archiveFile.good()) {
2344 throw DataException("extractData Error: problem reading from archive file");
2345 }
2346
2347 //
2348 // Close archive file
2349 archiveFile.close();
2350
2351 if (!archiveFile.good()) {
2352 throw DataException("extractData Error: problem closing archive file");
2353 }
2354
2355 //
2356 // Construct an appropriate Data object
2357 DataAbstract* tempData;
2358 switch (dataType) {
2359 case 0:
2360 // DataEmpty
2361 tempData=new DataEmpty();
2362 break;
2363 case 1:
2364 // DataConstant
2365 tempData=new DataConstant(fspace,dataPointShape,dataVec);
2366 break;
2367 case 2:
2368 // DataTagged
2369 tempData=new DataTagged(fspace,dataPointShape,tagNumbers,dataVec);
2370 break;
2371 case 3:
2372 // DataExpanded
2373 tempData=new DataExpanded(fspace,dataPointShape,dataVec);
2374 break;
2375 }
2376 shared_ptr<DataAbstract> temp_data(tempData);
2377 m_data=temp_data;
2378 }
2379
2380 ostream& escript::operator<<(ostream& o, const Data& data)
2381 {
2382 o << data.toString();
2383 return o;
2384 }

Properties

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

  ViewVC Help
Powered by ViewVC 1.1.26