unstable but might run ...
[physik/posic.git] / pair_correlation_calc.c
index 2f93ab7..efe81ad 100644 (file)
@@ -39,14 +39,24 @@ int main(int argc,char **argv) {
                return -1;
        }
 
+       memset(&moldyn,0,sizeof(t_moldyn));
+
+       printf("[pair corr calc] reading save file ...\n");
        ret=moldyn_read_save_file(&moldyn,argv[1]);
        if(ret) {
                printf("[pair corr calc] exit!\n");
                return ret;
        }
 
+       moldyn.cutoff*=2;
+       moldyn.cutoff_square*=4;
+
        dr=atof(argv[2]);
-       slots=(int)(moldyn.cutoff/dr);
+       slots=moldyn.cutoff/dr;
+       printf("[pair corr calc]\n");
+       printf("  slots: %d\n",slots);
+       printf("  cutoff: %f\n",moldyn.cutoff);
+       printf("  dr: %f\n",dr);
 
        stat=(double *)malloc(3*slots*sizeof(double));
        if(stat==NULL) {
@@ -54,30 +64,22 @@ int main(int argc,char **argv) {
                return -1;
        }
 
+       /* link cell init */
+       link_cell_init(&moldyn,VERBOSE);
+
        calculate_pair_correlation(&moldyn,dr,stat);
 
-       fd=open("pair_corr_func_ab.txt",
-               O_WRONLY|O_CREAT|O_TRUNC,S_IRUSR|S_IWUSR);
-       dprintf(fd,"# type a - type b bonds\n");
-       for(i=0;i<slots;i++)
-               dprintf(fd,"%f %f\n",i*dr,stat[i]);
-       close(fd);
-               
-       fd=open("pair_corr_func_aa.txt",
+       fd=open("pair_corr_func.txt",
                O_WRONLY|O_CREAT|O_TRUNC,S_IRUSR|S_IWUSR);
-       dprintf(fd,"# type a - type a bonds\n");
+       dprintf(fd,"#r #ab #aa #bb\n");
        for(i=0;i<slots;i++)
-               dprintf(fd,"%f %f\n",i*dr,stat[slots+i]);
+               dprintf(fd,"%f %f %f %f\n",
+                       i*dr,stat[i],stat[slots+i],stat[2*slots+i]);
        close(fd);
-               
-       fd=open("pair_corr_func_bb.txt",
-               O_WRONLY|O_CREAT|O_TRUNC,S_IRUSR|S_IWUSR);
-       dprintf(fd,"# type a - type b bonds\n");
-       for(i=0;i<slots;i++)
-               dprintf(fd,"%f %f\n",i*dr,stat[2*slots+i]);
-       close(fd);
-               
+
        free(stat);
 
+       moldyn_free_save_file(&moldyn);
+
        return 0;
 }