我试着通过改变很多变量来解决这个问题,但都没有奏效.以下是我的代码:

import numpy as np
import matplotlib.pyplot as plt
import math as mt
from numba import njit


N=1000
J=1
h=0.5

plt.rcParams['figure.dpi']=100

plt.xlabel('T', fontsize=14)

ns=100000000

s=np.empty([N,0])
for i in range(N):
    s[i]=np.random.choice([-1,1]) 
    


E=0
M=0
for i in range(-1,N-1):
    E=E-J*(s[i]*s[i+1]+s[i]*s[i-1]) 
    M=M+s[i]
    

Energy=np.empty([0,0])

C=np.empty([0,0])

Magne=np.empty([0,0])

chi=np.empty([0,0])
@njit
def average(k , E, M, ns, x, Energy, C, Magne, chi, s):
    Ev=0.
    Ev2=0.
    Mv=0.
    Mv2=0.
    
    for z in range(1,ns+1):
        
        i=np.random.randint(-1, high=N-1)
        
        dE=2*J*(s[i]*s[i+1]+s[i]*s[i-1])
        dM=2*h*s[i]
        
        pace=1./(1+mt.exp(k*(dE+dM)))
        #pace=min(1,mt.exp(-k*dE))
        
        
        if np.random.random()<pace:
            s[i]=-s[i]
            E=E+dE
            M=M-dM/h

        if z>x:
            Ev=Ev+(E-h*M)
            Ev2=Ev2+(E-h*M)**2
            Mv=Mv +M
            Mv2=Mv2 +M**2
        
        
    Ev=Ev/(ns-x)
    Ev2=Ev2/(ns-x)
    varE=(Ev2-Ev**2)
    Energy=np.append(Energy,Ev)
    
    C=np.append(C,varE*(k**2))
    
    Mv=Mv/(ns-x)
    Mv2=Mv2/(ns-x)
    varM=(Mv2-Mv**2)
    Magne=np.append(Magne,Mv)
    
    chi=np.append(chi,varM*k)
    
    return E, M, Energy, C, s, Magne, chi


b=np.arange(0.3,1.,0.01)

x=0.75*ns

for k in b:
    
    E, M, Energy, C, s, Magne, chi =average(k,E,M,ns, x,Energy,C, Magne, chi, s)

         
plt.plot(b,Energy, '.', color='r',linestyle='--')
plt.plot(b,C, '.', color='b',linestyle='--')
plt.title('Energy and heat capacity')
plt.legend(['E','C'])
plt.show()

plt.rcParams['figure.dpi']=100

plt.xlabel(r'$\beta$', fontsize=14)


plt.plot(b,Magne, '.',linestyle='--', color='r')
plt.plot(b,chi, '.', color='b',linestyle='--')
plt.title('Magnetization and susceptibility')
plt.legend(['M',r'$\chi$'])
plt.show()

我一直收到错误:

Cannot unify float64 and array(float64, 1d, C) for 'Mv2.3', defined at c:\users\usuario\documents\python scripts\ex13hw7.py (65)

File "ex13hw7.py", line 65: def average(k , E, M, ns, x, Energy, C, Magne, chi, s): <source elided> if z>x: Ev=Ev+(E-h*M) ^

期间:在c:\USERS\usuario\Documents\PYTHON SCRIPTS\ex13hw7.py中键入赋值(65)

有人知道问题出在哪里吗?

我try 更改变量的名称并在函数中包含/删除参数.

推荐答案

以下是使用Numba编译的代码的修复版本:

import math as mt

import matplotlib.pyplot as plt
import numpy as np
from numba import njit

N = 1000
J = 1
h = 0.5

plt.rcParams["figure.dpi"] = 100

plt.xlabel("T", fontsize=14)

ns = 100000000

s = np.empty(N)        # <-- don't use np.empty([N, 0])
for i in range(N):
    s[i] = np.random.choice([-1, 1])

E = 0
M = 0
for i in range(-1, N - 1):
    E = E - J * (s[i] * s[i + 1] + s[i] * s[i - 1])
    M = M + s[i]

Energy = np.array([])  # <-- don't use np.empty([0, 0])
C = np.array([])
Magne = np.array([])
chi = np.array([])


@njit
def average(k, E, M, ns, x, Energy, C, Magne, chi, s):
    Ev = 0.0
    Ev2 = 0.0
    Mv = 0.0
    Mv2 = 0.0

    for z in range(1, ns + 1):
        i = np.random.randint(-1, high=N - 1)

        dE = 2 * J * (s[i] * s[i + 1] + s[i] * s[i - 1])
        dM = 2 * h * s[i]

        pace = 1.0 / (1 + mt.e ** (k * (dE + dM)))  # <-- use math.e ** x

        if np.random.random() < pace:
            s[i] = -s[i]
            E = E + dE
            M = M - dM / h

        if z > x:
            Ev = Ev + (E - h * M)
            Ev2 = Ev2 + (E - h * M) ** 2
            Mv = Mv + M
            Mv2 = Mv2 + M**2

    Ev = Ev / (ns - x)
    Ev2 = Ev2 / (ns - x)
    varE = Ev2 - Ev**2
    Energy = np.append(Energy, Ev)

    C = np.append(C, varE * (k**2))

    Mv = Mv / (ns - x)
    Mv2 = Mv2 / (ns - x)
    varM = Mv2 - Mv**2
    Magne = np.append(Magne, Mv)

    chi = np.append(chi, varM * k)

    return E, M, Energy, C, s, Magne, chi


b = np.arange(0.3, 1.0, 0.1)
x = 0.75 * ns

for k in b:
    E, M, Energy, C, s, Magne, chi = average(k, E, M, ns, x, Energy, C, Magne, chi, s)


print("Finished!")

plt.plot(b, Energy, ".", color="r", linestyle="--")
plt.plot(b, C, ".", color="b", linestyle="--")
plt.title("Energy and heat capacity")
plt.legend(["E", "C"])
plt.show()

plt.rcParams["figure.dpi"] = 100

plt.xlabel(r"$\beta$", fontsize=14)
plt.plot(b, Magne, ".", linestyle="--", color="r")
plt.plot(b, chi, ".", color="b", linestyle="--")
plt.title("Magnetization and susceptibility")
plt.legend(["M", r"$\chi$"])
plt.show()

显示了这两个图表:

enter image description here

enter image description here

Python相关问答推荐

Python-Polars:如何用两个值的平均值填充NA?

如何对行使用分段/部分.diff()或.pct_change()?

从Python调用GMP C函数时的分段错误和内存泄漏

在Python中使用一行try

Django注释:将时差转换为小数或小数

如何使用PyTest根据self 模拟具有副作用的属性

Python:记录而不是在文件中写入询问在多文件项目中记录的最佳实践

请从Python访问kivy子部件的功能需要帮助

使用Ubuntu、Python和Weasyprint的Docker文件-venv的问题

Python plt.text中重叠,包adjust_text不起作用,如何修复?

优化在numpy数组中非零值周围创建缓冲区的函数的性能

三个给定的坐标可以是矩形的点吗

更改matplotlib彩色条的字体并勾选标签?

Excel图表-使用openpyxl更改水平轴与Y轴相交的位置(Python)

使用setuptools pyproject.toml和自定义目录树构建PyPi包

如果条件不满足,我如何获得掩码的第一个索引并获得None?

实现神经网络代码时的TypeError

将pandas导出到CSV数据,但在此之前,将日期按最小到最大排序

如何在达到end_time时自动将状态字段从1更改为0

在matplotlib中使用不同大小的标记顶部添加批注