X-Git-Url: https://hackdaworld.org/gitweb/?p=physik%2Fising.git;a=blobdiff_plain;f=ising.c;h=d3690247f9a94136609a8495d27ad637656ba65a;hp=e43b4be0a7df9ec7134f9cd5d668f651904e2dc7;hb=HEAD;hpb=78fcc4451906fd4fbc2830fc72fd511563d9bf98

diff --git a/ising.c b/ising.c
index e43b4be..d369024 100644
--- a/ising.c
+++ b/ising.c
@@ -7,32 +7,63 @@
  *
  */
 
+#define _GNU_SOURCE
 #include <stdio.h>
 #include <stdlib.h>
 #include <string.h>
 #include <math.h>
 #include <time.h>
+#include <sys/types.h>
+#include <sys/stat.h>
+#include <fcntl.h>
+#include <unistd.h>
+
 #include "dfbapi.h"
 
 #define X 200
 #define Y 200
-#define MAX_T 100
+#define I 20
+#define S 1
+
+int usage(void)
+{
+ puts("usage:");
+ puts("");
+ puts("-h \t show this help");
+ puts("-o <file> \t save T - M values to file");
+ puts("-i <value> \t specify itteration value");
+ puts("-x <value> \t # x lattice sites");
+ puts("-y <value> \t # y lattice sites");
+ puts("-s <value> \t spin interaction strength");
+ puts("-r \t run in interactive mode (still in work)");
+ puts("-d <value> \t refresh display rate");
+ puts("-t <value> maximal temperature");
+
+ return 1;
+}
 
 int main(int argc, char **argv)
 {
  unsigned char *atom;
- char *arg_v[4];
- char m_text[30];
- char t_text[30];
- char b_text[30];
+ char *arg_v[7];
+ char m_text[64];
+ char t_text[64];
+ char b_text[64];
+ char s_text[64];
+ char output_file[64];
+ int of_fd=0;
  int max_x,x_c,max_y,y_c;
  int i;
+ int itt;
  int count_p;
- unsigned int T;
- double S;
- int M;
+ double T;
+ double s;
+ int M=0;
  double beta;
  double delta_e;
+ int runfree=0;
+ double max_T;
+ int dr=1;
 
  /* random stuff*/
  srand(time(0));
@@ -40,59 +71,166 @@ int main(int argc, char **argv)
  /* display stuff */
  d2_lattice d2_l;
 
- /* we will parse argv later ... */
  max_x=X;
  max_y=Y;
-
+ itt=I;
+ s=S;
+ strcpy(output_file,"");
+ max_T=0;
+ /* parse argv */
+ for(i=1;i<argc;i++)
+ {
+  if(argv[i][0]=='-')
+  {
+   switch(argv[i][1])
+   {
+    case 'x':
+     max_x=atoi(argv[++i]);
+     break;
+    case 'y':
+     max_y=atoi(argv[++i]);
+     break;
+    case 'i':
+     itt=atoi(argv[++i]);
+     break;
+    case 's':
+     s=atof(argv[++i]);
+     break;
+    case 'o':
+     strcpy(output_file,argv[++i]);
+     break;
+    case 'r':
+     runfree=1;
+     break;
+    case 'd':
+     dr=atoi(argv[++i]);
+     break;
+    case 't':
+     max_T=atof(argv[++i]);
+     break;
+    default:
+     usage();
+     return -1;
+   }
+  } else usage();
+ }
+ 
+ /* prepare lattice */  	
  d2_lattice_init(&argc,argv,&d2_l,max_x,max_y);
-
  atom=(unsigned char *)(malloc(max_x*max_y*sizeof(unsigned char)));
- 
  d2_l.status=atom;
- 
+
+ if(strcmp(output_file,""))
+ {
+  of_fd=open(output_file,O_WRONLY|O_CREAT);
+  if(of_fd<1)
+  {
+   puts("can't open output file ...");
+   return -1;
+  }
+ }
+
  /* begin at T=0 M=1 situation */
  memset(atom,0,max_x*max_y*sizeof(unsigned char));
 
- S=1; /* i have no idea! */
+ if(max_T==0) max_T=3.0*s; 
+
+ for(T=.05;T<max_T;T+=.05)
+ {
+  beta=1.0/T; /* k_B = 1 */
+  for(i=0;i<itt;i++)
+  {
+
+ M=0;
+ for(x_c=0;x_c<max_x;x_c++)
+ {
+  for(y_c=0;y_c<max_y;y_c++)
+  {
+   count_p=0;
+   if((*(atom+x_c+((y_c+max_y+1)%max_y)*max_x))&1) ++count_p;
+   if((*(atom+x_c+((y_c+max_y-1)%max_y)*max_x))&1) ++count_p;
+   if((*(atom+((max_x+x_c+1)%max_x)+y_c*max_x))&1) ++count_p;
+   if((*(atom+((max_x+x_c-1)%max_x)+y_c*max_x))&1) ++count_p;
+   if(((*(atom+x_c+y_c*max_x))&1)==0) count_p=4-count_p;
+   delta_e=(2*count_p-4)*s;
+   if(delta_e<=0) *(atom+x_c+y_c*max_x)=(*(atom+x_c+y_c*max_x)+1)&1;
+   else
+   {
+    if(1.0*rand()/RAND_MAX<exp(-1.0*delta_e*beta))
+     *(atom+x_c+y_c*max_x)=(*(atom+x_c+y_c*max_x)+1)&1;
+   }
+   if((*(atom+x_c+y_c*max_x))&1) ++M;
+  }
+ }
  
- for(T=1;T<MAX_T;T++)
+ if(i%dr==0)
  {
-  beta=20.0/T; /* k_B = 1 */
-  /* do 10 itterations, we will need more */
-  for(i=0;i<10;i++)
+  sprintf(t_text," temp: %.3f",T);
+  arg_v[1]=t_text;
+  sprintf(b_text," beta: %.3f",beta);
+  arg_v[2]=b_text;
+  arg_v[3]=NULL;
+  sprintf(s_text," interaction strength: %.3f",s);
+  arg_v[4]=s_text;
+  arg_v[5]=NULL;
+  sprintf(m_text," magnetization: %.3f",1.0-2.0*M/(max_x*max_y));
+  arg_v[6]=m_text;
+  d2_lattice_draw(&d2_l,0,0,6,arg_v);
+ }
+
+  }
+  if(of_fd) dprintf(of_fd,"%f %f\n",T,1.0-2.0*M/(max_x*max_y));
+ }
+
+ for(T=max_T;T>0;T-=.05)
+ {
+  beta=1.0/T; /* k_B = 1 */
+  for(i=0;i<itt;i++)
   {
 
+ M=0;
  for(x_c=0;x_c<max_x;x_c++)
  {
   for(y_c=0;y_c<max_y;y_c++)
   {
    count_p=0;
-   M=0;
-   if((*(atom+x_c+((2*y_c+1)%max_y)*max_x))&1) ++count_p;
-   if((*(atom+x_c+((2*y_c-1)%max_y)*max_x))&1) ++count_p;
-   if((*(atom+((2*x_c+1)%max_x)+y_c*max_x))&1) ++count_p;
-   if((*(atom+((2*x_c-1)%max_x)+y_c*max_x))&1) ++count_p;
+   if((*(atom+x_c+((y_c+max_y+1)%max_y)*max_x))&1) ++count_p;
+   if((*(atom+x_c+((y_c+max_y-1)%max_y)*max_x))&1) ++count_p;
+   if((*(atom+((max_x+x_c+1)%max_x)+y_c*max_x))&1) ++count_p;
+   if((*(atom+((max_x+x_c-1)%max_x)+y_c*max_x))&1) ++count_p;
    if(((*(atom+x_c+y_c*max_x))&1)==0) count_p=4-count_p;
-   delta_e=(2*count_p-4)*S;
-   if(delta_e<0) *(atom+x_c+y_c*max_x)=(*(atom+x_c+y_c*max_x)+1)&1;
+   delta_e=(2*count_p-4)*s;
+   if(delta_e<=0) *(atom+x_c+y_c*max_x)=(*(atom+x_c+y_c*max_x)+1)&1;
    else
    {
     if(1.0*rand()/RAND_MAX<exp(-1.0*delta_e*beta))
      *(atom+x_c+y_c*max_x)=(*(atom+x_c+y_c*max_x)+1)&1;
    }
-   if((*(atom+x_c+((2*y_c+1)%max_y)*max_x))&1) ++M;
+   if((*(atom+x_c+y_c*max_x))&1) ++M;
   }
  }
-  sprintf(t_text,"T = %d",T);
+
+ if(i%dr==0)
+ {
+  sprintf(t_text," temp = %.3f",T);
   arg_v[1]=t_text;
-  sprintf(b_text,"b = %f",beta);
+  sprintf(b_text," beta = %.3f",beta);
   arg_v[2]=b_text;
-  sprintf(m_text,"M = %f",1.0-2.0*M/(max_x*max_y));
-  arg_v[3]=m_text;
-  d2_lattice_draw(&d2_l,0,0,3,arg_v);
+  arg_v[3]=NULL;
+  sprintf(s_text," interaction strength: %.3f",s);
+  arg_v[4]=s_text;
+  arg_v[5]=NULL;
+  sprintf(m_text," magnetization: %.3f",1.0-2.0*M/(max_x*max_y));
+  arg_v[6]=m_text;
+  d2_lattice_draw(&d2_l,0,0,6,arg_v);
+ }
+  
   }
+  if(of_fd) dprintf(of_fd,"%f %f\n",T,1.0-2.0*M/(max_x*max_y));
  }
 
+ if(of_fd) close(of_fd);
+
  getchar();
 
  return 1;