1. The problem statement, all variables and given/known data Lo,Im stuck on how to retrieve the specific heat capacity from an MC simulation, with the metropolis algorithm. I want my graph to look something like this: https://i.stack.imgur.com/NXeXs.png 2. Relevant equations C_v = ((<E^2>-<E>^2)/T^2 3. The attempt at a solution My code is similar to this guy: but without magnetization To retrieve C_v i wrote the following code to 4:19 from the video above. Cv = (E2/16384 - E.*E/16384^2)./(T(i).^2); plot(T,Cv); Where E2 is given by the code: function [E2] = ising_energy(u,J) [L1,L2] = size(u); E2=0; for i=1:L1 for j=1:L2 f1i=mod(i,L1)+1; f2j=mod(j,L2)+1; forward_neighbors=u(f1i,j)+u(i,f2j); E2=E2-(J*u(i,j)*forward_neighbors)*(J*u(i,j)*forward_neighbors); end;clear j end;clear i However my graph looks more like the plot of <E> vs T.