+ v3_sub(&dist_ij,btom,&(atom[i]));
+ d_ij=v3_norm(&dist_ij);
+ if(d_ij<=S) {
+ if(d_ij<=R) {
+ f_c=1.0;
+ df_c=0.0;
+ }
+ else {
+ s_r=S-R;
+ arg1=PI*(d_ij-R)/s_r;
+ f_c=0.5+0.5*cos(arg1);
+ df_c=-0.5*sin(arg1)*(PI/(s_r*d_ij));
+ }
+ }
+ else
+ continue;
+
+
+ /* end 2 body stuff */
+