%Demo_IsingHysteresis - Show the hysteresis
% of magnetization in simulated 3D IsingModel
clear all; help Demo_IsingHysteresis;
tic
% Lattice size
lsize = [20 20 20];
% Boltzmann constant; easiest to set to 1
kb = 1;
%Particle spin
spin = 1/2;
% Strength of interaction
J = 1;
% Temperature
fprintf('Simulating 20x20x20 spin 1/2 particles with J = 1.\n')
T = input('Input temperature (Suggested: 0.5): ');
% Field goes from + to - and back
Bz = cos(linspace(0,2*pi,50));
% Random Lattice
nspins = 2*spin + 1;
lattice = -spin + floor(rand(lsize)*nspins);
beta=1./(kb*T);
fprintf('Simulating...\n')
for tstep = 1:length(Bz)
% Status Update
if (mod(tstep,10) == 0)
fprintf(' [%g s] %g%% done\n',toc,100*tstep/length(Bz));
end
M_ct = 0; % Cumulative total for averaging
maxsteps = 1;
for istep=1:maxsteps
% Take a Monte Carlo step
lattice = MCStep(lattice, beta, J, Bz(tstep), spin);
% Update the magnetization and running totals
M = avgMagnetization(lattice);
M_ct = M_ct + M;
end
% Record the energy, magnetization, and fluctuations in magnetization
M_record(tstep) = M_ct/istep;
end
