Question

Я знаю, что это довольно распространенный вопрос, но поиск других ответов на этот вопрос не помог мне, так как я не совсем понимаю проблему.Я получаю саму ошибку, но эта проблема возникает только в одном из выпусков моего кода, где эта строка одинакова в обоих случаях.

Вот первая версия:

import numpy as np
import numpy.random as rnd
import matplotlib.pyplot as plt

#S[i, j] = 1 (spin up)
#S[i, j] = -1 (spin down)

def InitSpins(S, N):
    for i in range(N):
        for j in range(N):
            S[i, j] = 2*rnd.randint(2)-1
    return S

def GenerateMove(S, N):
    i = rnd.randint(N)
    j = rnd.randint(N)
    return i, j, S[i][j], - S[i][j]

def ComputeEnergy(S, N):
    E = 0.0
    for i in range(N):
        for j in range(N):
            E += -J*S[i,j]*S[(i+1)%N,j] - J*S[i,j]*S[i,(j+1)%N] - H*S[i,j]
    return E

def ComputeDiffEnergy(S, i, j, old, new, N):
    Eold = ComputeEnergy(S, N)
    S[i,j] = new
    Enew = ComputeEnergy(S, N)
    S[i,j] = old
    return Enew - Eold

def MCStep(S, N, E, nacc, T):
    i,j,old,new = GenerateMove(S, N)
    DE = ComputeDiffEnergy(S, i, j, old, new, N)
    if np.any(DE <= 0) or np.any(rnd.random() < np.exp(-DE/T)):
        AcceptMove(S, i, j, old, new, N)
        E += DE
        nacc += 1
#    elif 
#        AcceptMove(S, i, j, old, new, N)
#        E += DE
#        nacc += 1
    else:
        RejectMove(S, i, j, old, new, N)
    return E, nacc

def AcceptMove(S, i, j, old, new, N):
    S[i, j] = new

def RejectMove(S, i, j, old, new, N):
    S[i, j] = old

def ComputeX(S, N):
   X = 0.0
   for i in range(N):
       for j in range(N):
           X += S[i,j]
   X += 1/(N**2)*np.sum(X)
   return X

N = 10 #dimension of lattice
NIter = 10000 #iterations for production run
NEquil = NIter//10 #iterations in actual calculation
NT = 100 #number of time steps
T = 2.4
H = 0.0 #set outside magnetisation
J = 1.0 #set internal magnetisation
S = np.empty([N, N]) #set initlal spin array

S = InitSpins(S, N)
print(S)
print('energy',ComputeEnergy(S, N))

# Equilibration:
nacc = 0
E = ComputeEnergy(S, N)
for i in range(NEquil):
    E, nacc = MCStep(S, N, E, nacc, T)

# Production run
nacc = 0
sum_E = 0.0
sum_E2 = 0.0
E = ComputeEnergy(S, N)
for i in range(NIter):
    E, nacc = MCStep(S, N, E, nacc, T)
    sum_E += E
    sum_E2 += E**2

X = ComputeX(S, N)

def SumX(X, N):
    sum_X = 0.0
    sum_X2 = 0.0
    for i in range(NIter):
        sum_X += X
        sum_X2 += X**2
    return sum_X, sum_X2

sum_X, sum_X2 = SumX(X, N)

# Calculate averages
av_E = sum_E/float(NIter)
av_E2 = sum_E2/float(NIter)
av_X = sum_X/float(NIter)
av_X2 = sum_X2/float(NIter)
CV = 1/(1*(T**2))*(av_E2-av_E**2)
chi = 1/1*(T**2)*(av_X2-av_X**2)
#M = 

# print results
print("acceptance ratio",nacc/float(NIter))
print("average energy",av_E)
print("heat capacity",CV)
print("magnetic susceptibility", chi)
#print("magnetisation, )

Это, хотя и неполное, дает мне результаты.Затем я попытался заставить Python работать в диапазоне T (температура), чтобы получить график изменения теплоемкости.Это нарушает мою GenerateMove функцию с ошибкой, указанной в заголовке для return i,j, S[i][j], -S[i][j].

Обновленный код, который не работает:

import numpy as np
import numpy.random as rnd
import matplotlib.pyplot as plt

#S[i, j] = 1 (spin up)
#S[i, j] = -1 (spin down)

def InitSpins(S, N):
    for i in range(N):
        for j in range(N):
            S[i, j] = 2*rnd.randint(2)-1
    return S

def GenerateMove(S, N):
    i = rnd.randint(N)
    j = rnd.randint(N)
    return i, j, S[i][j], - S[i][j]


def ComputeEnergy(S, N):
    E = 0.0
    for i in range(N):
        for j in range(N):
            E += -J*S[i,j]*S[(i+1)%N,j] - J*S[i,j]*S[i,(j+1)%N] - H*S[i,j]
    return E

def ComputeDiffEnergy(S, i, j, old, new, N):
    Eold = ComputeEnergy(S, N)
    S[i,j] = new
    Enew = ComputeEnergy(S, N)
    S[i,j] = old
    return Enew - Eold

def MCStep(S, N, E, nacc, T):
    i,j,old,new = GenerateMove(S, N)
    DE = ComputeDiffEnergy(S, i, j, old, new, N)
    if np.any(DE <= 0) or np.any(rnd.random() < np.exp(-DE/T)):
        AcceptMove(S, i, j, old, new, N)
        E += DE
        nacc += 1
#    elif 
#        AcceptMove(S, i, j, old, new, N)
#        E += DE
#        nacc += 1
    else:
        RejectMove(S, i, j, old, new, N)
    return E, nacc

def AcceptMove(S, i, j, old, new, N):
    S[i, j] = new

def RejectMove(S, i, j, old, new, N):
    S[i, j] = old

def ComputeX(S, N):
   X = 0.0
   for i in range(N):
       for j in range(N):
           X += S[i,j]
   X += 1/(N**2)*np.sum(X)
   return X

N = 10 #dimension of lattice
NIter = 10000 #iterations for production run
NEquil = NIter//10 #iterations in actual calculation
NT = 100 #number of time steps
T = np.linspace(0.1,5.0,NT) #set temperaure range
C = np.zeros(NT) #inital heat capacity
H = 0.0 #set outside magnetisation
J = 1.0 #set internal magnetisation
S = np.empty([N, N]) #set initlal spin array

S = InitSpins(S, N)
print(S)
print('energy',ComputeEnergy(S, N))

## Equilibration:
#nacc = 0
#E = ComputeEnergy(S, N)
#for i in range(NEquil):
#    E, nacc = MCStep(S, N, E, nacc, T)
#
## Production run
#nacc = 0
#sum_E = 0.0
#sum_E2 = 0.0
#E = ComputeEnergy(S, N)
#for i in range(NIter):
#    E, nacc = MCStep(S, N, E, nacc, T)
#    sum_E += E
#    sum_E2 += E**2

X = ComputeX(S, N)

def SumX(X, N):
    sum_X = 0.0
    sum_X2 = 0.0
    for i in range(NIter):
        sum_X += X
        sum_X2 += X**2
    return sum_X, sum_X2

sum_X, sum_X2 = SumX(X, N)

for t in range(NT):
    nacc = 0
    sum_E  = 0.0
    sum_E2 = 0.0
    S = InitSpins(S, N)
    E = ComputeEnergy(S, N)
    for i in range(NIter):
        S = MCStep(S, N, E, nacc, T[t])
        for i in range(NEquil):
            S = MCStep(S, N, E, nacc, T[t])
            E, nacc = MCStep(S, N, E, nacc, T)
            sum_E += E
            sum_E2 += E**2
        av_E = sum_E/float(NIter)
        av_E2 = sum_E2/float(NIter)
        C[t] = (sum_E/(float(NIter*N*N)) - sum_E*sum_E/(NIter*NIter*N*N))/T[t]**2

# Calculate averages
av_X = sum_X/float(NIter)
av_X2 = sum_X2/float(NIter)
CV = 1/(1*(T**2))*(av_E2-av_E**2)
chi = 1/1*(T**2)*(av_X2-av_X**2)
#M = 

# print results
print("acceptance ratio",nacc/float(NIter))
print("average energy",av_E)
print("heat capacity",CV)
print("magnetic susceptibility", chi)
#print("magnetisation, )

#plotting resultts
fig = plt.figure()
ax = fig.add_subplot(111)
ax.plot(T, C)
plt.show()

РЕДАКТИРОВАТЬ: вот полныйtraceback

Traceback (most recent call last):

  File "<ipython-input-3-52258536ecbb>", line 1, in <module>
    runfile('F:/adv. numerical project graph.py', wdir='F:')

  File "C:\Users\User\Anaconda3\lib\site-packages\spyder_kernels\customize\spydercustomize.py", line 704, in runfile
    execfile(filename, namespace)

  File "C:\Users\User\Anaconda3\lib\site-packages\spyder_kernels\customize\spydercustomize.py", line 108, in execfile
    exec(compile(f.read(), filename, 'exec'), namespace)

  File "F:/adv. numerical project graph.py", line 121, in <module>
    S = MCStep(S, N, E, nacc, T[t])

  File "F:/adv. numerical project graph.py", line 42, in MCStep
    i,j,old,new = GenerateMove(S, N)

  File "F:/adv. numerical project graph.py", line 24, in GenerateMove
    return i, j, S[i][j], - S[i][j]

IndexError: tuple index out of range

Sheldore · Answer 1 · 30 декабря 2018

Вы излишне усложнили свой код.Расширение вашего первого кода для фиксированной температуры в диапазоне температур довольно просто.Все, что вам нужно было сделать, - это объединить код, выполняющий инициализацию решетки, уравновешивание MC и производство MC, в цикле for, работающем в заданном диапазоне температур.Я создал списки для хранения теплоемкости, ци и средней энергии при каждой температуре.

Ниже описано, как это делается.Я выбираю 30 значений температуры вместо 100, чтобы быстро дать вам рабочее решение.Я только показываю соответствующий код.Остальная часть кода с определениями функций остается прежней.

Мне не нравится поведение прогнозируемой теплоемкости вашего MC.Предполагается, что увеличивается с температурой.Вам необходимо перепроверить свою реализацию.У тебя все есть сейчас.

NT = 30 # number of temperature steps
T_list = np.linspace(0.1, 5.0, NT) # set temperaure range
av_E_list, CV_list, chi_list = [], [], []

for T in T_list: # <--- Loop over different temperatures
    print ("Performing MC for T = %.2f" %T)
    S = np.empty([N, N]) #set initlal spin array
    S = InitSpins(S, N)

    # Equilibration:
    nacc = 0
    E = ComputeEnergy(S, N)
    for i in range(NEquil):
        E, nacc = MCStep(S, N, E, nacc, T)

    # Production run
    nacc = 0
    sum_E = 0.0
    sum_E2 = 0.0
    E = ComputeEnergy(S, N)
    for i in range(NIter):
        E, nacc = MCStep(S, N, E, nacc, T)
        sum_E += E
        sum_E2 += E**2

    X = ComputeX(S, N)
    sum_X, sum_X2 = SumX(X, N)

    # Calculate averages
    av_E = sum_E/float(NIter)
    av_E2 = sum_E2/float(NIter)
    av_X = sum_X/float(NIter)
    av_X2 = sum_X2/float(NIter)
    CV = 1/(1*(T**2))*(av_E2-av_E**2)
    chi = 1/1*(T**2)*(av_X2-av_X**2)
    CV_list.append(CV)
    chi_list.append(chi)
    av_E_list.append(av_E)

# Plotting 
plt.plot(T_list, CV_list, '-b*')

Python: TypeError: объект 'int' не является подписным Монте-Карло

Пожалуйста, войдите или зарегистрируйтесь чтобы ответить на этот вопрос.

1 Ответ

Пожалуйста, войдите или зарегистрируйтесь что бы добавить комментарий.

Python: TypeError: объект 'int' не является подписным Монте-Карло

Пожалуйста, войдите или зарегистрируйтесь чтобы ответить на этот вопрос.

1 Ответ

Пожалуйста, войдите или зарегистрируйтесь что бы добавить комментарий.

Похожие темы