added still buggy post processing functions

diff --git a/moldyn.c b/moldyn.c
index e554a2f..aebd604 100644 (file)
--- a/moldyn.c
+++ b/moldyn.c
@@ -1805,3 +1805,102 @@ int moldyn_bc_check(t_moldyn *moldyn) {
 
        return 0;
 }
+
+/*
+ * postprocessing functions
+ */
+#define LINE_MAX 128
+int read_line(int fd,char *line) {
+
+       int count,ret;
+
+       count=0;
+
+       while(1) {
+               if(count==LINE_MAX) return count;
+               ret=read(fd,line+count,1);
+               if(ret<0) return ret;
+               if(line[count]=='\n') {
+                       line[count]='\0';
+                       return count+1;
+               }
+               count+=1;
+       }
+}
+
+int calc_fluctuations(double start,double end,char *file) {
+
+       int fd;
+       int count,ret;
+       double time,pot,kin,tot;
+       double p_m,k_m,t_m;
+       double p2_m,k2_m,t2_m;
+       double p_sum,k_sum,t_sum;
+       char buf[LINE_MAX];
+
+       fd=open(file,O_RDONLY);
+       if(fd<0) {
+               perror("[moldyn] post proc open");
+               return fd;
+       }
+
+       /* first calc the averages */
+       p_sum=0.0;
+       k_sum=0.0;
+       t_sum=0.0;
+       count=0;
+       while(1) {
+               ret=read_line(fd,buf);
+               if(ret<0) break;
+printf("%d\n",ret);
+               if(buf[0]=='#') continue;
+               sscanf(buf,"%lf %lf %lf %lf",&time,&kin,&pot,&tot);
+printf("%f %f %f %f\n",time,pot,kin,tot);
+               //if(time>end) break;
+               if((time>=start)&(time<=end)) {
+                       p_sum+=pot;
+                       k_sum+=kin;
+                       t_sum+=tot;
+                       count+=1;
+               }
+       }
+
+       p_m=p_sum/count;
+       k_m=k_sum/count;
+       t_m=t_sum/count;
+
+       /* mean square fluctuations */
+       if(lseek(fd,SEEK_SET,0)<0) {
+               perror("[moldyn] lseek");
+               return -1;
+       }
+       count=0;
+       p_sum=0.0;
+       k_sum=0.0;
+       t_sum=0.0;
+       while(1) {
+               ret=read_line(fd,buf);
+               if(ret<0) break;
+               if(buf[0]=='#') continue;
+               sscanf(buf,"%lf %lf %lf %lf",&time,&kin,&pot,&tot);
+               if(time>end) break;
+               if((time>=start)&(time<=end)) {
+                       p_sum+=((pot-p_m)*(pot-p_m));
+                       k_sum+=((kin-k_m)*(kin-k_m));
+                       t_sum+=((tot-t_m)*(tot-t_m));
+                       count+=1;
+               }
+       }
+
+       p2_m=p_sum/count;
+       k2_m=k_sum/count;
+       t2_m=t_sum/count;
+
+       printf("[moldyn] fluctuations (%f - %f):\n",start,end);
+       printf("  - averages   : %f %f %f\n",k_m,p_m,t_m);
+       printf("  - mean square: %f %f %f\n",k2_m,p2_m,t2_m);
+
+       close(fd);
+
+       return 0;
+}

diff --git a/moldyn.h b/moldyn.h
index 5160668..fc9bce7 100644 (file)
--- a/moldyn.h
+++ b/moldyn.h
@@ -427,4 +427,7 @@ int check_per_bound(t_moldyn *moldyn,t_3dvec *a);
 
 int moldyn_bc_check(t_moldyn *moldyn);
 
+int read_line(int fd,char *line);
+int calc_fluctuations(double start,double end,char *file);
+
 #endif