escript Revision_
DataVector.h
Go to the documentation of this file.
1
2/*****************************************************************************
3*
4* Copyright (c) 2003-2020 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-2017 by Centre for Geoscience Computing (GeoComp)
14* Development from 2019 by School of Earth and Environmental Sciences
15**
16*****************************************************************************/
17
18#ifndef __ESCRIPT_DATAVECTOR_H__
19#define __ESCRIPT_DATAVECTOR_H__
20
21#include "system_dep.h"
22#include "DataTypes.h"
23#include "Assert.h"
24#include "DataVectorAlt.h"
25#include "DataVectorTaipan.h"
26
27#ifdef ESYS_HAVE_BOOST_NUMPY
28#include <boost/python.hpp>
29#include <boost/python/numpy.hpp>
30#endif
31
32// ensure that nobody else tries to instantiate the complex version
34
35
36namespace escript {
37
38// Functions in DataTypes:: which manipulate DataVectors
39namespace DataTypes
40{
41
42 // This is the main version we had
43 //typedef DataVectorTaipan DataVector;
46
62 void
63 pointToStream(std::ostream& os, const RealVectorType::ElementType* data,const ShapeType& shape, int offset, bool needsep=true, const std::string& sep=",");
64
80 void
81 pointToStream(std::ostream& os, const CplxVectorType::ElementType* data,const ShapeType& shape, int offset, bool needsep=true, const std::string& sep=",");
82
83
84#ifdef ESYS_HAVE_BOOST_NUMPY
85
86 void
87 pointToNumpyArrayOld(boost::python::numpy::ndarray& dataArray, const RealVectorType::ElementType* data, const ShapeType& shape, long offset, long numsamples, long dpps, long numdata);
88
89
90 void
91 pointToNumpyArrayOld(boost::python::numpy::ndarray& dataArray, const CplxVectorType::ElementType* data, const ShapeType& shape, long offset, long numsamples, long dpps, long numdata);
92
93
105 void
106 pointToNumpyArray(boost::python::numpy::ndarray& dataArray, const RealVectorType::ElementType* data, const ShapeType& shape, int offset, int d, int index);
107
108
109 void
110 pointToNumpyArray(boost::python::numpy::ndarray& dataArray, const CplxVectorType::ElementType* data, const ShapeType& shape, int offset, int d, int index);
111#endif
112
121 std::string
122 pointToString(const RealVectorType& data,const ShapeType& shape, int offset, const std::string& prefix);
123
124
125 std::string
126 pointToString(const CplxVectorType& data,const ShapeType& shape, int offset, const std::string& prefix);
127
137 void copyPoint(RealVectorType& dest, vec_size_type doffset, vec_size_type nvals, const RealVectorType& src, vec_size_type soffset);
138
148 void copyPoint(CplxVectorType& dest, vec_size_type doffset, vec_size_type nvals, const CplxVectorType& src, vec_size_type soffset);
149
156
171 template <class VEC>
172 // ESCRIPT_DLL_API
173 void
174 copySlice(VEC& left,
175 const ShapeType& leftShape,
176 typename VEC::size_type leftOffset,
177 const VEC& other,
178 const ShapeType& otherShape,
179 typename VEC::size_type otherOffset,
180 const RegionLoopRangeType& region)
181 {
182 //
183 // Make sure vectors are not empty
184
185 ESYS_ASSERT(!left.size()==!1, "left data is empty."); //Note: !1=0, but clang returns an error if the rhs is 0 here
186 ESYS_ASSERT(!other.size()==!1, "other data is empty.");
187
188 //
189 // Check the vector to be sliced from is compatible with the region to be sliced,
190 // and that the region to be sliced is compatible with this vector:
191 ESYS_ASSERT(checkOffset(leftOffset,left.size(),noValues(leftShape)),
192 "offset incompatible with this vector.");
193 ESYS_ASSERT(otherOffset+noValues(leftShape)<=other.size(),
194 "offset incompatible with other vector.");
195
196 ESYS_ASSERT(getRank(otherShape)==region.size(),
197 "slice not same rank as vector to be sliced from.");
198
200 "slice shape not compatible shape for this vector.");
201
202 //
203 // copy the values in the specified region of the other vector into this vector
204
205 // the following loops cannot be parallelised due to the numCopy counter
206 int numCopy=0;
207
208 switch (region.size()) {
209 case 0:
210 /* this case should never be encountered,
211 as python will never pass us an empty region.
212 here for completeness only, allows slicing of a scalar */
213// (*m_data)[leftOffset+numCopy]=(*other.m_data)[otherOffset+other.relIndex()];
214
215 left[leftOffset+numCopy]=other[otherOffset];
216 numCopy++;
217 break;
218 case 1:
219 for (int i=region[0].first;i<region[0].second;i++) {
220 left[leftOffset+numCopy]=other[otherOffset+getRelIndex(otherShape,i)];
221 numCopy++;
222 }
223 break;
224 case 2:
225 for (int j=region[1].first;j<region[1].second;j++) {
226 for (int i=region[0].first;i<region[0].second;i++) {
227/* (*m_data)[leftOffset+numCopy]=(*other.m_data)[otherOffset+other.relIndex(i,j)];*/
228 left[leftOffset+numCopy]=other[otherOffset+getRelIndex(otherShape,i,j)];
229 numCopy++;
230 }
231 }
232 break;
233 case 3:
234 for (int k=region[2].first;k<region[2].second;k++) {
235 for (int j=region[1].first;j<region[1].second;j++) {
236 for (int i=region[0].first;i<region[0].second;i++) {
237// (*m_data)[leftOffset+numCopy]=(*other.m_data)[otherOffset+other.relIndex(i,j,k)];
238 left[leftOffset+numCopy]=other[otherOffset+getRelIndex(otherShape,i,j,k)];
239 numCopy++;
240 }
241 }
242 }
243 break;
244 case 4:
245 for (int l=region[3].first;l<region[3].second;l++) {
246 for (int k=region[2].first;k<region[2].second;k++) {
247 for (int j=region[1].first;j<region[1].second;j++) {
248 for (int i=region[0].first;i<region[0].second;i++) {
249/* (*m_data)[leftOffset+numCopy]=(*other.m_data)[otherOffset+other.relIndex(i,j,k,l)];*/
250 left[leftOffset+numCopy]=other[otherOffset+getRelIndex(otherShape,i,j,k,l)];
251 numCopy++;
252 }
253 }
254 }
255 }
256 break;
257 default:
258 std::stringstream mess;
259 mess << "Error - (copySlice) Invalid slice region rank: " << region.size();
260 throw DataException(mess.str());
261 }
262 }
263
278 template<typename VEC>
279 // ESCRIPT_DLL_API
280 void
281 copySliceFrom(VEC& left,
282 const ShapeType& leftShape,
283 typename VEC::size_type leftOffset,
284 const VEC& other,
285 const ShapeType& otherShape,
286 typename VEC::size_type otherOffset,
287 const RegionLoopRangeType& region)
288 {
289 //
290 // Make sure vectors are not empty
291
292 ESYS_ASSERT(left.size()!=0, "this vector is empty.");
293 ESYS_ASSERT(other.size()!=0, "other vector is empty.");
294
295 //
296 // Check this vector is compatible with the region to be sliced,
297 // and that the region to be sliced is compatible with the other vector:
298
299 ESYS_ASSERT(checkOffset(otherOffset,other.size(),noValues(otherShape)),
300 "offset incompatible with other vector.");
301 ESYS_ASSERT(leftOffset+noValues(otherShape)<=left.size(),
302 "offset incompatible with this vector.");
303
304 ESYS_ASSERT(getRank(leftShape)==region.size(),
305 "slice not same rank as this vector.");
306
307 ESYS_ASSERT(getRank(otherShape)==0 || noValues(otherShape)==noValues(getResultSliceShape(region)),
308 "slice shape not compatible shape for other vector.");
309
310 //
311 // copy the values in the other vector into the specified region of this vector
312
313 // allow for case where other vector is a scalar
314 if (getRank(otherShape)==0) {
315
316 // the following loops cannot be parallelised due to the numCopy counter
317 int numCopy=0;
318
319 switch (region.size()) {
320 case 0:
321 /* this case should never be encountered,
322 as python will never pass us an empty region.
323 here for completeness only, allows slicing of a scalar */
324 //(*m_data)[leftOffset+relIndex()]=(*other.m_data)[otherOffset];
325 left[leftOffset]=other[otherOffset];
326 numCopy++;
327 break;
328 case 1:
329 for (int i=region[0].first;i<region[0].second;i++) {
330 left[leftOffset+getRelIndex(leftShape,i)]=other[otherOffset];
331 numCopy++;
332 }
333 break;
334 case 2:
335 for (int j=region[1].first;j<region[1].second;j++) {
336 for (int i=region[0].first;i<region[0].second;i++) {
337 left[leftOffset+getRelIndex(leftShape,i,j)]=other[otherOffset];
338 numCopy++;
339 }
340 }
341 break;
342 case 3:
343 for (int k=region[2].first;k<region[2].second;k++) {
344 for (int j=region[1].first;j<region[1].second;j++) {
345 for (int i=region[0].first;i<region[0].second;i++) {
346 left[leftOffset+getRelIndex(leftShape,i,j,k)]=other[otherOffset];
347 numCopy++;
348 }
349 }
350 }
351 break;
352 case 4:
353 for (int l=region[3].first;l<region[3].second;l++) {
354 for (int k=region[2].first;k<region[2].second;k++) {
355 for (int j=region[1].first;j<region[1].second;j++) {
356 for (int i=region[0].first;i<region[0].second;i++) {
357 left[leftOffset+getRelIndex(leftShape,i,j,k,l)]=other[otherOffset];
358 numCopy++;
359 }
360 }
361 }
362 }
363 break;
364 default:
365 std::stringstream mess;
366 mess << "Error - (copySliceFrom) Invalid slice region rank: " << region.size();
367 throw DataException(mess.str());
368 }
369
370 } else {
371
372 // the following loops cannot be parallelised due to the numCopy counter
373 int numCopy=0;
374
375 switch (region.size()) {
376 case 0:
377 /* this case should never be encountered,
378 as python will never pass us an empty region.
379 here for completeness only, allows slicing of a scalar */
380 //(*m_data)[leftOffset+relIndex()]=(*other.m_data)[otherOffset+numCopy];
381 left[leftOffset]=other[otherOffset+numCopy];
382 numCopy++;
383 break;
384 case 1:
385 for (int i=region[0].first;i<region[0].second;i++) {
386 left[leftOffset+getRelIndex(leftShape,i)]=other[otherOffset+numCopy];
387 numCopy++;
388 }
389 break;
390 case 2:
391 for (int j=region[1].first;j<region[1].second;j++) {
392 for (int i=region[0].first;i<region[0].second;i++) {
393 left[leftOffset+getRelIndex(leftShape,i,j)]=other[otherOffset+numCopy];
394 numCopy++;
395 }
396 }
397 break;
398 case 3:
399 for (int k=region[2].first;k<region[2].second;k++) {
400 for (int j=region[1].first;j<region[1].second;j++) {
401 for (int i=region[0].first;i<region[0].second;i++) {
402 left[leftOffset+getRelIndex(leftShape,i,j,k)]=other[otherOffset+numCopy];
403 numCopy++;
404 }
405 }
406 }
407 break;
408 case 4:
409 for (int l=region[3].first;l<region[3].second;l++) {
410 for (int k=region[2].first;k<region[2].second;k++) {
411 for (int j=region[1].first;j<region[1].second;j++) {
412 for (int i=region[0].first;i<region[0].second;i++) {
413 left[leftOffset+getRelIndex(leftShape,i,j,k,l)]=other[otherOffset+numCopy];
414 numCopy++;
415 }
416 }
417 }
418 }
419 break;
420 default:
421 std::stringstream mess;
422 mess << "Error - (copySliceFrom) Invalid slice region rank: " << region.size();
423 throw DataException(mess.str());
424 }
425
426 }
427
428 }
429}
430
431} // end of namespace
432
433#endif // __ESCRIPT_DATAVECTOR_H__
434
#define ESYS_ASSERT(a, b)
EsysAssert is a MACRO that will throw an exception if the boolean condition specified is false.
Definition Assert.h:79
Definition DataException.h:28
Definition DataVectorAlt.h:37
T ElementType
Definition DataVectorAlt.h:43
void pointToStream(std::ostream &os, const RealVectorType::ElementType *data, const ShapeType &shape, int offset, bool needsep=true, const std::string &sep=",")
Display a single value (with the specified shape) from the data.
void fillComplexFromReal(const RealVectorType &r, CplxVectorType &c)
copy data from a real vector to a complex vector The complex vector will be resized as needed and any...
Definition DataVector.cpp:585
int noValues(const ShapeType &shape)
Calculate the number of values in a datapoint with the given shape.
Definition DataTypes.cpp:91
bool checkOffset(vec_size_type offset, int size, int noval)
Definition DataTypes.h:331
std::string pointToString(const RealVectorType &data, const ShapeType &shape, int offset, const std::string &prefix)
Display a single value (with the specified shape) from the data.
Definition DataVector.cpp:492
int getRank(const DataTypes::ShapeType &shape)
Return the rank (number of dimensions) of the given shape.
Definition DataTypes.h:225
std::vector< std::pair< int, int > > RegionLoopRangeType
Definition DataTypes.h:46
std::vector< int > ShapeType
The shape of a single datapoint.
Definition DataTypes.h:44
void copySlice(VEC &left, const ShapeType &leftShape, typename VEC::size_type leftOffset, const VEC &other, const ShapeType &otherShape, typename VEC::size_type otherOffset, const RegionLoopRangeType &region)
Copy a data slice specified by the given region and offset from the "other" vector into the "left" ve...
Definition DataVector.h:174
DataTypes::ShapeType getResultSliceShape(const RegionType &region)
Determine the shape of the specified slice region.
Definition DataTypes.cpp:173
void copySliceFrom(VEC &left, const ShapeType &leftShape, typename VEC::size_type leftOffset, const VEC &other, const ShapeType &otherShape, typename VEC::size_type otherOffset, const RegionLoopRangeType &region)
Copy data into a slice specified by the given region and offset in the left vector from the other vec...
Definition DataVector.h:281
vec_size_type getRelIndex(const DataTypes::ShapeType &shape, vec_size_type i)
Compute the offset (in 1D vector) of a given subscript with a shape.
Definition DataTypes.h:240
void copyPoint(RealVectorType &dest, vec_size_type doffset, vec_size_type nvals, const RealVectorType &src, vec_size_type soffset)
Copy a point from one vector to another. Note: This version does not check to see if shapes are the s...
long vec_size_type
Definition DataTypes.h:49
Definition AbstractContinuousDomain.cpp:23