/[escript]/trunk/escriptcore/src/DataExpanded.cpp
ViewVC logotype

Contents of /trunk/escriptcore/src/DataExpanded.cpp

Parent Directory Parent Directory | Revision Log Revision Log


Revision 6144 - (show annotations)
Wed Apr 6 05:25:13 2016 UTC (2 years, 7 months ago) by caltinay
File size: 43556 byte(s)
last round of namespacing defines.

1
2 /*****************************************************************************
3 *
4 * Copyright (c) 2003-2016 by The University of Queensland
5 * http://www.uq.edu.au
6 *
7 * Primary Business: Queensland, Australia
8 * Licensed under the Apache License, version 2.0
9 * http://www.apache.org/licenses/LICENSE-2.0
10 *
11 * Development until 2012 by Earth Systems Science Computational Center (ESSCC)
12 * Development 2012-2013 by School of Earth Sciences
13 * Development from 2014 by Centre for Geoscience Computing (GeoComp)
14 *
15 *****************************************************************************/
16
17 #include "Data.h"
18 #include "DataConstant.h"
19 #include "DataException.h"
20 #include "DataExpanded.h"
21 #include "DataVectorOps.h"
22 #include "DataTagged.h"
23
24 #include <limits>
25
26 #ifdef ESYS_HAVE_NETCDF
27 #include <netcdfcpp.h>
28 #endif
29
30 using namespace std;
31 using namespace escript::DataTypes;
32
33 #ifdef SLOWSHARECHECK
34 #define CHECK_FOR_EX_WRITE do {\
35 if (isShared()) {\
36 std::ostringstream ss;\
37 ss << "Attempt to modify shared object. Line " << __LINE__ << " in "\
38 << __FILE__;\
39 abort();\
40 throw DataException(ss.str());\
41 }\
42 } while(0);
43 #else
44 #define CHECK_FOR_EX_WRITE
45 #endif
46
47 namespace escript {
48
49 DataExpanded::DataExpanded(const WrappedArray& value,
50 const FunctionSpace& what)
51 : parent(what,value.getShape())
52 {
53 // initialise the data array for this object
54 initialise(what.getNumSamples(),what.getNumDPPSample(), value.isComplex());
55 // copy the given value to every data point
56 copy(value);
57 }
58
59 DataExpanded::DataExpanded(const DataExpanded& other)
60 : parent(other.getFunctionSpace(), other.getShape()),
61 m_data_r(other.m_data_r), m_data_c(other.m_data_c)
62 {
63 m_iscompl=other.m_iscompl;
64 }
65
66 DataExpanded::DataExpanded(const DataConstant& other)
67 : parent(other.getFunctionSpace(), other.getShape())
68 {
69 // initialise the data array for this object
70 initialise(other.getNumSamples(),other.getNumDPPSample(), other.isComplex());
71 // DataConstant only has one value, copy this to every data point
72 copy(other);
73 }
74
75 DataExpanded::DataExpanded(const DataTagged& other)
76 : parent(other.getFunctionSpace(), other.getShape())
77 {
78 // initialise the data array for this object
79 initialise(other.getNumSamples(),other.getNumDPPSample(), other.isComplex());
80 // for each data point in this object, extract and copy the corresponding
81 // data value from the given DataTagged object
82 if (isComplex())
83 {
84 DataTypes::cplx_t dummy=0;
85 #pragma omp parallel for
86 for (int i=0; i<m_noSamples; i++) {
87 for (int j=0; j<m_noDataPointsPerSample; j++) {
88 try {
89 DataTypes::copyPoint(getTypedVectorRW(dummy), getPointOffset(i,j),
90 getNoValues(), other.getTypedVectorRO(dummy),
91 other.getPointOffset(i,j));
92 } catch (std::exception& e) {
93 cerr << e.what() << endl;
94 }
95 }
96 }
97
98
99 }
100 else
101 {
102 DataTypes::real_t dummy=0;
103 #pragma omp parallel for
104 for (int i=0; i<m_noSamples; i++) {
105 for (int j=0; j<m_noDataPointsPerSample; j++) {
106 try {
107 DataTypes::copyPoint(getTypedVectorRW(dummy), getPointOffset(i,j),
108 getNoValues(), other.getTypedVectorRO(dummy),
109 other.getPointOffset(i,j));
110 } catch (std::exception& e) {
111 cerr << e.what() << endl;
112 }
113 }
114 }
115 }
116 }
117
118 DataExpanded::DataExpanded(const DataExpanded& other,
119 const DataTypes::RegionType& region)
120 : parent(other.getFunctionSpace(),DataTypes::getResultSliceShape(region))
121 {
122 // initialise this Data object to the shape of the slice
123 initialise(other.getNumSamples(),other.getNumDPPSample(), other.isComplex());
124 // copy the data
125 DataTypes::RegionLoopRangeType region_loop_range =
126 DataTypes::getSliceRegionLoopRange(region);
127 if (isComplex())
128 {
129 DataTypes::cplx_t dummy=0;
130 #pragma omp parallel for
131 for (int i=0; i<m_noSamples; i++) {
132 for (int j=0; j<m_noDataPointsPerSample; j++) {
133 try {
134 DataTypes::copySlice(getTypedVectorRW(dummy), getShape(),
135 getPointOffset(i,j), other.getTypedVectorRO(dummy),
136 other.getShape(),
137 other.getPointOffset(i,j),
138 region_loop_range);
139 } catch (std::exception& e) {
140 cerr << e.what() << endl;
141 }
142 }
143 }
144 }
145 else
146 {
147 DataTypes::real_t dummy=0;
148 #pragma omp parallel for
149 for (int i=0; i<m_noSamples; i++) {
150 for (int j=0; j<m_noDataPointsPerSample; j++) {
151 try {
152 DataTypes::copySlice(getTypedVectorRW(dummy), getShape(),
153 getPointOffset(i,j), other.getTypedVectorRO(dummy),
154 other.getShape(),
155 other.getPointOffset(i,j),
156 region_loop_range);
157 } catch (std::exception& e) {
158 cerr << e.what() << endl;
159 }
160 }
161 }
162 }
163 }
164
165 DataExpanded::DataExpanded(const FunctionSpace& what,
166 const DataTypes::ShapeType &shape,
167 const DataTypes::RealVectorType &data)
168 : parent(what,shape)
169 {
170 ESYS_ASSERT(data.size()%getNoValues()==0,
171 "DataExpanded Constructor - size of supplied data is not a multiple of shape size.");
172
173 if (data.size() == getNoValues()) {
174 RealVectorType& vec=m_data_r;
175 // create the view of the data
176 initialise(what.getNumSamples(),what.getNumDPPSample(), false);
177 // now we copy this value to all elements
178 for (int i=0; i<getLength();) {
179 for (unsigned int j=0;j<getNoValues();++j,++i) {
180 vec[i]=data[j];
181 }
182 }
183 } else {
184 // copy the data in the correct format
185 m_data_r = data;
186 }
187 }
188
189 DataExpanded::DataExpanded(const FunctionSpace& what,
190 const DataTypes::ShapeType &shape,
191 const DataTypes::CplxVectorType &data)
192 : parent(what,shape)
193 {
194 ESYS_ASSERT(data.size()%getNoValues()==0,
195 "DataExpanded Constructor - size of supplied data is not a multiple of shape size.");
196
197 if (data.size() == getNoValues()) {
198 CplxVectorType& vec=m_data_c;
199 // create the view of the data
200 initialise(what.getNumSamples(),what.getNumDPPSample(), true);
201 // now we copy this value to all elements
202 for (int i=0; i<getLength();) {
203 for (unsigned int j=0;j<getNoValues();++j,++i) {
204 vec[i]=data[j];
205 }
206 }
207 } else {
208 // copy the data in the correct format
209 m_data_c = data;
210 }
211 }
212
213
214 DataExpanded::DataExpanded(const FunctionSpace& what,
215 const DataTypes::ShapeType &shape,
216 const DataTypes::real_t v)
217 : parent(what,shape)
218 {
219 initialise(what.getNumSamples(),what.getNumDPPSample(), false);
220 DataTypes::RealVectorType& vec=m_data_r;
221 // now we copy this value to all elements
222 const int L=getLength();
223 #pragma omp parallel for
224 for (int i=0; i<L; ++i) {
225 vec[i]=v;
226 }
227 }
228
229 DataExpanded::DataExpanded(const FunctionSpace& what,
230 const DataTypes::ShapeType &shape,
231 const DataTypes::cplx_t v)
232 : parent(what,shape)
233 {
234 initialise(what.getNumSamples(),what.getNumDPPSample(), true);
235 DataTypes::CplxVectorType& vec=m_data_c;
236
237 // now we copy this value to all elements
238 const int L=getLength();
239 #pragma omp parallel for
240 for (int i=0; i<L; ++i) {
241 vec[i]=v;
242 }
243 }
244
245
246 DataExpanded::~DataExpanded()
247 {
248 }
249
250 DataAbstract* DataExpanded::deepCopy() const
251 {
252 return new DataExpanded(*this);
253 }
254
255 DataAbstract* DataExpanded::getSlice(const DataTypes::RegionType& region) const
256 {
257 return new DataExpanded(*this,region);
258 }
259
260 void DataExpanded::setSlice(const DataAbstract* value,
261 const DataTypes::RegionType& region)
262 {
263 const DataExpanded* tempDataExp=dynamic_cast<const DataExpanded*>(value);
264 if (!tempDataExp)
265 throw DataException("Programming error - casting to DataExpanded.");
266
267 CHECK_FOR_EX_WRITE;
268
269 // get shape of slice
270 DataTypes::ShapeType shape(DataTypes::getResultSliceShape(region));
271 DataTypes::RegionLoopRangeType region_loop_range =
272 DataTypes::getSliceRegionLoopRange(region);
273 // check shape
274 if (getRank() != region.size())
275 throw DataException("Error - Invalid slice region.");
276
277 if (tempDataExp->getRank() > 0 && !DataTypes::checkShape(value->getShape(), shape))
278 throw DataException(DataTypes::createShapeErrorMessage(
279 "Error - Couldn't copy slice due to shape mismatch.",
280 shape, value->getShape()));
281
282 // copy the data from the slice into this object
283 DataTypes::RealVectorType& vec=getVectorRW();
284 const ShapeType& mshape=getShape();
285 const DataTypes::RealVectorType& tVec=tempDataExp->getVectorRO();
286 const ShapeType& tShape=tempDataExp->getShape();
287 #pragma omp parallel for
288 for (int i=0; i<m_noSamples; i++) {
289 for (int j=0; j<m_noDataPointsPerSample; j++) {
290 DataTypes::copySliceFrom(vec, mshape, getPointOffset(i,j), tVec,
291 tShape, tempDataExp->getPointOffset(i,j),
292 region_loop_range);
293 }
294 }
295 }
296
297 void DataExpanded::copy(const DataConstant& value)
298 {
299 ESYS_ASSERT(checkShape(getShape(), value.getShape()),
300 createShapeErrorMessage("Error - Couldn't copy due to shape mismatch.", value.getShape(), getShape()));
301 if (isComplex())
302 {
303 if (value.isComplex())
304 {
305 // copy a single value to every data point in this object
306 #pragma omp parallel for
307 for (int i=0; i<m_noSamples; i++) {
308 for (int j=0; j<m_noDataPointsPerSample; j++) {
309 DataTypes::copyPoint(getTypedVectorRW((cplx_t)(0)), getPointOffset(i,j),
310 getNoValues(), value.getTypedVectorRO((cplx_t)(0)), 0);
311 }
312 }
313 }
314 else // value is real
315 {
316 throw DataException("Programming error - DataExpanded::copy source and target must be the same complexity.");
317 }
318 }
319 else
320 {
321 if (value.isComplex())
322 {
323 throw DataException("Programming error - DataExpanded::copy source and target must be the same complexity.");
324 }
325 else
326 {
327 real_t dummy=0;
328 // copy a single value to every data point in this object
329 #pragma omp parallel for
330 for (int i=0; i<m_noSamples; i++) {
331 for (int j=0; j<m_noDataPointsPerSample; j++) {
332 DataTypes::copyPoint(getTypedVectorRW(dummy), getPointOffset(i,j),
333 getNoValues(), value.getTypedVectorRO(dummy), 0);
334 }
335 }
336 }
337 }
338
339 }
340
341 void DataExpanded::copy(const WrappedArray& value)
342 {
343 // check the input shape matches this shape
344 if (!DataTypes::checkShape(getShape(),value.getShape()))
345 throw DataException(DataTypes::createShapeErrorMessage(
346 "Error - (DataExpanded) Cannot copy due to shape mismatch.",
347 value.getShape(),getShape()));
348 getVectorRW().copyFromArray(value, getNumDPPSample()*getNumSamples());
349 }
350
351 void DataExpanded::initialise(int noSamples, int noDataPointsPerSample, bool cplx)
352 {
353 this->m_iscompl=cplx;
354 if (noSamples==0) //retain the default empty object
355 return;
356
357 if (cplx)
358 {
359 // resize data array to the required size
360 m_data_c.resize(noSamples*noDataPointsPerSample*getNoValues(), 0.0, noDataPointsPerSample*getNoValues());
361 }
362 else
363 {
364 // resize data array to the required size
365 m_data_r.resize(noSamples*noDataPointsPerSample*getNoValues(), 0.0, noDataPointsPerSample*getNoValues());
366 }
367 }
368
369 bool
370 DataExpanded::hasNaN() const
371 {
372 bool haveNaN=false;
373 if (isComplex())
374 {
375 #pragma omp parallel for
376 for (DataTypes::CplxVectorType::size_type i=0;i<m_data_c.size();++i)
377 {
378 if (std::isnan(m_data_c[i].real()) || std::isnan(m_data_c[i].imag()))
379 {
380 #pragma omp critical
381 {
382 haveNaN=true;
383 }
384 }
385 }
386 }
387 else
388 {
389 #pragma omp parallel for
390 for (DataTypes::RealVectorType::size_type i=0;i<m_data_r.size();++i)
391 {
392 if (std::isnan(m_data_r[i]))
393 {
394 #pragma omp critical
395 {
396 haveNaN=true;
397 }
398 }
399 }
400 }
401 return haveNaN;
402 }
403
404
405 void
406 DataExpanded::replaceNaN(DataTypes::real_t value) {
407 CHECK_FOR_EX_WRITE
408 if (isComplex())
409 {
410 #pragma omp parallel for
411 for (DataTypes::CplxVectorType::size_type i=0;i<m_data_c.size();++i)
412 {
413 if (std::isnan(m_data_c[i].real()) || std::isnan(m_data_c[i].imag()))
414 {
415 m_data_c[i] = value;
416 }
417 }
418 }
419 else
420 {
421 #pragma omp parallel for
422 for (DataTypes::RealVectorType::size_type i=0;i<m_data_r.size();++i)
423 {
424 if (std::isnan(m_data_r[i]))
425 {
426 m_data_r[i] = value;
427 }
428 }
429 }
430 }
431
432 void
433 DataExpanded::replaceNaN(DataTypes::cplx_t value) {
434 CHECK_FOR_EX_WRITE
435 if (isComplex())
436 {
437 #pragma omp parallel for
438 for (DataTypes::CplxVectorType::size_type i=0;i<m_data_c.size();++i)
439 {
440 if (std::isnan(m_data_c[i].real()) || std::isnan(m_data_c[i].imag()))
441 {
442 m_data_c[i] = value;
443 }
444 }
445 }
446 else
447 {
448 complicate();
449 replaceNaN(value);
450 }
451 }
452
453 string DataExpanded::toString() const
454 {
455 stringstream ss;
456 FunctionSpace fs=getFunctionSpace();
457
458 int offset=0;
459 for (int i=0; i<m_noSamples; i++) {
460 for (int j=0; j<m_noDataPointsPerSample; j++) {
461 offset = getPointOffset(i,j);
462 stringstream suffix;
463 suffix << "( id: " << i << ", ref: "
464 << fs.getReferenceIDOfSample(i) << ", pnt: " << j << ")";
465 ss << (isComplex()?
466 DataTypes::pointToString(getTypedVectorRO((cplx_t)0), getShape(), offset, suffix.str())
467 :
468 DataTypes::pointToString(getTypedVectorRO((real_t)0), getShape(), offset, suffix.str()));
469 if (!(i==(m_noSamples-1) && j==(m_noDataPointsPerSample-1))) {
470 ss << endl;
471 }
472 }
473 }
474 string result(ss.str());
475 if (result.empty())
476 return "(data contains no samples)\n";
477
478 return result;
479 }
480
481 DataTypes::RealVectorType::size_type DataExpanded::getPointOffset(int sampleNo,
482 int dataPointNo) const
483 {
484 DataTypes::RealVectorType::size_type blockSize=getNoValues();
485 ESYS_ASSERT((isComplex()?
486 ((sampleNo >= 0) && (dataPointNo >= 0) && (m_data_c.size() > 0))
487 :
488 ((sampleNo >= 0) && (dataPointNo >= 0) && (m_data_r.size() > 0))),
489 "(DataBlocks2D) Index value out of range.");
490 DataTypes::RealVectorType::size_type temp=(sampleNo*m_noDataPointsPerSample+dataPointNo)*blockSize;
491 ESYS_ASSERT((isComplex()?
492 (temp <= (m_data_c.size()-blockSize))
493 :
494 (temp <= (m_data_r.size()-blockSize))), "Index value out of range.");
495
496 return temp;
497 }
498
499
500 void DataExpanded::complicate()
501 {
502 if (!isComplex())
503 {
504 fillComplexFromReal(m_data_r, m_data_c);
505 this->m_iscompl=true;
506 m_data_r.resize(0,0,1);
507 }
508 }
509
510
511 DataTypes::RealVectorType::size_type DataExpanded::getLength() const
512 {
513 return m_data_r.size();
514 }
515
516 void DataExpanded::copyToDataPoint(int sampleNo, int dataPointNo, const DataTypes::cplx_t value)
517 {
518 if (!isComplex())
519 {
520 throw DataException("Programming error - attempt to set complex value on real data.");
521 }
522 CHECK_FOR_EX_WRITE;
523 // Get the number of samples and data-points per sample.
524 int numSamples = getNumSamples();
525 int numDataPointsPerSample = getNumDPPSample();
526 int dataPointRank = getRank();
527 ShapeType dataPointShape = getShape();
528 if (numSamples*numDataPointsPerSample > 0) {
529 //TODO: global error handling
530 if (sampleNo >= numSamples || sampleNo < 0)
531 throw DataException("DataExpanded::copyDataPoint: invalid sampleNo.");
532 if (dataPointNo >= numDataPointsPerSample || dataPointNo < 0)
533 throw DataException("DataExpanded::copyDataPoint: invalid dataPointNo.");
534
535 DataTypes::CplxVectorType::size_type offset = getPointOffset(sampleNo, dataPointNo);
536 DataTypes::CplxVectorType& vec = getTypedVectorRW(cplx_t(0));
537 if (dataPointRank==0) {
538 vec[offset] = value;
539 } else if (dataPointRank==1) {
540 for (int i=0; i<dataPointShape[0]; i++) {
541 vec[offset+i] = value;
542 }
543 } else if (dataPointRank==2) {
544 for (int i=0; i<dataPointShape[0]; i++) {
545 for (int j=0; j<dataPointShape[1]; j++) {
546 vec[offset+getRelIndex(dataPointShape,i,j)] = value;
547 }
548 }
549 } else if (dataPointRank==3) {
550 for (int i=0; i<dataPointShape[0]; i++) {
551 for (int j=0; j<dataPointShape[1]; j++) {
552 for (int k=0; k<dataPointShape[2]; k++) {
553 vec[offset+getRelIndex(dataPointShape,i,j,k)] = value;
554 }
555 }
556 }
557 } else if (dataPointRank==4) {
558 for (int i=0; i<dataPointShape[0]; i++) {
559 for (int j=0; j<dataPointShape[1]; j++) {
560 for (int k=0; k<dataPointShape[2]; k++) {
561 for (int l=0; l<dataPointShape[3]; l++) {
562 vec[offset+getRelIndex(dataPointShape,i,j,k,l)] = value;
563 }
564 }
565 }
566 }
567 }
568 }
569 }
570
571
572 void DataExpanded::copyToDataPoint(int sampleNo, int dataPointNo, const DataTypes::real_t value)
573 {
574 if (isComplex())
575 {
576 copyToDataPoint(sampleNo, dataPointNo, cplx_t(value));
577 return;
578 }
579 CHECK_FOR_EX_WRITE;
580 // Get the number of samples and data-points per sample.
581 int numSamples = getNumSamples();
582 int numDataPointsPerSample = getNumDPPSample();
583 int dataPointRank = getRank();
584 ShapeType dataPointShape = getShape();
585 if (numSamples*numDataPointsPerSample > 0) {
586 //TODO: global error handling
587 if (sampleNo >= numSamples || sampleNo < 0)
588 throw DataException("DataExpanded::copyDataPoint: invalid sampleNo.");
589 if (dataPointNo >= numDataPointsPerSample || dataPointNo < 0)
590 throw DataException("DataExpanded::copyDataPoint: invalid dataPointNo.");
591
592 DataTypes::RealVectorType::size_type offset = getPointOffset(sampleNo, dataPointNo);
593 DataTypes::RealVectorType& vec = getVectorRW();
594 if (dataPointRank==0) {
595 vec[offset] = value;
596 } else if (dataPointRank==1) {
597 for (int i=0; i<dataPointShape[0]; i++) {
598 vec[offset+i] = value;
599 }
600 } else if (dataPointRank==2) {
601 for (int i=0; i<dataPointShape[0]; i++) {
602 for (int j=0; j<dataPointShape[1]; j++) {
603 vec[offset+getRelIndex(dataPointShape,i,j)] = value;
604 }
605 }
606 } else if (dataPointRank==3) {
607 for (int i=0; i<dataPointShape[0]; i++) {
608 for (int j=0; j<dataPointShape[1]; j++) {
609 for (int k=0; k<dataPointShape[2]; k++) {
610 vec[offset+getRelIndex(dataPointShape,i,j,k)] = value;
611 }
612 }
613 }
614 } else if (dataPointRank==4) {
615 for (int i=0; i<dataPointShape[0]; i++) {
616 for (int j=0; j<dataPointShape[1]; j++) {
617 for (int k=0; k<dataPointShape[2]; k++) {
618 for (int l=0; l<dataPointShape[3]; l++) {
619 vec[offset+getRelIndex(dataPointShape,i,j,k,l)] = value;
620 }
621 }
622 }
623 }
624 }
625 }
626 }
627
628 void DataExpanded::copyToDataPoint(int sampleNo, int dataPointNo,
629 const WrappedArray& value)
630 {
631 CHECK_FOR_EX_WRITE;
632 // Get the number of samples and data-points per sample
633 int numSamples = getNumSamples();
634 int numDataPointsPerSample = getNumDPPSample();
635 // check rank
636 if (value.getRank() != getRank())
637 throw DataException("Rank of value does not match Data object rank");
638
639 if (numSamples*numDataPointsPerSample > 0) {
640 //TODO: global error handling
641 if (sampleNo >= numSamples || sampleNo < 0)
642 throw DataException("DataExpanded::copyDataPoint: invalid sampleNo.");
643 if (dataPointNo >= numDataPointsPerSample || dataPointNo < 0)
644 throw DataException("DataExpanded::copyDataPoint: invalid dataPointNoInSample.");
645
646 if (isComplex())
647 {
648 DataTypes::CplxVectorType::size_type offset = getPointOffset(sampleNo, dataPointNo);
649 DataTypes::CplxVectorType& vec = getTypedVectorRW(cplx_t(0));
650 vec.copyFromArrayToOffset(value, offset, 1);
651 }
652 else
653 {
654 DataTypes::RealVectorType::size_type offset = getPointOffset(sampleNo, dataPointNo);
655 DataTypes::RealVectorType& vec = getTypedVectorRW(real_t(0));
656 vec.copyFromArrayToOffset(value, offset, 1);
657 }
658 }
659 }
660
661 void DataExpanded::symmetric(DataAbstract* ev)
662 {
663 const int numSamples = getNumSamples();
664 const int numDataPointsPerSample = getNumDPPSample();
665 DataExpanded* temp_ev = dynamic_cast<DataExpanded*>(ev);
666 if (!temp_ev)
667 throw DataException("DataExpanded::symmetric: casting to DataExpanded failed (probably a programming error).");
668
669 const ShapeType& shape = getShape();
670 const ShapeType& evShape = temp_ev->getShape();
671 if (isComplex())
672 {
673 const DataTypes::CplxVectorType& vec = getTypedVectorRO((DataTypes::cplx_t)0);
674 DataTypes::CplxVectorType& evVec = temp_ev->getTypedVectorRW((DataTypes::cplx_t)0);
675 #pragma omp parallel for
676 for (int sampleNo = 0; sampleNo < numSamples; sampleNo++) {
677 for (int dataPointNo = 0; dataPointNo < numDataPointsPerSample; dataPointNo++) {
678 escript::symmetric(vec, shape,
679 getPointOffset(sampleNo,dataPointNo), evVec, evShape,
680 ev->getPointOffset(sampleNo,dataPointNo));
681 }
682 }
683 }
684 else
685 {
686 const DataTypes::RealVectorType& vec = getTypedVectorRO(0.0);
687 DataTypes::RealVectorType& evVec = temp_ev->getTypedVectorRW(0.0);
688 #pragma omp parallel for
689 for (int sampleNo = 0; sampleNo < numSamples; sampleNo++) {
690 for (int dataPointNo = 0; dataPointNo < numDataPointsPerSample; dataPointNo++) {
691 escript::symmetric(vec, shape,
692 getPointOffset(sampleNo,dataPointNo), evVec, evShape,
693 ev->getPointOffset(sampleNo,dataPointNo));
694 }
695 }
696 }
697 }
698
699 void DataExpanded::antisymmetric(DataAbstract* ev)
700 {
701 const int numSamples = getNumSamples();
702 const int numDataPointsPerSample = getNumDPPSample();
703 DataExpanded* temp_ev=dynamic_cast<DataExpanded*>(ev);
704 if (!temp_ev)
705 throw DataException("DataExpanded::antisymmetric: casting to DataExpanded failed (probably a programming error).");
706
707 const ShapeType& shape = getShape();
708 const ShapeType& evShape = temp_ev->getShape();
709 if (isComplex())
710 {
711 const DataTypes::CplxVectorType& vec = getTypedVectorRO((DataTypes::cplx_t)0);
712 DataTypes::CplxVectorType& evVec = temp_ev->getTypedVectorRW((DataTypes::cplx_t)0);
713 #pragma omp parallel for
714 for (int sampleNo = 0; sampleNo < numSamples; sampleNo++)
715 {
716 for (int dataPointNo = 0; dataPointNo < numDataPointsPerSample; dataPointNo++)
717 {
718 escript::antisymmetric(vec, shape,
719 getPointOffset(sampleNo,dataPointNo), evVec, evShape,
720 ev->getPointOffset(sampleNo,dataPointNo));
721 }
722 }
723 }
724 else
725 {
726 const DataTypes::RealVectorType& vec = getTypedVectorRO(0.0);
727 DataTypes::RealVectorType& evVec = temp_ev->getTypedVectorRW(0.0);
728 #pragma omp parallel for
729 for (int sampleNo = 0; sampleNo < numSamples; sampleNo++) {
730 for (int dataPointNo = 0; dataPointNo < numDataPointsPerSample; dataPointNo++)
731 {
732 escript::antisymmetric(vec, shape,
733 getPointOffset(sampleNo,dataPointNo), evVec, evShape,
734 ev->getPointOffset(sampleNo,dataPointNo));
735 }
736 }
737 }
738 }
739
740
741 void DataExpanded::hermitian(DataAbstract* ev)
742 {
743 const int numSamples = getNumSamples();
744 const int numDataPointsPerSample = getNumDPPSample();
745 DataExpanded* temp_ev = dynamic_cast<DataExpanded*>(ev);
746 if (!temp_ev)
747 throw DataException("DataExpanded::hermitian: casting to DataExpanded failed (probably a programming error).");
748 if (!isComplex() || !temp_ev->isComplex())
749 {
750 throw DataException("DataExpanded::hermitian: do not call this method with real data");
751 }
752 const ShapeType& shape = getShape();
753 const ShapeType& evShape = temp_ev->getShape();
754 const DataTypes::CplxVectorType& vec = getTypedVectorRO((DataTypes::cplx_t)0);
755 DataTypes::CplxVectorType& evVec = temp_ev->getTypedVectorRW((DataTypes::cplx_t)0);
756 #pragma omp parallel for
757 for (int sampleNo = 0; sampleNo < numSamples; sampleNo++) {
758 for (int dataPointNo = 0; dataPointNo < numDataPointsPerSample; dataPointNo++) {
759 escript::hermitian(vec, shape,
760 getPointOffset(sampleNo,dataPointNo), evVec, evShape,
761 ev->getPointOffset(sampleNo,dataPointNo));
762 }
763 }
764 }
765
766 void DataExpanded::antihermitian(DataAbstract* ev)
767 {
768 const int numSamples = getNumSamples();
769 const int numDataPointsPerSample = getNumDPPSample();
770 DataExpanded* temp_ev=dynamic_cast<DataExpanded*>(ev);
771 if (!temp_ev)
772 throw DataException("DataExpanded::antihermitian: casting to DataExpanded failed (probably a programming error).");
773 if (!isComplex() || !temp_ev->isComplex())
774 {
775 throw DataException("DataExpanded::antihermitian: do not call this method with real data");
776 }
777 const ShapeType& shape = getShape();
778 const ShapeType& evShape = temp_ev->getShape();
779 const DataTypes::CplxVectorType& vec = getTypedVectorRO((DataTypes::cplx_t)0);
780 DataTypes::CplxVectorType& evVec = temp_ev->getTypedVectorRW((DataTypes::cplx_t)0);
781 #pragma omp parallel for
782 for (int sampleNo = 0; sampleNo < numSamples; sampleNo++)
783 {
784 for (int dataPointNo = 0; dataPointNo < numDataPointsPerSample; dataPointNo++)
785 {
786 escript::antihermitian(vec, shape,
787 getPointOffset(sampleNo,dataPointNo), evVec, evShape,
788 ev->getPointOffset(sampleNo,dataPointNo));
789 }
790 }
791 }
792
793
794
795 void DataExpanded::trace(DataAbstract* ev, int axis_offset)
796 {
797 const int numSamples = getNumSamples();
798 const int numDataPointsPerSample = getNumDPPSample();
799 DataExpanded* temp_ev = dynamic_cast<DataExpanded*>(ev);
800 if (!temp_ev)
801 throw DataException("DataExpanded::trace: casting to DataExpanded failed (probably a programming error).");
802 const ShapeType& shape=getShape();
803 const ShapeType& evShape=temp_ev->getShape();
804 if (isComplex())
805 {
806 const DataTypes::CplxVectorType& vec=getVectorROC();
807 DataTypes::CplxVectorType& evVec=temp_ev->getVectorRWC();
808
809 #pragma omp parallel for
810 for (int sampleNo = 0; sampleNo < numSamples; sampleNo++) {
811 for (int dataPointNo = 0; dataPointNo < numDataPointsPerSample; dataPointNo++) {
812 escript::trace(vec, shape, getPointOffset(sampleNo,dataPointNo),
813 evVec, evShape, ev->getPointOffset(sampleNo,dataPointNo),
814 axis_offset);
815 }
816 }
817 }
818 else
819 {
820 const DataTypes::RealVectorType& vec=getVectorRO();
821 DataTypes::RealVectorType& evVec=temp_ev->getVectorRW();
822
823 #pragma omp parallel for
824 for (int sampleNo = 0; sampleNo < numSamples; sampleNo++) {
825 for (int dataPointNo = 0; dataPointNo < numDataPointsPerSample; dataPointNo++) {
826 escript::trace(vec, shape, getPointOffset(sampleNo,dataPointNo),
827 evVec, evShape, ev->getPointOffset(sampleNo,dataPointNo),
828 axis_offset);
829 }
830 }
831 }
832 }
833
834 void DataExpanded::transpose(DataAbstract* ev, int axis_offset)
835 {
836 const int numSamples = getNumSamples();
837 const int numDataPointsPerSample = getNumDPPSample();
838 DataExpanded* temp_ev=dynamic_cast<DataExpanded*>(ev);
839 if (!temp_ev)
840 throw DataException("DataExpanded::transpose: casting to DataExpanded failed (probably a programming error).");
841 const ShapeType& shape = getShape();
842 if (isComplex())
843 {
844 const DataTypes::CplxVectorType& vec = getVectorROC();
845 DataTypes::CplxVectorType& evVec = temp_ev->getVectorRWC();
846 const ShapeType& evShape = temp_ev->getShape();
847 #pragma omp parallel for
848 for (int sampleNo = 0; sampleNo < numSamples; sampleNo++) {
849 for (int dataPointNo = 0; dataPointNo < numDataPointsPerSample; dataPointNo++) {
850 escript::transpose(vec, shape,
851 getPointOffset(sampleNo,dataPointNo), evVec, evShape,
852 ev->getPointOffset(sampleNo,dataPointNo), axis_offset);
853 }
854 }
855 }
856 else
857 {
858 const DataTypes::RealVectorType& vec = getVectorRO();
859 DataTypes::RealVectorType& evVec = temp_ev->getVectorRW();
860 const ShapeType& evShape = temp_ev->getShape();
861 #pragma omp parallel for
862 for (int sampleNo = 0; sampleNo < numSamples; sampleNo++) {
863 for (int dataPointNo = 0; dataPointNo < numDataPointsPerSample; dataPointNo++) {
864 escript::transpose(vec, shape,
865 getPointOffset(sampleNo,dataPointNo), evVec, evShape,
866 ev->getPointOffset(sampleNo,dataPointNo), axis_offset);
867 }
868 }
869 }
870 }
871
872 void DataExpanded::swapaxes(DataAbstract* ev, int axis0, int axis1)
873 {
874 const int numSamples = getNumSamples();
875 const int numDataPointsPerSample = getNumDPPSample();
876 DataExpanded* temp_ev=dynamic_cast<DataExpanded*>(ev);
877 if (!temp_ev)
878 throw DataException("Error - DataExpanded::swapaxes: casting to DataExpanded failed (probably a programming error).");
879 const ShapeType& shape=getShape();
880 const ShapeType& evShape=temp_ev->getShape();
881 if (isComplex())
882 {
883 const DataTypes::CplxVectorType& vec=getVectorROC();
884 DataTypes::CplxVectorType& evVec=temp_ev->getVectorRWC();
885
886 #pragma omp parallel for
887 for (int sampleNo = 0; sampleNo < numSamples; sampleNo++) {
888 for (int dataPointNo = 0; dataPointNo < numDataPointsPerSample; dataPointNo++) {
889 escript::swapaxes(vec, shape,
890 getPointOffset(sampleNo,dataPointNo), evVec, evShape,
891 ev->getPointOffset(sampleNo,dataPointNo), axis0, axis1);
892 }
893 }
894 }
895 else
896 {
897 const DataTypes::RealVectorType& vec=getVectorRO();
898 DataTypes::RealVectorType& evVec=temp_ev->getVectorRW();
899
900 #pragma omp parallel for
901 for (int sampleNo = 0; sampleNo < numSamples; sampleNo++) {
902 for (int dataPointNo = 0; dataPointNo < numDataPointsPerSample; dataPointNo++) {
903 escript::swapaxes(vec, shape,
904 getPointOffset(sampleNo,dataPointNo), evVec, evShape,
905 ev->getPointOffset(sampleNo,dataPointNo), axis0, axis1);
906 }
907 }
908 }
909 }
910
911 void DataExpanded::eigenvalues(DataAbstract* ev)
912 {
913 const int numSamples = getNumSamples();
914 const int numDataPointsPerSample = getNumDPPSample();
915 DataExpanded* temp_ev=dynamic_cast<DataExpanded*>(ev);
916 if (!temp_ev)
917 throw DataException("DataExpanded::eigenvalues: casting to DataExpanded failed (probably a programming error).");
918 const ShapeType& evShape=temp_ev->getShape();
919 const ShapeType& shape=getShape();
920 if (isComplex())
921 {
922 const DataTypes::CplxVectorType& vec=getVectorROC();
923 DataTypes::CplxVectorType& evVec=temp_ev->getVectorRWC();
924
925 #pragma omp parallel for
926 for (int sampleNo = 0; sampleNo < numSamples; sampleNo++) {
927 for (int dataPointNo = 0; dataPointNo < numDataPointsPerSample; dataPointNo++) {
928 escript::eigenvalues(vec, shape,
929 getPointOffset(sampleNo,dataPointNo), evVec, evShape,
930 ev->getPointOffset(sampleNo,dataPointNo));
931 }
932 }
933 }
934 else
935 {
936 const DataTypes::RealVectorType& vec=getVectorRO();
937 DataTypes::RealVectorType& evVec=temp_ev->getVectorRW();
938
939 #pragma omp parallel for
940 for (int sampleNo = 0; sampleNo < numSamples; sampleNo++) {
941 for (int dataPointNo = 0; dataPointNo < numDataPointsPerSample; dataPointNo++) {
942 escript::eigenvalues(vec, shape,
943 getPointOffset(sampleNo,dataPointNo), evVec, evShape,
944 ev->getPointOffset(sampleNo,dataPointNo));
945 }
946 }
947 }
948 }
949
950 void DataExpanded::eigenvalues_and_eigenvectors(DataAbstract* ev,
951 DataAbstract* V, double tol)
952 {
953 const int numSamples = getNumSamples();
954 const int numDataPointsPerSample = getNumDPPSample();
955 DataExpanded* temp_ev=dynamic_cast<DataExpanded*>(ev);
956 if (!temp_ev)
957 throw DataException("DataExpanded::eigenvalues_and_eigenvectors: casting to DataExpanded failed (probably a programming error).");
958
959 DataExpanded* temp_V=dynamic_cast<DataExpanded*>(V);
960 if (!temp_V)
961 throw DataException("DataExpanded::eigenvalues_and_eigenvectors: casting to DataExpanded failed (probably a programming error).");
962
963 const DataTypes::RealVectorType& vec = getVectorRO();
964 const ShapeType& shape = getShape();
965 DataTypes::RealVectorType& evVec = temp_ev->getVectorRW();
966 const ShapeType& evShape = temp_ev->getShape();
967 DataTypes::RealVectorType& VVec = temp_V->getVectorRW();
968 const ShapeType& VShape = temp_V->getShape();
969 #pragma omp parallel for
970 for (int sampleNo = 0; sampleNo < numSamples; sampleNo++) {
971 for (int dataPointNo = 0; dataPointNo < numDataPointsPerSample; dataPointNo++) {
972 escript::eigenvalues_and_eigenvectors(vec, shape,
973 getPointOffset(sampleNo,dataPointNo), evVec, evShape,
974 ev->getPointOffset(sampleNo,dataPointNo), VVec, VShape,
975 V->getPointOffset(sampleNo,dataPointNo), tol);
976 }
977 }
978 }
979
980
981 int DataExpanded::matrixInverse(DataAbstract* out) const
982 {
983 DataExpanded* temp=dynamic_cast<DataExpanded*>(out);
984 if (!temp)
985 throw DataException("DataExpanded::matrixInverse: casting to DataExpanded failed (probably a programming error).");
986
987 if (getRank() != 2)
988 throw DataException("DataExpanded::matrixInverse: input must be rank 2.");
989
990 const int numdpps=getNumDPPSample();
991 const int numSamples = getNumSamples();
992 const DataTypes::RealVectorType& vec=m_data_r;
993 int errcode=0;
994 #pragma omp parallel
995 {
996 int errorcode=0;
997 LapackInverseHelper h(getShape()[0]);
998 #pragma omp for
999 for (int sampleNo = 0; sampleNo < numSamples; sampleNo++)
1000 {
1001 // not sure I like all those virtual calls to getPointOffset
1002 DataTypes::RealVectorType::size_type offset=getPointOffset(sampleNo,0);
1003 int res=escript::matrix_inverse(vec, getShape(), offset,
1004 temp->getVectorRW(), temp->getShape(), offset, numdpps, h);
1005 if (res > errorcode) {
1006 errorcode=res;
1007 #pragma omp critical
1008 {
1009 // I'm not especially concerned which error gets reported
1010 // as long as one is
1011 errcode=errorcode;
1012 }
1013 }
1014 }
1015 }
1016 return errcode;
1017 }
1018
1019 void DataExpanded::setToZero()
1020 {
1021 // Q: Surely there is a more efficient way to do this????
1022 // Why is there no memset here? Parallel issues?
1023 // A: This ensures that memory is touched by the correct thread.
1024 CHECK_FOR_EX_WRITE;
1025 const int numSamples = getNumSamples();
1026 const int numDataPointsPerSample = getNumDPPSample();
1027 const DataTypes::RealVectorType::size_type n = getNoValues();
1028 #pragma omp parallel for
1029 for (int sampleNo = 0; sampleNo < numSamples; sampleNo++) {
1030 for (int dataPointNo = 0; dataPointNo < numDataPointsPerSample; dataPointNo++) {
1031 double* p = &m_data_r[getPointOffset(sampleNo,dataPointNo)];
1032 for (int i=0; i<n; ++i)
1033 p[i] = 0.;
1034 }
1035 }
1036 }
1037
1038 void DataExpanded::dump(const std::string fileName) const
1039 {
1040 #ifdef ESYS_HAVE_NETCDF
1041 const int ldims=2+DataTypes::maxRank;
1042 const NcDim* ncdims[ldims];
1043 NcVar *var, *ids;
1044 int rank = getRank();
1045 int type= getFunctionSpace().getTypeCode();
1046 int ndims =0;
1047 long dims[ldims];
1048 const double* d_ptr=&(m_data_r[0]);
1049 const DataTypes::ShapeType& shape = getShape();
1050 JMPI mpiInfo(getFunctionSpace().getDomain()->getMPI());
1051 const std::string newFileName(mpiInfo->appendRankToFileName(fileName));
1052 // netCDF error handler
1053 NcError err(NcError::verbose_nonfatal);
1054 NcFile dataFile(newFileName.c_str(), NcFile::Replace);
1055 if (!dataFile.is_valid())
1056 throw DataException("DataExpanded::dump: opening of netCDF file for output failed.");
1057 if (!dataFile.add_att("type_id",2))
1058 throw DataException("DataExpanded::dump: appending data type to netCDF file failed.");
1059 if (!dataFile.add_att("rank",rank))
1060 throw DataException("DataExpanded::dump: appending rank attribute to netCDF file failed.");
1061 if (!dataFile.add_att("function_space_type",type))
1062 throw DataException("DataExpanded::dump: appending function space attribute to netCDF file failed.");
1063 ndims=rank+2;
1064 if ( rank >0 ) {
1065 dims[0]=shape[0];
1066 if (! (ncdims[0] = dataFile.add_dim("d0",shape[0])) )
1067 throw DataException("DataExpanded::dump: appending ncdim 0 to netCDF file failed.");
1068 }
1069 if ( rank >1 ) {
1070 dims[1]=shape[1];
1071 if (! (ncdims[1] = dataFile.add_dim("d1",shape[1])) )
1072 throw DataException("DataExpanded::dump: appending ncdim 1 to netCDF file failed.");
1073 }
1074 if ( rank >2 ) {
1075 dims[2]=shape[2];
1076 if (! (ncdims[2] = dataFile.add_dim("d2", shape[2])) )
1077 throw DataException("DataExpanded::dump: appending ncdim 2 to netCDF file failed.");
1078 }
1079 if ( rank >3 ) {
1080 dims[3]=shape[3];
1081 if (! (ncdims[3] = dataFile.add_dim("d3", shape[3])) )
1082 throw DataException("DataExpanded::dump: appending ncdim 3 to netCDF file failed.");
1083 }
1084 dims[rank]=getFunctionSpace().getNumDataPointsPerSample();
1085 if (! (ncdims[rank] = dataFile.add_dim("num_data_points_per_sample", dims[rank])) )
1086 throw DataException("DataExpanded::dump: appending num_data_points_per_sample to netCDF file failed.");
1087 dims[rank+1]=getFunctionSpace().getNumSamples();
1088 if (! (ncdims[rank+1] = dataFile.add_dim("num_samples", dims[rank+1])) )
1089 throw DataException("DataExpanded::dump: appending num_sample to netCDF file failed.");
1090
1091 if (getFunctionSpace().getNumSamples() > 0) {
1092 if (! ( ids = dataFile.add_var("id", ncInt, ncdims[rank+1])) )
1093 throw DataException("DataExpanded::dump: appending reference id to netCDF file failed.");
1094 const dim_t* ids_p=getFunctionSpace().borrowSampleReferenceIDs();
1095 if (! (ids->put(ids_p,dims[rank+1])) )
1096 throw DataException("DataExpanded::dump: copy reference id to netCDF buffer failed.");
1097 if (! ( var = dataFile.add_var("data", ncDouble, ndims, ncdims)) )
1098 throw DataException("DataExpanded::dump: appending variable to netCDF file failed.");
1099 if (! (var->put(d_ptr,dims)) )
1100 throw DataException("DataExpanded::dump: copy data to netCDF buffer failed.");
1101 }
1102 #else
1103 throw DataException("DataExpanded::dump: not configured with netCDF. Please contact your installation manager.");
1104 #endif // ESYS_HAVE_NETCDF
1105 }
1106
1107 void DataExpanded::setTaggedValue(int tagKey,
1108 const DataTypes::ShapeType& pointshape,
1109 const DataTypes::RealVectorType& value,
1110 int dataOffset)
1111 {
1112 CHECK_FOR_EX_WRITE;
1113 if (isComplex())
1114 {
1115 CplxVectorType tv;
1116 fillComplexFromReal(value, tv);
1117 setTaggedValue(tagKey, pointshape, tv, dataOffset);
1118 return;
1119 }
1120 const int numSamples = getNumSamples();
1121 const int numDataPointsPerSample = getNumDPPSample();
1122 const DataTypes::RealVectorType::size_type n = getNoValues();
1123 const real_t* in = &value[0+dataOffset];
1124
1125 if (value.size() != n)
1126 throw DataException("DataExpanded::setTaggedValue: number of input values does not match number of values per data points.");
1127
1128 #pragma omp parallel for
1129 for (int sampleNo = 0; sampleNo < numSamples; sampleNo++) {
1130 if (getFunctionSpace().getTagFromSampleNo(sampleNo) == tagKey) {
1131 for (int dataPointNo = 0; dataPointNo < numDataPointsPerSample; dataPointNo++) {
1132 real_t* p = &m_data_r[getPointOffset(sampleNo,dataPointNo)];
1133 for (int i=0; i<n; ++i)
1134 p[i] = in[i];
1135 }
1136 }
1137 }
1138 }
1139
1140
1141 void DataExpanded::setTaggedValue(int tagKey,
1142 const DataTypes::ShapeType& pointshape,
1143 const DataTypes::CplxVectorType& value,
1144 int dataOffset)
1145 {
1146 CHECK_FOR_EX_WRITE;
1147 if (!isComplex())
1148 {
1149 throw DataException("Programming Error - Attempt to set a complex value on a real object.");
1150 }
1151 const int numSamples = getNumSamples();
1152 const int numDataPointsPerSample = getNumDPPSample();
1153 const DataTypes::CplxVectorType::size_type n = getNoValues();
1154 const DataTypes::cplx_t* in = &value[0+dataOffset];
1155
1156 if (value.size() != n)
1157 throw DataException("DataExpanded::setTaggedValue: number of input values does not match number of values per data points.");
1158
1159 #pragma omp parallel for
1160 for (int sampleNo = 0; sampleNo < numSamples; sampleNo++) {
1161 if (getFunctionSpace().getTagFromSampleNo(sampleNo) == tagKey) {
1162 for (int dataPointNo = 0; dataPointNo < numDataPointsPerSample; dataPointNo++) {
1163 cplx_t* p = &m_data_c[getPointOffset(sampleNo,dataPointNo)];
1164 for (int i=0; i<n; ++i)
1165 p[i] = in[i];
1166 }
1167 }
1168 }
1169 }
1170
1171
1172 void DataExpanded::reorderByReferenceIDs(dim_t *reference_ids)
1173 {
1174 CHECK_FOR_EX_WRITE;
1175 const int numSamples = getNumSamples();
1176 const DataTypes::RealVectorType::size_type n = getNoValues() * getNumDPPSample();
1177 FunctionSpace fs=getFunctionSpace();
1178
1179 for (int sampleNo = 0; sampleNo < numSamples; sampleNo++) {
1180 const index_t id_in = reference_ids[sampleNo];
1181 const index_t id = fs.getReferenceIDOfSample(sampleNo);
1182 if (id != id_in) {
1183 bool matched=false;
1184 for (int sampleNo2 = sampleNo+1; sampleNo2 < numSamples; sampleNo2++) {
1185 if (id == reference_ids[sampleNo2]) {
1186 double* p = &m_data_r[getPointOffset(sampleNo,0)];
1187 double* p2 = &m_data_r[getPointOffset(sampleNo2,0)];
1188 for (int i=0; i<n; i++) {
1189 const double rtmp=p[i];
1190 p[i] = p2[i];
1191 p2[i] = rtmp;
1192 }
1193 reference_ids[sampleNo]=id;
1194 reference_ids[sampleNo2]=id_in;
1195 matched=true;
1196 break;
1197 }
1198 }
1199 if (!matched)
1200 throw DataException("DataExpanded::reorderByReferenceIDs: unable to reorder sample data by reference ids");
1201 }
1202 }
1203 }
1204
1205 DataTypes::RealVectorType& DataExpanded::getVectorRW()
1206 {
1207 CHECK_FOR_EX_WRITE;
1208 return m_data_r;
1209 }
1210
1211 const DataTypes::RealVectorType& DataExpanded::getVectorRO() const
1212 {
1213 return m_data_r;
1214 }
1215
1216 DataTypes::CplxVectorType& DataExpanded::getVectorRWC()
1217 {
1218 CHECK_FOR_EX_WRITE;
1219 return m_data_c;
1220 }
1221
1222 const DataTypes::CplxVectorType& DataExpanded::getVectorROC() const
1223 {
1224 return m_data_c;
1225 }
1226
1227 DataTypes::RealVectorType& DataExpanded::getTypedVectorRW(DataTypes::real_t dummypar)
1228 {
1229 CHECK_FOR_EX_WRITE;
1230 return m_data_r;
1231 }
1232
1233 const DataTypes::RealVectorType& DataExpanded::getTypedVectorRO(DataTypes::real_t dummypar) const
1234 {
1235 return m_data_r;
1236 }
1237
1238 DataTypes::CplxVectorType& DataExpanded::getTypedVectorRW(DataTypes::cplx_t dummypar)
1239 {
1240 CHECK_FOR_EX_WRITE;
1241 return m_data_c;
1242 }
1243
1244 const DataTypes::CplxVectorType& DataExpanded::getTypedVectorRO(DataTypes::cplx_t dummypar) const
1245 {
1246 return m_data_c;
1247 }
1248
1249
1250 //void DataExpanded::randomFill(long seed)
1251 //{
1252 // CHECK_FOR_EX_WRITE;
1253 //
1254 // DataVector& dv=getVectorRW();
1255 // const size_t dvsize=dv.size();
1256 // randomFillArray(seed, &(dv[0]), dvsize);
1257 //}
1258
1259 } // end of namespace
1260

Properties

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

  ViewVC Help
Powered by ViewVC 1.1.26