X-Git-Url: https://hackdaworld.org/gitweb/?p=physik%2Fposic.git;a=blobdiff_plain;f=posic.c;h=6b92a5e51b0428de417e34392189a6aa246abffb;hp=7c1c41848aff58df3068542f0e6b7b9541bce733;hb=153ffa76a126130faca5be28aefad32f4c8b623d;hpb=792f14f95b47989f7f12df0ea70b54619be016ee diff --git a/posic.c b/posic.c index 7c1c418..6b92a5e 100644 --- a/posic.c +++ b/posic.c @@ -5,144 +5,38 @@ * */ -#include - -#include "moldyn.h" -#include "math/math.h" -#include "init/init.h" -#include "visual/visual.h" +/* main include file */ #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/lj_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)); +/* functions */ - visual_init(&vis,fb); - a=LEN_X; - b=LEN_Y; - c=LEN_Z; - /* set for 'bounding atoms' */ - //vis.dim.x=a*LC_SI; - //vis.dim.y=b*LC_SI; - //vis.dim.z=c*LC_SI; +/* main code */ - t=TEMPERATURE; +int parse_config_file() { - printf("placing silicon atoms ... "); - count=create_lattice(DIAMOND,SI,M_SI,LC_SI,a,b,c,&si); - printf("(%d) ok!\n",count); - - /* testing purpose - count=2; - si=malloc(2*sizeof(t_atom)); - si[0].r.x=0.16e-9; - si[0].r.y=0; - si[0].r.z=0; - si[0].element=SI; - si[0].mass=M_SI; - si[1].r.x=-0.16e-9; - si[1].r.y=0; - si[1].r.z=0; - si[1].element=SI; - si[1].mass=M_SI; - */ - - 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: %.40f [J]\n",e); - printf("3/2 N k T = %.40f [J]\n",1.5*count*K_BOLTZMANN*t); + return 0; +} - /* check total momentum */ - p=get_total_p(si,count); - printf("total momentum: %.30f [Ns]\n",v3_norm(&p)); +int main(int argc,char **argv) { - /* check potential energy */ - md.count=count; - md.atom=si; - md.potential=potential_lennard_jones; - md.force=force_lennard_jones; - md.cutoff=R_CUTOFF; - md.cutoff_square=(R_CUTOFF*R_CUTOFF); - md.pot_params=&lj; - md.integrate=velocity_verlet; - md.time_steps=RUNS; - md.tau=TAU; - md.status=0; - md.visual=&vis; - md.write=WRITE_FILE; + t_moldyn md; - lj.sigma6=LJ_SIGMA_SI*LJ_SIGMA_SI; - help=lj.sigma6*lj.sigma6; - lj.sigma6*=help; - lj.sigma12=lj.sigma6*lj.sigma6; - lj.epsilon=LJ_EPSILON_SI; + t_lj_params *lj; + t_ho_params *ho; + t_tersoff_mult_params *tp; + t_albe_mult_params *ap; - u=get_e_pot(&md); + lj=NULL; + ho=NULL; + tp=NULL; + ap=NULL; - printf("potential energy: %.40f [J]\n",u); - printf("total energy (1): %.40f [J]\n",e+u); - printf("total energy (2): %.40f [J]\n",get_total_energy(&md)); + memset(&md,0,sizeof(t_moldyn)); - md.dim.x=a*LC_SI; - md.dim.y=b*LC_SI; - md.dim.z=c*LC_SI; - printf("estimated accurate time step: %.30f [s]\n", - estimate_time_step(&md,3.0,t)); - - - /* - * let's do the actual md algorithm now - * - * integration of newtons equations - */ - - moldyn_integrate(&md); - - printf("total energy (after integration): %.40f [J]\n", - get_total_energy(&md)); - - /* close */ - - visual_tini(&vis); - - rand_close(&random); - return 0; }