Time Evolution of wave function when Hamiltonian is time dependent

ITensor Julia Questions
, , ,
1

Consider I have a Hamiltonian(H(t)) which is time-dependent. I want to find the time evolution of a wave function which has the initial condition that is the ground state(\psi(0)) of the Hamiltonian(H(0)). How do I find \psi(t)?
For this, I need to solve the first-order differential equation, which is the Schrodinger equation. Most of the examples in using ITensor I found for time-independent Hamiltonian. I found a way in ITensor discourse to do this using trotter gate, but it is not helping me as in my case my Hamiltonian is for bosonic sites with particle number conserved has onsite interaction term like a_i^{\dagger}a_i^{\dagger}a_ia_i or n_i(n_i-1).

I just want to know if there exists any way that I can solve this kind of problem. Like in ED it can be very easily done using Odeint in Python, where I just pass Hamiltonian as a function of time and initial condition, and it will do all the jobs. I will be very much thankful as well as helpful if you can suggest and guide me on a way to find the time evolution that I want. For information ultimately I need to find an expectation value of an observable \langle\psi|J|\psi\rangle(t), which I can do I hope if I am able to find the time evolution of \psi(t) according to Schrodinger equation of time-dependent Hamiltonian.

Note: I have already asked this question with other details in the CPP version of ITensor(Time Evolution of wave function when Hamiltonian is time dependent.), but I am not getting any reply back. I generally prefer to write my code in CPP, but I will be more than happy to write my code in Julia if I can solve the problem as I mentioned in Julia version of ITensor.

2

What have you tried? And why do you think TEBD won’t work?

3

Thank you for your quick response.
I have tried using the Cpp code, which I mentioned in the NOTE with a link to the first question.

I have tried using Trotter decomposition, but whenever I use op for the four operators, an error occurs as follows

((1.0 / 2.0) * U_ONSITE) * op(sites, “Adag”, i) * op(sites, “Adag”, i) * op(sites, “A”, i) * op(sites, “A”, i);

or

((1.0 / 2.0) * U_ONSITE) * op(sites, “N”, i) * op(sites, “N”, i);

problem arises like

From line 912, file itensor.cc

Mismatched QN Index arrows

Mismatched QN Index arrows
Aborted (core dumped)

I also tried my code with TDVP, but it also did not work for me as in each iteration I need to update my H(t), and there appears to be an index mismatch in global sub-space expansion.

I just want to know if there exists any way both in CPP or Julia (I would prefer CPP, but I can write my code in Julia if I can solve my problem there) that I can solve this kind of problem. Like in ED solving Schrodinger equation can be very easily done using Odeint in Python, where I just pass Hamiltonian as a function of time and initial condition, and it will do all the jobs. TDVP looks similar, but I failed to do it for time-dependent Hamiltonian. I will be very much thankful as well as helpful if you can suggest and guide me on a way to find the time evolution that I want. For information ultimately I need to find an expectation value of an observable \langle\psi|J|\psi\rangle(t), which I can do I hope if I am able to find the time evolution of \psi(t) according to Schrodinger equation of time-dependent Hamiltonian.

4

To perform a matrix multiplication of two on-site operators you should use either of the following:

julia> using ITensors

julia> s = siteinds("S=1/2", 10; conserve_qns=true);

julia> @show op("S+ * S-", s, 1);
op("S+ * S-", s, 1) = ITensor ord=2
Dim 1: (dim=2|id=406|"S=1/2,Site,n=1")' <Out>
 1: QN("Sz",1) => 1
 2: QN("Sz",-1) => 1
Dim 2: (dim=2|id=406|"S=1/2,Site,n=1") <In>
 1: QN("Sz",1) => 1
 2: QN("Sz",-1) => 1
NDTensors.BlockSparse{Float64, Vector{Float64}, 2}
 2×2
Block(1, 1)
 [1:1, 1:1]
 1.0

julia> @show product(op("S+", s, 1), op("S-", s, 1));
product(op("S+", s, 1), op("S-", s, 1)) = ITensor ord=2
Dim 1: (dim=2|id=406|"S=1/2,Site,n=1")' <Out>
 1: QN("Sz",1) => 1
 2: QN("Sz",-1) => 1
Dim 2: (dim=2|id=406|"S=1/2,Site,n=1") <In>
 1: QN("Sz",1) => 1
 2: QN("Sz",-1) => 1
NDTensors.BlockSparse{Float64, Vector{Float64}, 2}
 2×2
Block(1, 1)
 [1:1, 1:1]
 1.0

op("S+", s, 1) and op("S-", s, 1) make ITensors with indices (s[1]', dag(s[1])), and * performs a tensor contraction over the shared indices. They both have the same sets of indices so * tries to contract down to a scalar over both indices, but the indices don’t have the proper arrow directions (only opposite arrows can contract) so you are getting an error. The alternatives I show above perform the matrix multiplication you were hoping to perform.

5

Thank you very much!
I have successfully done the time evolution. In the CPP version, the onsite interaction term for my case can be written as

op(sites, "Adag*Adag*A*A", i) * op(sites, "Id", i + 1);

and then can be pushed back to the gate.

Just for information, is there anything available to do the time evolution which is easier to use and more generalized? Like in TDVP for time-independent Hamiltonian, I can just pass the Hamiltonian and the initial state, and it can do its job. In ED in Python passing the the Hamiltonian as a function of time and initial state in ODEINT also works well. I am just finding a similar inbuilt package like this for MPS kind of things.