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}$

Improve this question
$\endgroup$
12
  • $\begingroup$ Check vars beceause is nonsense ? $\endgroup$
    – Mariusz Iwaniuk
    Commented Jul 27, 2018 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
    Commented Jul 27, 2018 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
    Commented Jul 27, 2018 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
    Commented Jul 27, 2018 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
    Commented Jul 27, 2018 at 13:03

1 Answer 1

Reset to default
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

Improve this answer
$\endgroup$
11
  • $\begingroup$ Thanks. Yes this is very helpful. Is it possible to plot vars as function of P, instead of vars(t)? $\endgroup$
    – AtoZ
    Commented Jul 28, 2018 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
    Commented Jul 28, 2018 at 12:34
  • 1
    $\begingroup$ Ok, P= {0, 0, 0, 0, 1, 1, 0, 0, 0, 0}. Where is the parameter? $\endgroup$
    – Alex Trounev
    Commented Jul 28, 2018 at 12:53
  • 1
    $\begingroup$ Plot3D[vars[p0,t],{p0,10},{t,0,25}] ? $\endgroup$
    – Alex Trounev
    Commented Jul 28, 2018 at 13:14
  • 1
    $\begingroup$ I updated the solution, added 3D figure. It's probably what you wanted to see. $\endgroup$
    – Alex Trounev
    Commented Jul 29, 2018 at 4:58

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.

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