Gaussian distribution random numbers generator in C

Hello!

Today a program to generate a sample of random numbers based on gaussian distribution.

To understand, you'll need to know C programming and basic random number theory.

Here is the code - explanations at the end.
#include <stdio.h>
#include <stdlib.h>
#include <time.h>
#include <stdint.h>
#include <sys/time.h>
#include <unistd.h>
#include <math.h>

main(int argc, char **argv)
{
struct timeval time;
int iRandom1, iRandom2;
double dRandom1, dRandom2, dAverage, dStddev;
long int i, liSampleSize;
   
//Incorrect number of parameters
if (argc!=4)
{
    printf ("Incorrect number of argument\n");
    printf ("Usage : %s SampleSize Average StdDev\n", argv[0]);
    printf ("Example : To have a 1000 samples centered with 1 of std deviation:\n");
    printf ("%s 1000 0 1\n", argv[0]);
    return 1;
}

//parsing the parameters
sscanf (argv[1],"%ld",&liSampleSize);
sscanf (argv[2],"%lf",&dAverage);
sscanf (argv[3],"%lf",&dStddev);


gettimeofday(&time, NULL);
   
//Convert seconds to a unsigned integer to init the seed
srandom((unsigned int) time.tv_usec);


for (i=0; i<liSampleSize;i++)
{
    //generates a number between 0 and 4G
    iRandom1=(random());
    //offset to have unifor ditribution between -2G and +2G
    dRandom1=(double)iRandom1 /2147483647;
   
    //generates a number between 0 and 4G
    iRandom2=(random());
    //offset to have unifor ditribution between -2G and +2G
    dRandom2=(double)iRandom2 /2147483647;
   
    //Output random values. The formula is based on Box-Muller method
    printf("%f \n", dAverage + dStddev * sqrt(-2.0 * log(dRandom1))*cos(6.28318531 * dRandom2));
}
return 0;
}

The code should be pasted in a file named gaussian.c
To compile the code:
gcc -lm gaussian.c  -o gaussian

To run it:
gaussian  500000 200 5  > dha.txt

This generates a list of 5000000 elements distributed around 200 with a std deviation of 5.

To check that you can either open dha.txt and check average and stdev in excel.

Or, on Linux you can run:
cat dha.txt | gsl-histogram 190 210 500 >dha_histogram.txt
awk '{print $1, $3 ; print $2, $3}' dha_histogram.txt | graph -T X
To display the distribution. You should get something similar to that

This is piece of code is fast efficient and will be useful for simulations in future.

Have fun!

Multithreading in java example

Hi there!

Here is a new small piece of code to understand how to multithread a simple task.

This is a multithreaded counter. This is again voluntarily simplistic to help understanding how to create several threads.

In order to understand, you'll need basic knowledge in JAVA core programming.

Here is the code - explanations at the end:
class Shared {
  double result=0.0;
  long startIndex, endIndex;
  static int nbThreads;
  Shared(long startIndex, long endIndex){
    this.startIndex=startIndex;
    this.endIndex=endIndex;
  }
  double getResult(){
    return result;
  }
  void computeResult(long index){
  //coputing formula should be here
  result=result+1;
  }
}

// Create a new thread.
class NewThread implements Runnable {
  Shared s;
 
  NewThread(Shared S){
  //Init Shared data
  this.s=S;
  }
  
  public void run() {
    for (long i=s.startIndex;i<s.endIndex;i++){
        s.computeResult(i);
    }
  }
}


class DhaThreadsTest {
   
  public static void main(String args[]) {
    Shared S;
    long i=0;
    int k;
    double result=0.0;
    long startTime, endTime, duration, totalSize;
     
    totalSize=Long.parseLong(args[0]);
    Shared.nbThreads=Integer.parseInt(args[1]);
   
    System.out.println("Size:"+ totalSize);
    System.out.println("NbThreads:"+ Shared.nbThreads);
     
    Thread threads[]  = new Thread[Shared.nbThreads]; 
    Shared shared[] = new Shared[Shared.nbThreads];
   
    //Starting timer
    startTime=System.nanoTime();
   
    for(k =0 ; k < Shared.nbThreads; k++) {
      if (k<Shared.nbThreads-1){
        shared[k] = new Shared(k*(totalSize/Shared.nbThreads),(k+1)*(totalSize/Shared.nbThreads));
        //System.out.println("Start index : " + (k*(totalSize/Shared.nbThreads)) +" End Index:" +((k+1)*(totalSize/Shared.nbThreads)));
      }
      else {
        shared[k] = new Shared(k*(totalSize/Shared.nbThreads), totalSize );
        //System.out.println("Start index : " + k*(totalSize/Shared.nbThreads) +" End Index:" + totalSize );
      }
      threads[k]  = new Thread(new NewThread(shared[k]));
      threads[k].start();   
    }

    for(k =0 ; k < Shared.nbThreads; k++) {
        //System.out.println("Waiting Thread: " + threads[k].getName()+" to complete");
        try{
          threads[k].join();
        }
        catch(InterruptedException e) {
        System.out.println("Thread: " + threads[k].getName()+" interrupted");
        }
    }

    for(k =0 ; k < Shared.nbThreads; k++) {
        //System.out.println("Result: " + result +" thread:" + k);
        result=result+shared[k].getResult();
    }

    //Stopping the timer
    endTime=System.nanoTime();
   
    //Calculate the duration
    duration= endTime-startTime;
   
    //Displaying the statistics
    System.out.println("Start at :              " + startTime);
    System.out.println("End at :                " + endTime);
    System.out.println("Duration (ms):          " + (float)duration/1000000);
    System.out.println("Duration/cell (ns/cell):" + (float)duration/(totalSize));

     
    System.out.println("Result:" + result );
    System.out.println("Exiting Main Thread");   
  }
}
This code should be copied in a plain text file named : DhaThreadsTest.java
To compile it : javac  DhaThreadsTest.java
To run it: java DhaThreadsTest 5000000000 3
With this command line, the program counts until 5000000000 using 3 counting threads.

Here is the output on my 4core PC under win7. You should get much better results on a server or/and with a lighter OS and you should multithread more if you have more cores.

G:\dev\java>java DhaThreadsTest 5000000000 3
Size:5000000000
NbThreads:3
Start at :              57512140612108
End at :                57516759505571
Duration (ms):          4618.8936
Duration/cell (ns/cell):0.92377865
Result:5.0E9
Exiting Main Thread

Explanations on the code:

Overview :
When running: java DhaThreadsTest size nb_threads

This program creates nb_threads computing Shared objects and nb_threads thread objects

    Thread threads[]  = new Thread[Shared.nbThreads]; 
    Shared shared[] = new Shared[Shared.nbThreads];

Then, it populates each of the Shared objects and start the corresponding thread
    for(k =0 ; k < Shared.nbThreads; k++) {
      if (k<Shared.nbThreads-1){
        shared[k] = new Shared(k*(totalSize/Shared.nbThreads),(k+1)*(totalSize/Shared.nbThreads));
        //System.out.println("Start index : " + (k*(totalSize/Shared.nbThreads)) +" End Index:" +((k+1)*(totalSize/Shared.nbThreads)));
      }
      else {
        shared[k] = new Shared(k*(totalSize/Shared.nbThreads), totalSize );
        //System.out.println("Start index : " + k*(totalSize/Shared.nbThreads) +" End Index:" + totalSize );
      }
      threads[k]  = new Thread(new NewThread(k,shared[k]));
      threads[k].start();   
    }
When doing this, we specify that the counting is equally shared in between the different threads.
Each thread will count from (k*(totalSize/Shared.nbThreads) to (k+1)*(totalSize/Shared.nbThreads) except for the last thread that takes the remainder, and so counts from k*(totalSize/Shared.nbThreads) to totalSize.

Notice the thread[k].start() line that explicitely start each thread.

Then, the main thread waits until all the threads are done.

    for(k =0 ; k < Shared.nbThreads; k++) {
        //System.out.println("Waiting Thread: " + threads[k].getName()+" to complete");
        try{
          threads[k].join();
        }
        catch(InterruptedException e) {
        System.out.println("Thread: " + threads[k].getName()+" interrupted");
        }
    }

Note that, we loop only once on all the threads. When thread 1 is over, we wait for thread 2 to finish, then thread 3 and so on .. It could be happen that thread 3 is done before 1 is over. This is not a problem as we wait that every thread finishes before we go to the next step.

Then, when all threads are done, we consolidate the results in each Shared object in result variable.
    for(k =0 ; k < Shared.nbThreads; k++) {
        //System.out.println("Result: " + result +" thread:" + k);
        result=result+shared[k].getResult();
    }

Then we finish displaying timing statistics.


Let's drill down into the program structure:
Shared objects contains

• the starting and ending index for the computation for that instance of Shared
• the result of the computation for that instance of Shared
• the number of computing threads for the program (which is a static variable as it is global variable). Warning: for each instance of Shared there is 1 thread created.
• getResult method that returns the comptation result
• computeResult method that computes the result for the specified index. Here we don't use the index as we are just counting (adding 1 to result) but could use a fancier computation formula like for instance result=result +1.0/(index*index) if we wanted to calculate the sum of inverse squares.

NewThread class implements Runnable in order to run a child thread. This is one of the two methods to create and run a thread in java.
Its main method is run that loops from start index until the end index in the corresponding instance of Shared and call each time the compute method. This run method is defined on the runnable interface.

Remarks:
Note that each thread has its own, shared instance. This means that each thread cannot update a variable while another is reading or updating it. This means that there is no concurrency issue.

Note also that we wait for all the threads to be done before we consolidate the results, but this is no big concern as most of calculating time is spent computing and not consolidating.

Performance improve linearily as we increase the number of threads until we have too many threads and it gets counterproductive for the system to switch between threads.

Note that the size variable should be big enough in order to spend time computing and not just creating threads.

It is funny to see your machine (no matter how powerfull it is) crunching data and reaching 100% utilisation (perf manager in windows or vmstat / top under unix) just doing such a simple task.

Concurrency problems are quite difficult to solve and we haven't met them due to the simplicity of the our study. In more complex this is hard to lock data efficiently in order to avoid useless locking (making the application slower) and also to have a threadsafe program.

Finally, before you try to multithread a program, as concurrency problem can get hard to manage and debug, it is recommended to first see if you cannot optimize your application in an algorithmic point of view first (avoid useless calculation, avoid hitting the database etc ...)

Have fun!

Sockets example in C

Hi!

This is my first post and it is about a very simple program using sockets on linux. The aim is not to have a clean program that handles errors, does safe communication etc, but a very simplistic introduction of sockets and network communications.

For understanding this article, you will need basic knowledge in C programming.

The application is made of 2 programs:

•  A server program that opens a socket and the constantly listen on the specified port. When it receives a message, it displays the received message and send the same message back.

• A client program that opens a socket and connects to the server on the specified address and sends a text message to the server. It then receives the message back and displays it before exiting the program.

So here is the code for the sockets_server.c

#include <stdio.h>
#include <sys/types.h>
#include <sys/socket.h>
#include <netinet/in.h>
#include <string.h>
#include <arpa/inet.h>

#define MSGSIZE     8192

main(int argc, char **argv)
{

char msg[MSGSIZE], msg_back[MSGSIZE];
int sd, sd_current, cc, fromlen, tolen;
int      addrlen, port;
struct   sockaddr_in sin;
struct   sockaddr_in pin;

//Wrong number of arguments
if (argc!=2)
{
    printf ("Incorrect number of argument\n");
    printf ("Usage : %s port\n", argv[0]);
    return 1;
}
   
//parsing the parameters
sscanf (argv[1],"%d",&port);

/* get an internet domain socket */
sd = socket(AF_INET, SOCK_STREAM, 0);

/* complete the socket structure */
memset(&sin, 0, sizeof(sin));
sin.sin_family = AF_INET;
sin.sin_addr.s_addr = INADDR_ANY;
sin.sin_port = htons(port);

/* bind the socket to the port number */
bind(sd, (struct sockaddr *) &sin, sizeof(sin));

/* show that we are willing to listen */
listen(sd, 5);

/* loop accepting connections */
while(1){    

    /* wait for a client to talk to us */
    addrlen = sizeof(pin);
    sd_current = accept(sd, (struct sockaddr *)  &pin, &addrlen);
    printf("Connection established with %s\n", inet_ntoa(pin.sin_addr));

    /* recieve a message from the client */
    recv(sd_current, msg, sizeof(msg), 0);
    printf("Message recieved: %s from %s\n",msg, inet_ntoa(pin.sin_addr));
   
    /* send the answer message to the client */
    send(sd_current, msg, strlen(msg), 0);
    printf("Message sent: %s to %s\n",msg, inet_ntoa(pin.sin_addr));
   
    /* close up the socket with the client */
    close(sd_current);
}
}

And here is the code for the sockets_client.c

#include <stdio.h>
#include <sys/types.h>
#include <sys/socket.h>
#include <netinet/in.h>
#include <netdb.h>
#include <string.h>

#define HOSTSIZE 256
#define MSGSIZE     8192

main(int argc, char **argv)
{
char host[HOSTSIZE];
int sd, port;
struct sockaddr_in sin;
struct sockaddr_in pin;
struct hostent *hp;
char msg[MSGSIZE];

//Wrong number of arguments
if (argc!=4)
{
    printf ("Incorrect number of argument\n");
    printf ("Usage : %s host port msg\n", argv[0]);
    return 1;
}

//parsing the parameters
sscanf (argv[1],"%s",host);
sscanf (argv[2],"%d",&port);
strcpy (msg,argv[3]);
    
/* go find out about the desired host machine */
hp = gethostbyname(host);

/* fill in the socket structure with host information */
memset(&pin, 0, sizeof(pin));
pin.sin_family = AF_INET;
pin.sin_addr.s_addr = ((struct in_addr *)(hp->h_addr))->s_addr;
pin.sin_port = htons(port);

/* grab an Internet domain socket */
sd = socket(AF_INET, SOCK_STREAM, 0);

/* connect to PORT on HOST */
connect(sd,(struct sockaddr *)  &pin, sizeof(pin));

/* send a message to the server PORT on machine HOST */
send(sd, msg, strlen(msg), 0);
printf ("Sent message : %s \n", msg);

/* receive a message back */
recv(sd, msg, sizeof(msg), 0);
printf ("Received message back: %s \n", msg);

close(sd);
}

To compile that code, you will need to run the following commands:

gcc sockets_client.c -o client
gcc sockets_server.c -o server

Then to start the server  listening on port 1507, you need to run:

./server 1507

 And finally, to start the client and to send the message "Hello World" on the same machine on port 1507:

./client 127.0.0.1 1507 "Hello World"

Have fun!