0
$\begingroup$

This is an ODE, I want to plot P on x-axis and vars on y-axis. My code works fine but gives an empty plot

Subscript[V, 0] = -2.5;
Subscript[k, 0] = Pi/2;
\[Epsilon] = 0.05;
n = 10;
B = Sqrt[3];
Subscript[t, fin] = 25;

A = Table[Subscript[V, i], {i, n}];
For[j = 1, j < n + 1, j++, A[[j]] = 0];
A[[n/2]] = Subscript[V, 0]*(1 + \[Epsilon]); 
A[[n/2 + 1]] = Subscript[V, 0]*(1 - \[Epsilon]);
P = Table[Subscript[\[Alpha], i], {i, n}];
For[j = 1, j < n + 1, j++, P[[j]] = 0];
P[[n/2]] = 1; P[[n/2 + 1]] = 1;
Subscript[\[Phi], 0][t_] := Subscript[\[Phi], 1][t];
Subscript[\[Phi], n + 1][t_] := Subscript[\[Phi], n][t];

eqns = Table[{Sqrt[-1]*Subscript[\[Phi], i]'[t] == 
 A[[i]]*Subscript[\[Phi], i][t] - Subscript[\[Phi], i + 1][t] - 
  Subscript[\[Phi], i - 1][
   t] + (P[[i]]*Abs[Subscript[\[Phi], i][t]]^2*
    Subscript[\[Phi], i][t]), Subscript[\[Phi], i][0] == B}, {i, 
n}];
vars = Table[
Sum[Abs[\[Phi][i, Subscript[t, fin]]]^2, {i, n/2 + 2, n}], {i, 
25}]/Table[Sum[Abs[\[Phi][i, 0]]^2, {i, n/2 - 1}], {i, 25}];

sol = NDSolve[eqns, vars, {t, 0, 25}, 
Method -> {"ExplicitRungeKutta", "StiffnessTest" -> False}, 
MaxSteps -> \[Infinity], AccuracyGoal -> 8, PrecisionGoal -> 8];
Plot=ListPlot[Table[vars /. First[%], {t, 0, 250}], 
PlotRange -> All, ImageSize -> 400]

with vars explicitly given by vars$=\frac{\sum_{n>\frac{M}{2}+1}|\phi_n(t_{fin})|^2}{\sum_{n<\frac{M}{2}}|\phi_n(0)|^2}$

$\endgroup$
  • $\begingroup$ Check vars beceause is nonsense ? $\endgroup$ – Mariusz Iwaniuk Jul 27 '18 at 10:33
  • $\begingroup$ No, its the transmission coefficient. So essentially it will be a plot of transmission coefficient as a function of P. did I make a mistake in writing it appropriately you mean? $\endgroup$ – AtoZ Jul 27 '18 at 10:37
  • 1
    $\begingroup$ Replace by my code: vars = Table[Subscript[\[Phi], i][t], {i, 0, n}];Plot[Evaluate[vars /. sol // ReIm], {t, 0, 25}, PlotRange -> All, ImageSize -> 400, PlotLegends -> Table[\[Phi][i][t], {i, 0, n}]]? $\endgroup$ – Mariusz Iwaniuk Jul 27 '18 at 10:49
  • $\begingroup$ vars=$\frac{\sum_{n>\frac{M}{2}+1}|\phi_n(t_{fin})|^2}{\sum_{n<\frac{M}{2}}|\phi_n(0)|^2}$ $\endgroup$ – AtoZ Jul 27 '18 at 10:54
  • 1
    $\begingroup$ If I check your equation a see a lot of subscripted functions like Subscript[\[Phi], 6][t] were as your variables look like (Abs[\[Phi][7, 25]]^2) +.... Your main problem is due to the use of subscripted variables! $\endgroup$ – Ulrich Neumann Jul 27 '18 at 13:03
1
$\begingroup$

It is possible that this code will be useful

Subscript[V, 0] = -2.5;
Subscript[k, 0] = Pi/2;
\[Epsilon] = 0.05;
n = 10;
B = Sqrt[3];
(*Subscript[t,fin]=25;*)

A = Table[Subscript[V, i], {i, n}];
For[j = 1, j < n + 1, j++, A[[j]] = 0];
A[[n/2]] = Subscript[V, 0]*(1 + \[Epsilon]);
A[[n/2 + 1]] = Subscript[V, 0]*(1 - \[Epsilon]);
P = Table[Subscript[\[Alpha], i], {i, n}];
For[j = 1, j < n + 1, j++, P[[j]] = 0];
P[[n/2]] = 1; P[[n/2 + 1]] = 1;
Subscript[\[Phi], 0][t_] := Subscript[\[Phi], 1][t];
Subscript[\[Phi], n + 1][t_] := Subscript[\[Phi], n][t];

eqns = Table[{Sqrt[-1]*Subscript[\[Phi], i]'[t] == 
     A[[i]]*Subscript[\[Phi], i][t] - Subscript[\[Phi], i + 1][t] - 
      Subscript[\[Phi], i - 1][
       t] + (P[[i]]*Abs[Subscript[\[Phi], i][t]]^2*
        Subscript[\[Phi], i][t]), Subscript[\[Phi], i][0] == B}, {i, 
    n}];


sol = NDSolve[eqns, 
   Table[Subscript[\[Phi], i][t], {i, 1, n}], {t, 0, 25}, 
   Method -> {"ExplicitRungeKutta", "StiffnessTest" -> False}, 
   MaxSteps -> \[Infinity], AccuracyGoal -> 8, PrecisionGoal -> 8];
f[i_, t_] := Abs[Subscript[\[Phi], i][t]]^2 /. sol
A = Sum[f[i, t] /. t -> 0, {i, 1, n/2 - 1}];
Plot[Sum[f[i, t], {i, n/2 + 2, n}]/A, {t, 0, 25}, 
 PlotLabel -> "vars(t)"]

Table[Plot[
  Evaluate[Abs[Subscript[\[Phi], i][t]]^2 /. sol], {t, 0, 25}, 
  PlotLabel -> i, AxesLabel -> Automatic], {i, 1, n}]

fig1

If there is a dependence on the parameter, then use the following code

Subscript[V, 0] = -2.5;
Subscript[k, 0] = Pi/2;
\[Epsilon] = 0.05;
n = 10;
B = Sqrt[3];
Subscript[t, fin] = 25;

A = Table[Subscript[V, i], {i, n}];
For[j = 1, j < n + 1, j++, A[[j]] = 0];
A[[n/2]] = Subscript[V, 0]*(1 + \[Epsilon]);
A[[n/2 + 1]] = Subscript[V, 0]*(1 - \[Epsilon]);
P = Table[Subscript[\[Alpha], i], {i, n}];
For[j = 1, j < n + 1, j++, P[[j]] = 0];
P[[n/2]] = p0; P[[n/2 + 1]] = p0;
Subscript[\[Phi], 0][t_] := Subscript[\[Phi], 1][t];
Subscript[\[Phi], n + 1][t_] := Subscript[\[Phi], n][t];

eqns = Table[{Sqrt[-1]*Subscript[\[Phi], i]'[t] == 
     A[[i]]*Subscript[\[Phi], i][t] - Subscript[\[Phi], i + 1][t] - 
      Subscript[\[Phi], i - 1][
       t] + (P[[i]]*Abs[Subscript[\[Phi], i][t]]^2*
        Subscript[\[Phi], i][t]), Subscript[\[Phi], i][0] == B}, {i, 
    n}];


F = Table[
   ParametricNDSolveValue[eqns, 
    Subscript[\[Phi], i], {t, 0, 25}, {p0}], {i, 1, n}];



{Plot3D[Evaluate[Abs[F[[n/2]][p0][t]]^2], {p0, 0, 10}, {t, 0, 25}, 
  Mesh -> None, ColorFunction -> Hue, 
  PlotLabel -> "\!\(\*SubscriptBox[\(\[Phi]\), \(5\)]\)(p0,t)", 
  AxesLabel -> {"p0", "t", ""}, PlotRange -> All], 
 Plot3D[Evaluate[
   Sum[Abs[F[[i]][p0][t]]^2, {i, n/2 + 2, n}]/
    Sum[Abs[F[[i]][p0][0]]^2, {i, 1, n/2 - 1}]], {p0, 0, 10}, {t, 0, 
   25}, Mesh -> None, ColorFunction -> Hue, PlotLabel -> "vars(p0,t)",
   AxesLabel -> {"p0", "t", ""}, PlotRange -> All]}

fig2

$\endgroup$
  • $\begingroup$ Thanks. Yes this is very helpful. Is it possible to plot vars as function of P, instead of vars(t)? $\endgroup$ – AtoZ Jul 28 '18 at 12:08
  • 1
    $\begingroup$ This is possible, but by definition P = 0. There must be parameters or at least one parameter, eg For[j = 1, j < n + 1, j++, P[[j]] = p0]; $\endgroup$ – Alex Trounev Jul 28 '18 at 12:34
  • 1
    $\begingroup$ Ok, P= {0, 0, 0, 0, 1, 1, 0, 0, 0, 0}. Where is the parameter? $\endgroup$ – Alex Trounev Jul 28 '18 at 12:53
  • 1
    $\begingroup$ Plot3D[vars[p0,t],{p0,10},{t,0,25}] ? $\endgroup$ – Alex Trounev Jul 28 '18 at 13:14
  • 1
    $\begingroup$ I updated the solution, added 3D figure. It's probably what you wanted to see. $\endgroup$ – Alex Trounev Jul 29 '18 at 4:58

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.