Vacancy Concentration in Monte Carlo Model

Introduction to Monte Carlo Methods

In Section 7.4 we introduced the random walk model as an extremely simplified model of atomic movement compared to molecular dynamics (MD) models (if you don't remember our justification for the random walk model, feel free to review that page). Now, we will extend the random walk model to take into account energy differences at different locations.

In the random walk model, we assume atoms are located on a lattice and they randomly "jump" to neighboring sites due collisions constantly pushing them in random directions. However, in reality, not all jumps are equally likely. For an atom to jump to a neighboring location, the bonds with atoms at the current location must be broken, which requires energy, specifically kinetic energy. If an atom has more kinetic energy than the potential energy (PE) of its bonds with neighboring atoms, it will break those bonds and move. The more bonds that must be broken, the less likely that the atom will have enough energy to break them and make the jump. The energy of the new bonds being formed need to be taken into account as well. As a simplification, we assume that the atom only needs enough KE to make up the difference between the PE of the two locations. In reality, the atom would likely need more KE than this, because at intermediate positions the PE would likely be higher than either the starting or ending location, but we still get qualitatively correct behavior with this simplification. We then model atoms as "attempting" to jump and succeeding if they randomly have enough energy to overcome the PE difference between the locations.

The atom making the jump on the left ends up with 3 fewer bonds. So, it would need enough KE to break three bonds. The atom making the jump on the right has the same number of bonds before and after the jump. So, in our model, it will always make that jump.

Figure 8.5.1 The atom making the jump on the left ends up with 3 fewer bonds. So, it would need enough KE to break three bonds. The atom making the jump on the right has the same number of bonds before and after the jump. So, in our model, it will always make that jump.

The actual algorithm, which is a Monte Carlo method, works like this:

  • The atom calculates the PE of its current location, $U_0$, which is equal to the number of neighbors it has times the PE associated with each bond.
  • The atom "jumps" to a new neighboring location
  • The atom calculates the PE of the new location, $U_1$.
  • If $U_1 < U_0$, the atom always "accepts" the jump. This is because even if the atom had zero KE, it would feel forces accelerating it towards a neighboring location with lower PE.
  • If $U_1 >= U_0$ and the atom randomly has KE greater than or equal to the to overcome the PE difference, it "accepts" the jump. If not, it "rejects" the jump and goes back to the prior location.

How likely is it for an atom to have enough energy to make a jump? On the previous page, Section 8.4, we saw that the distribution of KE among atoms in a MD simulation is an exponential distribution with a variance (how far the tail extends) dependent on temperature. So, we use an exponential distribution to determine how likely the atom is to have enough KE to overcome the PE difference between the starting and ending locations. Figure 8.5.2 shows two KE distributions at different temperatures, illustrating that at higher temperatures, atoms are much more likely to have high KE. Mathematically, the probability of the atom having enough energy to make a jump is:

$$P=\exp(\frac{-\Delta U} {k_b T}) \tag{8.5.1}$$

Where:

  • $P$ is the probability of making the jump
  • $\Delta U$ is the energy needed to make the jump (in our case the difference between the PE before/after a jump). This is often called the activation energy
  • $k_b$ is the Boltzman's constant, which has units of $\pu{J*K-1}$
  • $T$ is the temperature in units of $\pu{K}$
The two histograms show the KE probability distributions for atoms at two different temperature. These were generated using a slightly modified version of @ref(CB10579). The distribution on the left is for a lower temperature. The large majority of the probability mass is for KE < 1.5 (blue).  A small amount is for 1.5 < KE < 3 (green) and almost none for KE > 3 (red). In contrast, in the distribution at higher temperature, atoms are much more likely to have high KE.

Figure 8.5.2 The two histograms show the KE probability distributions for atoms at two different temperature. These were generated using a slightly modified version of NetLogo model 8.4.2. The distribution on the left is for a lower temperature. The large majority of the probability mass is for KE < 1.5 (blue). A small amount is for 1.5 < KE < 3 (green) and almost none for KE > 3 (red). In contrast, in the distribution at higher temperature, atoms are much more likely to have high KE.

Equilibrium Vacancy Concentration

NetLogo model 8.5.1 below uses the algorithm described above. Explore what happens as temperature and bond-energy are varied and use it to answer the question we raised at the beginning of the chapter Section 8.3.

Exercise 8.5.1: Vacancy Concentration
Not Currently Assigned

  1. At a given bond-energy, how does the vacancy concentration change as temperature is increased/decreased? Why?

  2. At a given temperature, how does the vacancy concentration change as bond-energy is increased/decreased? Why?

  3. Try to explain in your own words why the equilibrium vacancy concentration is zero at low temperatures but non-zero at higher temperatures, even though vacancies increase the PE of the system.