0
$\begingroup$

I am trying to plot the solution given in the code with respect to "delc". Now the problem is that it can be plot for a particular value of "t" like t=10,20,50,60,.. upto 100 but what i need is mean of all these plots. Is it possible to plot that by adjusting the y-axis so that i can get one single graph in which time also varies upto 100.

 w1 = 1;
 gma1 = 0.005;
 n1 = 1;
 gma2 = 0.005;
 G1 = 0.005;
 k1 = .1;
 k2 = 0.1;
 a1 = 0.07;
 a2 = 0.58;
 k0 = 0.1;
 Q1 = 1.268;
 del0 = 1;
 N1 = 1;
 ome = 1;
 M1 = del0*(1 - Cos[ome*t]);
 s = ParametricNDSolveValue[{V11'[t] - V21[t]*w1 - V12[t]*w1 == 0, 
   V12'[t] - V22[t]*w1 + w1*V11[t] + gma1*V12[t] - 
  Sqrt[2]*G1*a1*V13[t] - Sqrt[2]*G1*a2*V14[t] == 0, 
V13'[t] - V23[t]*w1 + k1*V13[t] + 
  Sqrt[2]*G1*a2*V11[t] - (-G1*Q1 + delc)*V14[t] == 0, 
V14'[t] - V24[t]*w1 + k1*V14[t] - 
  Sqrt[2]*G1*a1*V11[t] - (G1*Q1 - delc)*V13[t] == 0, 
V21'[t] + V11[t]*w1 + gma1*V21[t] - Sqrt[2]*G1*a1*V31[t] - 
  Sqrt[2]*G1*a2*V41[t] - w1*V22[t] == 0, 
V22'[t] + V12[t]*w1 + gma1*V22[t] - Sqrt[2]*G1*a1*V32[t] - 
  Sqrt[2]*G1*a2*V42[t] + w1*V21[t] + gma1*V22[t] - 
  Sqrt[2]*G1*a1*V23[t] - Sqrt[2]*G1*a2*V24[t] - gma2*(2*n1 + 1) ==
  0, V23'[t] + w1*V13[t] + gma1*V23[t] - Sqrt[2]*G1*a1*V33[t] - 
  Sqrt[2]*G1*a2*V43[t] + k1*V23[t] + 
  Sqrt[2]*G1*a2*V21[t] - (-G1*Q1 + delc)*V24[t] == 0, 
V24'[t] + V14[t]*w1 + gma1*V24[t] - Sqrt[2]*G1*a1*V34[t] - 
  Sqrt[2]*G1*a2*V44[t] + k1*V24[t] - 
  Sqrt[2]*G1*a1*V21[t] - (G1*Q1 - delc)*V23[t] == 0, 
V31'[t] + k1*V31[t] + 
  Sqrt[2]*G1*a2*V11[t] - (-G1*Q1 + delc)*V41[t] - w1*V32[t] == 0, 
V32'[t] + k1*V32[t] + 
  Sqrt[2]*G1*a2*V12[t] - (-G1*Q1 + delc)*V42[t] + w1*V31[t] - 
  Sqrt[2]*G1*a2*V34[t] - Sqrt[2]*G1*a1*V33[t] + gma1*V32[t] == 0, 
V33'[t] + k1*V33[t] + 
  Sqrt[2]*G1*a2*V13[t] - (-G1*Q1 + delc)*V43[t] + k1*V33[t] + 
  Sqrt[2]*G1*a2*V31[t] - (-G1*Q1 + delc)*V34[t] - k0 == 0, 
V34'[t] + k1*V34[t] + 
  Sqrt[2]*G1*a2*V14[t] - (-G1*Q1 + delc)*V44[t] + k1*V34[t] - 
  Sqrt[2]*G1*a1*V31[t] - (G1*Q1 - delc)*V33[t] == 0, 
V41'[t] + k1*V41[t] - 
  Sqrt[2]*G1*a1*V11[t] - (G1*Q1 - delc)*V31[t] - w1*V42[t] == 0, 
V42'[t] + k1*V42[t] + 
  Sqrt[2]*G1*a1*V12[t] - (G1*Q1 - delc)*V32[t] + w1*V41[t] - 
  Sqrt[2]*G1*a2*V44[t] - Sqrt[2]*G1*a1*V43[t] + gma1*V42[t] == 0, 
V43'[t] + k1*V43[t] - 
  Sqrt[2]*G1*a1*V13[t] - (G1*Q1 - delc)*V33[t] + k1*V43[t] + 
  Sqrt[2]*G1*a2*V41[t] - (-G1*Q1 + delc)*V44[t] == 0, 
V44'[t] + k1*V44[t] - 
  Sqrt[2]*G1*a1*V14[t] - (G1*Q1 - delc)*V34[t] + k1*V44[t] - 
  Sqrt[2]*G1*a1*V41[t] - (G1*Q1 - delc)*V43[t] - k0 == 0, 
V11[0] == 1, V12[0] == 1, V13[0] == 0, V14[0] == 0, V21[0] == 0, 
V22[0] == 1, V23[0] == 0, V24[0] == 0, V31[0] == 0, V32[0] == 0, 
V33[0] == 0, V34[0] == 0, V41[0] == 0, V42[0] == 0, V43[0] == 0, 
V44[0] == 0}, Function[t,1/2*(V11[t] + V22[t] - 2*V12[t])^(-1)], {t, 0, 
100}, delc];
 Plot[s[delc][60], {delc, 0, 2}, Frame -> True, 
 PlotRange -> {All, {0.545, .5475}}, 
  FrameLabel -> {Style["\[CapitalDelta]c", Bold, 20], 
   Style["\!\(\*SubscriptBox[\(S\), \(q\)]\)", Bold, 20], 
    Style["t = 60", Bold, 20]}, 
  FrameTicksStyle -> Directive[FontSize -> 20], GridLines -> Automatic]
$\endgroup$
6
  • $\begingroup$ Moreover, the return value of ParametricNDSolveValue is a ParametricFunction not a Rule; so 1/2*(V11[t] + V22[t] - 2*V12[t])^(-1) /. s does not make sense at all. Instead, you can call f = ParametricNDSolveValue[eq, t \[Function] 1/2*(V11[t] + V22[t] - 2*V12[t])^(-1), {t, 0, 100}, delc] and plot f[delc][t]. $\endgroup$ – Henrik Schumacher Feb 3 '20 at 7:07
  • $\begingroup$ yes, it is delc. Now i have change it. $\endgroup$ – vini Feb 3 '20 at 7:08
  • $\begingroup$ @Henrik now please see the code again. $\endgroup$ – vini Feb 3 '20 at 7:22
  • $\begingroup$ I took the freedom to make another edit to your post. Please have a look. $\endgroup$ – Henrik Schumacher Feb 3 '20 at 7:29
  • $\begingroup$ change s = ParametricNDSolveValue[...] to s = ParametricNDSolve[...], set a value for the parameter delc (say param=.1) and use Evaluate[1/2*(V11[param][t] + V22[param][t] - 2*V12[delc][t])^(-1) /. s] in the first argument of Plot? $\endgroup$ – kglr Feb 3 '20 at 7:30
1
$\begingroup$
functions = {V11, V12, V13, V14, V21, V22, V23, V24, V31, V32, V33, 
   V34, V41, V42, V43, V44};

ClearAll[foo1]
foo1[d_, t_] := Through[s[d][t]]


delc = .1;
Plot[Evaluate@foo1[delc, t], {t, 0, 60}, Frame -> True, 
 FrameLabel -> {Style["Time", Bold, 20], Style[Subscript["S", "q"], Italic, Bold, 20]}, 
 PlotLegends -> Through[Through[functions[delc]][t]],
 ImageSize -> Large, FrameTicksStyle -> 20, PlotStyle -> Thick]

enter image description here

Plot[Mean@foo1[delc, t], {t, 0, 60}, Frame -> True, 
 FrameLabel -> {Style["Time", Bold, 20], Style[Subscript["S", "q"], Italic, Bold, 20]}, 
 PlotLegends -> {"Mean"}, ImageSize -> Large, FrameTicksStyle -> 20, 
 PlotStyle -> Thick]

enter image description here

Plot[Evaluate[foo1[delc, 50]], {delc, 0, 1}, Frame -> True, 
 FrameLabel -> {Style["Time", Bold, 20], Style[Subscript["S", "q"], Italic, Bold, 20]}, 
 PlotLegends -> LineLegend[Automatic, Through[Through[functions[ \[Delta]]][50]], 
   LegendLabel -> "t=50"], 
 ImageSize -> Large, FrameTicksStyle -> 20, PlotStyle -> Thick]

enter image description here

To plot 1/2*(V11[delc][t] + V22[delc][t] - 2*V12[delc][t])^(-1) you can do:

ClearAll[foo2]
foo2[d_, t_] := 1/2 /{1, 1, -2}.Through[s[d][[{1, 6, 2}]]@t];

Plot[foo2[delc,t], {t, 0, 60}, Frame -> True, 
 FrameLabel -> {Style["Time", Bold, 20], Style[Subscript["S", "q"], Italic, Bold, 20]}, 
 PlotLegends -> {Style[ToString[1/2*(V11[delc][t] + V22[delc][t] - 2*V12[delc][t])^(-1), 
     TraditionalForm], 20]}, 
 ImageSize -> Large, FrameTicksStyle -> 20, PlotStyle -> Thick]

enter image description here

$\endgroup$
2
  • $\begingroup$ i need the graph of "Sq" with respect to "delc" not with time. $\endgroup$ – vini Feb 3 '20 at 8:13
  • $\begingroup$ @vini, pls see the updated version. $\endgroup$ – kglr Feb 3 '20 at 8:37

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.