Sustainable Scientific Software Development

Europython 2017

Alice Harpole

soton logo soton logo

import unittest

def squared(x):
	return x*x

class test_units(unittest.TestCase):

    def test_squared(self):
        self.assertTrue(squared(-5) == 25)
        self.assertTrue(squared(1e5) == 1e10)
        self.assertRaises(TypeError, squared, "A string")

data = rand(80,80)                # declare some random data

def func(a):                       # function to apply to data
    return a**2 * numpy.sin(a)

output = func(data)               # calculate & plot some function of random data
plt.imshow(output);  plt.colorbar();   plt.show()
Input is $0 \leq x \leq 1$, so output must be $0 \leq f(x) \leq \sin(1) \simeq 0.841$ $$\overline{f(x)} = \int_0^1 f(x) \,dx \simeq 0.223$$

def test_limits(a):
    if numpy.all(a >= 0.) and numpy.all(a <= 0.842): return True
    return False

def test_average(a):
    if numpy.isclose(numpy.average(a), 0.223, rtol=5.e-2): return True
    return False

if test_limits(output):
	print('Function output within correct limits')
else:
	print('Function output is not within correct limits')
if test_average(output):
	print('Function output has correct average')
else:
	print('Function output does not have correct average')
Function output within correct limits
Function output has correct average

	# use trapezium rule to find integral of sin x between 0,1
	hs = numpy.array([1. / (4. * 2.**n) for n in range(8)])
	errors = numpy.zeros_like(hs)

	for i, h in enumerate(hs):
	    xs = numpy.arange(0., 1.+h, h)
	    ys = numpy.sin(xs)

	    # use trapezium rule to approximate integral of sin(x)
	    integral_approx = sum((xs[1:] - xs[:-1]) *
				  0.5 * (ys[1:] + ys[:-1]))
	    errors[i] = -numpy.cos(1) + numpy.cos(0) - integral_approx

	plt.loglog(hs, errors, 'x', label='Error')
	plt.plot(hs, 0.1*hs**2, label=r'$h^2$')
	plt.xlabel(r'$h$'); plt.ylabel('error')