X-Git-Url: https://hackdaworld.org/gitweb/?p=physik%2Fposic.git;a=blobdiff_plain;f=posic.c;h=6b92a5e51b0428de417e34392189a6aa246abffb;hp=997127e75a0a6490a31396fe3abb1c6f01b6b99c;hb=HEAD;hpb=737bbb619a53036e1b2f4eddeecaecb5ae088b38 diff --git a/posic.c b/posic.c index 997127e..6b92a5e 100644 --- a/posic.c +++ b/posic.c @@ -1,176 +1,41 @@ /* * posic.c - precipitation process of silicon carbide in silicon * - * author: Frank Zirkelbach + * author: Frank Zirkelbach * */ - + +/* main include file */ + #include "posic.h" -#define RAND(max) (max*(0.5-(1.0*rand()/RAND_MAX+1))); +/* functions */ + + + +/* main code */ + +int parse_config_file() { + + return 0; +} int main(int argc,char **argv) { - t_atom *si; - //t_si *c; - int i,j,runs,amount_si; - double time; - int fd1,fd2; - unsigned char xyz[128]; - unsigned char scr[128]; - unsigned char ppm[128]; - - double tau,tau2,m,m2; - double deltax,deltay,deltaz,distance; - double deltax2,deltay2,deltaz2,tmp; - double lj1,lj2,lj; - - /* silicon */ - amount_si=AMOUNT_SI; - printf("simulating %d silicon atoms\n",amount_si); - si=malloc(amount_si*sizeof(t_atom)); - if(!si) { - perror("silicon malloc"); - return -1; - } - memset(si,0,amount_si*sizeof(t_atom)); - m=SI_M; m2=2.0*m; - - /* init */ - printf("placing silicon atoms\n"); - for(i=0;iLX) si[i].x-=LEN_X; - else if(si[i].x<-LX) si[i].x+=LEN_X; - si[i].y+=(tau2*si[i].fy/m2); - if(si[i].y>LY) si[i].y-=LEN_Y; - else if(si[i].y<-LY) si[i].y+=LEN_Y; - si[i].z+=(tau2*si[i].fz/m2); - if(si[i].z>LZ) si[i].z-=LEN_Z; - else if(si[i].z<-LZ) si[i].z+=LEN_Z; - /* calculation of velocities v(t+h/2) */ - si[i].vx+=(tau*si[i].fx/m2); - si[i].vy+=(tau*si[i].fy/m2); - si[i].vz+=(tau*si[i].fz/m2); - } - for(i=0;iLX) deltax-=LEN_X; - else if(-deltax>LX) deltax+=LEN_X; - deltax2=deltax*deltax; - deltay=si[i].y-si[j].y; - if(deltay>LY) deltay-=LEN_Y; - else if(-deltay>LY) deltay+=LEN_Y; - deltay2=deltay*deltay; - deltaz=si[i].z-si[j].z; - if(deltaz>LZ) deltaz-=LEN_Z; - else if(-deltaz>LZ) deltaz+=LEN_Z; - deltaz2=deltaz*deltaz; - distance=deltax2+deltay2+deltaz2; - if(distance<=R2_CUTOFF) { - tmp=1.0/distance; // 1/r^2 - lj1=tmp; // 1/r^2 - tmp*=tmp; // 1/r^4 - lj1*=tmp; // 1/r^6 - tmp*=tmp; // 1/r^8 - lj2=tmp; // 1/r^8 - lj1*=tmp; // 1/r^14 - lj1*=LJ_SIGMA_12; - lj2*=LJ_SIGMA_06; - lj=-2*lj1+lj2; - si[i].fx=lj*deltax; - si[i].fy=lj*deltay; - si[i].fz=lj*deltaz; - si[i].fx=-lj*deltax; - si[i].fy=-lj*deltay; - si[i].fz=-lj*deltaz; - } - } - /* calculation of new velocities v(t+h) */ - si[i].vx+=(tau*si[i].fx/m2); - si[i].vy+=(tau*si[i].fy/m2); - si[i].vz+=(tau*si[i].fz/m2); - } - - if(!(runs%10)) { - - /* rasmol script & xyz file */ - sprintf(xyz,"./saves/si-%.15f.xyz",time); - sprintf(ppm,"./video/si-%.15f.ppm",time); - fd1=open(xyz,O_WRONLY|O_CREAT|O_TRUNC); - if(fd1<0) { - perror("rasmol xyz file open"); - return -1; - } - dprintf(fd2,"load xyz %s\n",xyz); - dprintf(fd2,"spacefill 200\n"); - dprintf(fd2,"rotate x 11\n"); - dprintf(fd2,"rotate y 13\n"); - dprintf(fd2,"set ambient 20\n"); - dprintf(fd2,"set specular on\n"); - dprintf(fd2,"write ppm %s\n",ppm); - dprintf(fd2,"zap\n"); - dprintf(fd1,"%d\nsilicon\n",amount_si); - for(i=0;i