1
$\begingroup$

Good morning,

I have an ODE linear system $\dot{\vec f} = A \vec f$ depending parametrically to two real parameters, s and p. I'd like to get a density map of the solutions of $f_A/(f_A+f_B)$ at $t=20$, in the space of parameters ${s, p} \in [0; 1]$.

I post my try:

rA = 2.1
rB = 1.0
fA := ParametricNDSolveValue[{A'[ t] == (rA - 1 + p) A[ t] + s B[t], B'[t] == (1 - p) A[t] + (rB - s) B[t], A[0] == 2, B[0] == 2}, A, {t, 0, 20}, {{s, p} \[Element] Reals}]
fB := ParametricNDSolveValue[{A'[ t] == (rA - 1 + p) A[ t] + s B[t], B'[t] == (1 - p) A[t] + (rB - s) B[t], A[0] == 2, B[0] == 2}, B, {t, 0, 20}, {{s, p} \[Element] Reals}]
DensityPlot[fA[s][p]/(fA[s][p] + fB[s][p]), {s, 0, 1}, {p, 0, 1}]

Here there's what I get

The plot I obtain can't be what I need, as I never ask Mathematica to evaluate the functions at $t=20$. I'm really new to Mathematica and I've the impression to be missing something really basic about its semantics. I hope to learn from your help, thank you in advance!

|improve this question
$\endgroup$
  • $\begingroup$ I know that nowhere in my code there is $t=20$, but when I tried few things to do that I got errors like "ParametricNDSolve::dsvar: 20 cannot be used as a variable" $\endgroup$ – Orso Jul 8 '14 at 9:06
  • 1
    $\begingroup$ I cut and pasted your code and the plot works? Try Clear[t] first, it looks like you may have executed t=20 before. $\endgroup$ – Ymareth Jul 8 '14 at 9:16
  • $\begingroup$ Hi @Ymareth and thank you for your help. Indeed, I already do Clear[t] before the lines I posted but it doesn't change anything in my case... $\endgroup$ – Orso Jul 8 '14 at 9:21
  • $\begingroup$ You're right, I pasted it into a new notebook and I obtain a plot. The problem is that nobody told the machine to plot at $t=20$, so I guess the plot is not at all what I need... $\endgroup$ – Orso Jul 8 '14 at 9:38
1
$\begingroup$

I may have misunderstood the aim. I post this as motivating:

I have re-written the system of equations (no guarantee error free):

mat = {{1.1 + p, s}, {1 - p, 1 - s}}; pf = 
 ParametricNDSolveValue[{{a'[t], b'[t]} == 
    mat.{a[t], b[t]}, {a[0], b[0]} == {2, 2}}, {a[t], b[t]}, {t, 0, 
   20}, {s, p}];
foi[par1_, par2_, u_] := (#1/(#2 + #1) & @@ pf[par1, par2]) /. t -> u

Looking at the solutions $a(t),b(t)$:

Manipulate[
 Column[{ParametricPlot[pf[v1, v2] /. t -> w, {w, 0, 20}, 
    PlotRange -> {{0, 10^12}, {0, 10^12}}, Frame -> True, 
    FrameLabel -> {"a", "b"}], 
   LogPlot[Evaluate[pf[v1, v2] /. t -> w], {w, 0, 20}, 
    PlotLegends -> {"a(t)", "b(t)"}
    , PlotRange -> {0, 10^12}, Frame -> True, 
    FrameLabel -> {"t", "a|b"}]}], {{v1, 0, "s"}, 0, 
  1}, {{v2, 0, "p"}, 0, 1}]

enter image description here

As expected with these linear system with constant positive coefficients->exponentials

The desired density plot for function (foi evaluated at t=20) above:

DensityPlot[foi[x, y, 20], {x, 0, 1}, {y, 0, 1}, 
 ColorFunction -> "SunsetColors", PlotLegends -> Automatic, 
 MeshFunctions -> (#3 &), Mesh -> {Range[0, 1, 0.1]}, 
 FrameLabel -> {"s", "p"}]

enter image description here

or in 3D:

Plot3D[foi[x, y, 20], {x, 0, 1}, {y, 0, 1}, 
 ColorFunction -> (ColorData["SunsetColors"][#3] &), 
 MeshFunctions -> (#3 &), Mesh -> {Range[0, 1, 0.1]}, 
 AxesLabel -> {"s", "p", 
   "\!\(\*FractionBox[\(a \((20)\)\), \(a \((20)\) + b \
\((20)\)\)]\)"}]

enter image description here

Apologies if I have misinterpreted your aim.

|improve this answer
$\endgroup$
  • $\begingroup$ Exactly, @ubpqdn, that's what I wanted! Thank you so much, I've learnt a lot from your solution! $\endgroup$ – Orso Jul 8 '14 at 13:48

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.