On 1/27/07, Keith Goodman <[EMAIL PROTECTED]> wrote:
I get slightly different results when I repeat a calculation. In a
long simulation the differences snowball and swamp the effects I am
trying to measure.

Here's a unit test for the problem. I am distributing it in hopes of
raising awareness of the problem. (What color should I make the
Repeatability Wristbands?)

I am sure others are having this problem without even knowing it.
"""Unit tests for repeatability.

The problem and possible causes are discussed in the following thread:
http://projects.scipy.org/pipermail/numpy-discussion/2007-January/025724.html

If the tests do no pass, then try removing ATLAS-SSE2 from your system. But
similar problems (not tested) may occur unless you remove all versions of
ATLAS (SSE and BASE, for example).
"""

import unittest
import numpy.matlib as M
M.seterr(divide='ignore')
M.seterr(invalid='ignore')

def calc(x, y): 
    return x * (x.T * y)
    
def load():

    # x data
    x = M.zeros((3,3))
    x[0,0] =  0.00301404794991108
    x[0,1] =  0.00264742266711118
    x[0,2] = -0.00112705028731085
    x[1,0] =  0.0228605377994491
    x[1,1] =  0.00337153112741583
    x[1,2] = -0.00823674912992519
    x[2,0] =  0.00447839875836716
    x[2,1] =  0.00274880280576514
    x[2,2] = -0.00161133933606597
         
    # y data
    y = M.zeros((3,1))
    y[0,0] = 0.000885398
    y[1,0] = 0.00667193
    y[2,0] = 0.000324727
    
    return x, y    

class Test_repeat(unittest.TestCase):
    "Test repeatability"
    
    def setUp(self):
        self.nsim = 100
        self.x, self.y = load()
        
    def test_repeat_1(self):
        "repeatability #1"
        z0 = calc(self.x, self.y)
        msg = ''
        result = True
        for i in xrange(self.nsim):
            z = calc(self.x, self.y)
            if (z != z0).any():
                msg =  'Max difference = %g' % abs(z - z0).max()       
                result = False 
                break
        self.assert_(result, msg)  

    def test_repeat_2(self):
        "repeatability #2" 
        z0 = calc(self.x, self.y)
        msg = ''
        result = True
        for i in xrange(self.nsim):
            x, y = load() 
            z = calc(x, y)
            if (z != z0).any():
                msg =  'Max difference = %g' % abs(z - z0).max()       
                result = False 
                break
        self.assert_(result, msg)      

    def test_repeat_3(self):
        "repeatability #3"  
        z0 = calc(100*self.x, self.y) / (100 * 100)
        msg = ''
        result = True
        for i in xrange(self.nsim):
            x, y = load() 
            z = calc(100*x, y) / (100 * 100)
            if (z != z0).any():
                msg =  'Max difference = %g' % abs(z - z0).max()       
                result = False 
                break
        self.assert_(result, msg) 

 
def testsuite():
    unit = unittest.TestLoader().loadTestsFromTestCase
    s = []
    s.append(unit(Test_repeat))    
    return unittest.TestSuite(s)

def run():   
    suite = testsuite()
    unittest.TextTestRunner(verbosity=2).run(suite)
_______________________________________________
Numpy-discussion mailing list
[email protected]
http://projects.scipy.org/mailman/listinfo/numpy-discussion

Reply via email to