for(i=0;i<lc->cells;i++)
list_destroy(&(moldyn->lc.subcell[i]));
- for(count=0;count<moldyn->count;count++) {
+ for(count=0;count<moedyn->count;count++) {
i=(atom[count].r.x+(moldyn->dim.x/2))/lc->x;
j=(atom[count].r.y+(moldyn->dim.y/2))/lc->y;
k=(atom[count].r.z+(moldyn->dim.z/2))/lc->z;
/* 2 body stuff */
+ /* we need: f_c, df_c, f_r, df_r */
+
v3_sub(&dist_ij,btom,&(atom[i]));
d_ij=v3_norm(&dist_ij);
if(d_ij<=S) {
- S=;
- R=;
+ /* determine the tersoff parameters */
+ if(atom[i].element!=btom->element) {
+ S=sqrt(TERSOFF_S[e1]*TERSOFF_S[e2]);
+ R=R_m;
A=;
lambda=;
B=;
betaN=;
if(d_ij<=R) {
- f_c=1.0;
- df_c=0.0;
+ df_r=-lambda*A*exp(-lambda*d_ij)/d_ij;
+ v3_scale(&force,&dist_ij,df_r);
+ v3_add(&(atom[i].f),&(atom[i].f),
+ &force);
}
else {
s_r=S-R;
f_c=0.5+0.5*cos(arg1);
df_c=-0.5*sin(arg1)*(PI/(s_r*d_ij));
f_r=A*exp(-lambda*d_ij);
- f_a=-B*exp(-mu*d_ij);
+ df_r=-lambda*f_r/d_ij;
+ scale=df_c*f_r+df_r*f_c;
+ v3_scale(&force,&dist_ij,scale);
+ v3_add(&(atom[i].f),&(atom[i].f),
+ &force);
}
}
else
v3_sub(&distance_jk,ktom,btom);
cos_theta=(d_ij2+d_ik*d_ik-d_jk*d_jk)/\
(2*d_ij*d_ik);
+ sin_theta=sqrt(1.0/\
+ (cos_theta*cos_theta));
theta=arccos(cos_theta);