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

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

Parent Directory Parent Directory | Revision Log Revision Log


Revision 113 - (hide annotations)
Mon Feb 28 07:06:33 2005 UTC (14 years, 9 months ago) by jgs
Original Path: trunk/esys2/escript/src/Data/Data.cpp
File size: 25505 byte(s)
*** empty log message ***

1 jgs 94 // $Id$
2     /*=============================================================================
3    
4     ******************************************************************************
5     * *
6     * COPYRIGHT ACcESS 2004 - All Rights Reserved *
7     * *
8     * This software is the property of ACcESS. No part of this code *
9     * may be copied in any form or by any means without the expressed written *
10     * consent of ACcESS. Copying, use or modification of this software *
11     * by any unauthorised person is illegal unless that *
12     * person has a software license agreement with ACcESS. *
13     * *
14     ******************************************************************************
15    
16     ******************************************************************************/
17    
18     #include "escript/Data/Data.h"
19    
20     #include <iostream>
21     #include <algorithm>
22     #include <vector>
23     #include <exception>
24     #include <functional>
25     #include <math.h>
26    
27     #include <boost/python/str.hpp>
28     #include <boost/python/extract.hpp>
29     #include <boost/python/long.hpp>
30    
31     #include "escript/Data/DataException.h"
32    
33     #include "escript/Data/DataExpanded.h"
34     #include "escript/Data/DataConstant.h"
35     #include "escript/Data/DataTagged.h"
36     #include "escript/Data/DataEmpty.h"
37     #include "escript/Data/DataArray.h"
38     #include "escript/Data/DataAlgorithm.h"
39     #include "escript/Data/FunctionSpaceFactory.h"
40     #include "escript/Data/AbstractContinuousDomain.h"
41 jgs 102 #include "escript/Data/UnaryFuncs.h"
42 jgs 94
43     using namespace std;
44     using namespace boost::python;
45     using namespace boost;
46     using namespace escript;
47    
48     Data::Data()
49     {
50     //
51     // Default data is type DataEmpty
52     DataAbstract* temp=new DataEmpty();
53 jgs 102 shared_ptr<DataAbstract> temp_data(temp);
54     m_data=temp_data;
55 jgs 94 }
56    
57     Data::Data(double value,
58     const tuple& shape,
59     const FunctionSpace& what,
60     bool expanded)
61     {
62     DataArrayView::ShapeType dataPointShape;
63     for (int i = 0; i < shape.attr("__len__")(); ++i) {
64     dataPointShape.push_back(extract<const int>(shape[i]));
65     }
66     DataArray temp(dataPointShape,value);
67     initialise(temp.getView(),what,expanded);
68     }
69    
70     Data::Data(double value,
71     const DataArrayView::ShapeType& dataPointShape,
72     const FunctionSpace& what,
73     bool expanded)
74     {
75     DataArray temp(dataPointShape,value);
76     pair<int,int> dataShape=what.getDataShape();
77     initialise(temp.getView(),what,expanded);
78     }
79    
80 jgs 102 Data::Data(const Data& inData)
81 jgs 94 {
82 jgs 102 m_data=inData.m_data;
83 jgs 94 }
84    
85     Data::Data(const Data& inData,
86     const DataArrayView::RegionType& region)
87     {
88     //
89 jgs 102 // Create Data which is a slice of another Data
90     DataAbstract* tmp = inData.m_data->getSlice(region);
91     shared_ptr<DataAbstract> temp_data(tmp);
92     m_data=temp_data;
93 jgs 94 }
94    
95     Data::Data(const Data& inData,
96     const FunctionSpace& functionspace)
97     {
98     if (inData.getFunctionSpace()==functionspace) {
99     m_data=inData.m_data;
100     } else {
101     Data tmp(0,inData.getPointDataView().getShape(),functionspace,true);
102     // Note for Lutz, Must use a reference or pointer to a derived object
103     // in order to get polymorphic behaviour. Shouldn't really
104     // be able to create an instance of AbstractDomain but that was done
105     // as a boost python work around which may no longer be required.
106     const AbstractDomain& inDataDomain=inData.getDomain();
107     if (inDataDomain==functionspace.getDomain()) {
108     inDataDomain.interpolateOnDomain(tmp,inData);
109     } else {
110     inDataDomain.interpolateACross(tmp,inData);
111     }
112     m_data=tmp.m_data;
113     }
114     }
115    
116     Data::Data(const DataTagged::TagListType& tagKeys,
117     const DataTagged::ValueListType & values,
118     const DataArrayView& defaultValue,
119     const FunctionSpace& what,
120     bool expanded)
121     {
122     DataAbstract* temp=new DataTagged(tagKeys,values,defaultValue,what);
123 jgs 102 shared_ptr<DataAbstract> temp_data(temp);
124     m_data=temp_data;
125 jgs 94 if (expanded) {
126     expand();
127     }
128     }
129    
130     Data::Data(const numeric::array& value,
131     const FunctionSpace& what,
132     bool expanded)
133     {
134     initialise(value,what,expanded);
135     }
136    
137     Data::Data(const DataArrayView& value,
138     const FunctionSpace& what,
139     bool expanded)
140     {
141     initialise(value,what,expanded);
142     }
143    
144     Data::Data(const object& value,
145     const FunctionSpace& what,
146     bool expanded)
147     {
148     numeric::array asNumArray(value);
149     initialise(asNumArray,what,expanded);
150     }
151    
152     Data::Data(const object& value,
153     const Data& other)
154     {
155     //
156     // Create DataConstant using the given value and all other parameters
157     // copied from other. If value is a rank 0 object this Data
158     // will assume the point data shape of other.
159     DataArray temp(value);
160     if (temp.getView().getRank()==0) {
161     //
162     // Create a DataArray with the scalar value for all elements
163     DataArray temp2(other.getPointDataView().getShape(),temp.getView()());
164     initialise(temp2.getView(),other.getFunctionSpace(),false);
165     } else {
166     //
167     // Create a DataConstant with the same sample shape as other
168     initialise(temp.getView(),other.getFunctionSpace(),false);
169     }
170     }
171    
172     escriptDataC
173     Data::getDataC()
174     {
175     escriptDataC temp;
176     temp.m_dataPtr=(void*)this;
177     return temp;
178     }
179    
180     escriptDataC
181     Data::getDataC() const
182     {
183     escriptDataC temp;
184     temp.m_dataPtr=(void*)this;
185     return temp;
186     }
187    
188     tuple
189     Data::getShapeTuple() const
190     {
191     const DataArrayView::ShapeType& shape=getDataPointShape();
192     switch(getDataPointRank()) {
193     case 0:
194     return make_tuple();
195     case 1:
196     return make_tuple(long_(shape[0]));
197     case 2:
198     return make_tuple(long_(shape[0]),long_(shape[1]));
199     case 3:
200     return make_tuple(long_(shape[0]),long_(shape[1]),long_(shape[2]));
201     case 4:
202     return make_tuple(long_(shape[0]),long_(shape[1]),long_(shape[2]),long_(shape[3]));
203     default:
204     throw DataException("Error - illegal Data rank.");
205     }
206     }
207    
208     void
209     Data::copy(const Data& other)
210     {
211     //
212     // Perform a deep copy
213     {
214     DataExpanded* temp=dynamic_cast<DataExpanded*>(other.m_data.get());
215     if (temp!=0) {
216     //
217     // Construct a DataExpanded copy
218     DataAbstract* newData=new DataExpanded(*temp);
219 jgs 102 shared_ptr<DataAbstract> temp_data(newData);
220     m_data=temp_data;
221 jgs 94 return;
222     }
223     }
224     {
225     DataTagged* temp=dynamic_cast<DataTagged*>(other.m_data.get());
226     if (temp!=0) {
227     //
228 jgs 102 // Construct a DataTagged copy
229 jgs 94 DataAbstract* newData=new DataTagged(*temp);
230 jgs 102 shared_ptr<DataAbstract> temp_data(newData);
231     m_data=temp_data;
232 jgs 94 return;
233     }
234     }
235     {
236     DataConstant* temp=dynamic_cast<DataConstant*>(other.m_data.get());
237     if (temp!=0) {
238     //
239     // Construct a DataConstant copy
240     DataAbstract* newData=new DataConstant(*temp);
241 jgs 102 shared_ptr<DataAbstract> temp_data(newData);
242     m_data=temp_data;
243 jgs 94 return;
244     }
245     }
246 jgs 102 {
247     DataEmpty* temp=dynamic_cast<DataEmpty*>(other.m_data.get());
248     if (temp!=0) {
249     //
250     // Construct a DataEmpty copy
251     DataAbstract* newData=new DataEmpty();
252     shared_ptr<DataAbstract> temp_data(newData);
253     m_data=temp_data;
254     return;
255     }
256     }
257 jgs 94 throw DataException("Error - Copy not implemented for this Data type.");
258     }
259    
260     void
261     Data::copyWithMask(const Data& other,
262     const Data& mask)
263     {
264     Data mask1;
265     Data mask2;
266    
267     mask1 = mask.wherePositive();
268     mask2.copy(mask1);
269    
270     mask1 *= other;
271     mask2 *= *this;
272     mask2 = *this - mask2;
273    
274     *this = mask1 + mask2;
275     }
276    
277     bool
278     Data::isExpanded() const
279     {
280     DataExpanded* temp=dynamic_cast<DataExpanded*>(m_data.get());
281     return (temp!=0);
282     }
283    
284     bool
285     Data::isTagged() const
286     {
287     DataTagged* temp=dynamic_cast<DataTagged*>(m_data.get());
288     return (temp!=0);
289     }
290    
291     bool
292     Data::isEmpty() const
293     {
294     DataEmpty* temp=dynamic_cast<DataEmpty*>(m_data.get());
295     return (temp!=0);
296     }
297    
298     bool
299     Data::isConstant() const
300     {
301     DataConstant* temp=dynamic_cast<DataConstant*>(m_data.get());
302     return (temp!=0);
303     }
304    
305     void
306     Data::expand()
307     {
308     if (isConstant()) {
309     DataConstant* tempDataConst=dynamic_cast<DataConstant*>(m_data.get());
310     DataAbstract* temp=new DataExpanded(*tempDataConst);
311 jgs 102 shared_ptr<DataAbstract> temp_data(temp);
312     m_data=temp_data;
313 jgs 94 } else if (isTagged()) {
314     DataTagged* tempDataTag=dynamic_cast<DataTagged*>(m_data.get());
315     DataAbstract* temp=new DataExpanded(*tempDataTag);
316 jgs 102 shared_ptr<DataAbstract> temp_data(temp);
317     m_data=temp_data;
318 jgs 94 } else if (isExpanded()) {
319     //
320     // do nothing
321     } else if (isEmpty()) {
322     throw DataException("Error - Expansion of DataEmpty not possible.");
323     } else {
324     throw DataException("Error - Expansion not implemented for this Data type.");
325     }
326     }
327    
328     void
329     Data::tag()
330     {
331     if (isConstant()) {
332     DataConstant* tempDataConst=dynamic_cast<DataConstant*>(m_data.get());
333     DataAbstract* temp=new DataTagged(*tempDataConst);
334 jgs 102 shared_ptr<DataAbstract> temp_data(temp);
335     m_data=temp_data;
336 jgs 94 } else if (isTagged()) {
337     // do nothing
338     } else if (isExpanded()) {
339     throw DataException("Error - Creating tag data from DataExpanded not possible.");
340     } else if (isEmpty()) {
341     throw DataException("Error - Creating tag data from DataEmpty not possible.");
342     } else {
343     throw DataException("Error - Tagging not implemented for this Data type.");
344     }
345     }
346    
347 jgs 102 void
348     Data::reshapeDataPoint(const DataArrayView::ShapeType& shape)
349     {
350     m_data->reshapeDataPoint(shape);
351     }
352    
353 jgs 94 Data
354     Data::wherePositive() const
355     {
356     return escript::unaryOp(*this,bind2nd(greater<double>(),0.0));
357     }
358    
359     Data
360 jgs 102 Data::whereNegative() const
361     {
362     return escript::unaryOp(*this,bind2nd(less<double>(),0.0));
363     }
364    
365     Data
366 jgs 94 Data::whereNonNegative() const
367     {
368     return escript::unaryOp(*this,bind2nd(greater_equal<double>(),0.0));
369     }
370    
371     Data
372 jgs 102 Data::whereNonPositive() const
373 jgs 94 {
374 jgs 102 return escript::unaryOp(*this,bind2nd(less_equal<double>(),0.0));
375 jgs 94 }
376    
377     Data
378     Data::whereZero() const
379     {
380     return escript::unaryOp(*this,bind2nd(equal_to<double>(),0.0));
381     }
382    
383     Data
384 jgs 102 Data::whereNonZero() const
385     {
386     return escript::unaryOp(*this,bind2nd(not_equal_to<double>(),0.0));
387     }
388    
389     Data
390 jgs 94 Data::interpolate(const FunctionSpace& functionspace) const
391     {
392     return Data(*this,functionspace);
393     }
394    
395     bool
396     Data::probeInterpolation(const FunctionSpace& functionspace) const
397     {
398     if (getFunctionSpace()==functionspace) {
399     return true;
400     } else {
401     const AbstractDomain& domain=getDomain();
402     if (domain==functionspace.getDomain()) {
403     return domain.probeInterpolationOnDomain(getFunctionSpace().getTypeCode(),functionspace.getTypeCode());
404     } else {
405     return domain.probeInterpolationACross(getFunctionSpace().getTypeCode(),functionspace.getDomain(),functionspace.getTypeCode());
406     }
407     }
408     }
409    
410     Data
411     Data::gradOn(const FunctionSpace& functionspace) const
412     {
413     if (functionspace.getDomain()!=getDomain())
414     throw DataException("Error - gradient cannot be calculated on different domains.");
415     DataArrayView::ShapeType grad_shape=getPointDataView().getShape();
416     grad_shape.push_back(functionspace.getDim());
417     Data out(0.0,grad_shape,functionspace,true);
418     getDomain().setToGradient(out,*this);
419     return out;
420     }
421    
422     Data
423     Data::grad() const
424     {
425     return gradOn(escript::function(getDomain()));
426     }
427    
428     int
429     Data::getDataPointSize() const
430     {
431     return getPointDataView().noValues();
432     }
433    
434     DataArrayView::ValueType::size_type
435     Data::getLength() const
436     {
437     return m_data->getLength();
438     }
439    
440     const DataArrayView::ShapeType&
441     Data::getDataPointShape() const
442     {
443     return getPointDataView().getShape();
444     }
445    
446     boost::python::numeric::array
447     Data::integrate() const
448     {
449     int index;
450     int rank = getDataPointRank();
451     DataArrayView::ShapeType shape = getDataPointShape();
452    
453     //
454     // calculate the integral values
455     vector<double> integrals(getDataPointSize());
456     AbstractContinuousDomain::asAbstractContinuousDomain(getDomain()).setToIntegrals(integrals,*this);
457    
458     //
459     // create the numeric array to be returned
460     // and load the array with the integral values
461     boost::python::numeric::array bp_array(1.0);
462     if (rank==0) {
463 jgs 108 bp_array.resize(1);
464 jgs 94 index = 0;
465     bp_array[0] = integrals[index];
466     }
467     if (rank==1) {
468     bp_array.resize(shape[0]);
469     for (int i=0; i<shape[0]; i++) {
470     index = i;
471     bp_array[i] = integrals[index];
472     }
473     }
474     if (rank==2) {
475     bp_array.resize(shape[0],shape[1]);
476     for (int i=0; i<shape[0]; i++) {
477     for (int j=0; j<shape[1]; j++) {
478     index = i + shape[0] * j;
479     bp_array[i,j] = integrals[index];
480     }
481     }
482     }
483     if (rank==3) {
484     bp_array.resize(shape[0],shape[1],shape[2]);
485     for (int i=0; i<shape[0]; i++) {
486     for (int j=0; j<shape[1]; j++) {
487     for (int k=0; k<shape[2]; k++) {
488     index = i + shape[0] * ( j + shape[1] * k );
489     bp_array[i,j,k] = integrals[index];
490     }
491     }
492     }
493     }
494     if (rank==4) {
495     bp_array.resize(shape[0],shape[1],shape[2],shape[3]);
496     for (int i=0; i<shape[0]; i++) {
497     for (int j=0; j<shape[1]; j++) {
498     for (int k=0; k<shape[2]; k++) {
499     for (int l=0; l<shape[3]; l++) {
500     index = i + shape[0] * ( j + shape[1] * ( k + shape[2] * l ) );
501     bp_array[i,j,k,l] = integrals[index];
502     }
503     }
504     }
505     }
506     }
507    
508     //
509     // return the loaded array
510     return bp_array;
511     }
512    
513     Data
514     Data::sin() const
515     {
516     return escript::unaryOp(*this,(Data::UnaryDFunPtr)::sin);
517     }
518    
519     Data
520     Data::cos() const
521     {
522     return escript::unaryOp(*this,(Data::UnaryDFunPtr)::cos);
523     }
524    
525     Data
526     Data::tan() const
527     {
528     return escript::unaryOp(*this,(Data::UnaryDFunPtr)::tan);
529     }
530    
531     Data
532     Data::log() const
533     {
534     return escript::unaryOp(*this,(Data::UnaryDFunPtr)::log10);
535     }
536    
537     Data
538     Data::ln() const
539     {
540     return escript::unaryOp(*this,(Data::UnaryDFunPtr)::log);
541     }
542    
543 jgs 106 Data
544     Data::sign() const
545 jgs 94 {
546 jgs 106 return escript::unaryOp(*this,escript::fsign);
547 jgs 94 }
548    
549 jgs 106 Data
550     Data::abs() const
551 jgs 94 {
552 jgs 106 return escript::unaryOp(*this,(Data::UnaryDFunPtr)::fabs);
553 jgs 94 }
554    
555 jgs 106 Data
556     Data::neg() const
557 jgs 94 {
558 jgs 106 return escript::unaryOp(*this,negate<double>());
559 jgs 94 }
560    
561 jgs 102 Data
562 jgs 106 Data::pos() const
563 jgs 94 {
564 jgs 102 return (*this);
565     }
566    
567     Data
568 jgs 106 Data::exp() const
569 jgs 102 {
570 jgs 106 return escript::unaryOp(*this,(Data::UnaryDFunPtr)::exp);
571 jgs 102 }
572    
573     Data
574 jgs 106 Data::sqrt() const
575 jgs 102 {
576 jgs 106 return escript::unaryOp(*this,(Data::UnaryDFunPtr)::sqrt);
577 jgs 102 }
578    
579 jgs 106 double
580     Data::Lsup() const
581 jgs 102 {
582 jgs 106 //
583     // set the initial absolute maximum value to zero
584     return algorithm(DataAlgorithmAdapter<AbsMax>(0));
585 jgs 102 }
586    
587 jgs 106 double
588     Data::sup() const
589 jgs 102 {
590 jgs 106 //
591     // set the initial maximum value to min possible double
592 jgs 108 return algorithm(DataAlgorithmAdapter<FMax>(numeric_limits<double>::max()*-1));
593 jgs 102 }
594    
595 jgs 106 double
596     Data::inf() const
597 jgs 102 {
598 jgs 106 //
599     // set the initial minimum value to max possible double
600     return algorithm(DataAlgorithmAdapter<FMin>(numeric_limits<double>::max()));
601 jgs 102 }
602    
603     Data
604 jgs 106 Data::maxval() const
605 jgs 102 {
606 jgs 113 //
607     // set the initial maximum value to min possible double
608 jgs 108 return dp_algorithm(DataAlgorithmAdapter<FMax>(numeric_limits<double>::max()*-1));
609 jgs 102 }
610    
611     Data
612 jgs 106 Data::minval() const
613 jgs 102 {
614 jgs 113 //
615     // set the initial minimum value to max possible double
616 jgs 106 return dp_algorithm(DataAlgorithmAdapter<FMin>(numeric_limits<double>::max()));
617 jgs 102 }
618    
619     Data
620 jgs 106 Data::length() const
621 jgs 102 {
622 jgs 106 return dp_algorithm(DataAlgorithmAdapter<Length>(0));
623 jgs 102 }
624    
625     Data
626 jgs 106 Data::trace() const
627 jgs 102 {
628 jgs 106 return dp_algorithm(DataAlgorithmAdapter<Trace>(0));
629 jgs 102 }
630    
631     Data
632 jgs 106 Data::transpose(int axis) const
633 jgs 102 {
634 jgs 106 // not implemented
635     throw DataException("Error - Data::transpose not implemented yet.");
636     return Data();
637 jgs 102 }
638    
639 jgs 104 void
640     Data::saveDX(std::string fileName) const
641     {
642     getDomain().saveDX(fileName,*this);
643     return;
644     }
645    
646 jgs 110 void
647     Data::saveVTK(std::string fileName) const
648     {
649     getDomain().saveVTK(fileName,*this);
650     return;
651     }
652    
653 jgs 102 Data&
654     Data::operator+=(const Data& right)
655     {
656 jgs 94 binaryOp(right,plus<double>());
657     return (*this);
658     }
659    
660 jgs 102 Data&
661     Data::operator+=(const boost::python::object& right)
662 jgs 94 {
663     binaryOp(right,plus<double>());
664     return (*this);
665     }
666    
667 jgs 102 Data&
668     Data::operator-=(const Data& right)
669 jgs 94 {
670     binaryOp(right,minus<double>());
671     return (*this);
672     }
673    
674 jgs 102 Data&
675     Data::operator-=(const boost::python::object& right)
676 jgs 94 {
677     binaryOp(right,minus<double>());
678     return (*this);
679     }
680    
681 jgs 102 Data&
682     Data::operator*=(const Data& right)
683 jgs 94 {
684     binaryOp(right,multiplies<double>());
685     return (*this);
686     }
687    
688 jgs 102 Data&
689     Data::operator*=(const boost::python::object& right)
690 jgs 94 {
691     binaryOp(right,multiplies<double>());
692     return (*this);
693     }
694    
695 jgs 102 Data&
696     Data::operator/=(const Data& right)
697 jgs 94 {
698     binaryOp(right,divides<double>());
699     return (*this);
700     }
701    
702 jgs 102 Data&
703     Data::operator/=(const boost::python::object& right)
704 jgs 94 {
705     binaryOp(right,divides<double>());
706     return (*this);
707     }
708    
709 jgs 102 Data
710     Data::powO(const boost::python::object& right) const
711 jgs 94 {
712     Data result;
713     result.copy(*this);
714     result.binaryOp(right,(Data::BinaryDFunPtr)::pow);
715     return result;
716     }
717    
718 jgs 102 Data
719     Data::powD(const Data& right) const
720 jgs 94 {
721     Data result;
722     result.copy(*this);
723     result.binaryOp(right,(Data::BinaryDFunPtr)::pow);
724     return result;
725     }
726    
727     //
728     // NOTE: It is essential to specify the namepsace this operator belongs to
729 jgs 102 Data
730     escript::operator+(const Data& left, const Data& right)
731 jgs 94 {
732     Data result;
733     //
734     // perform a deep copy
735     result.copy(left);
736     result+=right;
737     return result;
738     }
739    
740     //
741     // NOTE: It is essential to specify the namepsace this operator belongs to
742 jgs 102 Data
743     escript::operator-(const Data& left, const Data& right)
744 jgs 94 {
745     Data result;
746     //
747     // perform a deep copy
748     result.copy(left);
749     result-=right;
750     return result;
751     }
752    
753     //
754     // NOTE: It is essential to specify the namepsace this operator belongs to
755 jgs 102 Data
756     escript::operator*(const Data& left, const Data& right)
757 jgs 94 {
758     Data result;
759     //
760     // perform a deep copy
761     result.copy(left);
762     result*=right;
763     return result;
764     }
765    
766     //
767     // NOTE: It is essential to specify the namepsace this operator belongs to
768 jgs 102 Data
769     escript::operator/(const Data& left, const Data& right)
770 jgs 94 {
771     Data result;
772     //
773     // perform a deep copy
774     result.copy(left);
775     result/=right;
776     return result;
777     }
778    
779     //
780     // NOTE: It is essential to specify the namepsace this operator belongs to
781 jgs 102 Data
782     escript::operator+(const Data& left, const boost::python::object& right)
783 jgs 94 {
784     //
785     // Convert to DataArray format if possible
786     DataArray temp(right);
787     Data result;
788     //
789     // perform a deep copy
790     result.copy(left);
791     result+=right;
792     return result;
793     }
794    
795     //
796     // NOTE: It is essential to specify the namepsace this operator belongs to
797 jgs 102 Data
798     escript::operator-(const Data& left, const boost::python::object& right)
799 jgs 94 {
800     //
801     // Convert to DataArray format if possible
802     DataArray temp(right);
803     Data result;
804     //
805     // perform a deep copy
806     result.copy(left);
807     result-=right;
808     return result;
809     }
810    
811     //
812     // NOTE: It is essential to specify the namepsace this operator belongs to
813 jgs 102 Data
814     escript::operator*(const Data& left, const boost::python::object& right)
815 jgs 94 {
816     //
817     // Convert to DataArray format if possible
818     DataArray temp(right);
819     Data result;
820     //
821     // perform a deep copy
822     result.copy(left);
823     result*=right;
824     return result;
825     }
826    
827     //
828     // NOTE: It is essential to specify the namepsace this operator belongs to
829 jgs 102 Data
830     escript::operator/(const Data& left, const boost::python::object& right)
831 jgs 94 {
832     //
833     // Convert to DataArray format if possible
834     DataArray temp(right);
835     Data result;
836     //
837     // perform a deep copy
838     result.copy(left);
839     result/=right;
840     return result;
841     }
842    
843     //
844     // NOTE: It is essential to specify the namepsace this operator belongs to
845 jgs 102 Data
846     escript::operator+(const boost::python::object& left, const Data& right)
847 jgs 94 {
848     //
849     // Construct the result using the given value and the other parameters
850     // from right
851     Data result(left,right);
852     result+=right;
853     return result;
854     }
855    
856     //
857     // NOTE: It is essential to specify the namepsace this operator belongs to
858 jgs 102 Data
859     escript::operator-(const boost::python::object& left, const Data& right)
860 jgs 94 {
861     //
862     // Construct the result using the given value and the other parameters
863     // from right
864     Data result(left,right);
865     result-=right;
866     return result;
867     }
868    
869     //
870     // NOTE: It is essential to specify the namepsace this operator belongs to
871 jgs 102 Data
872     escript::operator*(const boost::python::object& left, const Data& right)
873 jgs 94 {
874     //
875     // Construct the result using the given value and the other parameters
876     // from right
877     Data result(left,right);
878     result*=right;
879     return result;
880     }
881    
882     //
883     // NOTE: It is essential to specify the namepsace this operator belongs to
884 jgs 102 Data
885     escript::operator/(const boost::python::object& left, const Data& right)
886 jgs 94 {
887     //
888     // Construct the result using the given value and the other parameters
889     // from right
890     Data result(left,right);
891     result/=right;
892     return result;
893     }
894    
895     //
896     // NOTE: It is essential to specify the namepsace this operator belongs to
897 jgs 102 //bool escript::operator==(const Data& left, const Data& right)
898     //{
899     // /*
900     // NB: this operator does very little at this point, and isn't to
901     // be relied on. Requires further implementation.
902     // */
903     //
904     // bool ret;
905     //
906     // if (left.isEmpty()) {
907     // if(!right.isEmpty()) {
908     // ret = false;
909     // } else {
910     // ret = true;
911     // }
912     // }
913     //
914     // if (left.isConstant()) {
915     // if(!right.isConstant()) {
916     // ret = false;
917     // } else {
918     // ret = true;
919     // }
920     // }
921     //
922     // if (left.isTagged()) {
923     // if(!right.isTagged()) {
924     // ret = false;
925     // } else {
926     // ret = true;
927     // }
928     // }
929     //
930     // if (left.isExpanded()) {
931     // if(!right.isExpanded()) {
932     // ret = false;
933     // } else {
934     // ret = true;
935     // }
936     // }
937     //
938     // return ret;
939     //}
940    
941     Data
942     Data::getItem(const boost::python::object& key) const
943 jgs 94 {
944 jgs 102 const DataArrayView& view=getPointDataView();
945 jgs 94
946 jgs 102 DataArrayView::RegionType slice_region=view.getSliceRegion(key);
947 jgs 94
948 jgs 102 if (slice_region.size()!=view.getRank()) {
949     throw DataException("Error - slice size does not match Data rank.");
950 jgs 94 }
951    
952 jgs 102 return getSlice(slice_region);
953 jgs 94 }
954    
955     Data
956 jgs 102 Data::getSlice(const DataArrayView::RegionType& region) const
957 jgs 94 {
958 jgs 102 return Data(*this,region);
959 jgs 94 }
960    
961     void
962 jgs 102 Data::setItemO(const boost::python::object& key,
963     const boost::python::object& value)
964 jgs 94 {
965 jgs 102 Data tempData(value,getFunctionSpace());
966     setItemD(key,tempData);
967     }
968    
969     void
970     Data::setItemD(const boost::python::object& key,
971     const Data& value)
972     {
973 jgs 94 const DataArrayView& view=getPointDataView();
974     DataArrayView::RegionType slice_region=view.getSliceRegion(key);
975     if (slice_region.size()!=view.getRank()) {
976     throw DataException("Error - slice size does not match Data rank.");
977     }
978 jgs 108 if (getFunctionSpace()!=value.getFunctionSpace()) {
979     setSlice(Data(value,getFunctionSpace()),slice_region);
980     } else {
981     setSlice(value,slice_region);
982     }
983 jgs 94 }
984    
985     void
986 jgs 102 Data::setSlice(const Data& value,
987     const DataArrayView::RegionType& region)
988 jgs 94 {
989 jgs 102 Data tempValue(value);
990     typeMatchLeft(tempValue);
991     typeMatchRight(tempValue);
992     m_data->setSlice(tempValue.m_data.get(),region);
993     }
994    
995     void
996     Data::typeMatchLeft(Data& right) const
997     {
998     if (isExpanded()){
999     right.expand();
1000     } else if (isTagged()) {
1001     if (right.isConstant()) {
1002     right.tag();
1003     }
1004     }
1005     }
1006    
1007     void
1008     Data::typeMatchRight(const Data& right)
1009     {
1010 jgs 94 if (isTagged()) {
1011     if (right.isExpanded()) {
1012     expand();
1013     }
1014     } else if (isConstant()) {
1015     if (right.isExpanded()) {
1016     expand();
1017     } else if (right.isTagged()) {
1018     tag();
1019     }
1020     }
1021     }
1022    
1023     void
1024     Data::setTaggedValue(int tagKey,
1025     const boost::python::object& value)
1026     {
1027     //
1028     // Ensure underlying data object is of type DataTagged
1029     tag();
1030    
1031     if (!isTagged()) {
1032     throw DataException("Error - DataTagged conversion failed!!");
1033     }
1034    
1035     //
1036     // Construct DataArray from boost::python::object input value
1037     DataArray valueDataArray(value);
1038    
1039     //
1040     // Call DataAbstract::setTaggedValue
1041     m_data->setTaggedValue(tagKey,valueDataArray.getView());
1042     }
1043    
1044 jgs 110 void
1045     Data::setRefValue(int ref,
1046     const boost::python::numeric::array& value)
1047     {
1048     //
1049     // Construct DataArray from boost::python::object input value
1050     DataArray valueDataArray(value);
1051    
1052     //
1053     // Call DataAbstract::setRefValue
1054     m_data->setRefValue(ref,valueDataArray);
1055     }
1056    
1057     void
1058     Data::getRefValue(int ref,
1059     boost::python::numeric::array& value)
1060     {
1061     //
1062     // Construct DataArray for boost::python::object return value
1063     DataArray valueDataArray(value);
1064    
1065     //
1066     // Load DataArray with values from data-points specified by ref
1067     m_data->getRefValue(ref,valueDataArray);
1068    
1069     //
1070     // Load values from valueDataArray into return numarray
1071    
1072     // extract the shape of the numarray
1073     int rank = value.getrank();
1074     DataArrayView::ShapeType shape;
1075     for (int i=0; i < rank; i++) {
1076     shape.push_back(extract<int>(value.getshape()[i]));
1077     }
1078    
1079     // and load the numarray with the data from the DataArray
1080     DataArrayView valueView = valueDataArray.getView();
1081    
1082     if (rank==0) {
1083     throw DataException("Data::getRefValue error: only rank 1 data handled for now.");
1084     }
1085     if (rank==1) {
1086     for (int i=0; i < shape[0]; i++) {
1087     value[i] = valueView(i);
1088     }
1089     }
1090     if (rank==2) {
1091     throw DataException("Data::getRefValue error: only rank 1 data handled for now.");
1092     }
1093     if (rank==3) {
1094     throw DataException("Data::getRefValue error: only rank 1 data handled for now.");
1095     }
1096     if (rank==4) {
1097     throw DataException("Data::getRefValue error: only rank 1 data handled for now.");
1098     }
1099    
1100     }
1101    
1102 jgs 94 /*
1103     Note: this version removed for now. Not needed, and breaks escript.cpp
1104     void
1105     Data::setTaggedValue(int tagKey,
1106     const DataArrayView& value)
1107     {
1108     //
1109     // Ensure underlying data object is of type DataTagged
1110     tag();
1111    
1112     if (!isTagged()) {
1113     throw DataException("Error - DataTagged conversion failed!!");
1114     }
1115    
1116     //
1117     // Call DataAbstract::setTaggedValue
1118     m_data->setTaggedValue(tagKey,value);
1119     }
1120     */
1121    
1122     ostream& escript::operator<<(ostream& o, const Data& data)
1123     {
1124     o << data.toString();
1125     return o;
1126     }

Properties

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

  ViewVC Help
Powered by ViewVC 1.1.26