X-Git-Url: https://hackdaworld.org/gitweb/?p=physik%2Fposic.git;a=blobdiff_plain;f=posic.c;h=6b92a5e51b0428de417e34392189a6aa246abffb;hp=926e20c4e43ee38965e0bb62e01b87b09f7b858d;hb=HEAD;hpb=83775c491117faa149281d0302fc8b8064d6b080 diff --git a/posic.c b/posic.c index 926e20c..6b92a5e 100644 --- a/posic.c +++ b/posic.c @@ -1,262 +1,41 @@ /* * posic.c - precipitation process of silicon carbide in silicon * - * author: Frank Zirkelbach + * author: Frank Zirkelbach * */ - -#include "moldyn.h" -#include "math/math.h" -#include "init/init.h" -#include "visual/visual.h" -#include "posic.h" - -int main(int argc,char **argv) { - - t_moldyn md; - - t_atom *si; - - t_visual vis; - - t_random random; - - int a,b,c; - double t,e,u; - double help; - t_3dvec p; - int count; - - t_lj_params lj; - - char fb[32]="saves/fcc_test"; - - /* init */ - - rand_init(&random,NULL,1); - random.status|=RAND_STAT_VERBOSE; - - /* testing random numbers */ - //for(a=0;a<1000000;a++) - // printf("%f %f\n",rand_get_gauss(&random), - // rand_get_gauss(&random)); - - visual_init(&vis,fb); +/* main include file */ - a=LEN_X; - b=LEN_Y; - c=LEN_Z; - - t=TEMPERATURE; - - printf("placing silicon atoms ... "); - //count=create_lattice(DIAMOND,Si,M_SI,LC_SI,a,b,c,&si); - //printf("(%d) ok!\n",count); - count=2; - si=malloc(2*sizeof(t_atom)); - si[0].r.x=2.0; - si[0].r.y=0; - si[0].r.z=0; - si[0].element=Si; - si[0].mass=14.0; - si[1].r.x=-2.0; - si[1].r.y=0; - si[1].r.z=0; - si[1].element=Si; - si[1].mass=14.0; - - printf("setting thermal fluctuations\n"); - //thermal_init(si,&random,count,t); - v3_zero(&(si[0].v)); - v3_zero(&(si[1].v)); - - /* check kinetic energy */ - - e=get_e_kin(si,count); - printf("kinetic energy: %f\n",e); - printf("3/2 N k T = %f\n",1.5*count*K_BOLTZMANN*t); - - /* check total momentum */ - p=get_total_p(si,count); - printf("total momentum: %f\n",v3_norm(&p)); - - /* check potential energy */ - md.count=count; - md.atom=si; - md.potential=potential_lennard_jones; - md.force=force_lennard_jones; - //md.cutoff_square=((LC_SI/4.0)*(LC_SI/4.0)); - md.cutoff_square=36.0; - md.pot_params=&lj; - md.integrate=velocity_verlet; - md.time_steps=RUNS; - md.tau=TAU; - md.status=0; - md.visual=&vis; - - lj.sigma6=3.0/16.0*LC_SI*LC_SI; - help=lj.sigma6*lj.sigma6; - lj.sigma6*=help; - lj.sigma12=lj.sigma6*lj.sigma6; - lj.epsilon=10000; - - u=get_e_pot(&md); +#include "posic.h" - printf("potential energy: %f\n",u); - printf("total energy (1): %f\n",e+u); - printf("total energy (2): %f\n",get_total_energy(&md)); +/* functions */ - md.dim.x=a*LC_SI; - md.dim.y=b*LC_SI; - md.dim.z=c*LC_SI; - /* - * let's do the actual md algorithm now - * - * integration of newtons equations - */ - /* visualize */ - //visual_atoms(&vis,0.0,si,count); - +/* main code */ - moldyn_integrate(&md); +int parse_config_file() { - printf("total energy (after integration): %f\n",get_total_energy(&md)); + return 0; +} - /* close */ +int main(int argc,char **argv) { - visual_tini(&vis); + t_moldyn md; - rand_close(&random); - + t_lj_params *lj; + t_ho_params *ho; + t_tersoff_mult_params *tp; + t_albe_mult_params *ap; - //printf("starting velocity verlet: "); - //fflush(stdout); + lj=NULL; + ho=NULL; + tp=NULL; + ap=NULL; - //for(runs=0;runsLX) 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); -// /* reset of forces */ -// si[i].fx=.0; -// si[i].fy=.0; -// si[i].fz=.0; -// } -// 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[j].fx+=lj*deltax; -// si[j].fy+=lj*deltay; -// si[j].fz+=lj*deltaz; -// } -// } -// } -// for(i=0;i