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

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

Parent Directory Parent Directory | Revision Log Revision Log


Revision 950 - (hide annotations)
Tue Feb 6 07:01:11 2007 UTC (12 years, 8 months ago) by gross
File size: 11247 byte(s)
escript data objects can now be saved to netCDF files, see http://www.unidata.ucar.edu/software/netcdf/.
Currently only constant data are implemented with expanded and tagged data to follow.
There are two new functions to dump a data object

   s=Data(...)
   s.dump(<filename>)

and to recover it

   s=load(<filename>, domain)

Notice that the function space of s is recovered but domain is still need. 

dump and load will replace archive and extract.

The installation needs now the netCDF installed. 


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