This was quite an interesting puzzle. It can be solved very quickly using an “idiot level solution” but at the same time it’s open to almost infinite optimisation. I chose to take the idiot level solution, since it still ran in 8 thousands of a second of CPU time, which is fast even for C. Very simple idea, we write a bang function to calculate factorials, we write an ncr function (that now I think of it should probably just have been a macro) to calculate the NCR value, which is just 3 factorials and a division. For factorials only going up to 100 we can afford to just calculate it every time. Then a pair of for loops try every value of n and r, incrementing a counter if the result is over 1000000. Print out the increment at the end. Didn’t even have to recompile.

Code can be found here: https://gist.github.com/2340044

Problem description can be found here: http://projecteuler.net/problem=53

To make up for the tiny amount of time I spent actually coding on this puzzle, let’s think about some potential optimisations that could be used to speed this program up, if we were being required to work out a similar piece of information for values up to (for example) 1000000 (which would be about 100000000 times more work for the computer I think). Firstly, and perhaps least helpfully, we can afford to push up the initial n value to 23, as stated in the problem description, the reason I chose not to do that here is that I’m always  wary about being right on the upper or lower bounds when I don’t need to be. I only have to be out by one or get an equality test slightly wrong and I’m left with an almost untraceable bug. Secondly, the r value could start a little higher in our tests. Without doing any real maths it’s obvious that we could start at at least 2, because when r=1, nCr=n. Similarly the r loop could stop a little sooner as well, because the last few values are just as small as the first few.

Moving onto actual optimisation rather than just bounds tightening, I have 3 main ways that the code could be sped up. Firstly, nCr values are symmetrical, ie for any given n, nCr and nC(n-r) have the same value. This cuts in half the number of actual calculations we have to do whilst still retaining 100% accuracy. Secondly, the most system intensive part of these calculations is the factorial function, which is very difficult to optimise in itself, there isn’t really any fast way to calculate factorials. However, over the course of this program, we calculate the factorials of the same number again and again. We can quite easily  store the result every time we do a calculation, then on future runs we just pull up the cached version rather than recalculating it. I’m not sure, but given how fast this code ran, I think GCC might actually have spotted this optimisation for me and implemented it, which is both clever and helpful of it.

Finally, probably the greatest possible optimisation we could have, is to remember that the nCr values are just the numbers from the appropriate  rows and columns of Pascals Triangle. Though the given formula is the fastest possible way of calculating an nCr value for arbitrarily large values of n and r, when you’re calculating it for a lot of consecutive values at the same time, you can just draw the triangle in memory and use that. Drawing the triangle like this reduces the calculation of each value to a single addition. This is made even faster by the fact that once we get a value over 1000000, we know all the values that come below it in the triangle are also over 1000000, so we can just add all of those to our count, and move on.

I’m sure there’s also a lot of possible mirco-optimisation but that’s usually not worth bothering with, GCC is much better at that sort of thing than programmers are, and it usually results in messy hard to maintain code.

On a related note, this was the 50th problem I solved on Project Euler, making me officially Level 2! A great achievement I think!

The problem description can be found here: http://projecteuler.net/problem=52

A very short problem can be solved with very little code. Mine can be found here: https://gist.github.com/2273876

54 lines of very spacious code. The primary challenge with this one was finding a good fast way to get the constituent digits of a number. My initial algorithm used an idiom that I fall back on quite a lot when Project Euler demands that I do something with the digits of a number, which is something like

while(n>0)

{

do_something_with(n%10);

n/=10

}

While I was programming I got worried that this wouldn’t account for trailing zeros (which turned out not to be a problem since there aren’t any in the answer, but I felt it had to be dealt with none the less). This was me being stupid as it fails to deal with LEADING zeros, which are generally ignored in Project Euler problems (as they are usually in all maths, since they don’t change the value of the number).

I also had to make a little use of the string.h memory function that I don’t especially like, but they were the easiest way to clear and compare the array’s holding the digit counts.

With the exception of minor syntax errors that were 2 second fixes (read the compiler error, go to the indicated line in the code, facepalm, fix) like having typed \ when I meant / there weren’t really any bugs this time around which I’m pretty pleased with. It more or less ran first time.

Added textures to my OpenGL cube today. Not too much difficult programming this time, since all the 3D magic has already been dealt with in the last problem, and that carries most of the hard maths. All my refactoring in the last section paid off quite a lot this time as well as all the code was beautifully easy to edit, no real bugs except a few mistyped types and attempts to run before I’d made all the necessary edits. There is a slight problem with warping, which I actually think is mostly just a result of not having a precise enough depth buffer and the like. Not going to reupload the code this time because there’s so much of it and it’s largely the same as last time. The big changes were adding individual vertex information for each corner for each side it appears on because they can’t be reused as they have different texture coordinates depending on which face they’re on. Also had to get the image into code which was quite easy since gimp provides a facility to save an image as C code.

A video can be found here http://youtu.be/MPniGcCiMmQ

On a separate note I’m not sure I’m going to bother with more  OpenGL work. It’s not really a section of computer science that I have that much interest in, despite how cool it looks, and it’s a lot of work to get very simple things happening. Potential things I might do in the future include learning C++ (no matter how much I hate the language it is very popular), more project euler, perhaps a little bit of GUI work with Qt once I’ve done C++ and also trying to find a more involving project that Project Euler that is still challenging rather than grunt work. The trouble with Project Euler is that each of the individual problems is too short to really get into and enjoy. Hopefully I’ll be able to find something but I’m not sure!

Actually quite an easy one to solve this time, though with a couple of lynch pins to get past. The biggest problem is the size of the numbers being dealt with. Took my a couple of goes to realise  that I needed to use not just a 64 bit variable to store my numbers, but one that didn’t waste a bit with a sign. I did spend some time messing around with trying to get bigger variables but didn’t have much success (the lowest ceiling I could determine analytically was just under 100 bits, so I’d have needed a 128 bit variable). Some research showed that this fairly enormous variable could also be combined with the fact that although the problem recommends 50 iterations, only about 25 are needed. I used 30 to save myself time working out the exact number.

The main other problem I has was slipping in how consistent I was with variables, so occasionally things got converted without me noticing that messed things up (causing overflow etc.).

Speed isn’t really a problem, even most interpreted languages should be able to brute force this one in under a second, never mind the one minute rule quoted by Project Euler, and I’m using C which is several orders of magnitude faster than most of those languages. Besides a very quick perusal of the solvers forum for this problem shows that no one was really using a better algorithm, only micro optimisation, which is both boring and relatively ineffective so not worth time worrying about in an academic problem like this.

My favourite part is probably the reverse function. It’s very elegant for some reason I can’t put my finger on .

Code is available in this gist: https://gist.github.com/2053561

Problem description can be found here: http://projecteuler.net/problem=55

This is really the point of OpenGL. There are much better ways of doing 2D graphics that don’t have anywhere near this level of complexity, but using OpenGL makes 3D graphics trivial compared to how difficult they would be to do manually. Unfortunately not trivial compared to regular 2D graphics, which is why this has taken so long. It’s taken me ages to work out all of the bugs I had turning my triangle into a cube, and OpenGL bugs are immensely hard to deal with because so much of the hard work is happening where you can’t see it to debug it, and most bugs only become visible well after the actual mistake is made, where the data that’s been generated is used or something similar.

I made this more difficult for myself by insisting on using my own code to handle the somewhat complicated maths that needs to be done on top of OpenGL, calculating matrices for example to convert coordinates relative to a model defined in code to coordinates that give pixel positions and depth. The primary reason for this was that the tutorial from which I’m learning OpenGL uses a library called GLM to do this, which is a C++ library, while I want to keep using raw C. All of the alternatives to GLM have a very different architecture to it, so I decided to build my own library with a similar structure to GLM (stipped to only the functions I use) but using C structures and functions rather than C++ classes and methods. This created an enormous new place for bugs to hide where it wouldn’t be easily apparent. In the end I had to get GLM and write a number of C++ programs to cross reference whether my functions where working, which for the most part they were (though there were a few radian-degrees mix-ups, none of them actually broke the programs, just meant the code and I were measuring angles differently).

The two biggest bugs were that the matrices I generated using this code were transposed from the way they were supposed to be (they were row major, not column major), which was solved with help from some nice people on the OpenGL forums, and the fact that Allegro, the library I’m using to give me an OpenGL context (basically a window I can draw to using OpenGL functions), was giving me an OpenGL context with no depth buffer (or, a depth buffer with zero bits of precision, which is just as bad) because it’s primarily meant to be a 2D graphics library so doesn’t care about depth, a bug which I spotted by comparing my code to a known working implementation of a similar cube that was different in a huge number of ways, and fixed using the advice of elias from the #allegro IRC channel.

The final code is available in this gist and a video of the cube spinning around can be found here.

Next step is to add textures to the cube, allowing me to show any image on its side rather than just some interpolated colours. I’ll start that more or less straight away, but in the event I get caught up again having to deal with horrible bugs I’ll mix in some work on Project Euler so this doesn’t just die for several months.

My first post on this blog about Duke of Edinburgh activity was a link to my entry to Ludum Dare 21, which happened a few months ago. Just last weekend Ludum Dare rolled around again and I decided to enter again. I did fairly badly I think because I made a complete waste of the first day and had to make my entire build in half the normal time. I think had started from where I did on the second day I might even have got some art/sound assets into the build this time.

Also this time I chose to make a time lapse of my programming on both days, recorded in two separate videos because I effectively started a separate project the second day.

My entry page can be seen here, with links to download windows and linux binaries of the game, the source and read my entry journal.

I decided to do a bit more work on the javascript counter for NaNoWriMo. It now does a fair bit more maths and is vastly more interactive. It can be played with here. The actual code used to make it can be seen using most browser by pressing Ctrl-U (it’s something else on IE, but then there’ no good reason to be using IE).

Hopefully a return to OpenGL by next week.

This weeks project is actually ridiculously simple, but it did take me quite a while. The reason for that is that this program was written in javascript, a language I don’t usually use and have never written more than a hello world program in. In order to write this I was googling just about every function as I typed it. The code is very simple, it tells the view of a web page how many words they would have to have written to be on track to win NaNoWriMo (a speed writing competition where those who enter aim to write 50,000 words in November).

In order to give the number of decimal places of accuracy that the code does, it calculates the number of milliseconds since November started, divides that by the number of milliseconds in November and multiplies by 50,000. It then uses the DOM to drop that number and a little relevant text into a span tag on the page. If it’s not November we just point out that NaNoWriMo isn’t currently happening.

A running example can be found here and the source can be viewed by using a browser to look at the source of that page.

I’ve been learning a bit about OpenGL recently. After a false start with OpenGL 1 (out of date now, though still in use) I got going and have a few problems left but am loving it so far. The current big problem I’m having is that GLM (the matrix math library built with OpenGL in mind0 is C++ only and I use C. I may even end up writing my own matrix library to handle it. None the less it’s going well so far and this is the most recent fruit of my labors, it’s a simple triangle in 3 colours that fades in and out over time. Though this is a ridiculously simple concept it actually requires quite a lot of work in the background because of the way OpenGL works.

Main C file:

#include <stdio.h>
#include <allegro5/allegro.h>
#include <math.h>

#define GL_GLEXT_PROTOTYPES
#include <GL/gl.h>
#include <GL/glext.h>

#define FPS 60
#define SCREEN_W 640
#define SCREEN_H 640

GLuint program;
GLint attribute_coord2d;
GLint attribute_colour;
GLint uniform_fade;
GLuint vbo_triangle;
unsigned int fadeval=0;

GLint makeshader(char * filename, GLenum type);
void print_log(GLuint object);

int init_resources();
void draw();

int main()
{
	ALLEGRO_DISPLAY *display = NULL;
	ALLEGRO_EVENT_QUEUE *event_queue = NULL;
	ALLEGRO_TIMER *timer = NULL;
	bool redraw=true;

	if(!al_init()) {
		fprintf(stderr, "failed to initialize allegro!\n");
		return -1;
	}

	timer = al_create_timer(1.0 / FPS);
	if(!timer) {
		fprintf(stderr, "failed to create timer!\n");
		return -1;
	}

	al_set_new_display_flags(ALLEGRO_OPENGL);
	display = al_create_display(SCREEN_W, SCREEN_H);
	if(!display) {
		fprintf(stderr, "failed to create display!\n");
		al_destroy_timer(timer);
		return -1;
	}

	event_queue = al_create_event_queue();
	if(!event_queue) {
		fprintf(stderr, "failed to create event_queue!\n");
		al_destroy_display(display);
		al_destroy_timer(timer);
		return -1;
	}

	al_register_event_source(event_queue, al_get_display_event_source(display));

	al_register_event_source(event_queue, al_get_timer_event_source(timer));

	al_clear_to_color(al_map_rgb(0,0,0));

	al_flip_display();

	al_start_timer(timer);

	init_resources();

	while(1)
	{
		ALLEGRO_EVENT ev;
		al_wait_for_event(event_queue, &ev);

		if(ev.type==ALLEGRO_EVENT_TIMER)
		{
			redraw=true;
			fadeval++;
		}
		else if(ev.type==ALLEGRO_EVENT_DISPLAY_CLOSE)
		{
			break;
		}

		if(redraw && al_is_event_queue_empty(event_queue))
		{
			redraw=false;

			draw(fadeval);
			al_flip_display();
		}
	}

	return 0;

}

int init_resources()
{
	GLint link_ok = GL_FALSE;
	GLuint vs;
	GLuint fs;

	glEnable(GL_BLEND);
	glBlendFunc(GL_SRC_ALPHA, GL_ONE_MINUS_SRC_ALPHA);

	GLfloat triangle_attributes[] = {
		0.0 , 0.8, 1.0, 0.0, 0.0,
		-0.8, -0.8, 0.0, 1.0, 0.0,
		0.8, -0.8, 0.0, 0.0, 1.0
	};

	glGenBuffers(1, &vbo_triangle);
	glBindBuffer(GL_ARRAY_BUFFER, vbo_triangle);

	glBufferData(GL_ARRAY_BUFFER, sizeof(triangle_attributes), triangle_attributes, GL_STATIC_DRAW);

	glBindBuffer(GL_ARRAY_BUFFER, 0);

	if(( vs = makeshader("vertex.v.glsl", GL_VERTEX_SHADER)) == 0)
		return 0;
	if(( fs = makeshader("fragment.f.glsl", GL_FRAGMENT_SHADER)) == 0)
		return 0;
	program = glCreateProgram();
	glAttachShader(program, vs);
	glAttachShader(program, fs);
	glLinkProgram(program);
	glGetProgramiv(program, GL_LINK_STATUS, &link_ok);
	if (!link_ok) {
		fprintf(stderr, "glLinkProgram:");
		print_log(program);
		return 0;
	}

	const char* attribute_name = "coord2d";
	attribute_coord2d = glGetAttribLocation(program, attribute_name);
	if (attribute_coord2d == -1) {
		fprintf(stderr, "Could not bind attribute %s\n", attribute_name);
		return 0;
	}
	attribute_name = "v_colour";
	attribute_colour = glGetAttribLocation(program, attribute_name);
	if(attribute_colour == -1)
	{
		fprintf(stderr, "Could not bind attribute %s\n", attribute_name);
		return 1;
	}

	attribute_name = "alpha";
	uniform_fade = glGetUniformLocation(program, attribute_name);
	if(uniform_fade == -1)
	{
		fprintf(stderr, "Could not bind uniform %s\n", attribute_name);
		return 0;
	}

	return 1;
}

void draw(unsigned int fadeval)
{
	/* Clear the background as white */
	glClearColor(1.0, 1.0, 1.0, 1.0);
	glClear(GL_COLOR_BUFFER_BIT);

	glUseProgram(program);

	glBindBuffer(GL_ARRAY_BUFFER, vbo_triangle);

	glEnableVertexAttribArray(attribute_coord2d);

	/* Describe our vertices array to OpenGL (it can't guess its format automatically) */
	glVertexAttribPointer(
		attribute_coord2d,      // attribute
		2,                      // number of elements per vertex, here (x,y)
		GL_FLOAT,               // the type of each element
		GL_FALSE,  	        // take our values as-is
		5 * sizeof(GLfloat),   //  stride
		0		        // offset of first element
	);
	glEnableVertexAttribArray(attribute_colour);

	glVertexAttribPointer(
		attribute_colour,
		3,
		GL_FLOAT,
		GL_FALSE,
		5 * sizeof(GLfloat),
		(GLvoid *) (2 * sizeof(GLfloat))
	);

	glUniform1f(uniform_fade, sinf((float)fadeval/360 * 2 * 3.14)/2 + 0.5);
//	printf("%f\n", sinf((float)fadeval/360 * 2 * 3.13)/2 + 0.5);

	/* Push each element in buffer_vertices to the vertex shader */
	glDrawArrays(GL_TRIANGLES, 0, 3);
	glDisableVertexAttribArray(attribute_coord2d);
	glDisableVertexAttribArray(attribute_colour);

}	

#define BLOCKSIZE 512
#define BLOCK (BLOCKSIZE*sizeof(char))

char * readfile(char * filename)
{
	int size;
	int read=0; //total we've read
	char * retval = malloc(size=(BLOCK));
	char * temp;
	FILE * in = fopen(filename, "r");

	if(!retval)
	{
		printf("Memory allocation failure in readfile\n");
		return NULL;
	}

	if ( ! in )
	{
		printf("%s\n", strerror(errno));
	}

	read+=fread(retval, sizeof(char), BLOCKSIZE, in);

	while(!feof(in))
	{
		temp=realloc(retval, size+=BLOCK);
		if(temp)
		{
			retval=temp;
		}
		else
		{
			printf("Memory allocation failure in readfile\n");
			return NULL;
		}
		read+=fread(retval+read, sizeof(char), BLOCK, in);
	}

	retval=realloc(retval, read+1);
	retval[read]='\0';
	return retval;
}

/**
 *  * Display compilation errors from the OpenGL shader compiler
 *   */
void print_log(GLuint object)
{
	GLint log_length = 0;
	if (glIsShader(object))
		glGetShaderiv(object, GL_INFO_LOG_LENGTH, &log_length);
	else if (glIsProgram(object))
		glGetProgramiv(object, GL_INFO_LOG_LENGTH, &log_length);
	else {
		fprintf(stderr, "printlog: Not a shader or a program\n");
		return;
	}

	char* log = (char*)malloc(log_length);

	if (glIsShader(object))
		glGetShaderInfoLog(object, log_length, NULL, log);
	else if (glIsProgram(object))
		glGetProgramInfoLog(object, log_length, NULL, log);

	fprintf(stderr, "%s", log);
	free(log);
}

GLint makeshader(char * filename, GLenum type)
{
	const GLchar * source = readfile(filename);
	GLuint shader = glCreateShader(type);
	GLint compile_ok = GL_FALSE;

	if(source==NULL)
	{
		printf("Failure reading file\n");
		return 0;
	}

	glShaderSource(shader, 1, &source, NULL);

	free((void *) source);

	glCompileShader(shader);
	glGetShaderiv(shader, GL_COMPILE_STATUS, &compile_ok);
	if(!compile_ok)
	{
		fprintf(stderr, "%s", filename);
		print_log(shader);
		glDeleteShader(shader);
		return 0;
	}
	return shader;
}

vertex.v.glsl

#version 120

attribute vec2 coord2d;
attribute vec3 v_colour;
varying vec3 f_colour;
void main(void) {
	gl_Position = vec4(coord2d, 0.0, 1.0);
	f_colour=v_colour;
}

fragment.f.glsl

#version 120

varying vec3 f_colour;
uniform float alpha;

void main(void) {
	gl_FragColor = vec4(f_colour.x, f_colour.y, f_colour.z, alpha);
}

Enjoy.

This wasn’t a difficult problem. I worked out that it can’t be an 8 or 9 digit number, because all 8 and 9 pandigital numbers are divisible by 3 (testable with the old trick of adding up all the digits and seeing if that number is divisible by 3). Then I just use a simple algorithm to produce all the 7 digit pandigital numbers very quickly which are tested for primality using a dead simple function (divide by every number up to the square root). This quickly spits out all the 7 digit pandigital primes of which I took the last one. This was the solution most people who I@ve seen used.

#include <stdio.h>
#include <stdbool.h>

#define LENGTH 7
bool prime(int n)
{
	int i;

	for(i=2;i*i<=n;i++)
	{
		if(n%i==0)
			return false;
	}
	return true;
}

void recurse(int depth, int sofar, bool used[LENGTH+1])
{
	int i=0;

	if(depth==LENGTH)
	{
		if(prime(sofar))
			printf("%d\n", sofar);
		return;
	}

	while(used[++i]);

	while(i<=LENGTH)
	{
		used[i]=true;
		recurse(depth+1, sofar*10+i, used);
		used[i]=false;
		while(used[++i]);
	}
}	

int main()
{
	bool array[LENGTH+2]={false};
	recurse(0, 0, array);
	return 0;
}