Sunday, 31 March 2013

Building a Linux OS with Buildroot and Running it using Qemu

So it's been a very long time since I blogged, and I thought I should write another post. Since my last post I've finished my PhD and I'm now working as a software engineer working on embedded linux systems. I've learnt a lot during my time working as a software engineer, and I've found some of the things I've been googling are quite obscure, and some of the tutorials are a little out of date.

One task that I was given was to create an operating system using the buildroot tool, and run it virtually using qemu. If you've found your way here from Google you probably know about these tools, briefly, buildroot is a set of Makefiles which allows you to compile Linux for specific target architectures and add in a number of packages depending on your needs. Qemu is an open source virtualisation tool similar to virtualbox.

By the way, I'm assuming you're using a linux operating system to do this. I used Ubuntu 10.04 to do this.

Setting up Buildroot

  • Download buildroot from here and type tar xf buildroot-2013.02.tar.gz
  • Change into the buildroot directory and type make menuconfig. You will see a blue screen with lots of configuration options.
  • Set target architecture to i386, and target architecture variant to i686.
  • I used linux 2.6.37.2 for my version, I got this to work but I see no reason why leaving the default configuration wouldn't work, however to follow my exact steps go into toolchain and set kernel headers to manually specified, then enter 2.6.37.2 in the linux version section.
  • Next go into system configuration and set Port to run a getty (login prompt) on to /dev/ttyS0, and then set the baud rate to 38400.
  • Go into Filesystem images and select an ext2 filesystem. I also used an iso image, but I think you need to build a kernel first before you can select an iso image.
  • Go into the Kernel menu, and select Linux kernel, select "custom tarball" for kernel version, then set the kernel url to http://ftp.lug.ro/kernel/linux/kernel/v2.6/linux-2.6.37.2.tar.gz, also select "Install kernel to /boot in target".
  • Still inside the Kernel menu set Defconfig name to "i386"
  • Now exit and save your configuration, then type make. This will build a filesystem and kernel and will take a while (~30 minutes).
  • When the compiler has finished type make linux-menuconfig and go to File systems, then select "Second extended fs support" (press space twice).
  • Exit linux-menuconfig and save the changes.
  • Now type make menuconfig. Go to Kernel configuration, and select "using a custom config file", then set the "Configuration file path" to "output/build/linux-2.6.37.2/.config".
  • Now exit and save the changes.
  • Type make and your new file system will build.
Now you have a linux operating system! You probably want to be able to connect to it when it's running.

Setting up Netwoking


  • To set up networking on the new system first add a couple more packages in buildroot.
  • Go to Package selection, and select show packages provided by busybox.
  • Go through and select dhcp and openssh, also select thttpd if you want.
  • Change to the directory with rootfs.ext2 and bzImage in it.
  • Type mkdir mnt.
  • Now run sudo mount -o loop rootfs.ext2 mnt/
  • Then cd mnt/etc/network and add the following lines to interfaces: 
auto eth0
iface eth0 inet dhcp


  • Now type sudo umount mnt
Now you're completely finished with buildroot, now to get it running using qemu.

Setting up Qemu
  • DO NOT install qemu using apt-get on the Ubuntu. I had problems with this.
  • I installed qemu using git, as follows:
git clone git://git.qemu-project.org/qemu.git
git submodule update --init pixman

  • Then the usual configure, make, sudo make install builds qemu (this can take a good few minutes too!).
  • You might need to do sudo apt-get install libglib2.0-dev if the configure stage complains about missing glib2-12.
Running as a VM

  • Now you can run the buildroot image as a VM using qemu.
  • Change to the directory containing bzImage and rootfs.ext2.
  • Run this command: qemu-system-i386 -kernel bzImage -hda rootfs.ext2 -boot c -m 128M -append "root=/dev/sda rw console=ttyS0,38400n8" -nographic -net nic -net user -redir tcp:2222::22 -redir tcp:3333::33 -no-reboot -localtime
  • You should see the buildroot login prompt, log in as root.
  • Type adduser <yourname> and enter password details.
  • Open up a terminal and type ssh -p 2222 <yourname><at>localhost and you can connect to the VM.
  • You can also set up thttpd to host content on port 33 on the VM and access content in your browser by navigating to localhost:3333.
That's it, you're done. You now have a working linux operating system that is running virtually on your machine. You can SSH in to it, and it also has internet access, (ping will probably fail because apparently qemu doesn't like ICMP, but you should be able to wget to test your connection). I hope this helps out whoever is reading, it took me a long time to figure out how to do all this.

If you get into any problems going through this then please leave a comment and I'll try to see if I can help.

Wednesday, 28 November 2012

Project Euler Q3 & Q4

OK, so moving on with the project Euler questions, we're up to question 3, which is where things start to get challenging.

Question 3 asks:


The prime factors of 13195 are 5, 7, 13 and 29.

What is the largest prime factor of the number 600851475143 ?

So, to solve this problem we're going to need a good way of generating primes. One fairly efficient way of doing this is to use an algorithm called the Sieve of Eratosthenes. This algorithm works by removing multiples of primes that have already been identified, this reduces the number of candidates that are possible prime numbers.

To illustrate this, consider the numbers 2 - 10:
2  3  4  5  6  7  8  9  10

obviously 2 is prime, so no other multiple of 2 can be prime. That leaves the list looking as follows:
2  3  4  5  6  7  8  9  10 

now 3 is tested against 5, 7 and 9, removing 9. Then all the primes below 10 have been found.

The following code is a Python prime number generator that uses the Sieve of Eratosthenes to find prime numbers:



def gen_primes():
    D = {}  

    # The running integer that's checked for primeness
    q = 2  

    while True:
        if q not in D:
            yield q        
            D[q * q] = [q]
        else:
            for p in D[q]:
                D.setdefault(p + q, []).append(p)
            del D[q]

        q += 1

Using this code we can now test which primes divide the large number provided by PE.

number=600851475143
divisors=[]

for prime in gen_primes():
        
        if(number%prime == 0):
                divisors.append(prime)
                number = number/prime           

        if (prime >= number):
                break

print max(divisors)

The code adds all the primes that divide the large number exactly into a list. At the end of the loop the largest divisor is printed to the screen. For this question the answer is 6857.


Next up is question 4, which asks:

A palindromic number reads the same both ways. The largest palindrome made from the product of two 2-digit numbers is 9009 = 91 99.

Find the largest palindrome made from the product of two 3-digit numbers.

My approach to this question was somewhat brute-force. I defined a function to test if a number was a palindrome by converting it to a string and using some of Python's built-in string manipulation tools. The source code for this function is given below:

  def is_palindrome(number):
        forward=str(number)
        backward=forward[::-1]

        if forward == backward:
                return True
        else:
                return False
so the code is pretty simple, convert the number to a string (called forward). reverse the string with the command backward=forward[::-1], and compare the resulting two strings. If forward is the same as backward then obviously we have a palindrome.

Then I just tested every combination of three digit numbers, storing the results that were palindromes. The following code outlines my method:

pals=[]
for num1 in xrange(100,1000):
        for num2 in xrange(100, 1000):
                prod = num1*num2
                if is_palindrome(prod):
                        pals.append(prod)


print max(pals)

I need to have nested loops to make sure I test all the combinations, even if I will test some combinations twice. This brute-force algorithm still runs pretty quickly on my laptop, producing a result in a couple of seconds. The answer, incidentally, is 906609.

Thanks for reading. If you've got any comments/suggestions, or you want to share your solution, please leave a comment at the bottom of the page.
 


Tuesday, 27 November 2012

Project Euler Q1 and Q2

So it's been a while again, last time I blogged about bisection method, which is a simple way to solve equations numerically.

This week I thought I'd share my solutions to the first couple of project Euler questions. These aren't particularly challenging, in particular I think the first question is just so you get used to the format.

I decided to tackle the problems on there using Python, as I wanted to learn the language. I found it quite a useful way to learn the language, in particular I was a little confused about generator functions when I was doing the tutorial, however I used question 2 as an opportunity to get my head around them.

So, to begin, the first question on Project Euler is:
If we list all the natural numbers below 10 that are multiples of 3 or 5, we get 3, 5, 6 and 9. The sum of these multiples is 23.

Find the sum of all the multiples of 3 or 5 below 1000. 

I solved this question with a simple loop and an if statement to test if the numbers were multiples of 3 or 5, the code is below:


tot=0
for number in xrange(1000):
        if (number%3 == 0) or (number%5==0):
                tot = tot + number

print tot

I think the code is pretty straightforward, the answer is 233168.

Now, moving on to the slightly more challenging question 2:

Each new term in the Fibonacci sequence is generated by adding the previous two terms. By starting with 1 and 2, the first 10 terms will be:

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

By considering the terms in the Fibonacci sequence whose values do not exceed four million, find the sum of the even-valued terms.

I decided that the Fibonacci sequence would be an ideal generator function. Python allows you to code functions that can be iterated over, returning one item at a time using the yield command. I think the best way to demonstrate this is by showing you the code!


def fib():
        a,b = 1, 1
        while 1:
                a,b = b,a+b
                yield a

tot=0
for num in fib():
        print num
        if num%2==0 and (num < 4000000):                
                tot+=num
        if (num > 4000000):
                break

print tot

The code works by iterating over the fib() function which returns the next fibonacci number in the sequence, if the number is even and less than 4 million it is added to the running total. Once the fibonacci numbers exceed 4 million the loop is broken, and the total is printed.

The answer is 4613732.

I'm going to blog about my solutions to the first 10 questions. It would be good to get some kind of discussion going about the different methods people have used to solve these problems, particularly as we get to the trickier questions.

Anyway, thanks for reading.


Tuesday, 20 November 2012

Bisection Method

Football Fan? Why not check out www.playerrank.co.uk to vote for your favourite ever player?

I thought I'd start writing my blog again, it's been a while but I've had a few hits, particularly for my Moore-Penrose  post, and my Richardson extrapolation posts. I decided I'd post about bisection method. This is a commonly taught numerical method for solving equations.

The method is robust, as long as you are able to bracket the region in which your solution exists it should be possible to find the solution you are looking for. This is where bisection method has an advantage over say, Newton's method, which will fail if it encounters stationary points.

To solve an equation using bisection method you have to be able to specify the upper and lower boundary for the root that you are calculating. It is important that only one root lies in this region, otherwise you may get errors. In the example code below the equation we solve has roots x=-2, 2, therefore the upper and lower bounds should only include one of the two.

Bisection method iteratively tests whether there is a sign change in one half of the bracketing interval. If there is the interval is reduced in size by a half, and the procedure is carried out again, reducing the bracketing region until the solution has been found to within an appropriate tolerance.

The diagram above illustrates the technique, the points L and U represent the lower and upper limits of the initial bracketing region. Point M is the midpoint of the interval, defined as (U+L)/2. In the code presented below f(U) and f(M) are evaluated, if there is a root between M and U then f(U)*f(M) will always be negative. In this case there is a root, so the algorithm sets L=M, halving the size of the interval, and repeats the procedure.

When the difference between U and L is sufficiently small the root is returned. There is also a test that limits the number of iterations, in case a problem is encountered.

The code is presented below:


#include <iostream>

double f(double x); //function to be solved

using namespace std;


int main(){

double LowerLimit;
cout << "Enter the lower limit:\n";
cin >> LowerLimit;
double UpperLimit;
cout << "Enter the upper limit:\n" ;
cin >> UpperLimit;

double tol;   //specify the required accuracy
cout << "Required Tolerance:\n";
cin >> tol;

double MidPoint;
int iter = 0;
int ITERMAX = 50; //limits the number of iterations

while( ((UpperLimit-LowerLimit)/2.0 > tol) && (iter < ITERMAX)   ){

MidPoint =(UpperLimit + LowerLimit)/2.0;

if (f(UpperLimit) * f(MidPoint) > 0)
{
UpperLimit = MidPoint;
}
else
{
LowerLimit = MidPoint;
}

iter++;
}

if(iter < ITERMAX){
cout << "x = " << MidPoint << endl;
cout << iter << " Iterations required\n";
}
else{
cout << "Too Many Iterations\n";
}

return 0;
}


double f(double x){
return x*x -4.0;
}

Here is an example output  when I tested the code earlier:

Even with a huge bracketing region the method converged to the required root in 40 iterations. Anyway, that's all from me today, please leave a comment below if you have any queries, or even a blog request.

Friday, 9 March 2012

Moore-Penrose Pseudo-inverse Best Fit

Football Fan? Why not check out www.playerrank.co.uk to vote for your favourite ever player?

Before we get started I just want to say that I've signed this blog up to technorati, in order to get some more hits. I need to "claim" this blog by posting the following code: HMR6GJ2P49B6

Now that's done let's get to the meat of today's blog...

I was in the pub on Monday with a few of my colleagues, where we were discussing our day's work. Two colleagues had been doing least-squares best fitting and one was complaining that the code she used wasn't working.

She was trying to fit an exponential curve to a data set using some kind of Chi-Squared minimizing function in IDL and it was all going horribly wrong. At this point I asked why she wasn't using the Moore-Penrose Pseudoinverse matrix to do a linear fit on the logged data values. To my surprise none of my colleagues had heard of this technique, so I thought it would make a good blog.

Suppose you had the following system of equations:

There is no one point where the lines all intercept, since the three equations are linearly independent, so a least squares solution is sought.

Decompose this system into matrix form, Ax=b. You now have a system that looks like this:

A = [2  -1]  x = [x] b = [ 2]
      [1   1]        [y]       [ 5]
      [6  -1]                   [-5]

To find the least squares value of x and y all you need to do is the following:


in this case the result is  x = -0.094595 and y = 2.445946. 

In IDL the code to do this is really easy:

x = (invert(transpose(A)##A)##transpose(A))##B


Anyway, that's all for today's blog, if you've got any questions please feel free to comment below, it would be good to hear from you.


  

Tuesday, 6 March 2012

BattleShips

So it's been a while since I last posted anything, sorry about that, I've been really busy finishing off my paper. In case you're interested in it there's a link here.

Anyway, back to coding stuff. It doesn't count as a Maths/Physics problem strictly, but I recently wrote a model answer for a scientific computing module. The question was to write a two player game where each player positions a "submarine" in a 5x5x5 grid and then take it in turns to shoot at the other submarine, until one of the players hits the other's submarine. Also if too many turns are taken with nobody hitting the other submarine the code should quit.

It's a pretty simple question, and you don't need to use arrays to solve this, but the students are still learning about arrays, so I thought I'd show them how to solve this using a 4d array. The code is once again C++, but I will post code in a different language some time.

#include<iostream>


using namespace std;




void space();


int main(){


//I used a 4d array so that both players could have Submarines in the same location
int volume[2][5][5][5] = {0};


//arrays to hold the location of the players' Submarines
int x[2], y[2], z[2];
int player; //For looping input
int x_guess, y_guess, z_guess; //coordinates of the latest guess, these variables will hold both players' guesses


//First I need to get the locations of the Subs off both players


for(player=0; player<2; player++){


//do..while loop to ensure that the inputs are valid
do{
cout << "Player " << player+1 << " input x, y, z coords of your submarine [0,4]\n";
cin >> x[player] >> y[player] >> z[player] ;
}while(  ( (x[player] > 4 ) || (x[player] < 0) ) || ( (y[player] > 4) || (y[player] < 0) ) || ( (z[player] > 4) || (z[player] < 0) ) );  


volume[player][x[player]][y[player]][z[player]] = 1; //Set the location in the array = 1
space(); // 100 lines blank space
}




int turn_count=0; //This is the loop variable that will keep track of which players' turn it is


/* I use turn_count%2 to get the array index of volume for the correct player.
player 1's sub is in the array volume[0][x][y][z]
player 2's sub is in the array volume[1][x][y][z]

so when turn_count is a multiple of 2, turn_count%2 references the location of player 1's sub
when turn_count is not a multiple of 2, turn_count%2 references the location of player 2's sub


This is because turn_count%2 can only equal 0 or 1 */


do{


//Same logic as at input, ensure that guesses are valid
do{
cout << "Player " << turn_count%2 + 1 << " guess the location of the other Player's submarine:\n";
cin >> x_guess >> y_guess >> z_guess;
}while(  ( (x_guess > 4) || (x_guess < 0) ) || ( (y_guess > 4) || (y_guess < 0) ) || ( (z_guess > 4) || (z_guess < 0) ) );  


//if the guess coords don't coincide with the other player's sub then tell the current player they missed
if(volume[(turn_count+1)%2][x_guess][y_guess][z_guess] == 0){
cout << "Sorry Player " << turn_count%2 +1 << ", you missed!\n";
}


turn_count++; //step turn_count by 1


}while((volume[turn_count%2][x_guess][y_guess][z_guess] !=1) && (turn_count < 10) );  //if the player hits the other sub or the players have taken too long then stop looping


if(volume[turn_count%2][x_guess][y_guess][z_guess] != 1){ //If the player hasn't hit the other's sub then they must have taken too long
//Give them some abuse
cout << endl << endl;
cout << "You're both rubbish, I give up, I'm bored.\n";
cout << "Your Sub locations were:\n";
cout << "Player 1:\t" << x[0] << "\t" << y[0] << "\t" << z[0] << endl;
cout << "Player 2:\t" << x[1] << "\t" << y[1] << "\t" << z[1] << endl;
}
else{ //Otherwise congratulate the winner
cout << endl << endl;
cout << "Congratulations Player " << (turn_count-1)%2 +1 << " You Win!! \n";
cout << "Your Sub locations were:\n";
cout << "Player 1:\t" << x[0] << "\t" << y[0] << "\t" << z[0] << endl;
cout << "Player 2:\t" << x[1] << "\t" << y[1] << "\t" << z[1] << endl;


}




return 0;
}




void space(){
for(int i=0; i < 100; i++){
cout << endl;
}
}   


It's quite a long piece of code, and I can't really think of another instance where this might be useful, except perhaps if you're setting a question for a coding class. Anyway if you got this far down thanks for reading, why not post a comment? 

Friday, 24 February 2012

Circle on a Plane

OK, so I thought I'd keep up my blog, for the 17 people who have read it so far! I told you guys this week I'd generate a random point that sits on a circle on a plane.

You might think this is a weird problem to encounter, and it is, but I came across this after I'd been solving the motion of a particle's gyro-center, and needed to assign a position to the particle at random. The code works by first generating a point on the plane perpendicular to the "gyrocenter", which in this example is positioned at x=10, y=20, z=50, and then scales the vector from the generated point to the gyrocenter down to the required radius of the circle.

So here's the code:
#include<iostream>
#include<fstream>
#include<cmath>
#include<stdlib.h>


using namespace std;


int main(){
/* generates random circle of radius 2 on the plane
   (Rx)x + (Ry)y + (Rz)z = (Rx^2 + Ry^2 + Rz^2)
*/
double Rx, Ry, Rz, rx, ry, rz, tot, rg, size, x0, y0, z0;


ofstream prt("plane.dat");


ofstream crc("circ.dat");


rg = 2.0;


Rx=10.0; Ry =20.0; Rz=50.0;




prt << Rx << "\t" << Ry << "\t" << Rz << endl;


for(int i = 0; i < 1000; i++){
tot = Rx*Rx + Ry*Ry + Rz*Rz;


rx = double(rand()%100);


tot = tot - Rx*rx;


ry = double(rand()%100);

tot = tot - Ry * ry; 


rz = tot/ double(Rz);


prt << rx << "\t" << ry << "\t" << rz << endl;


rx = Rx - rx; ry = Ry - ry; rz = Rz - rz;


size = sqrt(rx*rx + ry*ry + rz*rz );


rx = (rx/size)*rg;
ry = (ry/size)*rg;
rz = (rz/size)*rg;


crc << 1 << " " <<  rx << " " << ry << " " << rz << endl;


}


return 0;


}
And here is a plot of the result:
 
Anyway, if anyone's interested there you go. I can't think what to do next, why not set me a challenge?