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

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

Parent Directory Parent Directory | Revision Log Revision Log


Revision 2669 - (hide annotations)
Thu Sep 17 06:01:30 2009 UTC (10 years, 6 months ago) by jfenwick
File size: 11450 byte(s)
Added a bunch of valgrind supressions. Hopefully this will make it more useful.

Asking for an offset on a DataConstant containing no samples no longer fails under debug builds.
I've disabled the sample and dataPoint bounds checks since they don't make sense on DataConstant.
1 jgs 82
2 ksteube 1312 /*******************************************************
3 ksteube 1811 *
4 jfenwick 2548 * Copyright (c) 2003-2009 by University of Queensland
5 ksteube 1811 * 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 ksteube 1312
14 ksteube 1811
15 ksteube 1827 #include "Data.h"
16 jgs 474 #include "DataConstant.h"
17     #include "DataException.h"
18 robwdcock 682 #include "esysUtils/EsysAssert.h"
19 jgs 82
20     #include <iostream>
21     #include <boost/python/extract.hpp>
22 gross 1023 #ifdef USE_NETCDF
23 ksteube 1312 #include <netcdfcpp.h>
24 gross 1023 #endif
25 ksteube 1827 #ifdef PASO_MPI
26     #include <mpi.h>
27     #endif
28 jgs 82
29 matt 1319 #include <boost/python/extract.hpp>
30 jfenwick 1796 #include "DataMaths.h"
31 matt 1319
32 jfenwick 2271 // #define CHECK_FOR_EX_WRITE if (!checkNoSharing()) {throw DataException("Attempt to modify shared object");}
33    
34     #define CHECK_FOR_EX_WRITE if (!checkNoSharing()) {std::ostringstream ss; ss << " Attempt to modify shared object. line " << __LINE__ << " of " << __FILE__; ss << m_owners.size(); cerr << ss << endl; /* *((int*)0)=17; */throw DataException(ss.str());}
35    
36 jgs 82 using namespace std;
37 matt 1319 using namespace boost::python;
38 jgs 82
39     namespace escript {
40    
41 jfenwick 2271 DataConstant::DataConstant(const WrappedArray& value,
42 jgs 102 const FunctionSpace& what)
43 jfenwick 2271 : parent(what,value.getShape())
44 jgs 102 {
45 jfenwick 2271 m_data.copyFromArray(value,1);
46 jgs 102 }
47 jgs 82
48 jgs 102 DataConstant::DataConstant(const DataConstant& other)
49 jfenwick 2005 : parent(other.getFunctionSpace(),other.getShape())
50 jfenwick 2271 {
51 jgs 102 m_data=other.m_data;
52     }
53 jgs 82
54 jgs 102 DataConstant::DataConstant(const DataConstant& other,
55 jfenwick 1796 const DataTypes::RegionType& region)
56 jfenwick 2005 : parent(other.getFunctionSpace(),DataTypes::getResultSliceShape(region))
57 jgs 102 {
58     //
59     // allocate space for this new DataConstant's data
60 jfenwick 1796 int len = getNoValues();
61 jgs 151 m_data.resize(len,0.,len);
62 jgs 102 //
63     // create a view of the data with the correct shape
64 jfenwick 1796 DataTypes::RegionLoopRangeType region_loop_range=DataTypes::getSliceRegionLoopRange(region);
65 jgs 119 //
66     // load the view with the data from the slice
67 jfenwick 2271 DataTypes::copySlice(m_data,getShape(),0,other.getVectorRO(),other.getShape(),0,region_loop_range);
68 jgs 119 }
69    
70     DataConstant::DataConstant(const FunctionSpace& what,
71 jfenwick 1796 const DataTypes::ShapeType &shape,
72     const DataTypes::ValueType &data)
73 jfenwick 2005 : parent(what,shape)
74 jgs 119 {
75 jgs 102 //
76 jgs 119 // copy the data in the correct format
77     m_data=data;
78     //
79     // create the view of the data
80 jfenwick 1796 // DataArrayView tempView(m_data,shape);
81     // setPointDataView(tempView);
82 jgs 102 }
83 jgs 82
84 jgs 102 string
85     DataConstant::toString() const
86     {
87 jfenwick 1796 return DataTypes::pointToString(m_data,getShape(),0,"");
88 jgs 102 }
89 jgs 82
90 jfenwick 1799
91     DataAbstract*
92     DataConstant::deepCopy()
93     {
94     return new DataConstant(*this);
95     }
96    
97    
98 jfenwick 1796 DataTypes::ValueType::size_type
99 jgs 102 DataConstant::getPointOffset(int sampleNo,
100     int dataPointNo) const
101     {
102 jfenwick 2669 // We avoid this check for constant data due to issues manipulating
103     // data with no samples.
104    
105     // EsysAssert((validSamplePointNo(dataPointNo) && validSampleNo(sampleNo)),
106     // "Invalid index, sampleNo: " << sampleNo << " dataPointNo: " << dataPointNo);
107 jgs 122 //
108     // Whatever the coord's always return the same value as this is constant data.
109 jgs 102 return 0;
110     }
111 jgs 82
112 jfenwick 1796 DataTypes::ValueType::size_type
113 jfenwick 2005 DataConstant::getPointOffset(int sampleNo,
114     int dataPointNo)
115     {
116 jfenwick 2669 // We avoid this check for constant data due to issues manipulating
117     // data with no samples.
118    
119     // EsysAssert((validSamplePointNo(dataPointNo) && validSampleNo(sampleNo)),
120     // "Invalid index, sampleNo: " << sampleNo << " dataPointNo: " << dataPointNo);
121 jfenwick 2005 //
122     // Whatever the coord's always return the same value as this is constant data.
123 jfenwick 2669 //
124    
125 jfenwick 2005 return 0;
126     }
127    
128     DataTypes::ValueType::size_type
129 jgs 102 DataConstant::getLength() const
130     {
131     return m_data.size();
132     }
133 jgs 82
134 jgs 102 DataAbstract*
135 jfenwick 1796 DataConstant::getSlice(const DataTypes::RegionType& region) const
136 jgs 102 {
137     return new DataConstant(*this,region);
138     }
139 jgs 82
140 jgs 102 void
141     DataConstant::setSlice(const DataAbstract* value,
142 jfenwick 1796 const DataTypes::RegionType& region)
143 jgs 102 {
144     const DataConstant* tempDataConst=dynamic_cast<const DataConstant*>(value);
145     if (tempDataConst==0) {
146     throw DataException("Programming error - casting to DataConstant.");
147 jgs 82 }
148 jfenwick 2271 CHECK_FOR_EX_WRITE
149 matt 1319 //
150 jfenwick 1796 DataTypes::ShapeType shape(DataTypes::getResultSliceShape(region));
151     DataTypes::RegionLoopRangeType region_loop_range=DataTypes::getSliceRegionLoopRange(region);
152 jgs 108 //
153     // check shape:
154 jfenwick 1796 if (getRank()!=region.size()) {
155 jgs 108 throw DataException("Error - Invalid slice region.");
156     }
157 jfenwick 1796 if (getRank()>0 && !DataTypes::checkShape(value->getShape(),shape)) {
158     throw DataException (DataTypes::createShapeErrorMessage(
159     "Error - Couldn't copy slice due to shape mismatch.",shape,value->getShape()));
160 jgs 108 }
161 jfenwick 1796 // getPointDataView().copySliceFrom(tempDataConst->getPointDataView(),region_loop_range);
162 jfenwick 2271 DataTypes::copySliceFrom(m_data,getShape(),0,tempDataConst->getVectorRO(), tempDataConst->getShape(),0,region_loop_range);
163 jgs 102 }
164 jgs 82
165 jgs 123
166    
167 gross 580 void
168 ksteube 775 DataConstant::symmetric(DataAbstract* ev)
169     {
170     DataConstant* temp_ev=dynamic_cast<DataConstant*>(ev);
171     if (temp_ev==0) {
172     throw DataException("Error - DataConstant::symmetric: casting to DataConstant failed (propably a programming error).");
173     }
174 jfenwick 2271 DataMaths::symmetric(m_data,getShape(),0,temp_ev->getVectorRW(), temp_ev->getShape(),0);
175 ksteube 775 }
176    
177     void
178     DataConstant::nonsymmetric(DataAbstract* ev)
179     {
180     DataConstant* temp_ev=dynamic_cast<DataConstant*>(ev);
181     if (temp_ev==0) {
182     throw DataException("Error - DataConstant::nonsymmetric: casting to DataConstant failed (propably a programming error).");
183     }
184 jfenwick 2271 DataMaths::nonsymmetric(m_data,getShape(),0,temp_ev->getVectorRW(), temp_ev->getShape(),0);
185 ksteube 775 }
186    
187     void
188 gross 800 DataConstant::trace(DataAbstract* ev, int axis_offset)
189 ksteube 775 {
190     DataConstant* temp_ev=dynamic_cast<DataConstant*>(ev);
191     if (temp_ev==0) {
192 gross 800 throw DataException("Error - DataConstant::trace: casting to DataConstant failed (propably a programming error).");
193 ksteube 775 }
194 jfenwick 2271 ValueType& evVec=temp_ev->getVectorRW();
195 jfenwick 1796 const ShapeType& evShape=temp_ev->getShape();
196     DataMaths::trace(m_data,getShape(),0,evVec,evShape,0,axis_offset);
197 ksteube 775 }
198    
199     void
200 gross 804 DataConstant::swapaxes(DataAbstract* ev, int axis0, int axis1)
201 gross 800 {
202     DataConstant* temp_ev=dynamic_cast<DataConstant*>(ev);
203     if (temp_ev==0) {
204 gross 804 throw DataException("Error - DataConstant::swapaxes: casting to DataConstant failed (propably a programming error).");
205 gross 800 }
206 jfenwick 2271 DataMaths::swapaxes(m_data,getShape(),0,temp_ev->getVectorRW(), temp_ev->getShape(),0,axis0,axis1);
207 gross 800 }
208    
209     void
210 ksteube 775 DataConstant::transpose(DataAbstract* ev, int axis_offset)
211     {
212     DataConstant* temp_ev=dynamic_cast<DataConstant*>(ev);
213     if (temp_ev==0) {
214     throw DataException("Error - DataConstant::transpose: casting to DataConstant failed (propably a programming error).");
215     }
216 jfenwick 2271 DataMaths::transpose(m_data, getShape(),0, temp_ev->getVectorRW(),temp_ev->getShape(),0,axis_offset);
217 ksteube 775 }
218    
219     void
220 gross 580 DataConstant::eigenvalues(DataAbstract* ev)
221     {
222     DataConstant* temp_ev=dynamic_cast<DataConstant*>(ev);
223     if (temp_ev==0) {
224     throw DataException("Error - DataConstant::eigenvalues: casting to DataConstant failed (propably a programming error).");
225     }
226 jfenwick 2271 DataMaths::eigenvalues(m_data,getShape(),0,temp_ev->getVectorRW(), temp_ev->getShape(),0);
227 gross 580 }
228     void
229     DataConstant::eigenvalues_and_eigenvectors(DataAbstract* ev,DataAbstract* V,const double tol)
230     {
231     DataConstant* temp_ev=dynamic_cast<DataConstant*>(ev);
232     if (temp_ev==0) {
233     throw DataException("Error - DataConstant::eigenvalues_and_eigenvectors: casting to DataConstant failed (propably a programming error).");
234     }
235     DataConstant* temp_V=dynamic_cast<DataConstant*>(V);
236     if (temp_V==0) {
237     throw DataException("Error - DataConstant::eigenvalues_and_eigenvectors: casting to DataConstant failed (propably a programming error).");
238     }
239 jfenwick 2271 DataMaths::eigenvalues_and_eigenvectors(m_data, getShape(),0,temp_ev->getVectorRW(), temp_ev->getShape(),0,temp_V->getVectorRW(), temp_V->getShape(),0,tol);
240 gross 580 }
241    
242 gross 950 void
243 gross 1118 DataConstant::setToZero()
244     {
245 jfenwick 2271 CHECK_FOR_EX_WRITE
246 jfenwick 1796 DataTypes::ValueType::size_type n=m_data.size();
247 gross 1118 for (int i=0; i<n ;++i) m_data[i]=0.;
248     }
249    
250     void
251 gross 950 DataConstant::dump(const std::string fileName) const
252     {
253 gross 1149 #ifdef USE_NETCDF
254 jfenwick 1796 const NcDim* ncdims[DataTypes::maxRank];
255 gross 950 NcVar* var;
256 jfenwick 1796 int rank = getRank();
257 gross 950 int type= getFunctionSpace().getTypeCode();
258     int ndims =0;
259 jfenwick 1796 long dims[DataTypes::maxRank];
260 gross 1141 const double* d_ptr=&(m_data[0]);
261 jfenwick 1796 DataTypes::ShapeType shape = getShape();
262 jfenwick 1872 int mpi_iam=getFunctionSpace().getDomain()->getMPIRank();
263     int mpi_num=getFunctionSpace().getDomain()->getMPISize();
264 ksteube 1827 #ifdef PASO_MPI
265     MPI_Status status;
266     #endif
267 matt 1319
268 ksteube 1827 #ifdef PASO_MPI
269     /* Serialize NetCDF I/O */
270     if (mpi_iam>0) MPI_Recv(&ndims, 0, MPI_INT, mpi_iam-1, 81802, MPI_COMM_WORLD, &status);
271     #endif
272    
273 gross 950 // netCDF error handler
274     NcError err(NcError::verbose_nonfatal);
275     // Create the file.
276 ksteube 1827 char *newFileName = Escript_MPI_appendRankToFileName(fileName.c_str(), mpi_num, mpi_iam);
277     NcFile dataFile(newFileName, NcFile::Replace);
278 gross 950 // check if writing was successful
279 matt 1319 if (!dataFile.is_valid())
280 gross 950 throw DataException("Error - DataConstant:: opening of netCDF file for output failed.");
281 matt 1319 if (!dataFile.add_att("type_id",0) )
282 gross 950 throw DataException("Error - DataConstant:: appending data type to netCDF file failed.");
283 matt 1319 if (!dataFile.add_att("rank",rank) )
284 gross 950 throw DataException("Error - DataConstant:: appending rank attribute to netCDF file failed.");
285 matt 1319 if (!dataFile.add_att("function_space_type",type))
286 gross 950 throw DataException("Error - DataConstant:: appending function space attribute to netCDF file failed.");
287 gross 580
288 gross 950 if (rank == 0) {
289 matt 1319 if( ! (ncdims[0] = dataFile.add_dim("l", 1)) )
290 ksteube 1800 throw DataException("Error - DataConstant:: appending ncdimension 0 to netCDF file failed.");
291 gross 950 dims[0]=1,
292     ndims=1;
293     } else {
294     ndims=rank;
295     dims[0]=shape[0];
296 matt 1319 if (! (ncdims[0] = dataFile.add_dim("d0",shape[0])) )
297 ksteube 1800 throw DataException("Error - DataConstant:: appending ncdimension 0 to netCDF file failed.");
298 gross 950 if ( rank >1 ) {
299     dims[1]=shape[1];
300 matt 1319 if (! (ncdims[1] = dataFile.add_dim("d1",shape[1])) )
301 ksteube 1800 throw DataException("Error - DataConstant:: appending ncdimension 1 to netCDF file failed.");
302 gross 950 }
303     if ( rank >2 ) {
304     dims[2]=shape[2];
305 matt 1319 if (! (ncdims[2] = dataFile.add_dim("d2", shape[2])) )
306 ksteube 1800 throw DataException("Error - DataConstant:: appending ncdimension 2 to netCDF file failed.");
307 gross 950 }
308     if ( rank >3 ) {
309     dims[3]=shape[3];
310 matt 1319 if (! (ncdims[3] = dataFile.add_dim("d3", shape[3])) )
311 ksteube 1800 throw DataException("Error - DataConstant:: appending ncdimension 3 to netCDF file failed.");
312 gross 950 }
313     }
314 matt 1319
315 gross 950 if (! ( var = dataFile.add_var("data", ncDouble, ndims, ncdims)) )
316     throw DataException("Error - DataConstant:: appending variable to netCDF file failed.");
317 gross 1141 if (! (var->put(d_ptr,dims)) )
318 matt 1319 throw DataException("Error - DataConstant:: copy data to netCDF buffer failed.");
319 ksteube 1827 #ifdef PASO_MPI
320     if (mpi_iam<mpi_num-1) MPI_Send(&ndims, 0, MPI_INT, mpi_iam+1, 81802, MPI_COMM_WORLD);
321     #endif
322 gross 1023 #else
323     throw DataException("Error - DataConstant:: dump is not configured with netCDF. Please contact your installation manager.");
324     #endif
325 gross 950 }
326    
327 jfenwick 2271 // These used to be marked as inline in DataConstant.
328     // But they are marked virtual in DataReady
329     DataTypes::ValueType&
330     DataConstant::getVectorRW()
331     {
332     CHECK_FOR_EX_WRITE
333     return m_data;
334     }
335    
336     const DataTypes::ValueType&
337     DataConstant::getVectorRO() const
338     {
339     return m_data;
340     }
341    
342 jgs 82 } // end of namespace

Properties

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

  ViewVC Help
Powered by ViewVC 1.1.26