Julia-在短时间内制作数组{Float64,3} [英] Julia- Make Array{Float64,3} in a short time

问题描述

我想从矩阵值函数H（Lx，Ly，Lz）获得Array {Float64,3}，其中Lx，Ly，Lz是参数，H是（Lx×Ly×Lz）×（Lx ×Ly×Lz）矩阵。

I want to get Array{Float64,3} from a matrix-valued function H(Lx,Ly,Lz) ,where Lx,Ly,Lz are parameters and H is a (Lx×Ly×Lz)×(Lx×Ly×Lz) matrix.

示例代码为

using LinearAlgebra
eye(T::Type,n) = Diagonal{T}(I, n) 
eye(n) = eye(Float64,n)

function H(Lx,Ly,Lz) #def of H
N = Lx*Ly*Lz　

　mat_Htb = zeros(Complex{Float64},N,N) 

　for iz = 1:Lz
    　for ix = 1:Lx
        　for iy=1:Ly             
            　for dz in -1:1             
            　jz = iz + dz 
                　for dx in -1:1
                　jx = ix + dx
                    　for dy in -1:1
                     　jy = iy + dy
                        
                            　ii = (iz-1)*Lx*Ly + (ix-1)*Ly + (iy-1) + 1
                            　jj = (jz-1)*Lx*Ly + (jx-1)*Ly + (jy-1) + 1

                                　if 1 <= jx <= Lx && 1 <= jy <= Ly && 1 <= jz <= Lz
                                     
                                    　if dx == +1 && dy == 0 && dz == 0
                                        　mat_Htb[ii,jj] += im 　
                                    　end
                            
                                    　if dx == -1 && dy == 0 && dz == 0 
                                        　mat_Htb[ii,jj] += im/4 　
                                    　end
                        
                                    　if dx == 0 && dy == +1 && dz == 0 
                                        　mat_Htb[ii,jj] += im/2
                                    　end
                            
                                    　if dx == 0 && dy == -1 && dz == 0 
                                        　mat_Htb[ii,jj] += im
                                    　end
                        
                                    　if dx == 0 && dy == 0 && dz == +1 
                                        　mat_Htb[ii,jj] += -im
                                    　end
                                      
                                    　if dx == 0 && dy == 0 && dz == -1 
                                        　mat_Htb[ii,jj] += im*(3/7)
                                    　end
                                    
                                    　if dx == 0 && dy == 0 && dz == 0 
                                        　mat_Htb[ii,jj] += im
                                    　end
                                     
                        　end
                    　end
                　end
            　end
        　end
    　end
　end
                
    
　return mat_Htb
end


Lx = 10 #systemsize-parameters
Ly = 10
Lz = 10


ψ0 = Complex{Float64}[] #def of \psi0 ,(Lx×Ly×Lz)×1 vector

for iz = 1:Lz
    for ix = 1:Lx
        for iy=1:Ly         
            gauss = exp(-((ix-5)^2 + (iy-5)^2 + (iz-5)^2))
            push!(ψ0,gauss)
        end
    end
end

ψ(t) = exp((-im*t).*H(Lx,Ly,Lz))*ψ0 #time-evolution
abs2ψ(t) = abs2.(ψ(t)./norm(ψ(t))) #normalized density

然后，我尝试制作一个数组{Float64,3}

Then, I tried to make an Array{Float64,3} like this.

x = 1:Lx  # our value range
y = 1:Ly
z = 1:Lz
t = 15 #time

ρ(ix,iy,iz) = abs2ψ(t)[(iz-1)*Lx*Ly + (ix-1)*Ly + (iy-1) + 1] 
density = Float64[ρ(ix,iy,iz) for ix in x, iy in y,iz in z]

H（Lx，Ly，Lz），ψ（t），abs2ψ（t），ρ（ix， iy，iz）计算顺利。
但是密度大约需要30分钟。

H(Lx,Ly,Lz),ψ(t),abs2ψ(t),ρ(ix,iy,iz) are calculated smoothly. But the density takes about 30 minutes.

最后，我将对t进行循环计算。
所以我想减少计算时间。
您能告诉我如何解决此问题吗？

Ultimately, I will do loop-calculations for t. So I want to reduce the calculation time. Could you tell me how to solve this problem?

推荐答案

仍然有很多事情可能需要进行了改进，但以下版本应该已经比您的版本快得多。

There are still numerous things that would probably need to be improved, but the following version should already be much faster than yours.

要记住的关键是尝试不要多次重复计算同一件事（特别是如果花一些时间来计算）

The key thing to remember is to try and not recompute the same thing several times (especially if it takes some time to compute it, and you're going to re-use the result a large number of times).

在您的示例中，这适用于：

In your example, this applies to :

• H ，仅取决于 Lx ， Ly 和 Lz ，因此可以一劳永逸地计算


• ψ和abs2ψ，它们取决于 H 和 t ，因此应在每个时间步进行更新-但可以重新用于所有（ix，iy，iz）三元组。

• H, which only depends on Lx, Ly and Lz, and as such can be computed once and for all
• ψ and abs2ψ, which depend on H and t, and should therefore get updated at each time step -- but can be re-used for all (ix, iy, iz) triplets.

using LinearAlgebra

function H(Lx,Ly,Lz)
    N = Lx*Ly*Lz　

    mat_Htb = zeros(Complex{Float64},N,N)

    for iz = 1:Lz
        for ix = 1:Lx
            for iy=1:Ly
                for dz in -1:1
                    jz = iz + dz
                    for dx in -1:1
                        jx = ix + dx
                        for dy in -1:1
                            jy = iy + dy

                            ii = (iz-1)*Lx*Ly + (ix-1)*Ly + (iy-1) + 1
                            jj = (jz-1)*Lx*Ly + (jx-1)*Ly + (jy-1) + 1

                            if 1 <= jx <= Lx && 1 <= jy <= Ly && 1 <= jz <= Lz

                                if dx == +1 && dy == 0 && dz == 0
                                    mat_Htb[ii,jj] += im
                                end

                                if dx == -1 && dy == 0 && dz == 0
                                    mat_Htb[ii,jj] += im/4
                                end

                                if dx == 0 && dy == +1 && dz == 0
                                    mat_Htb[ii,jj] += im/2
                                end

                                if dx == 0 && dy == -1 && dz == 0
                                    mat_Htb[ii,jj] += im
                                end

                                if dx == 0 && dy == 0 && dz == +1
                                    mat_Htb[ii,jj] += -im
                                end

                                if dx == 0 && dy == 0 && dz == -1
                                    mat_Htb[ii,jj] += im*(3/7)
                                end

                                if dx == 0 && dy == 0 && dz == 0
                                    mat_Htb[ii,jj] += im
                                end

                            end
                        end
                    end
                end
            end
        end
    end

    return mat_Htb
end


function run(Lx, Ly, Lz)
    ψ0 = Complex{Float64}[] #def of \psi0 ,(Lx×Ly×Lz)×1 vector

    for iz = 1:Lz
        for ix = 1:Lx
            for iy=1:Ly
                gauss = exp(-((ix-5)^2 + (iy-5)^2 + (iz-5)^2))
                push!(ψ0,gauss)
            end
        end
    end

    x = 1:Lx  # our value range
    y = 1:Ly
    z = 1:Lz
    t = 15 #time

    H_ = H(Lx,Ly,Lz)
    ψ = exp((-im*t).*H_)*ψ0 #time-evolution
    abs2ψ = abs2.(ψ./norm(ψ)) #normalized density

    ρ(ix,iy,iz) = abs2ψ[(iz-1)*Lx*Ly + (ix-1)*Ly + (iy-1) + 1]
    density = Float64[ρ(ix,iy,iz) for ix in x, iy in y,iz in z]
end

Lx = 10 #systemsize-parameters
Ly = 10
Lz = 10
run(Lx, Ly, Lz)

---

对于像这样的问题，这些问题非常特定于您要优化的代码，我倾向于认为在朱莉娅的演讲论坛会更合适。

这篇关于Julia-在短时间内制作数组{Float64,3}的文章就介绍到这了，希望我们推荐的答案对大家有所帮助，也希望大家多多支持IT屋！