/[escript]/trunk/ripley/test/python/run_readWriteOnRipley.py
ViewVC logotype

Annotation of /trunk/ripley/test/python/run_readWriteOnRipley.py

Parent Directory Parent Directory | Revision Log Revision Log


Revision 6141 - (hide annotations)
Wed Apr 6 03:51:30 2016 UTC (2 years, 11 months ago) by caltinay
File MIME type: text/x-python
File size: 15813 byte(s)
more namespacing of defines.

1 caltinay 5008
2     ##############################################################################
3     #
4 jfenwick 5863 # Copyright (c) 2003-2016 by The University of Queensland
5 caltinay 5008 # http://www.uq.edu.au
6     #
7     # Primary Business: Queensland, Australia
8 jfenwick 6112 # Licensed under the Apache License, version 2.0
9     # http://www.apache.org/licenses/LICENSE-2.0
10 caltinay 5008 #
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 sshaw 5706 from __future__ import print_function, division
18    
19 jfenwick 5863 __copyright__="""Copyright (c) 2003-2016 by The University of Queensland
20 caltinay 5008 http://www.uq.edu.au
21     Primary Business: Queensland, Australia"""
22 jfenwick 6112 __license__="""Licensed under the Apache License, version 2.0
23     http://www.apache.org/licenses/LICENSE-2.0"""
24 caltinay 5008 __url__="https://launchpad.net/escript-finley"
25    
26     import os
27 caltinay 5012 import numpy as np
28 caltinay 5008 import esys.escriptcore.utestselect as unittest
29     from esys.escriptcore.testing import *
30     from esys.escript import *
31 caltinay 5077 from esys.ripley import *
32 caltinay 5008
33     try:
34     RIPLEY_WORKDIR=os.environ['RIPLEY_WORKDIR']
35     except KeyError:
36 caltinay 5051 RIPLEY_WORKDIR='/tmp'
37 caltinay 5008
38 caltinay 6141 HAVE_UNZIP = getEscriptParamInt('HAVE_UNZIP')
39    
40 caltinay 5008 #NE=4 # number elements, must be even
41 caltinay 5051 #for x in [int(sqrt(mpiSize)),2,3,5,7,1]:
42 caltinay 5008 # NX=x
43 caltinay 5051 # NY=mpiSize//x
44     # if NX*NY == mpiSize:
45 caltinay 5008 # break
46     #
47 caltinay 5051 #for x in [(int(mpiSize**(1/3.)),int(mpiSize**(1/3.))),(2,3),(2,2),(1,2),(1,1)]:
48 caltinay 5008 # NXb=x[0]
49     # NYb=x[1]
50 caltinay 5051 # NZb=mpiSize//(x[0]*x[1])
51     # if NXb*NYb*NZb == mpiSize:
52 caltinay 5008 # break
53    
54 caltinay 5051 mpiSize = getMPISizeWorld()
55     mpiRank = getMPIRankWorld()
56 caltinay 5012
57     def adjust(NE, ftype):
58 caltinay 5172 if ftype in (ContinuousFunction, Solution):
59 caltinay 5012 return [i+1 for i in NE]
60     return NE
61    
62    
63 sshaw 5059 class WriteBinaryGridTestBase(unittest.TestCase): #subclassing required
64 caltinay 5051 NX = 10*mpiSize-1
65 caltinay 5012 NZ = 10
66    
67     def generateUniqueData(self, ftype):
68     dim = self.domain.getDim()
69     FSx=ftype(self.domain).getX()
70     NE = adjust(self.NE, ftype)
71     # normalise and scale range of values
72     x = [FSx[i]-inf(FSx[i]) for i in range(dim)]
73     x = [(NE[i]-1)*(x[i]/sup(x[i])) for i in range(dim)]
74     xMax = [int(sup(x[i]))+1 for i in range(dim)]
75     nvals=NE[0]*NE[1]
76     data = x[0] + xMax[0]*x[1]
77 caltinay 5077 if self.datatype == DATATYPE_INT32:
78 caltinay 5012 data += 0.05
79     if dim > 2:
80     data = data + xMax[0]*xMax[1]*x[2]
81     nvals*=NE[2]
82    
83     grid = np.array(range(nvals), dtype=self.dtype).reshape(tuple(reversed(NE)))
84     return data, grid
85    
86 caltinay 5051 def writeThenRead(self, data, ftype, fcode):
87     filename = os.path.join(RIPLEY_WORKDIR, "_wgrid%dd%s"%(self.domain.getDim(),fcode))
88     filename = filename + self.dtype.replace('<','L').replace('>','B')
89     self.domain.writeBinaryGrid(data, filename, self.byteorder, self.datatype)
90     MPIBarrierWorld()
91     result = np.fromfile(filename, dtype=self.dtype).reshape(
92     tuple(reversed(adjust(self.NE,ftype))))
93     return result
94    
95 caltinay 5012 def test_writeGrid2D(self):
96     self.NE = [self.NX, self.NZ]
97     self.domain = Rectangle(self.NE[0], self.NE[1], d1=0)
98 caltinay 5172 for ftype,fcode in [(ReducedFunction,'RF'), (ContinuousFunction,'CF'), (Solution, 'Sol')]:
99 caltinay 5051 data, ref = self.generateUniqueData(ftype)
100     result = self.writeThenRead(data, ftype, fcode)
101     self.assertAlmostEquals(Lsup(ref-result), 0, delta=1e-9,
102     msg="Data doesn't match for "+str(ftype(self.domain)))
103 caltinay 5012
104     def test_writeGrid3D(self):
105     self.NE = [self.NX, self.NX, self.NZ]
106     self.domain = Brick(self.NE[0], self.NE[1], self.NE[2], d2=0)
107 caltinay 5172 for ftype,fcode in [(ReducedFunction,'RF'), (ContinuousFunction,'CF'), (Solution, 'Sol')]:
108 caltinay 5051 data, ref = self.generateUniqueData(ftype)
109     result = self.writeThenRead(data, ftype, fcode)
110     self.assertAlmostEquals(Lsup(ref-result), 0, delta=1e-9,
111     msg="Data doesn't match for "+str(ftype(self.domain)))
112 caltinay 5012
113     class Test_writeBinaryGridRipley_LITTLE_FLOAT32(WriteBinaryGridTestBase):
114 caltinay 5008 def setUp(self):
115 caltinay 5077 self.byteorder = BYTEORDER_LITTLE_ENDIAN
116     self.datatype = DATATYPE_FLOAT32
117 caltinay 5012 self.dtype = "<f4"
118    
119     class Test_writeBinaryGridRipley_LITTLE_FLOAT64(WriteBinaryGridTestBase):
120     def setUp(self):
121 caltinay 5077 self.byteorder = BYTEORDER_LITTLE_ENDIAN
122     self.datatype = DATATYPE_FLOAT64
123 caltinay 5012 self.dtype = "<f8"
124 caltinay 5008
125 caltinay 5012 class Test_writeBinaryGridRipley_LITTLE_INT32(WriteBinaryGridTestBase):
126     def setUp(self):
127 caltinay 5077 self.byteorder = BYTEORDER_LITTLE_ENDIAN
128     self.datatype = DATATYPE_INT32
129 caltinay 5012 self.dtype = "<i4"
130    
131     class Test_writeBinaryGridRipley_BIG_FLOAT32(WriteBinaryGridTestBase):
132     def setUp(self):
133 caltinay 5077 self.byteorder = BYTEORDER_BIG_ENDIAN
134     self.datatype = DATATYPE_FLOAT32
135 caltinay 5012 self.dtype = ">f4"
136    
137     class Test_writeBinaryGridRipley_BIG_FLOAT64(WriteBinaryGridTestBase):
138     def setUp(self):
139 caltinay 5077 self.byteorder = BYTEORDER_BIG_ENDIAN
140     self.datatype = DATATYPE_FLOAT64
141 caltinay 5012 self.dtype = ">f8"
142    
143     class Test_writeBinaryGridRipley_BIG_INT32(WriteBinaryGridTestBase):
144     def setUp(self):
145 caltinay 5077 self.byteorder = BYTEORDER_BIG_ENDIAN
146     self.datatype = DATATYPE_INT32
147 caltinay 5012 self.dtype = ">i4"
148    
149    
150 sshaw 5059 class ReadBinaryGridTestBase(unittest.TestCase): #subclassing required
151 caltinay 5051 """
152     The reader tests work in several stages:
153     1) create numpy array and write to temporary file (ref)
154     2) call readBinaryGrid with that filename
155     3) write the resulting Data object using writeBinaryGrid (test)
156     4) read the result using numpy and compare (ref) and (test)
157     As such, it is important to note that a working writeBinaryGrid() method
158     is assumed!
159     """
160 caltinay 5096 # set defaults which may be overridden in subclasses
161     NX = 10
162     NZ = 8
163     fspaces = [(ReducedFunction,'RF'), (ContinuousFunction,'CF')]
164     byteorder = BYTEORDER_NATIVE
165     datatype = DATATYPE_FLOAT64
166     dtype = "f8"
167 caltinay 5051 shape = ()
168     fill = -42.57
169 caltinay 5096 first = [0,0,0]
170     multiplier = [1,1,1]
171     reverse = [0,0,0]
172 caltinay 5051
173     def generateUniqueData(self, ftype):
174     dim = self.domain.getDim()
175 caltinay 5096 NE = adjust(self.Ndata, ftype)
176 caltinay 5051 nvals=NE[0]*NE[1]
177     if dim > 2:
178     nvals*=NE[2]
179     grid = np.array(range(nvals), dtype=self.dtype).reshape(tuple(reversed(NE)))
180     return grid
181    
182     def write(self, data, filename):
183     self.domain.writeBinaryGrid(data, filename, self.byteorder, self.datatype)
184    
185     def read(self, filename, ftype):
186 caltinay 5096 first = self.first[:self.domain.getDim()]
187     multiplier = self.multiplier[:self.domain.getDim()]
188     reverse = self.reverse[:self.domain.getDim()]
189     numValues=adjust(self.Ndata, ftype)
190 caltinay 5077 return readBinaryGrid(filename, ftype(self.domain),
191 caltinay 5051 shape=self.shape, fill=self.fill, byteOrder=self.byteorder,
192 caltinay 5096 dataType=self.datatype, first=first, numValues=numValues,
193     multiplier=multiplier, reverse=reverse)
194 caltinay 5051
195     def numpy2Data2Numpy(self, ref, ftype, fcode):
196     filename = os.path.join(RIPLEY_WORKDIR, "_rgrid%dd%s"%(self.domain.getDim(),fcode))
197     filename = filename + self.dtype.replace('<','L').replace('>','B')
198     if mpiRank == 0:
199     ref.tofile(filename)
200     MPIBarrierWorld()
201     # step 2 - read
202     data = self.read(filename, ftype)
203     MPIBarrierWorld()
204     # step 3 - write
205     self.write(data, filename) # overwrite is ok
206     MPIBarrierWorld()
207     result = np.fromfile(filename, dtype=self.dtype).reshape(
208     tuple(reversed(adjust(self.NE,ftype))))
209     return result
210    
211     def test_readGrid2D(self):
212 caltinay 5098 if self.multiplier[0] == 1:
213     self.NE = [self.NX*mpiSize-1, self.NZ*self.multiplier[1]]
214     else:
215     self.NE = [self.NX*mpiSize*self.multiplier[0]-1, self.NZ*self.multiplier[1]]
216 caltinay 5096 self.domain = Rectangle(self.NE[0], self.NE[1], d0=mpiSize, d1=1)
217     for ftype,fcode in self.fspaces:
218 caltinay 5098 self.Ndata = [self.NX*mpiSize-1, self.NZ]
219     if ftype==ContinuousFunction:
220     self.Ndata[1] = self.NZ-1
221 caltinay 5051 # step 1 - generate
222     ref = self.generateUniqueData(ftype)
223 caltinay 5096 # step 2 & 3
224 caltinay 5051 result = self.numpy2Data2Numpy(ref, ftype, fcode)
225 caltinay 5096 # apply transformations to be able to compare
226     if self.reverse[0]:
227     result = result[...,::-1]
228     if self.reverse[1]:
229     result = result[::-1,:]
230     for i in range(2):
231     ref = np.repeat(ref, self.multiplier[i], axis=1-i)
232 caltinay 5098
233     # if domain larger than data: add column(s)/row(s) with fill value
234     fill=np.array(self.fill, dtype=ref.dtype)
235     realNE = adjust(self.NE,ftype)
236     for d in range(2):
237     excess = realNE[d]-ref.shape[1-d]
238     if excess > 0:
239     shape = list(ref.shape)
240     shape[1-d] = excess
241     extra = fill * np.ones(shape)
242     if self.reverse[d]:
243     ref = np.append(extra, ref, axis=1-d)
244     else:
245     ref = np.append(ref, extra, axis=1-d)
246    
247 caltinay 5096 # step 4 - compare
248 caltinay 5051 self.assertAlmostEquals(Lsup(ref-result), 0, delta=1e-9,
249     msg="Data doesn't match for "+str(ftype(self.domain)))
250    
251     def test_readGrid3D(self):
252 caltinay 5098 if self.multiplier[0] == 1:
253     self.NE = [self.NX*mpiSize-1, self.NX*self.multiplier[1], self.NZ*self.multiplier[2]]
254     else:
255     self.NE = [self.NX*mpiSize*self.multiplier[0]-1,
256     self.NX*self.multiplier[1], self.NZ*self.multiplier[2]]
257 caltinay 5096 self.domain = Brick(self.NE[0], self.NE[1], self.NE[2], d0=mpiSize, d1=1, d2=1)
258     for ftype,fcode in self.fspaces:
259 caltinay 5098 self.Ndata = [self.NX*mpiSize-1, self.NX, self.NZ]
260     if ftype==ContinuousFunction:
261     self.Ndata[1] = self.NX-1
262     self.Ndata[2] = self.NZ-1
263 caltinay 5051 # step 1 - generate
264     ref = self.generateUniqueData(ftype)
265 caltinay 5096 # step 2 & 3
266 caltinay 5051 result = self.numpy2Data2Numpy(ref, ftype, fcode)
267 caltinay 5096 # apply transformations to be able to compare
268     if self.reverse[0]:
269     result = result[...,::-1]
270     if self.reverse[1]:
271     result = result[...,::-1,:]
272     if self.reverse[2]:
273     result = result[::-1,:,:]
274     for i in range(3):
275     ref = np.repeat(ref, self.multiplier[i], axis=2-i)
276 caltinay 5098
277     # if domain larger than data: add column(s)/row(s) with fill value
278     fill=np.array(self.fill, dtype=ref.dtype)
279     realNE = adjust(self.NE,ftype)
280     for d in range(3):
281     excess = realNE[d]-ref.shape[2-d]
282     if excess > 0:
283     shape = list(ref.shape)
284     shape[2-d] = excess
285     extra = fill * np.ones(shape)
286     if self.reverse[d]:
287     ref = np.append(extra, ref, axis=2-d)
288     else:
289     ref = np.append(ref, extra, axis=2-d)
290    
291 caltinay 5096 # step 4 - compare
292 caltinay 5051 self.assertAlmostEquals(Lsup(ref-result), 0, delta=1e-9,
293     msg="Data doesn't match for "+str(ftype(self.domain)))
294    
295    
296 caltinay 5096 # The following block tests the reader for different byte orders and data
297     # types with domain-filling data (i.e. multiplier=1, reverse=0 and N=NE)
298    
299 caltinay 5051 class Test_readBinaryGridRipley_LITTLE_FLOAT32(ReadBinaryGridTestBase):
300     def setUp(self):
301 caltinay 5077 self.byteorder = BYTEORDER_LITTLE_ENDIAN
302     self.datatype = DATATYPE_FLOAT32
303 caltinay 5051 self.dtype = "<f4"
304    
305     class Test_readBinaryGridRipley_LITTLE_FLOAT64(ReadBinaryGridTestBase):
306     def setUp(self):
307 caltinay 5077 self.byteorder = BYTEORDER_LITTLE_ENDIAN
308     self.datatype = DATATYPE_FLOAT64
309 caltinay 5051 self.dtype = "<f8"
310    
311     class Test_readBinaryGridRipley_LITTLE_INT32(ReadBinaryGridTestBase):
312     def setUp(self):
313 caltinay 5077 self.byteorder = BYTEORDER_LITTLE_ENDIAN
314     self.datatype = DATATYPE_INT32
315 caltinay 5051 self.dtype = "<i4"
316    
317     class Test_readBinaryGridRipley_BIG_FLOAT32(ReadBinaryGridTestBase):
318     def setUp(self):
319 caltinay 5077 self.byteorder = BYTEORDER_BIG_ENDIAN
320     self.datatype = DATATYPE_FLOAT32
321 caltinay 5051 self.dtype = ">f4"
322    
323     class Test_readBinaryGridRipley_BIG_FLOAT64(ReadBinaryGridTestBase):
324     def setUp(self):
325 caltinay 5077 self.byteorder = BYTEORDER_BIG_ENDIAN
326     self.datatype = DATATYPE_FLOAT64
327 caltinay 5051 self.dtype = ">f8"
328    
329     class Test_readBinaryGridRipley_BIG_INT32(ReadBinaryGridTestBase):
330     def setUp(self):
331 caltinay 5077 self.byteorder = BYTEORDER_BIG_ENDIAN
332     self.datatype = DATATYPE_INT32
333 caltinay 5051 self.dtype = ">i4"
334    
335 caltinay 5111 @unittest.skip("reverseX not supported yet")
336 caltinay 5096 class Test_readBinaryGridRipley_reverseX(ReadBinaryGridTestBase):
337     def setUp(self):
338     self.reverse = [1,0,0]
339 caltinay 5051
340 caltinay 5111 @unittest.skip("reverseY not supported yet")
341 caltinay 5096 class Test_readBinaryGridRipley_reverseY(ReadBinaryGridTestBase):
342     def setUp(self):
343     self.reverse = [0,1,0]
344    
345     class Test_readBinaryGridRipley_reverseZ(ReadBinaryGridTestBase):
346     def setUp(self):
347     self.reverse = [0,0,1]
348    
349     class Test_readBinaryGridRipley_multiplierX(ReadBinaryGridTestBase):
350     def setUp(self):
351     self.multiplier = [2,1,1]
352    
353     class Test_readBinaryGridRipley_multiplierY(ReadBinaryGridTestBase):
354     def setUp(self):
355     self.multiplier = [1,2,1]
356    
357     class Test_readBinaryGridRipley_multiplierZ(ReadBinaryGridTestBase):
358     def setUp(self):
359     self.multiplier = [1,1,2]
360    
361     class Test_readBinaryGridRipley_multiplierXYZ(ReadBinaryGridTestBase):
362     def setUp(self):
363     self.multiplier = [2,3,4]
364    
365    
366 caltinay 5012 @unittest.skipIf(getMPISizeWorld() > 1,
367     "Skipping compressed binary grid tests due to element stretching")
368     class Test_readBinaryGridZippedRipley(unittest.TestCase):
369     # constants
370 caltinay 5077 byteorder = BYTEORDER_NATIVE
371     datatype = DATATYPE_FLOAT64
372 caltinay 5012
373 caltinay 5008 def read(self, filename, FS, expected, zipped = False):
374     first = [0 for i in expected]
375     reverse = [0 for i in expected]
376     scale = [1 for i in expected]
377 caltinay 5012
378 caltinay 5008 if not zipped:
379 caltinay 5077 return readBinaryGrid(filename, FS, (), 50000,
380 caltinay 5008 self.byteorder, self.datatype, first, expected, scale, reverse)
381 caltinay 5012
382 caltinay 6141 if not HAVE_UNZIP:
383 caltinay 5008 raise unittest.SkipTest("unzip library not available (boost_iostreams)")
384     return ripleycpp._readBinaryGridFromZipped(filename, FS, (), 50000,
385     self.byteorder, self.datatype, first, expected, scale, reverse)
386    
387 caltinay 5012 def test_readCompressed2D(self):
388     NE = [9, 10]
389 caltinay 5008 domain = Rectangle(NE[0], NE[1], d1=0)
390 caltinay 5012 for filename, ftype in [("RectRedF%s.grid.gz", ReducedFunction),
391     ("RectConF%s.grid.gz", ContinuousFunction)]:
392 caltinay 5008 FS = ftype(domain)
393 caltinay 5051 filename = os.path.join("ref_data", filename%mpiSize)
394 caltinay 5012 unzipped = self.read(filename[:-3], FS, adjust(NE, ftype))
395     zipped = self.read(filename, FS, adjust(NE, ftype), True)
396     self.assertEqual(Lsup(zipped - unzipped), 0, "Data objects don't match for "+str(FS))
397 caltinay 5008
398 caltinay 5012 def test_readCompressed3D(self):
399     NE = [9, 9, 10]
400 caltinay 5008 domain = Brick(NE[0], NE[1], NE[2], d1=0, d2=0)
401     for filename, ftype in [("BrickRedF%s.grid.gz", ReducedFunction),
402     ("BrickConF%s.grid.gz", ContinuousFunction)]:
403     FS = ftype(domain)
404 caltinay 5051 filename = os.path.join("ref_data", filename%mpiSize)
405 caltinay 5012 unzipped = self.read(filename[:-3], FS, adjust(NE, ftype))
406     zipped = self.read(filename, FS, adjust(NE, ftype), True)
407 caltinay 5008 self.assertEqual(Lsup(zipped - unzipped), 0, "Data objects don't match for "+str(FS))
408    
409    
410     if __name__ == '__main__':
411     run_tests(__name__, exit_on_failure=True)
412    

  ViewVC Help
Powered by ViewVC 1.1.26