Ising Model FAQ

In trying to set up a function to determine the heat capacity we determined that we needed a value for T squared, but were unsure how to calculate T because we've been using KT/J instead.

1. $J$ is merely your unit of energy…it is no different than using kcal/mol, eV, etc. We are just setting this to one.
2. Similarly, you can report the value of the heat capacity, $C$, in units of $N k_{B}$ (where $N$ is the total number of spins and $k_{B}$ is the Boltzmann constant.
3. By combining these two considerations, you now have no "free" constants.

---

We were a little confused about the difference between Ml and Ei.

Denoting each site by $i$ (i.e., each of the $i$ represents a particular value of $col, row$ )

(1)

\begin{align} M = \sum_{i} \sigma_{i} \end{align}

(2)

\begin{align} E = J \sum_{<i,j>_{n.n.}} \sigma_{i} \sigma_{j} \end{align}

(where $<i,j>_{n.n.}$ denotes that site $i$ and $j$ are nearest neighbors (north/south/east/west))

---

Is there a command for summation in fortran?

Consider the following:

(3)

\begin{align} M = \sum_{i} \sigma_{i} = \sigma_{1} + \sigma_{2} + \sigma_{3} + ... + \sigma_{N} \end{align}

Since the summation is nothing other than an addition, we could perform this by accumulating each term in the sum using a do loop, e.g.,

M = 0.0
do row=1,16
  do col = 1,16
    M = M+spin(col,row)
  end do
end do

NOTE: When computers run calculations, they read from right to left. This means for the state of M=M +spin(col,row), it is taking the value of M, adding the value of the spin, and storing that value into M. For example, if M were 2 and the spin were -1, the computer would take 2 add -1, and store 1 as a new value for M. This is a useful trick as a counter as well. For example, lets say you wanted to actually keep track of how many data points you've collected for your Ising model. You could define a real number, lets say Q. In the part of code where you collect data (under your statement of "if (mod(nMC,50) .eq. 0) then…." You could type Q = Q + 1. This may not make sense at first, but if you read from right to left, it does. You are taking the previous value for Q, adding one, then making that your new value for Q. This means that each time you collect data, Q will increase by one. So if at the end of the program you wrote "write(*,*) Q", you would see the number of times you collected data. Of course you could figure this out on your own, but this is a very useful trick in computer programing.

---

After performing montecarlo integration for the first time how do I define the spin array after 10000 trials equilibrium spin? Can I say if (mod(nMC, 10000) .eq. 0) spin=spinequilibrium ??? how do I specify over which array to do montecarlo integration?

There are two ways to solve this:

1. For each value of kT, run two sets of "do nMC=1,10000" loops, one after the other
2. As you mention, just have one "do nMC=1,20000" loop, inside of which is a statement "if (nMC .gt. 10000) then … ", within which you collect statistics.

However, it both cases, there is only one spin-array. The purpose of the equilibration period is to make sure that it has a chance to get to the right temperature as the surroundings, analogous to physical systems that you might be familiar with (e.g., taking a frozen steak out of the freezer and letting it equilibrate to refrigerator temperature).

---

The magnetization and energy correlate to the graphs, but the molar heat capacity is giving my group some trouble. There is some resemblance to the graph, but there is not as sharp of a peak.

A common source of error is to calculate the heat capacity, $C$ within the sampling routine. The problem is that while you are indeed computing "instantaneous" values of $E$ and $E^2$, you don't actually know the average values, $<E>$ and $<E^2>$ until after all of the sampling is done. If you use $E$ and $E^2$, you are not actually calculating $C$ (you can prove this to yourself with a little algebra)

A better way is to wait until your final values of $<E>$ and $<E^2>$ are obtained, and then just use these (along with the formula you derived in the last homework set) to obtain $C$. In other words, you don't have to keep a running sum for $C$.

---

The graphs seem a little off, but I cannot for the life of me figure out what is wrong with these equations

The code sent in contained the following snippet:

      Q=Q+J*2*Exp(-Ei/kT)
      M=(M+spin(col,row))
      Et=Et+J*Ei*(exp(-Ei/kT))
      E=(Et/Q)/(16*16)

The problem is the erroneous "Boltzmann-type" average on lines 1 & 2. HOWEVER, the MC sampling procedure has already given you the proper distribution of states in the course of your sampling run. In other words, the MC sampling scheme "automatically" gives you the probability of a state occurring (because it had to occur for you to be able to sample it, right?) Therefore, you need only collect a simple average (mean) of E, M, E^2 over the course of your run. This is described in more detail in one of the readings.