3
$\begingroup$

I am simulating Ising model for NiO. I have simulated for 2d,3d,triangular lattices, and have tried to do the same with NiO model.

There are papers which say that the ground state energy is around -36 mev per particle. Where as I am getting -21 mev per particle (should I not consider Oxygen atoms?, then it will be -42mev). Is there something wrong I am doing. Following is the code.

def initial_state_Nio(N): 
    state = np.random.choice([-1, 1], (N, N))
    state[::2, ::2] = 0
    state[1::2, 1::2] = 0
    return state
def diag_nbrs(i,j,N): 
    return [((i+1)%N,(j+1)%N),((i+1)%N,(j-1)%N),((i-1)%N,(j+1)%N),((i-1)%N,(j-1)%N)]

def lat_nbrs(i,j,N): 
    return [(i,(j+2)%N),(i,(j-2)%N),((i+2)%N, j),((i-2)%N, j)]

def Energy_Nio(state, J1, J2, H):
    J1   = 2.3*10**(-3) #diagonal
    J2   = -21*10**(-3) #lateral coupling
    E = 0
    N = state.shape[0]
    for x in range(N):
        for y in range(N):
            if (x-y)%2:
                nbrs = diag_nbrs(x,y,N)
                for nbr in nbrs:
                    E += -state[x,y]*J1*state[nbr[0],nbr[1]]
                nbrs = lat_nbrs(x,y,N)
                for nbr in nbrs:
                    E += -state[x,y]*J2*state[nbr[0],nbr[1]]
    E/=2
    E -= H*state.sum()
    return E

def calcMag(state):
    return np.sum(state)

def step_update_Nio(state, beta, J1, J2,H,energy,mag,N):
    J1   = 2.3*10**(-3) #diagonal
    J2   = -21*10**(-3) #lateral coupling
    for i in range(N**2): #1 step per state on average
        dE = 0
        x = random.randint(0,N-1)
        y = random.randint(0,N-1)
        if (x-y)%2:
            nbrs = diag_nbrs(x,y,N)
            for nbr in nbrs:
                dE += 2*state[x,y]*J1*state[nbr[0],nbr[1]]
            nbrs = lat_nbrs(x,y,N)
            for nbr in nbrs:
                dE += 2*state[x,y]*J2*state[nbr[0],nbr[1]]
            dE += 2*H*state[x,y]
            if (dE <= 0):
                if state[x,y] == 1:
                    mag-=2
                else:
                    mag+=2
                energy += dE
                state[x, y] *= -1
            else:
                r = random.uniform(0,1)
                tau = np.exp(-dE*beta)
                if (r < tau) :
                    if state[x,y] == 1:
                        mag-=2
                    else:
                        mag+=2
                    energy += dE
                    state[x, y] *= -1
    return state,energy,mag

def run_Nio(state, steps, N, beta, J1, J2,H):
    J1   = 2.3*10**(-3) 
    J2   = -21*10**(-3)
    E = np.zeros(steps)
    M = np.zeros(steps)
    energy = Energy_Nio(state, J1, J2,H)
    mag = calcMag(state)
    for i in range(steps):
        state,energy,mag = step_update_Nio(state, beta, J1, J2,H,energy,mag,N)
        E[i] = energy
        M[i]= mag
    plt.plot(E)
    plt.show()
    plt.plot(M)
    plt.show()
    return state,E,M

J1 and J2 are shown here enter image description here

Also, if you can comment on magnetisation, capacity heat, magnetic susceptibility, that will be helpful too

$\endgroup$
1
$\begingroup$

enter image description here

As, $J2 = -21 meV$ is more dominant than $J1 = 2.3meV$, the system is antiferromagnetic in nature. The expected ground state energy per moelcule(NiO) is $-42meV$. When the state reaches equilibirum, two of the nearest neighbors are alike and two are opposite, cancelling the energy contirbutions of each other. The energy contribution is from the second nearest negihbours giving $-4*21 = -84 meV$, but as each bond(energy) is counted twice, it should be halved. The above calculations are for Ni atom which has magnetic moment. As oxygen does not couple with other spins, the average energy per site(atom) is $\frac{-42}{2} = -21 meV$.

I have run through the code. It is correct.

| cite | improve this answer | |
$\endgroup$

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.