ViewVC Help
View File | Revision Log | Show Annotations | View Changeset | Root Listing
root/group/trunk/dp/dp.c
(Generate patch)

Comparing trunk/dp/dp.c (file contents):
Revision 587 by xsun, Thu Jul 10 15:42:41 2003 UTC vs.
Revision 692 by xsun, Wed Aug 13 16:43:53 2003 UTC

# Line 1 | Line 1
1   #include<stdio.h>
2 + #include<sys/time.h>
3 + #include<unistd.h>
4 + #include<time.h>
5   #include<math.h>
6   #include<stdlib.h>
7 + #include<string.h>
8 + #include "mkl_vsl.h"
9 +
10 + #define SEED    1
11 + #define BRNG    VSL_BRNG_MCG31
12 + #define METHOD  0
13 + #define N       1
14 +
15   #define max_sites 15000
16   #define MYSEED 8973247
17  
# Line 13 | Line 24 | int main()
24   int attemp, nacc;
25  
26  
27 < int main()
27 > int main(int argc, char *argv[])
28   {
29    void getmag();
30    double eneri(double xi, double yi, double zi, double phii, double thetai, int i, int jb);
31    double toterg();
32    void adjust();
33 <  void mcmove();
33 >  void mcmove(VSLStreamStatePtr stream);
34    void store(FILE *outFile);
35    
36    double twopi, nnd, myran, kb;
37    double ent;
38    int icycl, ncycl, imove, nmoves, nsamp;
39    int nx, ny, which, i, j;
40 +  int readfile, mySeed;
41    
42 +  FILE *inFile;
43 +  char *endptr, s[100000];
44 +  double ux, uy, uz;
45 +
46 +  struct        timeval         now_time_val;
47 +  struct        timezone        time_zone;
48 +  struct        tm              *now_tm;
49 +  time_t                        now;
50 +
51 +  VSLStreamStatePtr stream;
52 +    
53    FILE *outFile;
54 <  outFile=fopen("dp_file","w");
54 >  outFile=fopen("dp_file1","a");
55    
33 –  srand48(MYSEED);
56    
57 <  // These two parameters set the size of the system:
57 >  //srand48(MYSEED);
58 >  //
59  
60 +  gettimeofday(&now_time_val, &time_zone);  /* get the time now  */
61 +  now = now_time_val.tv_sec;                /*  convert to epoch time */
62 +
63 +  mySeed = (int) now;
64 +
65 +  vslNewStream(&stream, BRNG, mySeed);
66 +  
67 +  // These two parameters set the size of the system:
68 +  
69    nnd = 7.0;
70    ny = 100;
71 <
71 >  
72    // we want the domains to be roughly square:
73 <
73 >  
74    nx = (int) ((double)ny / sqrt(3.0));
75 <
76 <
75 >  
76 >  
77    // periodic box boundaries:
78  
79    domsizex = sqrt(3.0) * nx * nnd;
# Line 57 | Line 89 | int main()
89    beta = 1.0 / (kb * 300.0);
90    z0 = 0.0;
91    deltaz = 0.1;
92 <  deltaphi = 0.1;
93 <  kz = kb;
94 <  ktheta = kb;
92 >  deltaphi = 0.1;
93 >  kz = kb;
94 >  ktheta = kb;
95    theta0 = M_PI / 2.0;
96 <    
96 >  
97    which = 0;
98    xlat = 0;
99  
100 <  for(i=0; i < nx; i++) {
100 >  if (argc > 1) {
101 >    if( strcmp(argv[1], "-resume") == 0)
102 >      readfile = 1;
103 >    else
104 >      readfile = 0;
105 >  } else
106 >    readfile = 0;
107  
70 –    xlat = xlat + 2;
71 –    ylat = 0;
108      
109 <    for(j=0; j<ny; j++) {    
109 >  if (readfile == 1) {
110 >    
111 >    inFile = fopen("./l_con","r");
112 >    if (inFile == NULL){
113 >      printf("Error opening file\n");
114 >      exit(-1);
115 >    }
116        
117 <      which = which + 1;
118 <      x[which] = i * sqrt(3.0) * nnd;
119 <      y[which] = j * nnd;
120 <      myran = drand48();
121 <      phi[which] = myran * twopi;
122 <      myran = drand48();
123 <      theta[which] = acos(2.0*myran-1.0);
124 <      myran = drand48();
125 <      z[which]=z0 + (2.0*myran-1.0)*deltaz;
126 <      mu[which] = strength;
127 <      
128 <      which = which + 1;
129 <
130 <      printf("%d\n", which);
117 >    while(!feof(inFile))
118 >      {
119 >        fgets(s, 500, inFile);
120 >        nsites = atoi(s);
121 >        fscanf(inFile,"%s",s);
122 >        xlat= atoi(s);
123 >        fscanf(inFile,"%s",s);
124 >        ylat= atoi(s);
125 >        for(i=1; i<=nsites; i++)
126 >          {
127 >            fscanf(inFile,"%s",s);
128 >            fscanf(inFile,"%s",s);
129 >            x[i] = strtod(s, &endptr);
130 >            fscanf(inFile,"%s",s);
131 >            y[i] = strtod(s, &endptr);
132 >            fscanf(inFile,"%s",s);
133 >            z[i] = strtod(s, &endptr);
134 >            fscanf(inFile,"%s",s);
135 >            phi[i] = strtod(s, &endptr);
136 >            fscanf(inFile,"%s",s);
137 >            ux = strtod(s, &endptr);
138 >            fscanf(inFile,"%s",s);
139 >            uy = strtod(s, &endptr);
140 >            fscanf(inFile,"%s",s);
141 >            uz = strtod(s, &endptr);
142 >            theta[i] = acos(uz);
143 >            mu[i] = strength;
144 >          }
145 >      }
146 >    fclose(inFile);
147 >    printf("%f\t%f\t%f\t%f\t%f\n", x[nsites], y[nsites], z[nsites], phi[nsites], theta[nsites]);
148 >  } else {
149  
150 <      x[which] = nnd* sqrt(3.0) * (2.0*i + 1.0) / 2.0;
151 <      y[which] = nnd* ( 2.0*j+1.0 ) / 2.0;
152 <      myran = drand48();
153 <      phi[which] = myran * twopi;
154 <      myran = drand48();
155 <      theta[which] = acos(2.0*myran-1.0);
156 <      myran = drand48();
157 <      z[which] = z0 + (2.0*myran-1.0) * deltaz;
158 <      mu[which] = strength;
159 <      ylat = ylat + 2;
150 >    for(i=0; i < nx; i++) {
151 >      
152 >      xlat = xlat + 2;
153 >      ylat = 0;
154 >      
155 >      for(j=0; j<ny; j++) {    
156 >        
157 >        which = which + 1;
158 >        x[which] = i * sqrt(3.0) * nnd;
159 >        y[which] = j * nnd;
160 >        
161 >        vdRngUniform( METHOD, stream, N, &myran, 0.0, 1.0);
162 >        //myran = drand48();
163 >        phi[which] = myran * twopi;
164 >        //myran = drand48();
165 >        vdRngUniform( METHOD, stream, N, &myran, 0.0, 1.0);
166 >        theta[which] = acos(2.0*myran-1.0);
167 >        //myran = drand48();
168 >        vdRngUniform( METHOD, stream, N, &myran, 0.0, 1.0);
169 >        z[which]=z0 + (2.0*myran-1.0)*deltaz;
170 >        mu[which] = strength;
171 >        
172 >        which = which + 1;
173 >
174 >        x[which] = nnd* sqrt(3.0) * (2.0*i + 1.0) / 2.0;
175 >        y[which] = nnd* ( 2.0*j+1.0 ) / 2.0;
176 >        //myran = drand48();
177 >        vdRngUniform( METHOD, stream, N, &myran, 0.0, 1.0);
178 >        phi[which] = myran * twopi;
179 >        //myran = drand48();
180 >        vdRngUniform( METHOD, stream, N, &myran, 0.0, 1.0);
181 >        theta[which] = acos(2.0*myran-1.0);
182 >        //myran = drand48();
183 >        vdRngUniform( METHOD, stream, N, &myran, 0.0, 1.0);
184 >        z[which] = z0 + (2.0*myran-1.0) * deltaz;
185 >        mu[which] = strength;
186 >        ylat = ylat + 2;
187 >      }
188      }
189 +    
190 +    nsites= which;
191    }
102 –
103 –  nsites= which;
192    
193    en=toterg();
194    printf("\nthe initial energy is : \t%f\n\n",en);
# Line 114 | Line 202 | int main()
202          {
203                  for(imove=1;imove<=nmoves;imove++)
204                  {
205 <                        mcmove();
205 >                        mcmove(stream);
206                  }
207  
208                  if((icycl%nsamp) == 0)
# Line 142 | Line 230 | int main()
230                  printf("\nTotal energy end of simulation : %f\nrunning energy : %f\ndifference : %f\n\n", ent, en, ent-en);
231          }
232          fclose(outFile);
233 +        vslDeleteStream( &stream );
234          return 0;
235   }
236  
# Line 275 | Line 364 | void mcmove()
364          
365   }
366  
367 < void mcmove()
367 > void mcmove( VSLStreamStatePtr stream )
368   {
369          int o;
370          double eno, zn, phin, thetan, enn, myran;
# Line 283 | Line 372 | void mcmove()
372  
373          attemp = attemp + 1;
374          jb = 1;
375 <        myran = drand48();
375 >        //myran = drand48();
376 >        vdRngUniform( METHOD, stream, N, &myran, 0.0, 1.0);
377          o = (int)(nsites*myran) + 1;
378          eno=eneri(x[o], y[o], z[o], phi[o], theta[o], o, jb);
379 <        myran = drand48();
379 >        //myran = drand48();
380 >        vdRngUniform( METHOD, stream, N, &myran, 0.0, 1.0);
381          zn = z[o] + (2.0*myran-1.0) * deltaz;
382 <        myran = drand48();
382 >        //myran = drand48();
383 >        vdRngUniform( METHOD, stream, N, &myran, 0.0, 1.0);
384          phin = phi[o] + (2.0*myran-1.0) * deltaphi;
385 <        myran = drand48();
385 >        //myran = drand48();
386 >        vdRngUniform( METHOD, stream, N, &myran, 0.0, 1.0);
387          thetan = theta[o] + (2.0*myran-1.0)*dtheta;
388          enn=eneri(x[o], y[o], zn, phin, thetan, o, jb);
389 <        myran = drand48();
389 >        //myran = drand48();
390 >        vdRngUniform( METHOD, stream, N, &myran, 0.0, 1.0);
391          if(myran<=(exp(-beta*(enn-eno))))
392          {
393                  nacc = nacc + 1;

Diff Legend

– Removed lines
+ Added lines
< Changed lines
> Changed lines