fixed 2bp, todo: adjust create lattice, then go for 3bp

[physik/posic.git] / moldyn.c

diff --git a/moldyn.c b/moldyn.c
index 7ea5d50..aaab0de 100644 (file)
--- a/moldyn.c
+++ b/moldyn.c
@@ -99,6 +99,12 @@ int set_dim(t_moldyn *moldyn,double x,double y,double z,u8 visualize) {
                moldyn->vis.dim.z=z;
        }
 
+       printf("[moldyn] dimensions in A:\n");
+       printf("  x: %f\n",moldyn->dim.x);
+       printf("  y: %f\n",moldyn->dim.y);
+       printf("  z: %f\n",moldyn->dim.z);
+       printf("  visualize simulation box: %s\n",visualize?"on":"off");
+
        return 0;
 }
 
@@ -802,6 +808,7 @@ int velocity_verlet(t_moldyn *moldyn) {
                v3_add(&(atom[i].v),&(atom[i].v),&delta);
        }
 
+moldyn_bc_check(moldyn);
        /* neighbour list update */
        link_cell_update(moldyn);
 
@@ -885,11 +892,12 @@ int potential_force_calc(t_moldyn *moldyn) {
                                        continue;
 
                                if((jtom->attr&ATOM_ATTR_2BP)&
-                                  (itom[i].attr&ATOM_ATTR_2BP))
+                                  (itom[i].attr&ATOM_ATTR_2BP)) {
                                        moldyn->func2b(moldyn,
                                                       &(itom[i]),
                                                       jtom,
                                                       bc_ij);
+                               }
 
                                /* 3 body potential/force */
 
@@ -1212,7 +1220,7 @@ int tersoff_mult_2bp(t_moldyn *moldyn,t_atom *ai,t_atom *aj,u8 bc) {
 
        /* f_a, df_a calc (again, same for ij and ji) | save for later use! */
        exchange->f_a=-B*exp(-mu*d_ij);
-       exchange->df_a=-mu*exchange->f_a/d_ij;
+       exchange->df_a=mu*exchange->f_a/d_ij;
 
        /* f_c, df_c calc (again, same for ij and ji) */
        if(d_ij<R) {
@@ -1354,7 +1362,7 @@ int tersoff_mult_post_2bp(t_moldyn *moldyn,t_atom *ai,t_atom *aj,u8 bc) {
        v3_scale(&force,&force,-0.5);
 
        /* add force */
-       v3_sub(&(ai->f),&(ai->f),&force);
+       v3_add(&(ai->f),&(ai->f),&force);
 
        return 0;
 }