*
*/
-#include <math.h>
-
-#include "moldyn.h"
-#include "math/math.h"
-#include "init/init.h"
-#include "visual/visual.h"
+/* main include file */
#include "posic.h"
+/* functions */
+
+
+
+/* main code */
+
+int parse_config_file() {
+
+ return 0;
+}
+
int main(int argc,char **argv) {
t_moldyn md;
- t_atom *si;
- t_visual vis;
- t_random random;
-
- int a,b,c;
- double e,u;
- double help;
- t_3dvec p;
- int count;
-
- t_lj_params lj;
- t_ho_params ho;
-
- /* parse arguments */
- a=moldyn_parse_argv(&md,argc,argv);
- if(a<0) return -1;
-
- /* init */
- moldyn_log_init(&md,&vis);
- 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));
-
- 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;
-
- /* init lattice */
- 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.23*sqrt(3.0)*LC_SI/2.0;
- si[0].r.y=0;
- si[0].r.z=0;
- si[0].element=SI;
- si[0].mass=M_SI;
- si[1].r.x=-si[0].r.x;
- si[1].r.y=0;
- si[1].r.z=0;
- si[1].element=SI;
- si[1].mass=M_SI;
- */
-
- /* moldyn init (now si is a valid address) */
- md.count=count;
- md.atom=si;
- md.potential=potential_lennard_jones;
- md.force=force_lennard_jones;
- //md.potential=potential_harmonic_oscillator;
- //md.force=force_harmonic_oscillator;
- md.cutoff=R_CUTOFF;
- md.cutoff_square=(R_CUTOFF*R_CUTOFF);
- md.pot_params=&lj;
- //md.pot_params=&ho;
- md.integrate=velocity_verlet;
- //md.time_steps=RUNS;
- //md.tau=TAU;
- md.status=0;
- md.visual=&vis;
-
- printf("setting thermal fluctuations (T=%f K)\n",md.t);
- thermal_init(&md,&random,count);
- //for(a=0;a<count;a++) v3_zero(&(si[0].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*md.t);
-
- /* check total momentum */
- p=get_total_p(si,count);
- printf("total momentum: %.30f [Ns]\n",v3_norm(&p));
-
- /* check potential energy */
- lj.sigma6=LJ_SIGMA_SI*LJ_SIGMA_SI;
- help=lj.sigma6*lj.sigma6;
- lj.sigma6*=help;
- lj.sigma12=lj.sigma6*lj.sigma6;
- lj.epsilon4=4.0*LJ_EPSILON_SI;
-
- ho.equilibrium_distance=0.25*sqrt(3.0)*LC_SI;
- ho.spring_constant=4.0*LJ_EPSILON_SI;
-
- u=get_e_pot(&md);
-
- 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));
-
- 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,md.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 */
-
- rand_close(&random);
-
- moldyn_shutdown(&md);
-
+
+ t_lj_params *lj;
+ t_ho_params *ho;
+ t_tersoff_mult_params *tp;
+ t_albe_mult_params *ap;
+
+ lj=NULL;
+ ho=NULL;
+ tp=NULL;
+ ap=NULL;
+
+ memset(&md,0,sizeof(t_moldyn));
+
+
return 0;
}