Keymer Lab in silico contact process

From OpenWetWare

Jump to: navigation, search

Home        Lab Webiste        In Vitro        Research        Journal Club        Hacks       


We start by the creation of the 1d array that will represent the 2d lattice:

  1. // the int * syntax indicates that this function will return a pointer  
  2. int * createlandscape2d(int x_size, int y_size)
  3. {
  4. 	// this creates the pointer to the array landscape, calloc allocates the size of array on the memory for the whole programm until it is freed
  5. 	int *landscape = calloc(y_size*x_size,y_size*x_size);
  6. 	// this sets all the values in the array to zero
  7. 	memset(landscape, 0, sizeof(landscape));
  8. 	printf("Your landscape of x = %d  and y = %d was successfully created.\n It is still empty. \n", x_size, y_size);
  9. 	// this returns the pointer to the array 
  10. 	return landscape;
  11. }


We then initialize this lattice with a single 'cell' in the middle:

  1. int * inilandscapesingle(int x_size, int y_size, int *ptom)
  2. {
  3. 	// this makes landscape the pointer to the array pointed at by ptom
  4. 	int *landscape = ptom;
  5. 	// this line creates a single cell at the middle of the lattice 
  6. 	landscape[((y_size*x_size)/2)+ (x_size /2)]  = 1;
  7. 	printf("Your landscape was successfully initialize with a single cell in the middle of the landscape\n");
  8. 	return landscape;
  9. }


We can also initialise the landscape at random by a certain occupancy rate. This uses a random function based on the Mersenne twister. The implementation used here can be found at this adress

  1. int * inilandscaperandom(int x_size, int y_size,int *ptom, double occupancy, long long seed)
  2. {
  3. 	int *landscape = ptom;
  4. 	// this intializes the Mersenne twister random function 
  5. 	init_genrand64(seed);
  6. 	int j =  0;
  7. 	// this is the loop that goes through the whole lattice 
  8. 	while (j < x_size*y_size)
  9. 	{
  10. 		// this line generates a random real between 0 and 1
  11. 		double r = genrand64_real2();
  12. 		// Here if this condition is verified the position in the array is filled with a 1
  13. 		if (r < occupancy) 
  14. 		{
  15. 			landscape[j] = 1;
  16. 			j++;
  17. 		}
  18. 		else 
  19. 		{
  20. 			j++;
  21. 		} 
  22.  
  23. 	}
  24. 	printf("Your landscape was successfully inintialized randomly with occupancy rate equal to %f \n",occupancy);
  25. 	return landscape;
  26. }

This is the randomwalk through the array to update random sites inside the array depending on their neighbours:

  1. int * randomwalk2(int x_size, int y_size,int *ptom,  double survival, long long seed)
  2. {
  3. 	int *landscape = ptom;
  4. 	init_genrand64(seed);
  5. 	int sites = 0;
  6. 	int size  = y_size * x_size;
  7. 	int randomx = 0;
  8. 	int neighbour = 0;
  9. 	// a loop running for total nuber of sites
  10. 	for (sites  = 0; sites < size; sites ++)
  11. 	{
  12. 		// creating a random number that will be the site of interest
  13. 		randomx = (int) floor((genrand64_real2()*size));
  14. 		// switch is like an if except faster in c so here it returns the value of site randomx
  15. 		switch(landscape[randomx])
  16. 		{
  17. 			case 0:
  18. 				// choosing one of four neighbours of site randomx at random
  19. 				switch((int) ceil(genrand64_real3()*4))
  20. 				{
  21. 					//Picking the neighbour on the right 
  22. 					case 1:
  23. 						neighbour = landscape[(randomx+1)%x_size + (randomx/x_size)*x_size];
  24. 						break;
  25.  
  26. 					//Picking the neighbour on the left
  27. 					case 2:
  28. 						if (randomx-1 < 0)
  29. 						{
  30. 							neighbour = landscape[x_size-1];
  31. 							break;
  32. 						}
  33. 						else
  34. 						{
  35. 							neighbour = landscape[(randomx-1)%x_size + (randomx/x_size)*x_size];
  36. 							break;
  37. 						}
  38.  
  39. 					// Picking the southern neighbour 
  40. 					case 3:
  41. 						neighbour = landscape[(randomx + x_size)%(x_size*y_size)];
  42. 						break;
  43.  
  44. 					//Picking the northen neighbour 
  45. 					case 4:
  46. 						if (randomx-x_size < 0)
  47. 						{
  48. 							neighbour = landscape[randomx+x_size*(y_size-1)];
  49. 							break;
  50. 						}
  51. 						else
  52. 						{
  53. 							neighbour = landscape[randomx-x_size];	
  54. 							break;
  55. 						}
  56.  
  57. 				}
  58. 				if (neighbour ==1)
  59. 				{
  60. 					landscape[randomx] = 1;
  61. 				}
  62. 				break;
  63. 			case 1:
  64. 				// computing if an occupied site should die or not depending on death rate 
  65. 				if (genrand64_real2() < survival)
  66. 				{	
  67. 					landscape[randomx] = 0;
  68. 				}
  69. 				break;
  70. 		} 
  71.  
  72. 	}
  73. 	return landscape;
  74. }

We then run the simulation for a certain time t :

  1. int runsimulationbmp(int x_size, int y_size,int *ptom, int timeofsim, double survival, long long seed)
  2. {
  3. 	init_genrand64(seed);
  4. 	int i = 0;
  5. 	// creating a string to stock the name of the files we wish to create 
  6. 	char name[11] = "";
  7. 	// runs for timesteps = timeofsim
  8. 	while (i < timeofsim)
  9. 	{
  10. 		// updates the name string with iteration number or time t
  11. 		sprintf(name,"images/simu%d.bmp",(i+200));
  12. 		// run the randomwalk algorithm on the array 
  13. 		ptom = randomwalk2( x_size, y_size, ptom, survival, genrand64_int63());
  14. 		// creates a bmp image this function is descibed bellow
  15. 		createfrommatrix( x_size, y_size, 0, 255, 0, ptom, name);
  16. 		i++;
  17. 	}
  18. 	printf("Your simulation ran for time %d is finished and bmp files were created\n",timeofsim );
  19. 	// frees the memory taken up by the array Very Important !!! 
  20. 	void free(void * ptom); 
  21. 	return 0;
  22. }

This is the function that allows the user to transform the lattice in a bitmap image and save it in a file:

  1. int createfrommatrix(int x_size, int y_size, int red, int green, int blue,int *ptom,char name[])
  2. {
  3. 	int *landscape = ptom;
  4. 	// creating a bmp object and it's space in memory
  5. 	BMP* bmp = BMP_Create( x_size, y_size, 24 );
  6. 	int j = 0;
  7. 	// loop that runs throughout the whole array fed into this function
  8. 	while (j < x_size*y_size)
  9. 	{
  10. 		if (landscape[j] == 1)
  11. 		{
  12. 			// allows the user to set rgb values to pixels 
  13. 			BMP_SetPixelRGB(bmp, j%x_size, j/x_size , red , green , blue );
  14. 			j++;
  15. 		}
  16. 		else if(landscape[j] == 0)
  17. 		{
  18. 			BMP_SetPixelRGB(bmp, j%x_size, j/x_size , 0, 40 , 0 );
  19. 			j++;
  20. 		}
  21. 		else
  22. 		{
  23. 			BMP_SetPixelRGB(bmp, j%x_size, j/x_size , 255 , 0 , 0 );
  24. 			j++;
  25. 		}
  26. 	}
  27. 	// creates the file at the indicated name location 
  28. 	BMP_WriteFile( bmp, name );
  29. 	// frees the memory used by the bmp 
  30. 	BMP_Free(bmp);
  31. 	return 0 ;
  32.  
  33. }

The code to run this model can be found here. If you are using a debian based system:

tar zxvf Contact_model.tar.gz

Then:

  1. cd Contact_model/
  2. make
  3. ./test 200 200 500 0.05 0.2

The first argument is the width, the second the height, the third the number of time steps you want to run it for, the fourth the occupancy rate, and the last the death rate.

Personal tools