XZ method tested in HO basis (does not work)

XZ_technique.jl

@@ -0,0 +1,54 @@
+using LinearAlgebra, Plots
+include("ho_basis.jl")
+include("p_space.jl")
+
+μ = 0.5
+l = 0
+V1 = -5
+R1 = sqrt(3)
+V2 = 2
+R2 = sqrt(10)
+n_max = 15
+
+ns = collect(0:n_max)
+ls = fill(l, n_max + 1)
+
+T = 1/μ * sp_T_matrix(ns, ls)
+V = V1 .* V_Gaussian.(R1, l, ns, transpose(ns)) + V2 .* V_Gaussian.(R2, l, ns, transpose(ns))
+
+n_EC = 4
+train_cs = (0.7 .+ 0.05 * randn(n_EC)) - 1im * (0.2 .+ 0.05 * randn(n_EC))
+target_cs = range(0.77, 0.22, 6)
+
+train_E = zeros(ComplexF64, n_EC)
+EC_basis = zeros(ComplexF64, (n_max + 1, length(train_cs)))
+exact_E = zeros(ComplexF64, length(target_cs))
+extrapolate_E = similar(exact_E)
+
+near_E = 0.2 + 0.2im
+
+for (j, c) in enumerate(train_cs)
+    H = T + c .* V
+    evals, evecs = eigen(H)
+    i = argmin(abs.(evals .- near_E))
+    train_E[j] = evals[i]
+    EC_basis[:, j] = evecs[:, i]
+end
+
+N_EC = transpose(EC_basis) * EC_basis
+
+for (j, c) in enumerate(target_cs)
+    exact_E[j] = quick_pole_E((p, q) -> c*(V1*g0(R1, p, q) + V2*g0(R2, p, q)), μ; cs_angle=0.5)
+
+    H = T + c .* V
+    H_EC = transpose(EC_basis) * H * EC_basis
+    evals = eigvals(H_EC, N_EC)
+    i = argmin(abs.(evals .- exact_E[j]))
+    extrapolate_E[j] = evals[i]
+end
+
+scatter(real.(train_E), imag.(train_E), label="training")
+scatter!(real.(exact_E), imag.(exact_E), label="exact")
+scatter!(real.(extrapolate_E), imag.(extrapolate_E), label="extrapolated")
+xlims!(-0.2,0.3)
+ylims!(-0.3,0.3)