2 * harmonic_oscillator.c - harmonic oscillator potential
4 * author: Frank Zirkelbach <frank.zirkelbach@physik.uni-augsburg.de>
12 #include <sys/types.h>
18 #include "../moldyn.h"
19 #include "../math/math.h"
20 #include "harmonic_oscillator.h"
22 int harmonic_oscillator(t_moldyn *moldyn,t_atom *ai,t_atom *aj,u8 bc) {
25 t_3dvec force,distance;
30 params=moldyn->pot_params;
31 sc=params->spring_constant;
32 equi_dist=params->equilibrium_distance;
36 v3_sub(&distance,&(aj->r),&(ai->r));
38 if(bc) check_per_bound(moldyn,&distance);
40 if(d<=moldyn->cutoff) {
41 energy=(0.5*sc*(d-equi_dist)*(d-equi_dist));
42 moldyn->energy+=energy;
45 /* f = -grad E; grad r_ij = -1 1/r_ij distance */
46 f=sc*(1.0-equi_dist/d);
47 v3_scale(&force,&distance,f);
48 v3_add(&(ai->f),&(ai->f),&force);
49 virial_calc(ai,&force,&distance);
50 virial_calc(aj,&force,&distance); /* f and d signe switched */
51 v3_scale(&force,&distance,-f);
52 v3_add(&(aj->f),&(aj->f),&force);