From cb177e7c208a85b45d77b09fcada23b62d0248b5 Mon Sep 17 00:00:00 2001 From: hackbard Date: Wed, 31 Jan 2007 16:35:26 +0000 Subject: [PATCH] create_lattice fixes, still virial probs ... --- moldyn.c | 85 +++++++++++++++++++++++-------------------------- moldyn.h | 3 +- report/report.h | 2 +- sic.c | 12 +++---- visual/visual.c | 3 +- 5 files changed, 49 insertions(+), 56 deletions(-) diff --git a/moldyn.c b/moldyn.c index 0d5027c..f6c3c81 100644 --- a/moldyn.c +++ b/moldyn.c @@ -339,22 +339,24 @@ int create_lattice(t_moldyn *moldyn,u8 type,double lc,int element,double mass, } moldyn->atom=ptr; atom=&(moldyn->atom[count]); - - v3_zero(&origin); + + /* no atoms on the boundaries (only reason: it looks better!) */ + origin.x=0.5*lc; + origin.y=0.5*lc; + origin.z=0.5*lc; switch(type) { case CUBIC: set_nn_dist(moldyn,lc); - origin.x=0.5*lc; - origin.y=0.5*lc; - origin.z=0.5*lc; ret=cubic_init(a,b,c,lc,atom,&origin); break; case FCC: + v3_scale(&origin,&origin,0.5); set_nn_dist(moldyn,0.5*sqrt(2.0)*lc); - ret=fcc_init(a,b,c,lc,atom,NULL); + ret=fcc_init(a,b,c,lc,atom,&origin); break; case DIAMOND: + v3_scale(&origin,&origin,0.25); set_nn_dist(moldyn,0.25*sqrt(3.0)*lc); ret=diamond_init(a,b,c,lc,atom,&origin); break; @@ -429,65 +431,55 @@ int cubic_init(int a,int b,int c,double lc,t_atom *atom,t_3dvec *origin) { int fcc_init(int a,int b,int c,double lc,t_atom *atom,t_3dvec *origin) { int count; - int i,j; + int i,j,k,l; t_3dvec o,r,n; t_3dvec basis[3]; - double help[3]; - double x,y,z; - - x=a*lc; - y=b*lc; - z=c*lc; - if(origin) v3_copy(&o,origin); - else v3_zero(&o); + count=0; + if(origin) + v3_copy(&o,origin); + else + v3_zero(&o); /* construct the basis */ - for(i=0;i<3;i++) { - for(j=0;j<3;j++) { - if(i!=j) help[j]=0.5*lc; - else help[j]=.0; - } - v3_set(&basis[i],help); - } + memset(basis,0,3*sizeof(t_3dvec)); + basis[0].x=0.5*lc; + basis[0].y=0.5*lc; + basis[1].x=0.5*lc; + basis[1].z=0.5*lc; + basis[2].y=0.5*lc; + basis[2].z=0.5*lc; - v3_zero(&r); - count=0; - /* fill up the room */ r.x=o.x; - while(r.xf),&(ai->f),&force); virial_calc(ai,&force,&distance); +if(force.x*distance.x<=0) printf("virial xx: %.15f -> %f %f %f\n",force.x*distance.x,distance.x,distance.y,distance.z); virial_calc(aj,&force,&distance); /* f and d signe switched */ } diff --git a/moldyn.h b/moldyn.h index f2ccdb5..3c53a86 100644 --- a/moldyn.h +++ b/moldyn.h @@ -334,7 +334,8 @@ typedef struct s_tersoff_mult_params { #define M_SI 28.08553 /* amu */ //#define LJ_SIGMA_SI ((0.25*sqrt(3.0)*LC_SI)/1.122462) /* A */ -#define LJ_SIGMA_SI (LC_SI/1.122462) /* A */ +//#define LJ_SIGMA_SI (LC_SI/1.122462) /* A */ +#define LJ_SIGMA_SI (0.5*sqrt(2.0)*LC_SI/1.122462) /* A */ #define LJ_EPSILON_SI (2.1678*EV) /* NA */ #define TM_R_SI (2.7e-10*METER) /* A */ diff --git a/report/report.h b/report/report.h index 3868515..e46ff2a 100644 --- a/report/report.h +++ b/report/report.h @@ -53,7 +53,7 @@ set ytic auto \n\ set title 'Energy per atom vs. time' \n\ set xlabel 'Time [fs]' \n\ set ylabel 'Energy per atom [eV]' \n\ -set terminal postscript eps enhanced color dashed lw 1 'Helvetica' 14 \n\ +set terminal postscript eps enhanced color solid lw 1 'Helvetica' 14 \n\ set output 'energy.eps' \n\ plot \"energy\" using 1:2 title 'Kinetic energy' with lines , \"energy\" using 1:3 title 'Potential energy' with lines , \"energy\" using 1:4 title 'Total energy' with lines\ "; diff --git a/sic.c b/sic.c index 2a2c7bc..f97ffe8 100644 --- a/sic.c +++ b/sic.c @@ -121,8 +121,8 @@ int main(int argc,char **argv) { set_pbc(&md,TRUE,TRUE,TRUE); /* create the lattice / place atoms */ - create_lattice(&md,CUBIC,LC_SI,SI,M_SI, - //create_lattice(&md,FCC,LC_SI,SI,M_SI, + //create_lattice(&md,CUBIC,LC_SI,SI,M_SI, + create_lattice(&md,FCC,LC_SI,SI,M_SI, //create_lattice(&md,DIAMOND,LC_SI,SI,M_SI, // ATOM_ATTR_1BP|ATOM_ATTR_2BP|ATOM_ATTR_3BP|ATOM_ATTR_HB, ATOM_ATTR_2BP|ATOM_ATTR_HB, @@ -151,8 +151,8 @@ int main(int argc,char **argv) { set_pressure(&md,ATM); /* set p/t scaling */ - set_pt_scale(&md,P_SCALE_BERENDSEN,0.001, - T_SCALE_BERENDSEN,100.0); + //set_pt_scale(&md,P_SCALE_BERENDSEN,0.001, + // T_SCALE_BERENDSEN,100.0); //set_pt_scale(&md,0,0,T_SCALE_DIRECT,1.0); //set_pt_scale(&md,P_SCALE_BERENDSEN,0.001,0,0); @@ -160,7 +160,7 @@ int main(int argc,char **argv) { thermal_init(&md,TRUE); /* create the simulation schedule */ - moldyn_add_schedule(&md,100001,1.0); + moldyn_add_schedule(&md,101,1.0); //moldyn_add_schedule(&md,501,1.0); //moldyn_add_schedule(&md,501,1.0); @@ -171,7 +171,7 @@ int main(int argc,char **argv) { moldyn_set_log_dir(&md,argv[1]); moldyn_set_report(&md,"Frank Zirkelbach","Test 1"); moldyn_set_log(&md,LOG_TOTAL_ENERGY,1); - moldyn_set_log(&md,VISUAL_STEP,1000); + moldyn_set_log(&md,VISUAL_STEP,1); moldyn_set_log(&md,CREATE_REPORT,0); /* diff --git a/visual/visual.c b/visual/visual.c index d6ce369..573ac4f 100644 --- a/visual/visual.c +++ b/visual/visual.c @@ -87,8 +87,7 @@ int visual_atoms(t_visual *v,double time,t_atom *atom,int n) { /* script file update */ dprintf(v->fd,"load xyz %s\n",file); - dprintf(v->fd,"spacefill\n"); - //dprintf(v->fd,"spacefill 100\n"); + dprintf(v->fd,"spacefill 200\n"); dprintf(v->fd,"rotate x 100\n"); dprintf(v->fd,"rotate y 10\n"); dprintf(v->fd,"set ambient 20\n"); -- 2.20.1