Recently on the blog I've been looking at some neat mathematical algorithms that we can implement in Sniff. The code to calculate Pi and find square roots are probably beyond Sniff's target audience in terms of stuff you'd expect them to know, but the code is simple, and uses maths that's simple enough that it should be reasonable for a strong KS4 student to grasp whats going on if you choose to use them as exercises.

Today I'm going to do something that's a part of the A level further maths syllabus - Taylor approximations. These are definitely more advanced, as they're derived from calculus. To understand them you need to know about derivatives, and even to use them you need radians so they're not really going to be much use for younger kids. However they're something that older students should absolutely be encouraged to implement in code: they're tricky to calculate by hand, somewhat confusing, but yield great results when you actually run them on a computer.

I'm not going to explain the maths behind these, but will skip to the take home:

We can approximate complex functions like sin() with polynomials. In the case of sin() it works best when x is a small number (close to 0), but if we keep adding terms then the approximation just keeps getting better and better.

To implement that in Sniff there are a couple of gotcha's: there's no power function, and no factorial function! This makes the implementation more difficult , but actually forces us to think about the problem and produces a better final solution.

Lets start with a chat about factorials. Every programming class teaches you to implement factorials as:

def f(x)

if x=1

return 1

return x*f(x-1)

clever huh? we define factorial of X based on factorial of X-1. Eventually X gets down to 1, which we already know the answer to (hint: its 1! by which I mean 1=1! by which I mean 1!=1! I mean!!!!! you get the idea!). This is recursion and its a big thing in programming classes.

Well not so fast... Scratch and Sniff can't do recursion so lets look at another way to write that:

make f number

when factorial

.make i number

.set f to 1

.repeat x using i

..set f to f * i

Here we've written the factorial without using recursion. Now its not "clever" - its just the

*obvious*way of calculating 1x2x3x4x5x... and so on. If you've never been taught how to write a factorial, this is probably what you would probably write. Not only is it clearer, and simpler, its*better!*Each time the function f is called it creates a copy of x, and a bunch of other data. To find f(100) the recursive function calls*itself*100 times, making 100 copies. Making those copies is slow, and uses up a lot of memory.
So lets get this clear: the fancy version that I was taught in programming classes is slower, uses more memory, and is harder to understand than the regular and "bleedin' obvious" version. Yep, sure is! The recursive version has NO advantages - the regular version is better in every respect. Advocates of recursion will point out that a really good compiler can maybe claw back the performance penalties, but that's really not the point.

The thing is recursion is a powerful idea. It's important, and does have uses (as an undergraduate I remember making fun of the physics students struggling to implement dozens of lines of code in Fortran that could be written elegantly in 6 lines of Pascal using recursion), but its also tricky to understand. Teachers in programming classes use factorial as an example because its simple, and they can't think of a

*good*example of when you should use recursion, so they use a*bad*example. Kids: if you're ever in a programming class and get "taught" to write a factorial, ask why you can't use the simpler and faster version!
But lets get back to that sin function - in fact its going to be better not to break out a factorial function at all. Look again at that equation. For each term the the denomentator factorial increase by 2. The fastest way to calculate the denominator for a term is to perform two multiples on the previous denominator - much quicker and easier than than doing a full factorial every time.

Once we spot that we release that power of the numerator increases by two so we can get the next numerator by multiplying the current one by x*x - much cheaper than a power.

If we put that together we get:

make x number

make s number

make s number

when approxSin

.make numerator number

.make denominator number

.make count number

.set s to x

.set numerator to x

.set denominator to 1

.set count to 1

.repeat 2

..set numerator to -numerator*x*x

..set denominator to denominator*(count+1)*(count+2)

..change count by 2

..

..change s by numerator/denominator

As this is an "advanced" example - at least mathematically, I've used local variables. These are optional, but once you realise you need them they're a big help - by defining variables inside a script we don't have to worry about other code using them at the same time - each script gets its own version.

The numerator starts of as x and the denominator as 1. Each time around the loop we multiply the numerator by -x*x. This not only increases the power by two, but flips the sign.

For the denominator we use the variable count to tell us how far the sequence of factorials we are. count starts as 1, then to get to the second term we multiply by 2 and 3. We increment count by 2 to indicate that the denominator is now 3!

Having calculated the new numerator and denominator we just add the fraction to the running total

And that's it. We go around a few times, and get the right answer. I've set repeat to 2 here (one less than the equation) as this is enough to get useful answer, but leaves enough error to see that there's something interesting going on.

when start

.make count number

.repeat 10 using count

..set x to (pi/2)*0.1 *count

..broadcast approxSin and wait

..say join "x"[x]

..say join "approx:" [s]

..say join "sin(x):" [sin of (x/pi*180)]

..say join "error:" [s- sin of (x/pi*180)]

..say ""

This script just tests it. We're working in radians, so set x to count from 0 to pi/2 in to 10 steps, calculating the approximation, the correct answer and the error. Note we have to use degrees with the built in Sniff sin function.

x0.15708

approx:0.156434

sin(x):0.156434

error:-1.3411e-07

x0.314159

approx:0.309017

sin(x):0.309017

error:-1.78814e-07

x0.628318

approx:0.587792

sin(x):0.587785

error:7.21216e-06

x0.942477

approx:0.809146

sin(x):0.809017

error:0.000128925

x1.25664

approx:0.952017

sin(x):0.951057

error:0.000960231

x1.5708

approx:1.00452

sin(x):1

error:0.00452483

Pretty impressive! The results are more accurate with small values of x (due to the way the approximation is derived), but even at pi/2 (90 degrees) the error is only 0.004. If we increase the number of iterations by only 1 additional term we can reduce that to 0.0001 which is probably good enough for a lot of applications (in fact Sniff uses code very similar to this to implement trig functions on the Parallax Propeller - the code uses much less memory than the official sin functions).

It's worth noting there's nothing here that couldn't be done just as easily in Scratch. I've used repeat using, and local variables which Scratch 1.4 doesn't have but they're purely for convenience, and the code would work just as well without them. So next time someone says that Scratch is too simplistic to do real work, as them if they've done Taylor approximations (or Fourier transforms, or Monte Carlo Simulation)

## No comments:

## Post a comment