Julia中Ising模型数据无法写入.dat文件的问题求助

Julia新手编写Ising模型时，其余功能正常，但存储磁化强度和能量的数组无法写入.dat文件，在文件写入行抛出UndefRefError错误。

相关写入代码片段：

file = open("lattice_eqm.dat", "w") for i in 1:length(TIME_ARR)     write(file, string(TIME_ARR[i]), "\t", string(MAGNETIZATION_PER_SPIN_ARR[i]), "\t", string(ENERGY_PER_SPIN_ARR[i]), "\n") end close(file)

完整代码：

function Ising(J_ising, T, L, EQUI_STEPS, SIM_STEPS)
    setprecision(50)

# rand(Bool, n, n) creates a matrix with random 1 and 0
    LATTICE_SPINS = 2*rand(Bool, L, L).-1 #hence this will make a matrix with random 1 and -1
    E, M = 0, 0

# now finding the intial magnetization and the intial energy
    E = -J_ising * sum((LATTICE_SPINS[i, j] * (LATTICE_SPINS[i, mod1(j+1, L)]) + LATTICE_SPINS[i, mod1(j-1, L)] + LATTICE_SPINS[mod1(i+1, L), j] + LATTICE_SPINS[mod1(i-1, L), j]) for i in 1:L for j in 1:L)
    M = sum(LATTICE_SPINS[i, j] for i in 1:L for j in 1:L)

    E_per_spin = E/(2*(L*L))
    M_per_spin = M/(L*L)
    println("The initial energy per spin is $E_per_spin and the inital magnetization per spin is $M_per_spin")

# define arrays that will contain final values
    ENERGY_PER_SPIN_ARR = Array{BigFloat}(undef, 1)
    MAGNETIZATION_PER_SPIN_ARR = Array{BigFloat}(undef, 1)
    TIME_ARR = Array{Int32}(undef, 1)
    

# now to reach eqm
    for t in 1:EQUI_STEPS

        i, j = rand(1:L), rand(1:L)

        delE = J_ising * (((LATTICE_SPINS[i, j] * (LATTICE_SPINS[i, mod1(j+1, L)]) + LATTICE_SPINS[i, mod1(j-1, L)] + LATTICE_SPINS[mod1(i+1, L), j] + LATTICE_SPINS[mod1(i-1, L), j])) - (((-LATTICE_SPINS[i, j] * (LATTICE_SPINS[i, mod1(j+1, L)]) + LATTICE_SPINS[i, mod1(j-1, L)] + LATTICE_SPINS[mod1(i+1, L), j] + LATTICE_SPINS[mod1(i-1, L), j]))))

        if delE < 0
            E += delE
            M += 2*LATTICE_SPINS[i, j]
        
        else
            prob = exp(-delE/T)
            random = rand()
            
            if random < prob
                LATTICE_SPINS[i, j] = -LATTICE_SPINS[i, j]
                E += delE
                M += 2*LATTICE_SPINS[i, j]
            end
        end
        push!(TIME_ARR, t)
        push!(MAGNETIZATION_PER_SPIN_ARR, M/(L*L))
        push!(ENERGY_PER_SPIN_ARR, E/(L*L))
    end

    # writing to dat file
    file = open("lattice_eqm.dat", "w")
    for i in 1:length(TIME_ARR)    
        write(file, string(TIME_ARR[i]), "\t", string(MAGNETIZATION_PER_SPIN_ARR[i]), "\t", string(ENERGY_PER_SPIN_ARR[i]), "\n")
    end
    close(file)
end
Ising(big(0.5), big(0.2), 50, 100000, 500000)

UndefRefError的根源是数组初始化方式错误：

• 使用Array{BigFloat}(undef, 1)和Array{Int32}(undef, 1)创建的数组，包含1个未初始化的内存位置（值为未定义的垃圾数据）。
• 后续调用push!时，仅在数组末尾追加有效元素，原数组的第一个未初始化元素仍会保留。
• 循环遍历整个数组时，访问第一个未定义元素就会触发UndefRefError。

1. 初始化空数组：将数组初始化改为创建对应类型的空数组，避免引入未定义元素。
2. 使用do语法安全处理文件：Julia推荐用open(...) do file的语法，自动处理文件关闭，避免手动close可能出现的问题。

修改后的完整代码：

function Ising(J_ising, T, L, EQUI_STEPS, SIM_STEPS)
    setprecision(50)

    # 创建随机自旋晶格：1和-1的L×L矩阵
    LATTICE_SPINS = 2 .* rand(Bool, L, L) .- 1
    E, M = 0, 0

    # 计算初始能量和磁化强度
    E = -J_ising * sum(LATTICE_SPINS[i,j] * (LATTICE_SPINS[i, mod1(j+1, L)] + LATTICE_SPINS[i, mod1(j-1, L)] + LATTICE_SPINS[mod1(i+1, L), j] + LATTICE_SPINS[mod1(i-1, L), j]) for i in 1:L, j in 1:L)
    M = sum(LATTICE_SPINS)

    E_per_spin = E / (2 * L^2)
    M_per_spin = M / L^2
    println("初始单自旋能量：$E_per_spin，初始单自旋磁化强度：$M_per_spin")

    # 初始化空数组存储结果
    ENERGY_PER_SPIN_ARR = BigFloat[]
    MAGNETIZATION_PER_SPIN_ARR = BigFloat[]
    TIME_ARR = Int32[]

    # 弛豫到平衡态
    for t in 1:EQUI_STEPS
        i, j = rand(1:L), rand(1:L)
        
        # 简化自旋翻转的能量变化计算
        neighbor_sum = LATTICE_SPINS[i, mod1(j+1, L)] + LATTICE_SPINS[i, mod1(j-1, L)] + LATTICE_SPINS[mod1(i+1, L), j] + LATTICE_SPINS[mod1(i-1, L), j]
        delE = 2 * J_ising * LATTICE_SPINS[i,j] * neighbor_sum

        if delE < 0
            # 能量降低，直接翻转
            LATTICE_SPINS[i,j] = -LATTICE_SPINS[i,j]
            E += delE
            M += 2 * LATTICE_SPINS[i,j]
        else
            # 按Metropolis准则决定是否翻转
            prob = exp(-delE / T)
            if rand() < prob
                LATTICE_SPINS[i,j] = -LATTICE_SPINS[i,j]
                E += delE
                M += 2 * LATTICE_SPINS[i,j]
            end
        end

        push!(TIME_ARR, Int32(t))
        push!(MAGNETIZATION_PER_SPIN_ARR, M / L^2)
        push!(ENERGY_PER_SPIN_ARR, E / L^2)
    end

    # 写入数据到文件（使用do语法自动关闭文件）
    open("lattice_eqm.dat", "w") do file
        for i in eachindex(TIME_ARR)
            write(file, string(TIME_ARR[i]), "\t", string(MAGNETIZATION_PER_SPIN_ARR[i]), "\t", string(ENERGY_PER_SPIN_ARR[i]), "\n")
        end
    end
end

Ising(big(0.5), big(0.2), 50, 100000, 500000)

额外优化：简化了delE的计算逻辑（原表达式可等价简化），同时用eachindex替代1:length(...)，更符合Julia的编程习惯。

内容的提问来源于stack exchange，提问作者ludvary