/[escript]/branches/schroedinger/escript/src/DataLazy.cpp
ViewVC logotype

Contents of /branches/schroedinger/escript/src/DataLazy.cpp

Parent Directory Parent Directory | Revision Log Revision Log


Revision 1899 - (show annotations)
Mon Oct 20 05:13:24 2008 UTC (10 years, 6 months ago) by jfenwick
File size: 24053 byte(s)
Branch commit.
Added some doco to DataLazy.
Made Data::integrate aware of DataLazy.


1
2 /*******************************************************
3 *
4 * Copyright (c) 2003-2008 by University of Queensland
5 * Earth Systems Science Computational Center (ESSCC)
6 * http://www.uq.edu.au/esscc
7 *
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
15 #include "DataLazy.h"
16 #ifdef USE_NETCDF
17 #include <netcdfcpp.h>
18 #endif
19 #ifdef PASO_MPI
20 #include <mpi.h>
21 #endif
22 #ifdef _OPENMP
23 #include <omp.h>
24 #endif
25 #include "FunctionSpace.h"
26 #include "DataTypes.h"
27 #include "Data.h"
28 #include "UnaryFuncs.h" // for escript::fsign
29
30 /*
31 How does DataLazy work?
32 ~~~~~~~~~~~~~~~~~~~~~~~
33
34 Each instance represents a single operation on one or two other DataLazy instances. These arguments are normally
35 denoted left and right.
36
37 A special operation, IDENTITY, stores an instance of DataReady in the m_id member.
38 This means that all "internal" nodes in the structure are instances of DataLazy.
39
40 Each operation has a string representation as well as an opgroup - eg G_IDENTITY, G_BINARY, ...
41 Note that IDENITY is not considered a unary operation.
42
43 I am avoiding calling the structure formed a tree because it is not guaranteed to be one (eg c=a+a).
44 It must however form a DAG (directed acyclic graph).
45 I will refer to individual DataLazy objects with the structure as nodes.
46
47 Each node also stores:
48 - m_readytype \in {'E','T','C','?'} ~ indicates what sort of DataReady would be produced if the expression was
49 evaluated.
50 - m_length ~ how many values would be stored in the answer if the expression was evaluated.
51 - m_buffsrequired ~ the larged number of samples which would need to be kept simultaneously in order to
52 evaluate the expression.
53 - m_samplesize ~ the number of doubles stored in a sample.
54
55 When a new node is created, the above values are computed based on the values in the child nodes.
56 Eg: if left requires 4 samples and right requires 6 then left+right requires 7 samples.
57
58 The resolve method, which produces a DataReady from a DataLazy, does the following:
59 1) Create a DataReady to hold the new result.
60 2) Allocate a vector (v) big enough to hold m_buffsrequired samples.
61 3) For each sample, call resolveSample with v, to get its values and copy them into the result object.
62
63 (In the case of OMP, multiple samples are resolved in parallel so the vector needs to be larger.)
64
65 resolveSample returns a Vector* and an offset within that vector where the result is stored.
66 Normally, this would be v, but for identity nodes their internal vector is returned instead.
67
68 The convention that I use, is that the resolve methods should store their results starting at the offset they are passed.
69
70 For expressions which evaluate to Constant or Tagged, there is a different evaluation method.
71 The collapse method invokes the (non-lazy) operations on the Data class to evaluate the expression.
72 */
73
74
75 using namespace std;
76 using namespace boost;
77
78 namespace escript
79 {
80
81 const std::string&
82 opToString(ES_optype op);
83
84 namespace
85 {
86
87 enum ES_opgroup
88 {
89 G_UNKNOWN,
90 G_IDENTITY,
91 G_BINARY, // pointwise operations with two arguments
92 G_UNARY // pointwise operations with one argument
93 };
94
95
96
97
98 string ES_opstrings[]={"UNKNOWN","IDENTITY","+","-","*","/","sin","cos","tan",
99 "asin","acos","atan","sinh","cosh","tanh","erf",
100 "asinh","acosh","atanh",
101 "log10","log","sign","abs","neg","pos","exp","sqrt",
102 "1/","where>0","where<0","where>=0","where<=0"};
103 int ES_opcount=32;
104 ES_opgroup opgroups[]={G_UNKNOWN,G_IDENTITY,G_BINARY,G_BINARY,G_BINARY,G_BINARY,G_UNARY,G_UNARY,G_UNARY, //9
105 G_UNARY,G_UNARY,G_UNARY,G_UNARY,G_UNARY,G_UNARY,G_UNARY, // 16
106 G_UNARY,G_UNARY,G_UNARY, // 19
107 G_UNARY,G_UNARY,G_UNARY,G_UNARY,G_UNARY,G_UNARY,G_UNARY,G_UNARY, // 27
108 G_UNARY,G_UNARY,G_UNARY,G_UNARY,G_UNARY};
109 inline
110 ES_opgroup
111 getOpgroup(ES_optype op)
112 {
113 return opgroups[op];
114 }
115
116 // return the FunctionSpace of the result of "left op right"
117 FunctionSpace
118 resultFS(DataAbstract_ptr left, DataAbstract_ptr right, ES_optype op)
119 {
120 // perhaps this should call interpolate and throw or something?
121 // maybe we need an interpolate node -
122 // that way, if interpolate is required in any other op we can just throw a
123 // programming error exception.
124
125
126 if (left->getFunctionSpace()!=right->getFunctionSpace())
127 {
128 throw DataException("FunctionSpaces not equal - interpolation not supported on lazy data.");
129 }
130 return left->getFunctionSpace();
131 }
132
133 // return the shape of the result of "left op right"
134 DataTypes::ShapeType
135 resultShape(DataAbstract_ptr left, DataAbstract_ptr right, ES_optype op)
136 {
137 if (left->getShape()!=right->getShape())
138 {
139 throw DataException("Shapes not the same - shapes must match for lazy data.");
140 }
141 return left->getShape();
142 }
143
144 // determine the number of points in the result of "left op right"
145 size_t
146 resultLength(DataAbstract_ptr left, DataAbstract_ptr right, ES_optype op)
147 {
148 switch (getOpgroup(op))
149 {
150 case G_BINARY: return left->getLength();
151 case G_UNARY: return left->getLength();
152 default:
153 throw DataException("Programmer Error - attempt to getLength() for operator "+opToString(op)+".");
154 }
155 }
156
157 // determine the number of samples requires to evaluate an expression combining left and right
158 int
159 calcBuffs(const DataLazy_ptr& left, const DataLazy_ptr& right, ES_optype op)
160 {
161 switch(getOpgroup(op))
162 {
163 case G_IDENTITY: return 1;
164 case G_BINARY: return max(left->getBuffsRequired(),right->getBuffsRequired()+1);
165 case G_UNARY: return max(left->getBuffsRequired(),1);
166 default:
167 throw DataException("Programmer Error - attempt to calcBuffs() for operator "+opToString(op)+".");
168 }
169 }
170
171
172 } // end anonymous namespace
173
174
175
176 // Return a string representing the operation
177 const std::string&
178 opToString(ES_optype op)
179 {
180 if (op<0 || op>=ES_opcount)
181 {
182 op=UNKNOWNOP;
183 }
184 return ES_opstrings[op];
185 }
186
187
188 DataLazy::DataLazy(DataAbstract_ptr p)
189 : parent(p->getFunctionSpace(),p->getShape()),
190 m_op(IDENTITY)
191 {
192 if (p->isLazy())
193 {
194 // I don't want identity of Lazy.
195 // Question: Why would that be so bad?
196 // Answer: We assume that the child of ID is something we can call getVector on
197 throw DataException("Programmer error - attempt to create identity from a DataLazy.");
198 }
199 else
200 {
201 m_id=dynamic_pointer_cast<DataReady>(p);
202 if(p->isConstant()) {m_readytype='C';}
203 else if(p->isExpanded()) {m_readytype='E';}
204 else if (p->isTagged()) {m_readytype='T';}
205 else {throw DataException("Unknown DataReady instance in DataLazy constructor.");}
206 }
207 m_length=p->getLength();
208 m_buffsRequired=1;
209 m_samplesize=getNumDPPSample()*getNoValues();
210 cout << "(1)Lazy created with " << m_samplesize << endl;
211 }
212
213
214
215 DataLazy::DataLazy(DataAbstract_ptr left, ES_optype op)
216 : parent(left->getFunctionSpace(),left->getShape()),
217 m_op(op)
218 {
219 if (getOpgroup(op)!=G_UNARY)
220 {
221 throw DataException("Programmer error - constructor DataLazy(left, op) will only process UNARY operations.");
222 }
223 DataLazy_ptr lleft;
224 if (!left->isLazy())
225 {
226 lleft=DataLazy_ptr(new DataLazy(left));
227 }
228 else
229 {
230 lleft=dynamic_pointer_cast<DataLazy>(left);
231 }
232 m_readytype=lleft->m_readytype;
233 m_length=left->getLength();
234 m_left=lleft;
235 m_buffsRequired=1;
236 m_samplesize=getNumDPPSample()*getNoValues();
237 }
238
239
240 DataLazy::DataLazy(DataAbstract_ptr left, DataAbstract_ptr right, ES_optype op)
241 : parent(resultFS(left,right,op), resultShape(left,right,op)),
242 m_op(op)
243 {
244 if (getOpgroup(op)!=G_BINARY)
245 {
246 throw DataException("Programmer error - constructor DataLazy(left, right, op) will only process BINARY operations.");
247 }
248 if (left->isLazy()) // the children need to be DataLazy. Wrap them in IDENTITY if required
249 {
250 m_left=dynamic_pointer_cast<DataLazy>(left);
251 }
252 else
253 {
254 m_left=DataLazy_ptr(new DataLazy(left));
255 }
256 if (right->isLazy())
257 {
258 m_right=dynamic_pointer_cast<DataLazy>(right);
259 }
260 else
261 {
262 m_right=DataLazy_ptr(new DataLazy(right));
263 }
264 char lt=m_left->m_readytype;
265 char rt=m_right->m_readytype;
266 if (lt=='E' || rt=='E')
267 {
268 m_readytype='E';
269 }
270 else if (lt=='T' || rt=='T')
271 {
272 m_readytype='T';
273 }
274 else
275 {
276 m_readytype='C';
277 }
278 m_length=resultLength(m_left,m_right,m_op);
279 m_samplesize=getNumDPPSample()*getNoValues();
280 m_buffsRequired=calcBuffs(m_left, m_right,m_op);
281 cout << "(3)Lazy created with " << m_samplesize << endl;
282 }
283
284
285 DataLazy::~DataLazy()
286 {
287 }
288
289
290 int
291 DataLazy::getBuffsRequired() const
292 {
293 return m_buffsRequired;
294 }
295
296
297 /*
298 \brief Evaluates the expression using methods on Data.
299 This does the work for the collapse method.
300 For reasons of efficiency do not call this method on DataExpanded nodes.
301 */
302 DataReady_ptr
303 DataLazy::collapseToReady()
304 {
305 if (m_readytype=='E')
306 { // this is more an efficiency concern than anything else
307 throw DataException("Programmer Error - do not use collapse on Expanded data.");
308 }
309 if (m_op==IDENTITY)
310 {
311 return m_id;
312 }
313 DataReady_ptr pleft=m_left->collapseToReady();
314 Data left(pleft);
315 Data right;
316 if (getOpgroup(m_op)==G_BINARY)
317 {
318 right=Data(m_right->collapseToReady());
319 }
320 Data result;
321 switch(m_op)
322 {
323 case ADD:
324 result=left+right;
325 break;
326 case SUB:
327 result=left-right;
328 break;
329 case MUL:
330 result=left*right;
331 break;
332 case DIV:
333 result=left/right;
334 break;
335 case SIN:
336 result=left.sin();
337 break;
338 case COS:
339 result=left.cos();
340 break;
341 case TAN:
342 result=left.tan();
343 break;
344 case ASIN:
345 result=left.asin();
346 break;
347 case ACOS:
348 result=left.acos();
349 break;
350 case ATAN:
351 result=left.atan();
352 break;
353 case SINH:
354 result=left.sinh();
355 break;
356 case COSH:
357 result=left.cosh();
358 break;
359 case TANH:
360 result=left.tanh();
361 break;
362 case ERF:
363 result=left.erf();
364 break;
365 case ASINH:
366 result=left.asinh();
367 break;
368 case ACOSH:
369 result=left.acosh();
370 break;
371 case ATANH:
372 result=left.atanh();
373 break;
374 case LOG10:
375 result=left.log10();
376 break;
377 case LOG:
378 result=left.log();
379 break;
380 case SIGN:
381 result=left.sign();
382 break;
383 case ABS:
384 result=left.abs();
385 break;
386 case NEG:
387 result=left.neg();
388 break;
389 case POS:
390 // it doesn't mean anything for delayed.
391 // it will just trigger a deep copy of the lazy object
392 throw DataException("Programmer error - POS not supported for lazy data.");
393 break;
394 case EXP:
395 result=left.exp();
396 break;
397 case SQRT:
398 result=left.sqrt();
399 break;
400 case RECIP:
401 result=left.oneOver();
402 break;
403 case GZ:
404 result=left.wherePositive();
405 break;
406 case LZ:
407 result=left.whereNegative();
408 break;
409 case GEZ:
410 result=left.whereNonNegative();
411 break;
412 case LEZ:
413 result=left.whereNonPositive();
414 break;
415 default:
416 throw DataException("Programmer error - do not know how to resolve operator "+opToString(m_op)+".");
417 }
418 return result.borrowReadyPtr();
419 }
420
421 /*
422 \brief Converts the DataLazy into an IDENTITY storing the value of the expression.
423 This method uses the original methods on the Data class to evaluate the expressions.
424 For this reason, it should not be used on DataExpanded instances. (To do so would defeat
425 the purpose of using DataLazy in the first place).
426 */
427 void
428 DataLazy::collapse()
429 {
430 if (m_op==IDENTITY)
431 {
432 return;
433 }
434 if (m_readytype=='E')
435 { // this is more an efficiency concern than anything else
436 throw DataException("Programmer Error - do not use collapse on Expanded data.");
437 }
438 m_id=collapseToReady();
439 m_op=IDENTITY;
440 }
441
442 /*
443 \brief Compute the value of the expression (binary operation) for the given sample.
444 \return Vector which stores the value of the subexpression for the given sample.
445 \param v A vector to store intermediate results.
446 \param offset Index in v to begin storing results.
447 \param sampleNo Sample number to evaluate.
448 \param roffset (output parameter) the offset in the return vector where the result begins.
449
450 The return value will be an existing vector so do not deallocate it.
451 If the result is stored in v it should be stored at the offset given.
452 Everything from offset to the end of v should be considered available for this method to use.
453 */
454 DataTypes::ValueType*
455 DataLazy::resolveUnary(ValueType& v, size_t offset, int sampleNo, size_t& roffset) const
456 {
457 // we assume that any collapsing has been done before we get here
458 // since we only have one argument we don't need to think about only
459 // processing single points.
460 if (m_readytype!='E')
461 {
462 throw DataException("Programmer error - resolveUnary should only be called on expanded Data.");
463 }
464 const ValueType* vleft=m_left->resolveSample(v,offset,sampleNo,roffset);
465 const double* left=&((*vleft)[roffset]);
466 double* result=&(v[offset]);
467 roffset=offset;
468 switch (m_op)
469 {
470 case SIN:
471 tensor_unary_operation(m_samplesize, left, result, ::sin);
472 break;
473 case COS:
474 tensor_unary_operation(m_samplesize, left, result, ::cos);
475 break;
476 case TAN:
477 tensor_unary_operation(m_samplesize, left, result, ::tan);
478 break;
479 case ASIN:
480 tensor_unary_operation(m_samplesize, left, result, ::asin);
481 break;
482 case ACOS:
483 tensor_unary_operation(m_samplesize, left, result, ::acos);
484 break;
485 case ATAN:
486 tensor_unary_operation(m_samplesize, left, result, ::atan);
487 break;
488 case SINH:
489 tensor_unary_operation(m_samplesize, left, result, ::sinh);
490 break;
491 case COSH:
492 tensor_unary_operation(m_samplesize, left, result, ::cosh);
493 break;
494 case TANH:
495 tensor_unary_operation(m_samplesize, left, result, ::tanh);
496 break;
497 case ERF:
498 #ifdef _WIN32
499 throw DataException("Error - Data:: erf function is not supported on _WIN32 platforms.");
500 #else
501 tensor_unary_operation(m_samplesize, left, result, ::erf);
502 break;
503 #endif
504 case ASINH:
505 #ifdef _WIN32
506 tensor_unary_operation(m_samplesize, left, result, escript::asinh_substitute);
507 #else
508 tensor_unary_operation(m_samplesize, left, result, ::asinh);
509 #endif
510 break;
511 case ACOSH:
512 #ifdef _WIN32
513 tensor_unary_operation(m_samplesize, left, result, escript::acosh_substitute);
514 #else
515 tensor_unary_operation(m_samplesize, left, result, ::acosh);
516 #endif
517 break;
518 case ATANH:
519 #ifdef _WIN32
520 tensor_unary_operation(m_samplesize, left, result, escript::atanh_substitute);
521 #else
522 tensor_unary_operation(m_samplesize, left, result, ::atanh);
523 #endif
524 break;
525 case LOG10:
526 tensor_unary_operation(m_samplesize, left, result, ::log10);
527 break;
528 case LOG:
529 tensor_unary_operation(m_samplesize, left, result, ::log);
530 break;
531 case SIGN:
532 tensor_unary_operation(m_samplesize, left, result, escript::fsign);
533 break;
534 case ABS:
535 tensor_unary_operation(m_samplesize, left, result, ::fabs);
536 break;
537 case NEG:
538 tensor_unary_operation(m_samplesize, left, result, negate<double>());
539 break;
540 case POS:
541 // it doesn't mean anything for delayed.
542 // it will just trigger a deep copy of the lazy object
543 throw DataException("Programmer error - POS not supported for lazy data.");
544 break;
545 case EXP:
546 tensor_unary_operation(m_samplesize, left, result, ::exp);
547 break;
548 case SQRT:
549 tensor_unary_operation(m_samplesize, left, result, ::sqrt);
550 break;
551 case RECIP:
552 tensor_unary_operation(m_samplesize, left, result, bind1st(divides<double>(),1.));
553 break;
554 case GZ:
555 tensor_unary_operation(m_samplesize, left, result, bind2nd(greater<double>(),0.0));
556 break;
557 case LZ:
558 tensor_unary_operation(m_samplesize, left, result, bind2nd(less<double>(),0.0));
559 break;
560 case GEZ:
561 tensor_unary_operation(m_samplesize, left, result, bind2nd(greater_equal<double>(),0.0));
562 break;
563 case LEZ:
564 tensor_unary_operation(m_samplesize, left, result, bind2nd(less_equal<double>(),0.0));
565 break;
566
567 default:
568 throw DataException("Programmer error - resolveUnary can not resolve operator "+opToString(m_op)+".");
569 }
570 return &v;
571 }
572
573
574
575
576
577 #define PROC_OP(X) \
578 for (int i=0;i<steps;++i,resultp+=getNoValues()) \
579 { \
580 tensor_binary_operation(chunksize, &((*left)[lroffset]), &((*right)[rroffset]), resultp, X); \
581 lroffset+=leftStep; \
582 rroffset+=rightStep; \
583 }
584
585 /*
586 \brief Compute the value of the expression (binary operation) for the given sample.
587 \return Vector which stores the value of the subexpression for the given sample.
588 \param v A vector to store intermediate results.
589 \param offset Index in v to begin storing results.
590 \param sampleNo Sample number to evaluate.
591 \param roffset (output parameter) the offset in the return vector where the result begins.
592
593 The return value will be an existing vector so do not deallocate it.
594 If the result is stored in v it should be stored at the offset given.
595 Everything from offset to the end of v should be considered available for this method to use.
596 */
597 // This method assumes that any subexpressions which evaluate to Constant or Tagged Data
598 // have already been collapsed to IDENTITY. So we must have at least one expanded child.
599 // If both children are expanded, then we can process them in a single operation (we treat
600 // the whole sample as one big datapoint.
601 // If one of the children is not expanded, then we need to treat each point in the sample
602 // individually.
603 DataTypes::ValueType*
604 DataLazy::resolveBinary(ValueType& v, size_t offset, int sampleNo, size_t& roffset) const
605 {
606 cout << "Resolve binary: " << toString() << endl;
607
608 size_t lroffset=0, rroffset=0; // offsets in the left and right result vectors
609 // first work out which of the children are expanded
610 bool leftExp=(m_left->m_readytype=='E');
611 bool rightExp=(m_right->m_readytype=='E');
612 bool bigloops=((leftExp && rightExp) || (!leftExp && !rightExp)); // is processing in single step?
613 int steps=(bigloops?1:getNumDPPSample());
614 size_t chunksize=(bigloops? m_samplesize : getNoValues()); // if bigloops, pretend the whole sample is a datapoint
615 int leftStep=((leftExp && !rightExp)? getNoValues() : 0);
616 int rightStep=((rightExp && !leftExp)? getNoValues() : 0);
617 // Get the values of sub-expressions
618 const ValueType* left=m_left->resolveSample(v,offset,sampleNo,lroffset);
619 const ValueType* right=m_right->resolveSample(v,offset+m_samplesize,sampleNo,rroffset); // Note
620 // the right child starts further along.
621 double* resultp=&(v[offset]); // results are stored at the vector offset we recieved
622 switch(m_op)
623 {
624 case ADD:
625 PROC_OP(plus<double>());
626 break;
627 case SUB:
628 PROC_OP(minus<double>());
629 break;
630 case MUL:
631 PROC_OP(multiplies<double>());
632 break;
633 case DIV:
634 PROC_OP(divides<double>());
635 break;
636 default:
637 throw DataException("Programmer error - resolveBinary can not resolve operator "+opToString(m_op)+".");
638 }
639 roffset=offset;
640 return &v;
641 }
642
643
644
645 /*
646 \brief Compute the value of the expression for the given sample.
647 \return Vector which stores the value of the subexpression for the given sample.
648 \param v A vector to store intermediate results.
649 \param offset Index in v to begin storing results.
650 \param sampleNo Sample number to evaluate.
651 \param roffset (output parameter) the offset in the return vector where the result begins.
652
653 The return value will be an existing vector so do not deallocate it.
654 */
655 // the vector and the offset are a place where the method could write its data if it wishes
656 // it is not obligated to do so. For example, if it has its own storage already, it can use that.
657 // Hence the return value to indicate where the data is actually stored.
658 // Regardless, the storage should be assumed to be used, even if it isn't.
659
660 // the roffset is the offset within the returned vector where the data begins
661 const DataTypes::ValueType*
662 DataLazy::resolveSample(ValueType& v, size_t offset, int sampleNo, size_t& roffset)
663 {
664 cout << "Resolve sample " << toString() << endl;
665 // collapse so we have a 'E' node or an IDENTITY for some other type
666 if (m_readytype!='E' && m_op!=IDENTITY)
667 {
668 collapse();
669 }
670 if (m_op==IDENTITY)
671 {
672 const ValueType& vec=m_id->getVector();
673 if (m_readytype=='C')
674 {
675 roffset=0;
676 return &(vec);
677 }
678 roffset=m_id->getPointOffset(sampleNo, 0);
679 return &(vec);
680 }
681 if (m_readytype!='E')
682 {
683 throw DataException("Programmer Error - Collapse did not produce an expanded node.");
684 }
685 switch (getOpgroup(m_op))
686 {
687 case G_UNARY: return resolveUnary(v, offset,sampleNo,roffset);
688 case G_BINARY: return resolveBinary(v, offset,sampleNo,roffset);
689 default:
690 throw DataException("Programmer Error - resolveSample does not know how to process "+opToString(m_op)+".");
691 }
692 }
693
694
695 // To simplify the memory management, all threads operate on one large vector, rather than one each.
696 // Each sample is evaluated independently and copied into the result DataExpanded.
697 DataReady_ptr
698 DataLazy::resolve()
699 {
700
701 cout << "Sample size=" << m_samplesize << endl;
702 cout << "Buffers=" << m_buffsRequired << endl;
703
704 if (m_readytype!='E') // if the whole sub-expression is Constant or Tagged, then evaluate it normally
705 {
706 collapse();
707 }
708 if (m_op==IDENTITY) // So a lazy expression of Constant or Tagged data will be returned here.
709 {
710 return m_id;
711 }
712 // from this point on we must have m_op!=IDENTITY and m_readytype=='E'
713 size_t threadbuffersize=m_samplesize*(max(1,m_buffsRequired)); // Each thread needs to have enough
714 // storage to evaluate its expression
715 int numthreads=1;
716 #ifdef _OPENMP
717 numthreads=getNumberOfThreads();
718 int threadnum=0;
719 #endif
720 ValueType v(numthreads*threadbuffersize);
721 cout << "Buffer created with size=" << v.size() << endl;
722 DataExpanded* result=new DataExpanded(getFunctionSpace(),getShape(), ValueType(getNoValues()));
723 ValueType& resvec=result->getVector();
724 DataReady_ptr resptr=DataReady_ptr(result);
725 int sample;
726 size_t outoffset; // offset in the output data
727 int totalsamples=getNumSamples();
728 const ValueType* res=0; // Vector storing the answer
729 size_t resoffset=0; // where in the vector to find the answer
730 #pragma omp parallel for private(sample,resoffset,outoffset,threadnum,res) schedule(static)
731 for (sample=0;sample<totalsamples;++sample)
732 {
733 cout << "################################# " << sample << endl;
734 #ifdef _OPENMP
735 res=resolveSample(v,threadbuffersize*omp_get_thread_num(),sample,resoffset);
736 #else
737 res=resolveSample(v,0,sample,resoffset); // res would normally be v, but not if its a single IDENTITY op.
738 #endif
739 cerr << "-------------------------------- " << endl;
740 outoffset=result->getPointOffset(sample,0);
741 cerr << "offset=" << outoffset << endl;
742 for (unsigned int i=0;i<m_samplesize;++i,++outoffset,++resoffset) // copy values into the output vector
743 {
744 resvec[outoffset]=(*res)[resoffset];
745 }
746 cerr << "*********************************" << endl;
747 }
748 return resptr;
749 }
750
751 std::string
752 DataLazy::toString() const
753 {
754 ostringstream oss;
755 oss << "Lazy Data:";
756 intoString(oss);
757 return oss.str();
758 }
759
760
761 void
762 DataLazy::intoString(ostringstream& oss) const
763 {
764 switch (getOpgroup(m_op))
765 {
766 case G_IDENTITY:
767 if (m_id->isExpanded())
768 {
769 oss << "E";
770 }
771 else if (m_id->isTagged())
772 {
773 oss << "T";
774 }
775 else if (m_id->isConstant())
776 {
777 oss << "C";
778 }
779 else
780 {
781 oss << "?";
782 }
783 oss << '@' << m_id.get();
784 break;
785 case G_BINARY:
786 oss << '(';
787 m_left->intoString(oss);
788 oss << ' ' << opToString(m_op) << ' ';
789 m_right->intoString(oss);
790 oss << ')';
791 break;
792 case G_UNARY:
793 oss << opToString(m_op) << '(';
794 m_left->intoString(oss);
795 oss << ')';
796 break;
797 default:
798 oss << "UNKNOWN";
799 }
800 }
801
802 // Note that in this case, deepCopy does not make copies of the leaves.
803 // Hopefully copy on write (or whatever we end up using) will take care of this.
804 DataAbstract*
805 DataLazy::deepCopy()
806 {
807 if (m_op==IDENTITY)
808 {
809 return new DataLazy(m_left); // we don't need to copy the child here
810 }
811 return new DataLazy(m_left->deepCopy()->getPtr(),m_right->deepCopy()->getPtr(),m_op);
812 }
813
814
815 DataTypes::ValueType::size_type
816 DataLazy::getLength() const
817 {
818 return m_length;
819 }
820
821
822 DataAbstract*
823 DataLazy::getSlice(const DataTypes::RegionType& region) const
824 {
825 throw DataException("getSlice - not implemented for Lazy objects.");
826 }
827
828 DataTypes::ValueType::size_type
829 DataLazy::getPointOffset(int sampleNo,
830 int dataPointNo) const
831 {
832 throw DataException("getPointOffset - not implemented for Lazy objects - yet.");
833 }
834
835 } // end namespace

  ViewVC Help
Powered by ViewVC 1.1.26