# Colloidal networks

Authors: Yoldas, Frank, Laura

## Plans

### NVT Monte Carlo Simulation

We will explore the behavior of colloidal networks of patchy particles using simulations and Wertheim theory.

## To Do

:::success
green is for finished tasks
:::
### Simulations
:::success
* Get the simulation running
* Implement visualization in Python
* Plot behavior of the energy and pressure as a function of time
* Plot behavior of the energy and pressure as a function of T
* Plot behavior of the energy and pressure as a function of denstiy
:::

### 21/05: Notes
* vizualisation code is done. Some notes on the code
* python keeps the location in mind when the patch file is read.
* new command learned: **np.dot** is for the matrix multiplication of numpy arrays.
* An example image that seems to match the linux code:

### 24/05: Notes
* plotted potential and pressure with regards to time.
* it is clear high epsilon needs more time to reach equilibria --> skip first 200 rows
* run the simulation for different epsilon and plot
* as $\epsilon$ increases the pressure has higher fluctuations, even reaching negative values.
* 
* as $\epsilon$ increases pot energy goes down but fluctuation isn't as high.
* 
* numerically for epsilon 10
* potential= -13.535601753333333 std = 0.09570560566061856
pressure= -0.03288044333333333 std= 0.07562040064039237
* the std for pressure is so high

* for epsilon 0.1, 5.5 and 10
 

### 02/06: Notes
**To Do**
- [x] plot $\Delta F$ with different values of epsilon and density
* Derivation of Gamma
$d\Gamma_{12} = (1/2\pi)^2 d\vec{r}_{12} d\theta_1 d\theta_2$
$d\Gamma_{12} = (1/2\pi)^2 \int_{\sigma}^{r}\int_{0}^{2\pi}rsin\theta dr d\theta \int_{-\theta_{max}}^{\theta_{max}} d\theta_1\int_{-\theta_{max}}^{\theta_{max}} d\theta_2$
$d\Gamma_{12} = \frac{(r_{max}^2-\sigma^2)\theta_{max}^2}{\pi}$
* Seems like there is an issue? Because delta F is positive? I expected the difference to get larger as we have higher $\epsilon$ and $n$ but it is positive because the term inside log is too small.
- [x] from $\delta F$ calculate the bonding probability, the number of bonds in the system and so the potential energy.
- [x] compare with data from the simulation
- 
- theory underestimates the number of bonds, probably becuase we took g{r} to be 1. In reality the particles are more densely packed.
- the full code is in name theory_changinT

### 03/06 Notes
- [x] incorporate g_r from Yuste to the above case
* It is done. Gives more matching results than the first order approximation. 
* changed the code in a way to find $p_b$ analytically. `fsolve` function choose the wrong solution occasionally. Since it works depending on the initial guess.

### 04/06 Notes

- [x] calculate pressure using Bianchi(2008) paper as a function of density and compare with simulation
- [ ] Find expression for free energy as a funciton of density to investigate phase behaviour

The expression for the free energy as a function of density.
$$F(\rho)=F^\mathrm{id}(\rho) + N \int_0^\rho \mathrm{d}\rho' \frac{P(\rho') - \rho' k_B T}{\rho'^2}$$

$$F^\mathrm{id}(\rho) = k_B T N\left[\log(\rho \Lambda^2) - 1\right]$$
Here, $\rho$ is the number density, $P(\rho)$ is the pressure at density $\rho$, and $\Lambda$ is the thermal wavelength. In practice, $\Lambda$ does not affect the phase behavior, so you can simply set it equal to $\sigma$.

$P(\rho) = -\frac{d(F)}{dV} = \frac{d(F/N)}{d\rho} \rho^2 = \frac{1}{N}\frac{dF}{d\rho} \rho^2$
(Check my math)

$\frac{d(F/V)}{d\rho} = \frac{d(F)}{d\rho}/V + F(\rho)\frac{d(1/V)}{d\rho}$