- #1
kelly0303
- 573
- 33
I am trying to solve this two level (Schrödinger) equation as a function of time:$$i\begin{pmatrix}
\dot{x}\\
\dot{y}
\end{pmatrix} = \begin{pmatrix}
0 & iW+dE_0sin(\omega t)\\
-iW+dE_0sin(\omega t) & \Delta
\end{pmatrix}\begin{pmatrix}
x\\
y
\end{pmatrix}$$
(I can go into more details about the Hamiltonian if needed). The initial conditions are ##x(0)=1##, ##y(0)=0##. This is my Mathematica code for it:
W = 10;
OmegaRabi = 1000;
omega = 700000;
delta = 10000;
eqns = {I*x'[t] == (I*W + OmegaRabi*Sin[omega*t])*y[t],
I*y'[t] == (-I*W + OmegaRabi*Sin[omega*t])*x[t] + y[t]*delta,
x[0] == 1, y[0] == 0};
sol = NDSolve[eqns, {x, y}, {t, 0, 4}][[1]]
Plot[{Abs[x[t]]^2 + Abs[y[t]]^2 /. sol}, {t, 0, 1}]
All the variables are as in the original equation except for OmegaRabi which is equal to ##dE_0##. The output of this code (which would be the total probability i.e. ##|x(t)|^2+|y(t)^2|##) should be constant 1. However I get what is seen in the plot below. I assume this has to do with some rounding numerical errors (?). Is there a way to fix it and ideally have it deviate from 1 much slower as a function of time? I am interested in the ##|y(t)|^2## as a function of time, and if one looks at that, the values for it are around 0.00002, so currently the error after one second seems to be 5 times bigger than the value I am interested in.
\dot{x}\\
\dot{y}
\end{pmatrix} = \begin{pmatrix}
0 & iW+dE_0sin(\omega t)\\
-iW+dE_0sin(\omega t) & \Delta
\end{pmatrix}\begin{pmatrix}
x\\
y
\end{pmatrix}$$
(I can go into more details about the Hamiltonian if needed). The initial conditions are ##x(0)=1##, ##y(0)=0##. This is my Mathematica code for it:
W = 10;
OmegaRabi = 1000;
omega = 700000;
delta = 10000;
eqns = {I*x'[t] == (I*W + OmegaRabi*Sin[omega*t])*y[t],
I*y'[t] == (-I*W + OmegaRabi*Sin[omega*t])*x[t] + y[t]*delta,
x[0] == 1, y[0] == 0};
sol = NDSolve[eqns, {x, y}, {t, 0, 4}][[1]]
Plot[{Abs[x[t]]^2 + Abs[y[t]]^2 /. sol}, {t, 0, 1}]
All the variables are as in the original equation except for OmegaRabi which is equal to ##dE_0##. The output of this code (which would be the total probability i.e. ##|x(t)|^2+|y(t)^2|##) should be constant 1. However I get what is seen in the plot below. I assume this has to do with some rounding numerical errors (?). Is there a way to fix it and ideally have it deviate from 1 much slower as a function of time? I am interested in the ##|y(t)|^2## as a function of time, and if one looks at that, the values for it are around 0.00002, so currently the error after one second seems to be 5 times bigger than the value I am interested in.
