How do I solve Schrödinger equation with a piecewise periodic Hamiltonian? [duplicate]

Question

My Hamiltonian is

HI(t) = 
  Piecewise[
    {{H1, n τ <= t <= (n + 1/2) τ}, 
     {H2, (n + 1/2) τ <= t <= (n + 1) τ}}]

where

  H1 = {{x, y, 0}, {y, x, 0}, {0, 0, x}}; 
  H2 = {{x, y, z}, {y, x, 0}, {z, 0, x}};

Lets take values for numerical solution:

x = 2; y = 3; z = 4; z = 4; 

and

n = 20; τ = 0.5; 

initial condition:

{0,1,0};    

I want to solve the time-dependent Schrödinger equation for the piecewise periodic Hamiltonian:

$$i\frac{\partial\psi(t)}{\partial t}=HI(t)\psi(t)$$

We can also think this problem as

$$i \psi'(t) - HI(t) \psi(t) = 0,$$

where $HI(t)$ is piecewise continuous function.

Answer 1

I think possible issue here was the fact that Piecewise is used here to distinguish between different matrices inside NDSolve. Since Piecewise remains unevaluated until you put in specific time values, this hides the matrix structure from NDSolve at the initial step where it parses the differential equation.

Also, I think the periodicity you want to achieve is formulated in a way that can't be used in an actual computation, so I would suggest rewriting it with Mod instead.

To be on the safe side, I used MapThread to put the Piecewise inside the matrices (acting element-wise). Then I also defined the solution vector element-wise and used Thread in the vectorial differential equation to convert it into a system of equations. Here is the result:

H1 = {{x, y, 0}, {y, x, 0}, {0, 0, x}}; H2 = {{x, y, z}, {y, x, 
   0}, {z, 0, x}};

x = 2; y = 3; z = 4; z = 4; n = 20; τ = 0.5;

h[t_] = MapThread[
  Piecewise[{{#1, Mod[t, τ] <= τ/2}, {#2, True}}] &, {H1, 
   H2}, 2]

$$\left( \begin{array}{ccc} 2 & 3 & \begin{array}{cc} \{ & \begin{array}{cc} 0 & (t \bmod 0.5)\leq 0.25 \\ 4 & \text{True} \\ \end{array} \\ \end{array} \\ 3 & 2 & 0 \\ \begin{array}{cc} \{ & \begin{array}{cc} 0 & (t \bmod 0.5)\leq 0.25 \\ 4 & \text{True} \\ \end{array} \\ \end{array} & 0 & 2 \\ \end{array} \right)$$

ψ[t_] = Through[Array["ψ", 3][t]];

sol[t_] = ψ[t] /. 
   First@NDSolve[{Thread[
       I ψ'[t] == h[t].ψ[t]], ψ[0] == {0, 1, 
        0}}, ψ[t], {t, 0, n τ}];

ParametricPlot3D[Re@sol[t], {t, 0, n τ}, PlotStyle -> Tube[.01]]

Answer 2

It is simple like that:

H1 = {{x, y, 0}, {y, x, 0}, {0, 0, x}};
H2 = {{x, y, z}, {y, x, 0}, {z, 0, x}};

x = 2; y = 3; z = 4;
n = 20; \[Tau] = 0.5;
HI[t_?NumericQ] := 
 Piecewise[{{H1, 
    n \[Tau] <= t <= (n + 1/2) \[Tau]}, {H2, (n + 1/2) \[Tau] <= 
     t <= (n + 1) \[Tau]}}]

NDSolve[{\[Psi][n \[Tau]] == {0, 1, 0}, 
 I \[Psi]'[t] == HI[t].\[Psi][t]}, \[Psi][t]
, {t, n \[Tau], (n + 1) \[Tau]}]