I tried to run Julia-1.6.0 dev version in my laptop. When running the LinearMap function, I encountered an error, which I have attached along with this query. I would like any suggestions to improve.
Thanks in advance,
Sai Krishna
1 Like
Please try to avoid posting screenshots and quote your code instead, as this makes it much easier to help you. Also try to make sure, you provide a complete working example, as you didn’t define multiplyHPsi anywhere. This thread may be helpful to you:
That said, does your code work on a previous version of Julia? What version of LinearMaps.jl are you using? (You can find that out with ]st LinearMaps in the REPL.)
2 Likes
Thank you for your reply. I will quote the entire working example now. The LinearMaps.jl version is 2.6.1.
The code is as follows:
function multiplyHPsi!(nPsi::Vector, Psi::Vector) # requires H2 to have been defined
length(nPsi)==length(Psi) || throw(DimensionMismatch())
dimN = length(Psi) # Dimension of the vector space
N = convert(Int64,log2(dimN)) # Number of spins
dimPsi = tuple(fill(2::Int,1,N)…) # dimsPsi =(2,2,2, …, 2)
a = collect(2:N)
a = [a;1]
cyclperm = tuple(a…) #tuple for cyclic permutation cyclperm = (2, 3, … N, 1)
Psi = reshape(Psi,dimPsi)
nPsi = zeros(dimPsi)
for n=1:N
Psi = reshape(Psi,(4,2^(N-2)))
nPsi = reshape(nPsi,(4,2^(N-2)))
nPsi += H2*Psi
Psi = reshape(Psi,dimPsi)
nPsi = reshape(nPsi,dimPsi)
Psi = permutedims(Psi,cyclperm)
nPsi = permutedims(nPsi,cyclperm)
end
nPsi = reshape(nPsi,dimN)
return nPsi
end
Ising model
function buildIsing(h=1.0) # Ising model with transverse magnetic field h (critical h=1 by default)
X = [0 1; 1 0]
Z = [1 0; 0 -1]
I = diagm(0 => ones(2))
XX = kron(X,X)
ZI = kron(Z,I)
IZ = kron(I,Z)
H2 = -(XX + h/2*(ZI+IZ))
return H2
end
XX model
function buildXX(h=0.0) # XX model with magnetix field h (h=0 by default)
X = [0. 1; 1 0]
Y = [0. -1; 1 0] # we avoid using complex numbers
Z = [1 0; 0 -1]
I = diagm(0 => ones(2))
ZI = kron(Z,I)
IZ = kron(I,Z)
H2 = -(kron(X,X)-kron(Y,Y)) + h/2*(ZI+IZ) # sXsX+sYsY + hZ
return H2
end
display(buildIsing())
display(buildXX())
h = 1
H2 = buildIsing(h) # Hamiltonian
H2 = buildXX()
D,U = eigen(H2)
shiftE = D[end]
H2 = H2 - shiftE*diagm(0 => ones(4)) # remember to re-shift the energy later on!
N = 10 # number of spins
Psi = rand(2^N)
nPsi = rand(2^N)
Psi = Psi/norm(Psi)
init_steps = 1
energy = ones(0)
initial_step=1
Nsteps = 0
@time nPsi = multiplyHPsi!(nPsi, Psi) # N = 20 took 4 seconds;
initial_step += Nsteps
Nsteps = 10
final_step = initial_step + Nsteps-1
for n=initial_step:final_step
nPsi = multiplyHPsi!(nPsi,Psi)
newenergy = real(Psi’nPsi)[1] + shiftEN
energy = [energy; newenergy]
Psi = nPsi/norm(nPsi)
print(n-initial_step+1, “:”, Nsteps, " ")
end
figure(“Power_method”,figsize=(14,4))
subplot(121) # Create the 1st axis of a 2x2 arrax of axes
grid(“on”) # Create a grid on the axis
title(“Power method”)
ax = gca()
xlabel(“all iterations”)
ylabel(“Energy”)
plot(energy, marker = “o”)
subplot(122) # Create the 1st axis of a 2x2 arrax of axes
grid(“on”) # Create a grid on the axis
title(“Power method”)
ax = gca()
xlabel(“most recent iterations”)
ylabel(“Energy”)
plot(initial_step:final_step, energy[initial_step:final_step], marker = “o”)
using LinearMaps; using Arpack
N = 14 # N=14 took 4.5 sec; N=16 took 35 sec; N=18 took 122 sec
m = 1 # how many states?
this produces a linear map H that can be called from eigs
H = LinearMap{Float64}(multiplyHPsi!, 2^N; isreal = true, ismutating = true) # use this for real symmetric H
#H = LinearMap(multiplyHPsi!, 2^N, isreal=false, ismutating=true, ishermitian=true) # use this for complex hermitian H
@time u, v = eigs(H; nev=m, which=:SR) #:SR stands for smallest real part
Energy = real(u[1] + shiftE*N) # let us undo the shift in energy
Psi = v
E = Energy
Note that isreal is not a keyword argument of LinearMap, that’s also the complaint of the MethodError. You should specify the eltype of your map as you do with LinearMap{Float64}. Just leave out the isreal will probably resolve your problem.
Thank you for the suggestion but I ran into the following error after removing the isreal. And how do you specify the eltype of the map. Please let me know that.
The following is the error:
syntax: invalid keyword argument syntax “true” around In[66]:4
Stacktrace:
[1] top-level scope at none:1
I am guessing, that you left a semicolon before the positional argument isreal instead of a comma. Replacing that with a comma should fix this error. It’s a bit of a confusing error message though, because 1.5 has a new feature, where variables after a semicolon in a function call are interpreted as x=x, i.e. a keyword argument where the value is that variable in the current scope. Do you want to open an issue on GitHub, because this should probably throw a more helpful error for true and false?
I was assuming that you understand a bit the basic concepts in Julia. When I said, leave out the keyword argument isreal, this implies also its value:
H = LinearMap{Float64}(multiplyHPsi!, 2^N; ismutating = true)
The eltype is what you put in between the curly braces here, LinearMap{ eltype }.
In this case, as it seems you want to diagonalize quantum Hamiltonians, it might be useful to specify that this is a symmetric (or hermitian) operator, i.e.
H = LinearMap{Float64}(multiplyHPsi!, 2^N; ismutating = true, issymmetric = true)
1 Like
But once I removed the isreal argument, I got the same error.
Can you post the line which you entered just before the error was thrown, and the error.
This is the line which I entered before the error :
H = LinearMap{Float64}(multiplyHPsi!; 2^N, true) # use this for real symmetric H
This is the error:
MethodError: no method matching (LinearMaps.FunctionMap{Float64,F1,F2} where F2 where F1)(::typeof(multiplyHPsi!), ::Int64, ::Bool)
Closest candidates are:
(LinearMaps.FunctionMap{T,F1,F2} where F2 where F1)(::Any, ::Int64; kwargs…) where T at /home/sai/.julia/packages/LinearMaps/RjB5x/src/functionmap.jl:20
(LinearMaps.FunctionMap{T,F1,F2} where F2 where F1)(::Any, ::Int64, !Matched::Int64; kwargs…) where T at /home/sai/.julia/packages/LinearMaps/RjB5x/src/functionmap.jl:21
(LinearMaps.FunctionMap{T,F1,F2} where F2 where F1)(::Any, ::Any, !Matched::Int64; kwargs…) where T at /home/sai/.julia/packages/LinearMaps/RjB5x/src/functionmap.jl:22
…
Stacktrace:
[1] LinearMap{Float64}(::Function, ::Int64, ::Vararg{Any,N} where N; kwargs::Base.Iterators.Pairs{Union{},Union{},Tuple{},NamedTuple{(),Tuple{}}}) at /home/sai/.julia/packages/LinearMaps/RjB5x/src/LinearMaps.jl:164
[2] LinearMap{Float64}(::Function, ::Int64, ::Vararg{Any,N} where N) at /home/sai/.julia/packages/LinearMaps/RjB5x/src/LinearMaps.jl:164
[3] top-level scope at In[84]:4
I am sorry for my mistake about the error. The actual error is:
syntax: invalid keyword argument syntax “true” around In[88]:4
Stacktrace:
[1] top-level scope at none:1
Why did you not copy the code I posted?
H = LinearMap{Float64}(multiplyHPsi!, 2^N; ismutating = true)
I think you want to read a bit about basic julia syntax.
1 Like
Thanks for your reply. I have used the syntax and read about the LinearMap function. Now, I need to use the eigs function to decompose the resulting matrix. The code is as follows:
N = 14 # N=14 took 4.5 sec; N=16 took 35 sec; N=18 took 122 sec
m = 1 # how many states?
H = LinearMap(multiplyHPsi!, 2^N, ismutating = true, ishermitian = true)
#H = LinearMap(multiplyHPsi!, 2^N, isreal=false, ismutating=true, ishermitian=true) # use this for complex hermitian H
@time u, v = eigs(H; nev=m, which=:SR) #:SR stands for smallest real part
#Energy = real(u[1] + shiftE*N) # let us undo the shift in energy
#Psi = v
#E = Energy
The error thrown is:
MethodError: no method matching multiplyHPsi!(::SubArray{Float64,1,Array{Float64,1},Tuple{UnitRange{Int64}},true}, ::SubArray{Float64,1,Array{Float64,1},Tuple{UnitRange{Int64}},true})
Stacktrace:
[1] A_mul_B!(::SubArray{Float64,1,Array{Float64,1},Tuple{UnitRange{Int64}},true}, ::LinearMaps.FunctionMap{Float64,typeof(multiplyHPsi!),Nothing}, ::SubArray{Float64,1,Array{Float64,1},Tuple{UnitRange{Int64}},true}) at /home/sai/.julia/packages/LinearMaps/RjB5x/src/functionmap.jl:112
[2] mul!(::SubArray{Float64,1,Array{Float64,1},Tuple{UnitRange{Int64}},true}, ::LinearMaps.FunctionMap{Float64,typeof(multiplyHPsi!),Nothing}, ::SubArray{Float64,1,Array{Float64,1},Tuple{UnitRange{Int64}},true}, ::Bool, ::Bool) at /home/sai/.julia/packages/LinearMaps/RjB5x/src/LinearMaps.jl:78
[3] mul! at /home/sai/.julia/packages/LinearMaps/RjB5x/src/LinearMaps.jl:76 [inlined]
[4] (::Arpack.var"#matvecA!#24"{LinearMaps.FunctionMap{Float64,typeof(multiplyHPsi!),Nothing}})(::SubArray{Float64,1,Array{Float64,1},Tuple{UnitRange{Int64}},true}, ::SubArray{Float64,1,Array{Float64,1},Tuple{UnitRange{Int64}},true}) at /home/sai/.julia/packages/Arpack/o35I5/src/Arpack.jl:156
[5] aupd_wrapper(::Type{T} where T, ::Arpack.var"#matvecA!#24"{LinearMaps.FunctionMap{Float64,typeof(multiplyHPsi!),Nothing}}, ::Arpack.var"#18#25", ::Arpack.var"#19#26", ::Int64, ::Bool, ::Bool, ::String, ::Int64, ::Int64, ::String, ::Float64, ::Int64, ::Int64, ::Array{Float64,1}) at /home/sai/.julia/packages/Arpack/o35I5/src/libarpack.jl:83
[6] _eigs(::LinearMaps.FunctionMap{Float64,typeof(multiplyHPsi!),Nothing}, ::UniformScaling{Bool}; nev::Int64, ncv::Int64, which::Symbol, tol::Float64, maxiter::Int64, sigma::Nothing, v0::Array{Float64,1}, ritzvec::Bool) at /home/sai/.julia/packages/Arpack/o35I5/src/Arpack.jl:181
[7] #eigs#16 at /home/sai/.julia/packages/Arpack/o35I5/src/Arpack.jl:66 [inlined]
[8] #eigs#9 at /home/sai/.julia/packages/Arpack/o35I5/src/Arpack.jl:45 [inlined]
[9] macro expansion at ./timing.jl:174 [inlined]
[10] top-level scope at ./In[195]:4
Sorry my mistake. I forgot write the function to accept AbstractVector instead Vector. This solved the issue.
Thank you everyone for your valuable time and inputs.
1 Like
Hi! I’m following Vidal’s PSI online course too.
Do you succeed in this notebook?
I found the following strange behavior.
If I define LinearMap like this:
H = LinearMap{Float64}(multiplyHPsi!, 2^N; ismutating=true, issymmetric = true)
I need to check whether it does what it should, so I test
Psi = randn(2^N)
nPsi = randn(2^N)
Psi = Psi / norm(Psi)
nPsi1 = multiplyHPsi!(nPsi, Psi)
nPsi2 = H(Psi)
Now, the result nPsi1 is OK and looks correct; however, nPsi2 looks like random outcomes and looks nothing close to nPsi1, as it should.
I’m confused what goes wrong here.
That’s odd; could you post the full code that gives you the issue, including a working version of multiplyHPsi!. You can also put it in a gist and just paste the link here.
Hi Jutho!
Thank you very much for the reply.
I’ve create and .jl file to reproduce the above-mentioned issue.
Here is the link:
https://gist.github.com/brucelyu/b038114b507ab5ca241f78a3f93eb0f2