# Predicting Networks (Part 2)

###### tags: `masters projects` `machine learning` `patchy particles` `glassy dynamics`

---

### Owners (the only one with the permission to edit the main text, others can comment)

Duncan, Frank, Laura

---

# Energy minimalisation

Idea: Take a snapshot: minize free energy based on some potential

Steepest descent:

C-code:
1) read in snapshot
2) calculate energy
3) calculate forces
4) slightly adjust positions
5) goto 3 (until forces are small enough)

new position = old position + $\gamma \mathbf{f}$

Bound particles:

Sum over all bound pairs::
$U(r) = k (r-r_0)^2$
(r_0 = 1.06)

Sum over all bound triplets:
$U(\theta_{ijk}) = k (\cos \theta_{ijk} - \cos 2\pi /3)^2$

Sum over Non-bound pairs:
$U(r) = K (1-r)^2$ if $r < 1$, use a big $K$.

Cell list

# Bond type predictability

|||
|:---------------------------------:|:---------------------------------:|
|||
|||
|||
|| |

# Cutoff bond times

|||
|:---------------------------------:|:---------------------------------:|

# Energy Minimalisation code results (Energy problem is fixed.)

I wrote code that attempts to relax a system of particles that can each undergo three bonds according to the following energy contributions:

For all particles: $U(r)=k1(1-r)^2$ for $r<1$
For bonded particles: $U(r)=k2(r-1.06)^2$ for $r>1$
For bonded particles: $U(\theta)=k3(\cos(\theta)-\cos(2\pi/3))^2$

To relax the system I calculate the force due to each of the energy contributions on each particle using $F(r)=-\nabla U(r)$.

For particles $i$ and $j$ we have: $r_{ij}=\sqrt{(x_j-x_i)^2+(y_j-y_i)^2}$.
For bonded particles $i$, $j$ and $k$ we also have: $\theta_{ijk}=\arccos(\frac{(x_j-x_i)(x_j-x_k)+(y_j-y_i)(y_j-y_k)}{r_{ji}*r_{jk}})$

Using this we can calculate the forces:

For all particles:
$F_{x_i}(x_i)=-2k1*\frac{x_j-x_i}{r_{ij}}*(1-r_{ij})$ for $r_{ij}<1$
$F_{y_i}(y_i)=-2k1*\frac{y_j-y_i}{r_{ij}}*(1-r_{ij})$ for $r_{ij}<1$
For bonded particles:
$F_{x_i}(x_i)=2k2*\frac{x_j-x_i}{r_{ij}}*(r_{ij}-1.06)$ for $r_{ij}>1$
$F_{y_i}(y_i)=2k2*\frac{y_j-y_i}{r_{ij}}*(r_{ij}-1.06)$ for $r_{ij}>1$
For forces due to the angle of bonded particles $i, j$ and $k$ we have:
$F_{x_i}(x_i)=-2k3*(\frac{(x_j-x_i)(x_j-x_k)+(y_j-y_i)(y_j-y_k)}{\sqrt{(x_j-x_i)^2+(y_j-y_i)^2}*\sqrt{(x_j-x_k)^2+(y_j-y_k)^2}}+0.5)*$
$(\frac{-(x_j-x_k)}{\sqrt{(x_j-x_i)^2+(y_j-y_i)^2}*\sqrt{(x_j-x_k)^2+(y_j-y_k)^2}}+\frac{(x_j-x_i)*((x_j-x_i)(x_j-x_k)+(y_j-y_i)(y_j-y_k))}{((x_j-x_i)^2+(y_j-y_i)^2)^{1.5}*\sqrt{(x_j-x_k)^2+(y_j-y_k)^2}})$

$F_{y_i}(y_i)=-2k3*(\frac{(x_j-x_i)(x_j-x_k)+(y_j-y_i)(y_j-y_k)}{\sqrt{(x_j-x_i)^2+(y_j-y_i)^2}*\sqrt{(x_j-x_k)^2+(y_j-y_k)^2}}+0.5)*$
$(\frac{-(y_j-y_k)}{\sqrt{(x_j-x_i)^2+(y_j-y_i)^2}*\sqrt{(x_j-x_k)^2+(y_j-y_k)^2}}+\frac{(y_j-y_i)*((x_j-x_i)(x_j-x_k)+(y_j-y_i)(y_j-y_k))}{((x_j-x_i)^2+(y_j-y_i)^2)^{1.5}*\sqrt{(x_j-x_k)^2+(y_j-y_k)^2}})$

$F_{x_k}(x_k)=-2k3*(\frac{(x_j-x_i)(x_j-x_k)+(y_j-y_i)(y_j-y_k)}{\sqrt{(x_j-x_i)^2+(y_j-y_i)^2}*\sqrt{(x_j-x_k)^2+(y_j-y_k)^2}}+0.5)*$
$(\frac{-(x_j-x_i)}{\sqrt{(x_j-x_i)^2+(y_j-y_i)^2}*\sqrt{(x_j-x_k)^2+(y_j-y_k)^2}}+\frac{(x_j-x_k)*((x_j-x_i)(x_j-x_k)+(y_j-y_i)(y_j-y_k))}{\sqrt{(x_j-x_i)^2+(y_j-y_i)^2}*((x_j-x_k)^2+(y_j-y_k)^2)^{1.5}})$

$F_{y_k}(y_k)=-2k3*(\frac{(x_j-x_i)(x_j-x_k)+(y_j-y_i)(y_j-y_k)}{\sqrt{(x_j-x_i)^2+(y_j-y_i)^2}*\sqrt{(x_j-x_k)^2+(y_j-y_k)^2}}+0.5)*$
$(\frac{-(y_j-y_i)}{\sqrt{(x_j-x_i)^2+(y_j-y_i)^2}*\sqrt{(x_j-x_k)^2+(y_j-y_k)^2}}+\frac{(y_j-y_k)*((x_j-x_i)(x_j-x_k)+(y_j-y_i)(y_j-y_k))}{\sqrt{(x_j-x_i)^2+(y_j-y_i)^2}*((x_j-x_k)^2+(y_j-y_k)^2)^{1.5}})$

Once I calculate all forces I update the snapshot by changing all particle coordinates: $x_{i,new}=x_{i,old}+F_{x_i}*dt$, $y_{i,new}=y_{i,old}+F_{y_i}*dt$.

I did this once with all energy and forces, once with just energy and forces due to particle distance and once with only energy and forces due to angles between bonds. For the k values I choose: k1=64.0 k2=16.0 and k3=2.0. The results are shown below:

|All energy and force||
|:---------------------------------:|:---------------------------------:|
|||
|Only distance||
|||
|Only angle||
|||

During relaxation we want the force to go down to zero, which would mean we have found a minimum in the energy. We also want to see the energy go down to some value and then remain at this value.
We see that this happens nicely when only considering the distance between particles, however when we also look at the angles we see that at a certain point the energy starts to go up again. And for the case where we only look at energy due to angles we see that after some time the force blows up and the simulation stops.

To test whether my force and energy calculation were right I copied my c code into Mathematica and made a plot where I can easily change the angle of the bond and the orientation. I then plot the forces on the two outer particles as an arrow and the energy is represented by a red dot, where the bigger the red dot the bigger the energy. I do this for various angles and orientations to see if all cases work. They do seem to work. I also compare this to the plots I have of the energy and force vs bond angle.
|||
|:---------------------------------:|:---------------------------------:|
|||
|||

To try and make sure the step size is not to big to find the minimum I divide the time step deltat by 10 every time the energy of a step is larger then the energy of the previous step.
|||
|:---------------------------------:|:---------------------------------:|
When I do this the energy and force decrease until the force is at around 100. From that point the system stops changing. So something must be going wrong since it seems that the force is pointing in a direction that increases system energy.

So far I have not been able to find my mistake, but I still have some things I can try.

## Fixed energy
For the relaxation results above I kept applying forces on the two outer particles of a bond pair angle. Now I tried to fix the problem by applying forces on all three particles involved in the bond. (Since then the total force is zero...) This solved the issue of the energy growing.

|Forces on all three particles||
|:---------------------------------:|:---------------------------------:|
|||
Now the energy drops to a static value and stays there while the force still goes to zero. In the final snapshot you see that the angles look better then for the case without taking angles into account. (A good example of this is the red bond bottom right.)