25th of January: Complex numbers are neither complicated nor non-existent

I've been asked about quantum computing. I had planned to write a little on the subject but I don't intend to go into incredible depth. In any case, before I can do so, I must cover a preliminary issue - complex numbers. This is important because quantum mechanics, which is the mathematical framework needed to describe quantum computing, cannot be formulated without them.

I know some people find the concept of complex numbers difficult. The complaint is that they seem artificial, but this shouldn't matter because in fact all of the systems of numbers we use in mathematics are artificial...

The natural numbers model our intuitive notion of counting (e.g. 1, 2, 3...), and then we extend this into the integers by inventing negative numbers (...-2, -1, 0, 1, 2...). The reason we invent negatives is to solve sums which don't have solutions within the positive numbers. More specifically, by introducing zero and assigning each natural number a corresponding negative number - an additive inverse - we end up with a system which is closed with respect to subtraction. Note that integers are not closed under division.

We can extend the integers further in order to describe quantities which are a part, or fraction, of something. We define a fraction as the ratio between two integers and call all such numbers the rational numbers. The rationals inherit their subtractive closure from the integers, but what about for division? we can assign to every rational r a corresponding reciprocal 1/r which will serve as a multiplicative inverse, giving us closure under division.

The rationals are still not sufficient, because there are numbers which cannot be written as a fraction of integers. Examples include the square root of 2, π, and e. So once again we extend this system to what are called the 'real numbers', purely because it is useful to do so. And I must emphasize this point: numbers like π do not exist in the real world! π is defined to be the ratio between the diameter and circumference of a unit circle, but measure a real circle and you'll only measure a finite decimal approximation to π. The real numbers are very good for modelling quantities in the universe and that is why we use them.

Now real numbers are weird. They do not appeal to our intuition of counting, and if you believe otherwise you probably don't understand them well enough. So the next step up - complex numbers - shouldn't shock you any more than the real numbers do (you are, of course, welcome to be shocked and confused by both. I only ask for consistency).

Let's solve some equations...

1. x² - 4 = 0


2. x² + 4 = 0

The first is easy: x = sqrt(4) = 2. What about the second? well: x = sqrt(-4). But wait, I can't do that, can I? after all 'sqrt(a)' is really just the solution to ?² = a, and how could and real number multiplied by itself ever result in a negative number such as -4?

Still, let's imagine that sqrt(-4) somehow does have a solution. We can be sure that such a solution isn't a real number, because there's surely no real number satisfying the equation, so I'll invent a new number which isn't a member of the set of real numbers (beware: 'real' is not an adjective, but a noun). I'll call it i, and say that i has the property that i² = -1. Now I can solve the second equation: x = 2*i (or just 2i). Let's check that works: (2i)² = 2²*i² = 4 * (-1) = -4. Yay!

A complex number is defined as follows: c = a + b*i, where a and b are real numbers and i is the 'imaginary unit' defined above. A nice geometric perspective is found by plotting a and b on a plane with a along the x-axis and b on the y-axis.

So why is i called imaginary? this is an unfortunate historical accident which has caused many misconceptions ever since. Ultimately it's Rene Descartes' fault, as he used this word in his writing in order to mock the idea. Most mathematicians at the time considered it to be absurd. At the time they might have been justified in thinking so. After all, what possible use could it have?

Well it turns out that if we use complex numbers then every polynomial equation of degree n has exactly n solutions, and this is important enough that it's known as the 'fundamental theorem of algebra'. For similar reasons complex numbers allow us to solve any linear differential equation of any order, which consequently has in applications physics, engineering, chemistry, biology, economics, and so on.

I haven't explained how to do algebra with complex numbers but I will introduce further concepts when needed. The rules are a little different compared with the reals, but since the real numbers are actually a subset of the complex numbers, all the algebra you learnt in school is still 'true' - it can all be derived as a special case of complex algebra.

08:00:00 - 3 Comments

22nd of January: Miscellaneous things

Brief personal post. Skip if boring.

I've been told to update this site by somebody who I'd have a hard time saying no to, so here it is. I have a few more interesting posts in the pipeline and I promise I'll get on with them. In fact, did you (plural) know I do requests?

• Moggers87 wants me to write stuff on fundamental mathematics for CS. I suspect Mr Zini will appreciate some of that too, so I will


• Mr Zini also wants to know about natural computation, so: neural nets, chemical computers (this is where I embarrass myself infront of the >= 2 chemists who read this blog), etc


• I request something on writing systems. Writing systems are fun.


• You will all have requested something about chaos theory. I know you haven't yet, but once you read it, you'll know that you would have requested it.
.

• I also want to expand on an old essay I wrote about ethics on Facebook, and then to have a huge flamewar debate in the comments.


• Finally, music. I have a new piano thingy in the works that might be recorded eventually. Anyone have a video camera that captures audio?

I've been quite busy but also battling against my general laziness. I like to think I'm becoming more productive this year, but I like to think a lot of things of which are of dubious truth.

New years was excellent spent with various Milton Keynesers, ex-Milton Keynesers, and one Canadian. Shortly after that there were exams. Only one (I know, you all hate me) - quantum computing. I went in there ready to power through it with all the maths I'd drilled into myself, only to find that half the marks required me to actually remember the details of certain algorithms. I was royally screwed, yet somehow passed anyway. This week I've been working on 3 assessments that are due in quite soon.

02:30:00 - 7 Comments

31st of October: England's not so ancient ancestry

I know I'm coming a bit late to this party, but last night I watched Question Time featuring Nick Griffin of the BNP - a nationalist political party here in the UK, for those who don't know). He was shown to be a liar, a racist and a homophobe (but I won't dwell on what we already know). Instead, something that made me laugh in particular was when he said that he represents the 'indigenous' English people, who he claims have an ancestry going back 17,000 years. He says 'science' backs him up on this, yet as the historian on the panel pointed out, at that time there was an ice age (actually a glaciation period) which meant that no humans inhabited England 17000 years ago. Neanderthals did inhabit Southern Europe, and what we would consider modern humans emerged from Africa towards the end of the ice age.

It's likely that a few waves of people have entered Europe at various times, partly or completely replacing the existing populations. I'm not sure where we can point to on a time line as the point at which English, Scottish, Welsh, and Irish (the ethnic groups he mentioned) people appeared, but it would be far more recent than Nick would like.

The majority of languages spoken in Europe are closely related and fall into a family known as Indo-European. These languages are hypothesised to descend from a common ancestor, named proto-indo-european. "indo?" you say... yep! because that family also includes Kurdish, Iranian, Pastho, Hindi, Urdu*, Punjabi, and many others - but don't tell the BNP that. Now the question of how that language family spread over its present geographic extent is quite controversial and we must be careful what conclusions we try to draw from it (without extra evidence from archaeology, history, etc).

Anyway, my point is this: the English language has only existed for maybe 2000 years at most, and its form has changed dramatically over that time. The speakers of what soon became English (Anglo-Saxons) came from Germanic tribes, but England was already inhabited before they turned up - by Celts (Irish, Scottish, Welsh, etc), and by speakers of Roman (Latin). The Celts, like the Anglo-Saxons, migrated here even earlier from mainland Europe, and there were already people here before them too - the Picts! and I'm sure the Picts weren't the first either...

Not to mention the fact that England was invaded by the Normans 1000 years ago, and has over the last few thousand years had all sorts of people coming in. We are an island, after all. So, who are the 'true English'? I'd wager that all of us have ancestry outside of England going back only a few generations, so the concept is essentially meaningless. What Nick Griffin wants to do is invent a false history that underlies a made up sense of identity, because he can't quite get away with flat out saying "white, blonde, blue-eyed people only". As much as an argumentum ad hitlerum makes me cringe, I'd have to point out that the Nazis tried this too.

I welcome any corrections to my history as it is only rough. My main purpose here was to discredit the absurd idea that the English people can trace their roots back 17,000 years as a continuous culture inhabiting these islands.

Now watch this.

* Urdu is classified as Indo-European because it has mostly Sanskritic origins, as does Hindi. However, unlike Hindi it also borrows a large amount of vocabulary from Arabic and Persian. The latter is Indo-European, but the former is not; it is Semitic, a branch of Afro-Asiatic. Naturally these languages are far closer to each other than their distant relatives in Europe, and they fall under a major IE branch called Indo-Iranian.

19:17:08 - 2 Comments

17th of October: Programmer maths

Good evening. So I'm in York now, and now that the chaos has settled into a sort of quiescence - for a few days, at least - I blog again.

I came across an article recently: "Programmers should know math.. just not all of it". It's a reasonable list but I found a number of commenters objecting very strongly to it. It may have been more appropriate to have said "computer scientists" rather than programmers (there's quite a bit more to CS than programming, and it all involves maths in one way or another), but even plain old programming requires a maths understanding beyond school level.

Certain specialist applications involve a lot of mathematics, like graphics, ai, and finance. But these aren't what I'm talking about. I can list a few things I use regularly in any kind of programming - even mundane work:

• Modular arithmetic: for cycling through a set of values


• Boolean algebra: for simplying conditional statements


• Set theory: for any kind of data handling


• Complexity analysis: for deciding if what I wrote is going to take longer than the lifetime of the universe to run (which I like to avoid)

Off the top of my head (there are surely more). All of these are covered in the first part of the article's list. And then if you're wanting to do anything that's actually interesting, then look to the second part of that list.

It's quite possible for someone to be a good programmer without any formal mathematics training. Such a person would likely still have a good intuitive grasp of many deep mathematics without being aware of it.

20:18:35 - 0 Comments

21st of September: Holidays in North Korea

Random thought: North Korea might be one of the safest places in the world for a foreigner to visit. While I can get stabbed in my country, in NK I'd be assigned at least one heavily armed guy to follow me everywhere I go and I'd be locked up in a hotel for most of the time anyway.

So long as I don't do anything stupid like try to incite an uprising or insult the Great Leader... what could possibly go wrong?

12:36:03 - 3 Comments

14th of September: Jump-starting computers and FreeBSD

This weekend I helped build a computer out of spare parts - which included "jump starting" it by shorting the power switch pins using a pair of scissors - and then built another one, also from spare parts. I've put FreeBSD on the second one. Don't know what Tomzini is doing with the first, but probably something wimpy like Windows XP (it has a real power switch now, in case anyone wonders. Also, there are photos of the scissor switch in action).

I've never used any BSDs, and in fact the only 'proper UNIX' I have used is Sun Solaris occasionally. One of the biggest culture shocks for FreeBSD has been that unlike the Linuxes I've used, it isn't configured as a desktop system out of the box. This is probably intentional as FreeBSD certainly comes from a different culture to, say, Fedora. I should also add that here are more desktop-friendly BSDs around like DragonFly BSD and DesktopBSD.

What I mean by the above is that upon first booting it presents the user with a text terminal (as an aside, text terminals in BSD so much nicer-looking than in Fedora). I had opted to install Xorg for graphics, but it doesn't start by default. Nor is there a login manager set up. Now I'm not afraid to dive in and set it up the way I want it myself, but this does get tiresome after a while and by then I would prefer to just make a choice in the installer, like "configure this as a desktop computer", or "configure this as a server", and then have it do the rest for me.

Anyway it isn't much work to change the default behaviour and I have learnt a few new things about Xorg along the way. The biggest problem I encountered actually turned out to be a hardware issue and nothing to do with FreeBSD. I had problems with installation as it gave errors while trying to install packages. After I had installed, I went back to sysinstall and tried again to add packages like Emacs, Bash, and Gnome. But once again there were errors.

I eventually worked out that it was something to do with my ancient DVD rom drive not being able to work with files above a certain size - a "READ_BIG MEDIUM ERROR" occured when I tried to install a package from the terminal. My solution was to extract the /packages directory from the disk onto my laptop - which had no trouble reading the disk - and then transfer them to the BSD machine over the network. Then I could use pkg_add to install anything I wanted.

Despite some initial frustration - which wasn't the operation system's fault - I actually quite like FreeBSD so far. The only thing that really bothers me is the installer. It certainly needs some work. The error messages it gave me were nothing short of useless, and though it keeps telling me to go to a 'debug screen', this is only actually visible upon first installation. Installation is something that modern Linuxes really do well, and FreeBSD seems behind here.

The next step will be to set up XMonad which I've been curious about for a while, and then I'll be using this as my main machine for a while to see how things go.

19:25:29 - 4 Comments

8th of September: Gambling, martingale, and other sillyness

I keep seeing this message in my Gmail spam. Looked at one way, it presents a really dumb gambling strategy, but looked at another way it's a really good money-making strategy (just not for you). Here's an excerpt:

yo mate, ok I`ll give you my trick but if you give it someone else I`ll fuckin kill you :)
you know in roulette you can bet on blacks or reds. If you bet $1 on black and it goes black you win $1 but
if it goes red you loose your $1.
So I found a way you can win everytime:

bet $1 on black if it goes black you win $1

now again bet $1 on black, if it goes red bet $3 on black, if it goes red again bet $8 on black,
if red again bet $20 on black, red again bet $52 on black (always multiple you previous lost bet around 2.5)
if now is black you win $52 so you have $104 and you bet:

$1 + $3 + $8 + $20 + $52 = $84 So you just won $20 :)

now when you won you start with $1 on blacks again etc etc. its always bound to go black
eventually (it`s 50/50) so that way you eventually always win. But there`s a catch. If you
start winning too much (like $1000 a day) casino will finally notice something and can ban
you. I was banned once on royal casino. So don`t be too greedy and don`t win more then $200
a day and you can do it for years. I think bigger casinos know that trick so I play for real
money on smaller ones, right now I play on lucky jackpot casino: [Website link removed] for more
then 3 months, I win $50-$200 a day and my account still works.

So, three things:

1. It's a kind martingale strategy and it won't work (more on this in a moment).
2. The email is formatted to look 'legitimate'
3. There's an address for a gambling website which the sender 'recommends' (where 'recommends' actually means 'owns').

Now I don't know much about gambling but I know a few things about probability. The error in the above strategy is due to what's called the 'gambler's fallacy', which is the belief that certain outcomes of independent trials are 'due', if they haven't occurred for a while.

A classic example is flipping a fair coin: there are only two outcomes, and each flip is independent of the previous ones. The probability of getting 10 heads in a row is very small - 1/1024. But imagine that this, by chance, does happen. Then I flip the coin an 11th time: what's the chance of it also being heads? By the gambler's fallacy one might be inclined to say that it "must come up tails this time" or "tails is about due". But this is incorrect; the probability of heads is still 1/2, because each coin flip is independent. Our intuition about these things can completely mislead us.

For the gambling trick in this email the idea is the same: that the desired outcome will eventually come up. Now it is true, of course, that it will come up eventually, but when this happens is completely unaffected by the number of times that it hasn't happened so far. Meanwhile it is suggested that the gambler place an exponentially increasing (repeatedly doubling) bet upon each failure, promising a large payoff when 'eventually' comes around. The problem is that either: you'll run out of money, or the house will have a limit which you'll hit, before you even get to that 'eventual' win. In short, the strategy is more likely to bankrupt you than make you rich.

I should point out that martingale strategies like this do work if you have infinite funding and no betting limits (but if you have an infinite supply of money then why would you even bother?).

Of course the people running the game are quite happy to let you go ahead and use your strategy, and the tone of the email makes it seem as though you are privileged to have been let in on this secret out of the goodness of the sender's heart - even though you've never met or heard of this person before. And who is that person? oh yeah, probably the one who owns the website that he suggests you try this trick on.

18:11:09 - 2 Comments

6th of September: Version 4

I must say I'm pleasantly surprised at the number of people who have asked me over the past year or so whether and when InsideReality would return. I stopped posting sometime in 2008 in order to concentrate on my final year at university. Now that I've made it out with a degree I intend to revive InsideReality.

Some of you will have seen the new design a long time ago and might wonder what took so long? I'm also aware that I promised it would all be finished in June, which I failed to keep. In truth it didn't really take very long at all - most of the code was written last week - I've simply been busy with many other things.

So, since I last posted here I've graduated from Newcastle, moved back to Milton Keynes for a few months and have got an offer to study for an MSc in Natural Computation at York university. Naturally I've accepted that offer; It's freaking York! they have ducks and everything. Some interesting things that you've missed (some of which might make for future blog posts): random wanderings round Newcastle at sunrise last summer, with photos to go with it; random cow visitations in Newcastle; PI day; hats and giant coats.

I don't know how frequently I'll be posting now but I do know that I'll spam the various social networks to let y'all know when I do update. I'm about to start an MSc which will certainly be a lot of work, but I feel like my time management has improved so we'll see how things go. The subject area - 'natural computing' as it's called - is one which is of deep interest to me and I hope to be able to write quite a bit about it over the coming year, along with the usual mixture of maths, computing, and philosophy (and the occasional bit of linguistics).

You'll notice content is a bit sparse at the moment. I'll hopefully be adding quite a bit more code over the coming week. I've dropped a lot of old blog posts just because they were either: outdated, poor quality, or silly ephemera. In any case, enjoy...

Update

Thanks for all the comments, feedback and bug reports. I have:

• Filled in author names on some old comments where they hadn't been copied properly


• Fixed garbled comments, also due to copying problems


• Corrected a formatting problem on code pages


• Put RSS feed back online


• Fixed hyperlinks on tag list page

I plan soon to have it store your names in a cookie when you comment so you don't have to retype them. Along the same lines, I might experiment with requiring a CAPTCHA to only be entered once and subsequently remembering (with cookies) that you're a human once you've proved it.

20:10:04 - 16 Comments

16th of August: Decentralize your data

I'm not entirely comfortable with the trend towards centralising larger amounts of personal data. For example right now you can use Google's services to store all your email, contacts, documents, newsgroup subscriptions, photos, and your calender.

These 'free' services are provided because they allow Google to sell targeted advertising placements to its real customers. But should we really trust Google with all that personal data? what if it were abused, or leaked, or what happens when the courts are able to force certain data to be disclosed?

This is one of several reasons I have for advocating more decentralisation, in order to counter the trend. I firmly believe that my data should remain on my computer alone, 'just in case'. I'm even going to move most of my email back to my webhosting provider and use Gmail less. Now of course, that may be less secure there than at Google (easier for some sysadmin to abuse and not get noticed), but at least it isn't being held indefinitely in a massive database, and sold on as a product to advertisers. More to the point, the messages will only stay on the server until next time I retrieve them through POP3 and delete the remote copies, meaning that the permanent copies reside on my computer alone, and placing me in control of them.

Likewise the data given away to social networking services is worrying. What's interesting is that people have started using these as an alternative to email, which is barely usable now thanks to spam. And the problem with that is that unlike email, these services are closed and proprietary, which would mean giving up control once more to large companies.

Before we had these, we had similar tools, like finger on Unix. Of course finger isn't searchable, which is one problem with decentralising: people like the convenience of being able to search for data, emails, people, etc.

Another aspect of this comes from DRM and similar efforts to control what can and can't be done with certain files. I'm not going to get into a long debate regarding music piracy, mostly because all of my music is paid for anyway. Right now there exist 3 digital copies of my whole collection, and I can make as many more as I like thanks to avoiding these sort of restrictions. I think that using open source software helps here, as DRM isn't exactly popular in the open source community.

Indeed I see Linux and other open source software as the perfect means of seizing back what is rightfully mine - my computer, and my data.

11:29:22 - 1 Comments

28th of June: Binomial distribution fun

So earlier today I needed to look up cumulative probabilities over a binomial distribution with n = 30... as you do. Problem was I couldn't find any tables online which went that far. The obvious solution was to use a spreadsheet, which I did, and it even had the function in-built. But then I was curious as to how to actually write an algorithm for it. More specifically, how to write a fast algorithm to do it.

The binomial distribution is used in a situation in which there is some event which can have one of two outcomes. The event is usually termed a 'trial', and one of the outcomes is defined to be 'success' while the other is 'failure' (though the terminology shouldn't be taken too literally). A probability p is assigned to a successful outcome, and consequently, q = 1 - p is the probability of failure.

If a trial is repeated n times then we can ask the question "what is the probability of k successes out of those n trials?". To work this out, the following reasoning is used:

1. There is an event which has probability p of occurring. The chance of it happening k times in a row is p^k.


2. If that event does happen k times during n trials, then the probability that it doesn't happen in any of the remaining n-k trials is q^(n-k) where q = 1 - p.


3. During n trials there are many different ways of arranging the k successes. For example, they could all happen in a row, or they could happen alternately in-between failures. In total there are n! / (n-k)!k! ways/orderings in which k successes may be selected from n trials. We will write this as nCk (the binomial coefficient).

Putting that together we get:

Pr(K = k) = nCk * p^k * q^(n-k)

Now for a more interesting question: what is the probability that there are at most k successes?

The solution is to add together the probabilities of 0 successes, 1 success, 2 successes, and so on, up to k successes, and this is called the cumulative probability. Naturally it would be far too tedious to calculate by hand, and that is where the computer program comes in handy.

As for the algorithm, it turns out (as can be expected) that computing nCk for each probability causes some problems. The most elegant way in my opinion of finding nCk is to use the following recursive function, where binom(n, c) is the same as nCk.

binom(n, 0) = 1
binom(n, n) = 1
nCn = binom(n - 1, k - 1) + binom(n - 1, k)

Which can be found by looking for a formula to find numbers on Pascal's triangle. But though this looks very nice, it also becomes very slow for a computer as (n - k) gets larger. So instead the solution used here is non-recursive and works by repeatedly doing a multiplication followed by a division. This is not only much faster than a recursive solution but avoids overflow by keeping the intermediate result small.

Program source code.

02:16:02 - 0 Comments

22nd of February: Faster Fibonacci calculation

The Fibonacci sequence is a sequence of numbers where each number is the sum of the preceding two. We begin the sequence with [1, 1] and obtain:

1, 1, 2, 3, 5, 8, 13, 21, 34, 55, 89, 144, etc...

Looking at the numbers in this sequence you will see that each value is the sum of the previous two values. The exception is the first two values, [1, 1]. We define these values initially and the rest of the sequence grows from there.

We can describe the Fibonacci sequence nicely using the following recurrence relation:

fib(1) = 1;
fib(2) = 1;
fib(n) = fib(n - 1) + fib(n - 2);

So to find the nth Fibonacci number we add together the n-1th and the n-2th. Following that rule we will eventually hit both of the base cases, fib(1) and fib(2), which have been defined explicitly. The recursion ends at the base case.

Now programming this description into a computer is quick and easy, but it is also a horribly inefficient method of calculating these numbers. Aside from the extensive use of stack space, it is also calculating all the intermediate values twice. The number of operations performed for the first method grows exponentially in n – its exact complexity is O(φn) where φ = (1 + sqrt(5)) / 2 (the Golden Ratio).

A better method is an iterative one. We write an algorithm which looks something like this:

Integer fib(Integer n)
begin:
	Integer a := 1;
	Integer b := 0;
	for Integer i :=  1 to n
		//Compute next a, b values
		Integer next_a := a + b;
		Integer next_b := a;
		//Update: this really only needs one temporary variable but I used two for clarity
		a := next_a;
		b := next_b;
	next i
	return b;
end.

This is simple enough and significantly faster than the recursive algorithm. There are two assignments at the start, and each loop iteration consists of 4 assignments. So 2 + 4n operations, which comes to O(n). But can we improve on the linear algorithm? Well of course, otherwise I wouldn't be writing any of this in the first place. Some algebraic manipulation will provide us with a better algorithm:

The second algorithm applies the following mapping n times:

a <- a + b
b <- a

With an initial state (a, b) = (1, 0).

This mapping is a special case of a more general mapping, which we will call T_pq:

a <- bq + aq + ap
b <- bp + aq

Setting p = 0 and q = 1 we arrive at the above Fibonacci mapping. It can be shown that the result of T²_pq, that is, applying T_pq twice, is the same as applying T_p'q' once. This second mapping would give us every 2nd pair of Fibonacci numbers, rather than every consecutive pair. If we can skip over some numbers while looking for the nth Fibonacci number then we will certainly get there faster.

If p' and q' can be expressed in terms of p and q then we can build an algorithm which calculates a Fibonacci number by 'homing' in on the result. Jump some distance, then jump half that distance, and half again, until reaching the desired value.

And so here's the promised algebra... (transcribed from notes made yesterday)

1) Apply T_pq once

a = bq + aq + ap = q + p
b = bp + aq = q

(initial values a = 1, b = 0, were substituted in to produce the final expressions)

2) Apply T_pq again

a = q² + (q + p)q + (q + p)p = 2p² + 2pq + p²
b = qp + (q+p)q = pq + q² + pq = q² + 2pq

So these are the values of a and b after applying the mapping twice. Replace q with q' and p with p'

3) Solve for p' and q'

2p² + 2pq + p² = p' + q'
q² + 2pq = q'

Here I compare the result from mapping once using the values p' and q', with mapping twice using the values p and q.

4) Then with some rearranging

q' = q² + 2pq
p' = q² + p²

And that's the algebra over with. Now using p' and q', an algorithm is born...

Integer fib(n)
begin:
	//Set up initial state
	Integer a := 1;
	Integer b := 0;
	Integer p := 0;
	Integer q := 1;
	Integer nc := c;

	while(nc > 0)
		nc := nc / 2; //Discarding any fractional part
		Integer next_q := q*q + 2*p*q;
		Integer next_p := q*q + p*p;
		
		q := next_q;
		p := next_p;
		
		if nc is odd then
			Integer next_a := b*q + a*q + a*p;
			Integer next_b := b*p + q*q;

			a := next_a;
			b := next_b;
		end if
	repeat

	//Return result
	if(n is even)
		return b;
	else
		return a;
end.

A couple of notes: when n is divided by two only the integer part should be kept. Secondly, notice that the loop results in a pair of integers (a, b). Which one we return as the result depends on whether the input, n, is even or odd. We have 5 assignments before the loop. The loop itself will repeat until n is equal to zero, and n is repeatedly divided by two on each iteration of the loop. The number of times the loop repeats is thus given by log₂(n), or intuitively, the number of times we can divide n by two before it n becomes zero. We have 5 + log₂(n) operations, giving a time complexity of O(log₂(n)).

This algorithm is basically an imperative style solution to exercise 1.19 in The Structure and Interpretation of Computer Programs - a book every computer scientist should read. The exercise is to find p' and q' in terms of p and q, which is what I shown above. I have written an implementation in c (source code). The book has a nice Scheme implementation of a similar algorithm.

01:48:09 - 8 Comments

16th of February: Red-eyed sky

Time for another post... just to show off another fractal image. I've been experimenting with full-colour rendering, but my colouring algorithm is still giving odd results. Still, I like this:

22:46:23 - 0 Comments