0
$\begingroup$

Given the system of differential equations:

$i\frac{d}{dt}c_{1}(t)=H_{12} E(t)e^{(iw_{12}t)}c_2(t)+H_{13} E(t)e^{(iw_{13}t)}c_3(t)$, $i\frac{d}{dt}c_{2}(t)=H_{21} E(t)e^{(iw_{21}t)}c_2(t)+H_{13} E(t)e^{(iw_{13}t)}c_3(t)$, $i\frac{d}{dt}c_{3}(t)=H_{31} E(t)e^{(iw_{31}t)}c_2(t)+H_{32} E(t)e^{(iw_{32}t)}c_2(t)$,

where $H_{kn}=0.005, w_{12}=-0.001, w_{21}=0.001,w_{13}=1.000,w_{31}=-1.000,w_{32}=0.999,w_{23}=0.999$ $E(t)=E_0 sech(\frac{t-400}{150})^2 sin(t)$, where $E_0=9.0$,

I wish to plot $|c_1(t)|^2,|c_2(t)|^2,|c_3(t)|^2 $vs $t$,

I'm relatively new in using Mathematica. I tried to solve this differential equation by using NDSolve. Although I get the output as an interpolating function, it still does not generate any plot. Any help will be highly helpful!

My code is

Clear[x, y, z, s, c1, c2, c3, t, w12, w21, w23, w23, w32, w13, w31, e0, ham]
w12 = -0.001;
w21 = 0.001;
w13 = -1.000;
w31 = 1.000;
w23 = -0.999;
w32 = 0.999;
ham = 0.005;
e0 = 9.0;
paramND = {
  I*D[c1[t], t] == 
    ham*(e0*Sech[(t - 400)/150]^2*Sin[t])*((Exp[I*w12*t]*c2[t]) + (Exp[I*w13*t]*c3[t])),
  I*D[c2[t], t] == 
    ham*(e0*Sech[(t - 400)/150]^2Sin[t])*((Exp[I*w21*t]*c1[t])+ (Exp[I*w23*t]*c3[t])),
  I*D[c3[t], t] == 
    ham*(e0*Sech[(t - 400)/150]^2*Sin[t])*((Exp[I*w31*t]*c1[t]) + (Exp[I*w32*t]*c2[t])),
   c1[0.01] == 1.0, c2[0.01] == 0.0, c3[0.01] == 0.0};
numSol = 
  NDSolve[paramND, {c1[t], c2[t], c3[t]}, {t, 0.01, 1}, 
    MaxSteps -> ∞, Method -> "ExplicitRungeKutta"] ;
paramPlot = {c1[t]^2, c2[t]^2, c3[t]^2} //. numSol;
Plot[Evaluate @ paramPlot, {t, 0.01, 800},
  PlotRange -> {0.01,1.0}, PlotStyle -> {{Black, Thick}}]
$\endgroup$
  • $\begingroup$ ps. you also have typo above, you have Clear=[...], there should be no = $\endgroup$ – Nasser Mar 9 '15 at 13:57
1
$\begingroup$

it works, you just have complex functions, also the range you had was not helping. It is also better I think to use {c1, c2, c3} instead of {c1[t], c2[t], c3[t]}

Clear[x, y, z, s, c1, c2, c3, t, w12, w21, w23, w23, w32, w13, w31, e0, ham]
w12 = -0.001;
w21 = 0.001;
w13 = -1.000;
w31 = 1.000;
w23 = -0.999;
w32 = 0.999;
ham = 0.005;
e0 = 9.0;
paramND = {I*D[c1[t], t] == ham*(e0*Sech[(t - 400)/150]^2*Sin[t])*
   ((Exp[I*w12*t]*c2[t]) + (Exp[I*w13*t]*c3[t])), 
   I*D[c2[t], t] == ham*(e0*Sech[(t - 400)/150]^2 Sin[t])*((Exp[I*w21*t]*c1[t]) + 
  (Exp[I*w23*t]*c3[t])), 
   I*D[c3[t], t] == ham*(e0*Sech[(t - 400)/150]^2*Sin[t])*((Exp[I*w31*t]*c1[t]) + 
   (Exp[I*w32*t]*c2[t])), c1[0.01] == 1.0, 
   c2[0.01] == 0.0, c3[0.01] == 0.0};
numSol = First@NDSolve[paramND, {c1, c2, c3}, {t, 0.01, 1}, MaxSteps -> \[Infinity], 
    Method -> "ExplicitRungeKutta"];
paramPlot = c3 /. numSol;
Plot[Re@paramPlot[t], {t, 0.01, 800}, PlotRange -> All, PlotStyle -> {{Black, Thick}}]

Mathematica graphics

Plot[{Re@paramPlot[t], Im@paramPlot[t]}, {t, 0.01, 1}, PlotRange -> All, 
 PlotTheme -> "Detailed"]

Mathematica graphics

Update for comment:

Plot[{Abs[(c1 /. numSol)[t]]^2, Abs[(c2 /. numSol)[t]]^2, 
 Abs[(c3 /. numSol)[t]]^2}, {t, 0.01, 800}, PlotRange -> All, PlotTheme -> "Detailed"]

Mathematica graphics

$\endgroup$
  • $\begingroup$ just a quick question if you could answer, how can I plot {|c1|^2,|c2|^2,|c3|^2 va time (t)} in the same plot ? I changed the time range from {t,0.01,800}, PlotRange->{0,1}. $\endgroup$ – Debopriyo Chaudhuri Mar 9 '15 at 14:09

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.