db=-0.5*b/(1.0+exchange->zeta_ij);
}
- /* force contribution */
+ /* force contribution for atom i */
scale=-0.5*(f_c*(df_r-b*df_a)+df_c*(f_r-b*f_a)); // - in albe formalism
v3_scale(&force,&(exchange->dist_ij),scale);
v3_add(&(ai->f),&(ai->f),&force);
- v3_sub(&(aj->f),&(aj->f),&force); // dri rij = - drj rij
+
+ /* force contribution for atom j */
+ v3_scale(&force,&force,-1.0); // dri rij = - drj rij
+ v3_add(&(aj->f),&(aj->f),&force);
+
+ /* virial */
+ virial_calc(aj,&force,&(exchange->dist_ij));
#ifdef DEBUG
if(moldyn->time>DSTART&&moldyn->time<DEND) {
}
#endif
- /* virial */
- virial_calc(ai,&force,&(exchange->dist_ij));
-
/* dzeta prefactor = - f_c f_a db, (* -0.5 due to force calc) */
exchange->pre_dzeta=0.5*f_a*f_c*db;