Field based annealing#

This example shows how to set up a field-based annealing protocol. Below we anneal square spin ice using a rotating field, whose strength is gradually decreased. Disorder introduces small variations in the coercive fields, which actually helps the annealing process by creating nucleation points in the lattice.

from flatspin.model import SquareSpinIceClosed

model = SquareSpinIceClosed(size=(25,25), disorder=0.05, use_opencl=True)
model.plot_vertex_mag();
../_images/036910df4ad5fe63f5d070477913e4d08e8d019ad5028780840b195359e09e75.svg

Rotating field protocol#

We use the Rotate encoder to set up the external rotating field. The timesteps parameter denotes the resolution of one full rotation, while H0 and H sets the minimum and maximum field strength, respectively. The length of the input array defines the number of rotations (20), while the values define the field strength for each rotation, where a value of 1.0 is mapped to H, and a value of 0.0 maps to H0.

from flatspin.encoder import Rotate

timesteps = 64
enc = Rotate(H=0.09, H0=0.06, timesteps=timesteps)
input = np.linspace(1, 0, 20)

h_ext = enc(input)
H = norm(h_ext, axis=1)

plt.plot(norm(h_ext, axis=1), label="norm(h_ext)")
plt.plot(h_ext[:,0], label="h_ext[0]")
plt.plot(h_ext[:,1], label="h_ext[1]")
plt.xlabel("time step")
plt.ylabel("h_ext [T]")
plt.legend(loc='upper left', bbox_to_anchor=(1.0, 1.0));
../_images/a813778883711b9a13911b869bbf13196568cd923c246c0735f39e3fc84db957.svg

Run the field protocol#

Below we iterate over each h_ext value in the field protocol, and update the model accordingly. We record the number of spin flips (steps) and dipolar energy (E_dip) per field value. At the end of each rotation, we also take a snapshot of the spin array.

# Start in polarized state
model.polarize()

# Record spins, number of spin flips and dipolar energy over time
spins = []
flips = []
E_dip = []
for i, h in enumerate(h_ext):
    model.set_h_ext(h)
    s = model.relax()
    if (i+1) % timesteps == 0:
        # Record spin state at the end of each rotation
        spins.append(model.spin.copy())
    flips.append(s)
    E_dip.append(model.total_dipolar_energy())

print(f"Completed {sum(flips)} steps")
Completed 22081 steps

Spin flips over time#

Here we plot the total number of spin flips per field strength, i.e., per field rotation. As can be seen, the strongest fields saturates the array by flipping every spin twice. As the field strength decreases, so does the number of spin flips. Eventually the field becomes too weak to flip any spins.

H = norm(h_ext, axis=-1).round(10)
df = pd.DataFrame({'H': H, 'flips': flips})
df.groupby('H', sort=False).sum().plot(legend=False)
plt.gca().invert_xaxis()
plt.ylabel("Spin flips")
plt.xlabel("Field strength [T]");
../_images/10c847ff21d7dcbcc64adbe21a62d5029f495a27b8cd03aa208ef3de9444f3de.svg

Dipolar energy#

The total dipolar energy is the sum of the dipolar fields acting anti-parallel to the spin magnetization. It is a measure of the total frustration in the system, and a good measure of how well the system has annealed.

plt.plot(E_dip)
[<matplotlib.lines.Line2D at 0x7f24abbd2e30>]
../_images/35671c9898b5c723b536057bc3154c2704c5fec6e0743a1d900074af38195541.svg

Animation of the annealing process#

Finally we animate the annealing process by plotting vertex magnetization at the end of each rotation. In the animation below we see the emergence of antiferromagnetic domains (white regions), which correspond to low energy states. The domains are separated by domain walls with a net positive magnetic moment.

fig, ax = plt.subplots()

def animate_spin(i):
    H = norm(h_ext[(i+1) * timesteps - 1])
    ax.set_title(f"H={H:.3f} [T]")
    model.set_spin(spins[i])
    model.plot_vertex_mag(ax=ax, replace=True)
        
anim = FuncAnimation(fig, animate_spin, frames=len(spins), interval=200, blit=False)
plt.close() # Only show the animation
HTML(anim.to_jshtml())