Ka Wa Yip
Ka Wa Yip

Reputation: 3021

Solving time-dependent Schrodinger equation using MATLAB ode45

The Schrodinger equation for a time-dependent Hamiltonian is:

$$i\hbar\frac{d}{dt}\psi(t) = H(t)\psi(t) \, .$$

I try to implement a solver for the Schrodinger equation for a time-dependent Hamiltonian in ode45. However, because the Hamiltonian $H(t)$ is dependent on time. I do not know how to do interpolation in ode45. Can you give me some hints?

psi0 = [0 1];
H = [1 0;0 1]*cos(t); %this is wrong, I do not know how to implement this and pass it to ode45
hbar = 1;
t    = [0:1:100];
[T, psi] = ode45(dpsi, t, psi);
function dpsi = f(t, psi, H, psi0)
dpsi = (1/i)*H*psi;

I also try to come up with a solution of matrix interpolation in MATLAB: Interpolation that involve a matrix.

Upvotes: 0

Views: 3105

Answers (1)

edwinksl
edwinksl

Reputation: 7093

H is just an identity matrix in your case, so we can just multiply it with the psi vector to get back the psi vector itself. Then, we bring i*hbar to the right-hand-side of the equation so that the final equation is in a form that ode45 accepts. Finally, we use the following code to solve for psi:

function schrodinger_equation

  psi0 = [0;1];
  hbar = 1;
  t = [0 100];
  [T,psi] = ode45(@(t,psi)dpsi(t,psi,hbar),t,psi0);

  for i = 1:length(psi0)
    figure
    plot(T,real(psi(:,i)),T,imag(psi(:,i)))
    xlabel('t')
    ylabel('Re(\psi) or Im(\psi)')
    title(['\psi_0 = ' num2str(psi0(i))])
    legend('Re(\psi)','Im(\psi)','Location','best')
  end

end

function rhs = dpsi(t,psi,hbar)
  rhs = 1/(1i*hbar)*cos(t).*ones(2,1);
end

Note that I have plotted the two components of psi separately and for each such plot, I have also plotted the real and imaginary components separately. Here are the plots for two different values of psi0:

enter image description here

enter image description here

Upvotes: 2

Related Questions