This tutorial demonstrates the behaviour of a two-level system upon action of periodic perturbation. Just evaluate the notebook. Then use the last input cell and change the parameters for your enjoyment!

S. Rochester and D. Budker; October 6, 1999

Here we define the Schroedinger Equation (SE) and the wave function in the 2-state basis:

SE[H_, Ψ_] := I (∂_t) Ψ == H .  Ψ

MatrixForm[Ψ[t_] = {a[t], b[t]}]

Explicit form of the Hamiltonian and the SE;
note that the Hamiltonian is not quite Hermitian because of the i Γ/2 :

H=(({{0, v0 (E^(I ω t))}, {v0 (E^(-I ω t)), ω0 - I Γ/2}}));
MatrixForm[
      eqs = {LogicalExpand[
              SE[H, Ψ[t]]][[1]], LogicalExpand[
              SE[H, Ψ[t]]][[2]]}]

Initial condition: all atoms are in the ground state at t=0

MatrixForm[Ψ0 = {1, 0}]
inits = {LogicalExpand[Ψ[
                  0] == Ψ0][[1]], LogicalExpand[Ψ[0] == Ψ0][[2]]};

Now we solve the SE explicitly. Here is the place where we put numerical values for various parameters. Start with Γ=0, ω=ω0.

sols = NDSolve[
            Union[eqs /. {Γ -> 0, ω0 -> 
                    20, ω -> 20, v0 -> 1}, inits], {a, 
              b}, {t, 0, 100}, 
            MaxSteps -> 10000][[1]];

Now we plot the solution. Here we show the probability of finding the system in the excited state. We see the Rabi oscillations. Note that the oscillation (Rabi) frequency is (Subscript[Ω, R])=2*v0 (the wavefunction's amplitude oscillates with v0).

Plot[{Abs[b[t]]^2 /. sols}, {t, 0, 10}, PlotRange -> {0, 1}, 
      AxesLabel ->  {t, Probability}]

Now let us "turn on" relaxation. First let's try Γ<<v0:

sols = NDSolve[
            Union[eqs /. {Γ ->  .3, ω0 -> 
                    20, ω -> 20, v0 -> 1}, inits], {a, 
              b}, {t, 0, 100}, 
            MaxSteps -> 10000][[1]];
Plot[{Abs[b[t]]^2 /. sols}, {t, 0, 10}, PlotRange -> {0, 1}, 
      AxesLabel ->  {t, Probability}]

We see damped oscillation. Now try Γ>>v0:

sols = NDSolve[
            Union[eqs /. {Γ -> 1, ω0 -> 
                    20, ω -> 20, v0 ->  .1}, inits], {a, 
              b}, {t, 0, 100}, 
            MaxSteps -> 10000][[1]];
Plot[{Abs[b[t]]^2 /. sols}, {t, 0, 10}, PlotRange -> {0,  .1}, 
      AxesLabel ->  {t, Probability}]

Now we see that there are no more Rabi oscillations. This the "overdamped" regime. Note that at first the upper state population grows as if there was no relaxation, but then it "saturates" at a small level (Subscript[P, max])~(((((2  v0)/Γ)))^2), and then eventually decays away.  The maximum upper state population occurs at a time (Subscript[t, max])~2π/Γ.

More things to try and think about:
    1. What happens if the frequency ω  is different from the resonance frequency ω0 ?
    2. For ω=ω0, what is the critical value of v0/Γ at which the "overdamped" regime  turns into oscillatory regime?
    3. In the above, we assumed that the relaxation occurs due to the decay into unobserved states. What will change if the excited state decays back to the ground state?