/* albe 3 body potential function (first k loop) */
int albe_mult_3bp_k1(t_moldyn *moldyn,
- t_atom *ai,t_atom *aj,t_atom *ak,u8 bc) {
+ t_atom *ai,t_atom *aj,t_atom *ak,u8 bc) {
t_albe_mult_params *params;
t_albe_exchange *exchange;
}
/* force contribution */
- scale=-0.5*(f_c*(df_r+b*df_a)+df_c*(f_r+b*f_a));
+ scale=-0.5*(f_c*(df_r-b*df_a)+df_c*(f_r-b*f_a));
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
/* virial */
virial_calc(ai,&force,&(exchange->dist_ij));
- /* dzeta prefactor = - 0.5 f_c f_a db */
- exchange->pre_dzeta=-0.5*f_a*f_c*db;
+ /* dzeta prefactor = - f_c f_a db, (* -0.5 due to force calc) */
+ exchange->pre_dzeta=0.5*f_a*f_c*db;
/* energy contribution */
- energy=0.5*f_c*(f_r+b*f_a);
+ energy=0.5*f_c*(f_r-b*f_a);
moldyn->energy+=energy;
ai->e+=energy;
/* albe 3 body potential function (second k loop) */
int albe_mult_3bp_k2(t_moldyn *moldyn,
- t_atom *ai,t_atom *aj,t_atom *ak,u8 bc) {
+ t_atom *ai,t_atom *aj,t_atom *ak,u8 bc) {
t_albe_mult_params *params;
t_albe_exchange *exchange;