/[escript]/branches/arrayview_from_1695_trunk/escript/src/DataConstant.cpp
ViewVC logotype

Contents of /branches/arrayview_from_1695_trunk/escript/src/DataConstant.cpp

Parent Directory Parent Directory | Revision Log Revision Log


Revision 1773 - (show annotations)
Tue Sep 9 02:52:26 2008 UTC (11 years ago) by jfenwick
File size: 12650 byte(s)
Branch commit

Fixed some bugs.
Code now passes all_tests.
Next step is to test under OMP and MPI


1
2 /* $Id$ */
3
4 /*******************************************************
5 *
6 * Copyright 2003-2007 by ACceSS MNRF
7 * Copyright 2007 by University of Queensland
8 *
9 * http://esscc.uq.edu.au
10 * Primary Business: Queensland, Australia
11 * Licensed under the Open Software License version 3.0
12 * http://www.opensource.org/licenses/osl-3.0.php
13 *
14 *******************************************************/
15
16 #include "DataConstant.h"
17 #include "DataException.h"
18 #include "esysUtils/EsysAssert.h"
19
20 #include <iostream>
21 #include <boost/python/extract.hpp>
22 #ifdef USE_NETCDF
23 #include <netcdfcpp.h>
24 #endif
25
26 #include <boost/python/extract.hpp>
27 #include "DataMaths.h"
28
29 using namespace std;
30 using namespace boost::python;
31
32 namespace escript {
33
34 DataConstant::DataConstant(const boost::python::numeric::array& value,
35 const FunctionSpace& what)
36 : DataAbstract(what,DataTypes::shapeFromNumArray(value))
37 {
38 // // extract the shape of the numarray
39 // DataTypes::ShapeType tempShape;
40 // for (int i=0; i < value.getrank(); i++) {
41 // tempShape.push_back(extract<int>(value.getshape()[i]));
42 // }
43
44 // get the space for the data vector
45 // int len = getNoValues();
46 // DataVector temp_data(len, 0.0, len);
47 // DataArrayView temp_dataView(temp_data, tempShape);
48 // temp_dataView.copy(value);
49
50 m_data.copyFromNumArray(value);
51 //
52
53 // copy the data in the correct format
54 // m_data=temp_data;
55 //
56 // create the view of the data
57 // DataArrayView tempView(m_data,temp_dataView.getShape());
58 // setPointDataView(tempView);
59 }
60
61 // DataConstant::DataConstant(const DataArrayView& value,
62 // const FunctionSpace& what)
63 // : DataAbstract(what)
64 // {
65 // //
66 // // copy the data in the correct format
67 // m_data=value.getData();
68 // //
69 // // create the view of the data
70 // DataArrayView tempView(m_data,value.getShape());
71 // setPointDataView(tempView);
72 // }
73
74 DataConstant::DataConstant(const DataConstant& other)
75 : DataAbstract(other.getFunctionSpace(),other.getShape())
76 { //
77 // copy the data in the correct format
78 m_data=other.m_data;
79 //
80 // // create the view of the data
81 // DataArrayView tempView(m_data,other.getPointDataView().getShape());
82 // setPointDataView(tempView);
83 }
84
85 DataConstant::DataConstant(const DataConstant& other,
86 const DataTypes::RegionType& region)
87 : DataAbstract(other.getFunctionSpace(),DataTypes::getResultSliceShape(region))
88 {
89 //
90 // get the shape of the slice to copy from
91 // DataTypes::ShapeType shape(DataTypes::getResultSliceShape(region));
92 //
93 // allocate space for this new DataConstant's data
94 int len = getNoValues();
95 m_data.resize(len,0.,len);
96 //
97 // create a view of the data with the correct shape
98 // DataArrayView tempView(m_data,shape);
99 DataTypes::RegionLoopRangeType region_loop_range=DataTypes::getSliceRegionLoopRange(region);
100 //
101 // load the view with the data from the slice
102 // tempView.copySlice(other.getPointDataView(),region_loop_range);
103 DataTypes::copySlice(m_data,getShape(),0,other.getVector(),other.getShape(),0,region_loop_range);
104 // setPointDataView(tempView);
105 }
106
107 DataConstant::DataConstant(const FunctionSpace& what,
108 const DataTypes::ShapeType &shape,
109 const DataTypes::ValueType &data)
110 : DataAbstract(what,shape)
111 {
112 //
113 // copy the data in the correct format
114 m_data=data;
115 //
116 // create the view of the data
117 // DataArrayView tempView(m_data,shape);
118 // setPointDataView(tempView);
119 }
120
121 string
122 DataConstant::toString() const
123 {
124 return DataTypes::pointToString(m_data,getShape(),0,"");
125 }
126
127 DataTypes::ValueType::size_type
128 DataConstant::getPointOffset(int sampleNo,
129 int dataPointNo) const
130 {
131 EsysAssert((validSamplePointNo(dataPointNo) && validSampleNo(sampleNo)),
132 "Invalid index, sampleNo: " << sampleNo << " dataPointNo: " << dataPointNo);
133 //
134 // Whatever the coord's always return the same value as this is constant data.
135 return 0;
136 }
137
138 DataTypes::ValueType::size_type
139 DataConstant::getLength() const
140 {
141 return m_data.size();
142 }
143
144 // DataArrayView
145 // DataConstant::getDataPoint(int sampleNo,
146 // int dataPointNo)
147 // {
148 // EsysAssert((validSamplePointNo(dataPointNo) && validSampleNo(sampleNo)),
149 // "Invalid index, sampleNo: " << sampleNo << " dataPointNo: " << dataPointNo);
150 // //
151 // // Whatever the coord's always return the same value as this is constant data.
152 // return getPointDataView();
153 // }
154
155 DataAbstract*
156 DataConstant::getSlice(const DataTypes::RegionType& region) const
157 {
158 return new DataConstant(*this,region);
159 }
160
161 void
162 DataConstant::setSlice(const DataAbstract* value,
163 const DataTypes::RegionType& region)
164 {
165 const DataConstant* tempDataConst=dynamic_cast<const DataConstant*>(value);
166 if (tempDataConst==0) {
167 throw DataException("Programming error - casting to DataConstant.");
168 }
169 //
170 DataTypes::ShapeType shape(DataTypes::getResultSliceShape(region));
171 DataTypes::RegionLoopRangeType region_loop_range=DataTypes::getSliceRegionLoopRange(region);
172 //
173 // check shape:
174 if (getRank()!=region.size()) {
175 throw DataException("Error - Invalid slice region.");
176 }
177 if (getRank()>0 && !DataTypes::checkShape(value->getShape(),shape)) {
178 throw DataException (DataTypes::createShapeErrorMessage(
179 "Error - Couldn't copy slice due to shape mismatch.",shape,value->getShape()));
180 }
181 // getPointDataView().copySliceFrom(tempDataConst->getPointDataView(),region_loop_range);
182 DataTypes::copySliceFrom(m_data,getShape(),0,tempDataConst->getVector(), tempDataConst->getShape(),0,region_loop_range);
183 }
184
185 int
186 DataConstant::archiveData(ofstream& archiveFile,
187 const DataTypes::ValueType::size_type noValues) const
188 {
189 return(m_data.archiveData(archiveFile, noValues));
190 }
191
192 int
193 DataConstant::extractData(ifstream& archiveFile,
194 const DataTypes::ValueType::size_type noValues)
195 {
196 return(m_data.extractData(archiveFile, noValues));
197 }
198
199 void
200 DataConstant::symmetric(DataAbstract* ev)
201 {
202 DataConstant* temp_ev=dynamic_cast<DataConstant*>(ev);
203 if (temp_ev==0) {
204 throw DataException("Error - DataConstant::symmetric: casting to DataConstant failed (propably a programming error).");
205 }
206 /* DataArrayView& thisView=getPointDataView();
207 DataArrayView& evView=ev->getPointDataView();*/
208 DataMaths::symmetric(m_data,getShape(),0,temp_ev->getVector(), temp_ev->getShape(),0);
209 }
210
211 void
212 DataConstant::nonsymmetric(DataAbstract* ev)
213 {
214 DataConstant* temp_ev=dynamic_cast<DataConstant*>(ev);
215 if (temp_ev==0) {
216 throw DataException("Error - DataConstant::nonsymmetric: casting to DataConstant failed (propably a programming error).");
217 }
218 /* DataArrayView& thisView=getPointDataView();
219 DataArrayView& evView=ev->getPointDataView();*/
220 DataMaths::nonsymmetric(m_data,getShape(),0,temp_ev->getVector(), temp_ev->getShape(),0);
221 }
222
223 void
224 DataConstant::trace(DataAbstract* ev, int axis_offset)
225 {
226 DataConstant* temp_ev=dynamic_cast<DataConstant*>(ev);
227 if (temp_ev==0) {
228 throw DataException("Error - DataConstant::trace: casting to DataConstant failed (propably a programming error).");
229 }
230 /* DataArrayView& thisView=getPointDataView();
231 DataArrayView& evView=ev->getPointDataView();*/
232 ValueType& evVec=temp_ev->getVector();
233 const ShapeType& evShape=temp_ev->getShape();
234 DataMaths::trace(m_data,getShape(),0,evVec,evShape,0,axis_offset);
235 }
236
237 void
238 DataConstant::swapaxes(DataAbstract* ev, int axis0, int axis1)
239 {
240 DataConstant* temp_ev=dynamic_cast<DataConstant*>(ev);
241 if (temp_ev==0) {
242 throw DataException("Error - DataConstant::swapaxes: casting to DataConstant failed (propably a programming error).");
243 }
244 // DataArrayView& thisView=getPointDataView();
245 // DataArrayView& evView=ev->getPointDataView();
246 DataMaths::swapaxes(m_data,getShape(),0,temp_ev->getVector(), temp_ev->getShape(),0,axis0,axis1);
247 }
248
249 void
250 DataConstant::transpose(DataAbstract* ev, int axis_offset)
251 {
252 DataConstant* temp_ev=dynamic_cast<DataConstant*>(ev);
253 if (temp_ev==0) {
254 throw DataException("Error - DataConstant::transpose: casting to DataConstant failed (propably a programming error).");
255 }
256 /* DataArrayView& thisView=getPointDataView();
257 DataArrayView& evView=ev->getPointDataView();*/
258 DataMaths::transpose(m_data, getShape(),0, temp_ev->getVector(),temp_ev->getShape(),0,axis_offset);
259 }
260
261 void
262 DataConstant::eigenvalues(DataAbstract* ev)
263 {
264 DataConstant* temp_ev=dynamic_cast<DataConstant*>(ev);
265 if (temp_ev==0) {
266 throw DataException("Error - DataConstant::eigenvalues: casting to DataConstant failed (propably a programming error).");
267 }
268 /* DataArrayView& thisView=getPointDataView();
269 DataArrayView& evView=ev->getPointDataView();*/
270 DataMaths::eigenvalues(m_data,getShape(),0,temp_ev->getVector(), temp_ev->getShape(),0);
271 }
272 void
273 DataConstant::eigenvalues_and_eigenvectors(DataAbstract* ev,DataAbstract* V,const double tol)
274 {
275 DataConstant* temp_ev=dynamic_cast<DataConstant*>(ev);
276 if (temp_ev==0) {
277 throw DataException("Error - DataConstant::eigenvalues_and_eigenvectors: casting to DataConstant failed (propably a programming error).");
278 }
279 DataConstant* temp_V=dynamic_cast<DataConstant*>(V);
280 if (temp_V==0) {
281 throw DataException("Error - DataConstant::eigenvalues_and_eigenvectors: casting to DataConstant failed (propably a programming error).");
282 }
283 // DataArrayView thisView=getPointDataView();
284 // DataArrayView evView=ev->getPointDataView();
285 // DataArrayView VView=V->getPointDataView();
286
287 DataMaths::eigenvalues_and_eigenvectors(m_data, getShape(),0,temp_ev->getVector(), temp_ev->getShape(),0,temp_V->getVector(), temp_V->getShape(),0,tol);
288 }
289
290 void
291 DataConstant::setToZero()
292 {
293 DataTypes::ValueType::size_type n=m_data.size();
294 for (int i=0; i<n ;++i) m_data[i]=0.;
295 }
296
297 void
298 DataConstant::dump(const std::string fileName) const
299 {
300 #ifdef PASO_MPI
301 throw DataException("Error - DataConstant:: dump is not implemented for MPI yet.");
302 #endif
303 #ifdef USE_NETCDF
304 const NcDim* ncdims[DataTypes::maxRank];
305 NcVar* var;
306 int rank = getRank();
307 int type= getFunctionSpace().getTypeCode();
308 int ndims =0;
309 long dims[DataTypes::maxRank];
310 const double* d_ptr=&(m_data[0]);
311 DataTypes::ShapeType shape = getShape();
312
313 // netCDF error handler
314 NcError err(NcError::verbose_nonfatal);
315 // Create the file.
316 NcFile dataFile(fileName.c_str(), NcFile::Replace);
317 // check if writing was successful
318 if (!dataFile.is_valid())
319 throw DataException("Error - DataConstant:: opening of netCDF file for output failed.");
320 if (!dataFile.add_att("type_id",0) )
321 throw DataException("Error - DataConstant:: appending data type to netCDF file failed.");
322 if (!dataFile.add_att("rank",rank) )
323 throw DataException("Error - DataConstant:: appending rank attribute to netCDF file failed.");
324 if (!dataFile.add_att("function_space_type",type))
325 throw DataException("Error - DataConstant:: appending function space attribute to netCDF file failed.");
326
327 if (rank == 0) {
328 if( ! (ncdims[0] = dataFile.add_dim("l", 1)) )
329 throw DataException("Error - DataConstant:: appending ncdimsion 0 to netCDF file failed.");
330 dims[0]=1,
331 ndims=1;
332 } else {
333 ndims=rank;
334 dims[0]=shape[0];
335 if (! (ncdims[0] = dataFile.add_dim("d0",shape[0])) )
336 throw DataException("Error - DataConstant:: appending ncdimsion 0 to netCDF file failed.");
337 if ( rank >1 ) {
338 dims[1]=shape[1];
339 if (! (ncdims[1] = dataFile.add_dim("d1",shape[1])) )
340 throw DataException("Error - DataConstant:: appending ncdimsion 1 to netCDF file failed.");
341 }
342 if ( rank >2 ) {
343 dims[2]=shape[2];
344 if (! (ncdims[2] = dataFile.add_dim("d2", shape[2])) )
345 throw DataException("Error - DataConstant:: appending ncdimsion 2 to netCDF file failed.");
346 }
347 if ( rank >3 ) {
348 dims[3]=shape[3];
349 if (! (ncdims[3] = dataFile.add_dim("d3", shape[3])) )
350 throw DataException("Error - DataConstant:: appending ncdimsion 3 to netCDF file failed.");
351 }
352 }
353
354 if (! ( var = dataFile.add_var("data", ncDouble, ndims, ncdims)) )
355 throw DataException("Error - DataConstant:: appending variable to netCDF file failed.");
356 if (! (var->put(d_ptr,dims)) )
357 throw DataException("Error - DataConstant:: copy data to netCDF buffer failed.");
358 #else
359 throw DataException("Error - DataConstant:: dump is not configured with netCDF. Please contact your installation manager.");
360 #endif
361 }
362
363 } // end of namespace

Properties

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

  ViewVC Help
Powered by ViewVC 1.1.26