- 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);
+ 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;
+
+ ho.equilibrium_distance=0.3e-9;
+ ho.spring_constant=1.0e-9;
+
+ 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,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));