0
votes

I have an ODE and I solve it with NDSolve, then I plot the solution on a simplex in 2D.

Valid XHTML http://ompldr.org/vY2c5ag/simplex.jpg

Then I need to transform (align or just plot) this simplex in 3D at coordinates (1,0,0),(0,1,0),(0,0,1), so it looks like this scheme:

Valid XHTML http://ompldr.org/vY2dhMg/simps.png

I use ParametricPlot to do my plot so far. Maybe all I need is ParametricPlot3D, but I don't know how to call it properly.

Here is my code so far:

Remove["Global`*"];
phi[x_, y_] = (1*x*y)/(beta*x + (1 - beta)*y);
betam = 0.5;
betaf = 0.5;
betam = s;
betaf = 0.1;
sigma = 0.25;
beta = 0.3;
i = 1;
Which[i == 1, {betam = 0.40,  betaf = 0.60,  betam = 0.1,
   betaf = 0.1,  sigma = 0.25 , tmax = 10} ];
eta[x2_, y2_, p2_] = (betam + betaf + sigma)*p2 - betam*x2 -
   betaf*y2 - phi[x2, y2];
syshelp = {x2'[t] == (betam + betaf + sigma)*p2[t] - betam*x2[t] -
   phi[x2[t], y2[t]] - eta[x2[t], y2[t], p2[t]]*x2[t],
   y2'[t] == (betaf + betam + sigma)*p2[t] - betaf*y2[t] -
   phi[x2[t], y2[t]] - eta[x2[t], y2[t], p2[t]]*y2[t],
   p2'[t] == -(betam + betaf + sigma)*p2[t] + phi[x2[t], y2[t]] -
   eta[x2[t], y2[t], p2[t]]*p2[t]};
initialcond = {x2[0] == a, y2[0] == b, p2[0] == 1 - a - b};
tmax = 50;

solhelp =
   Table[
      NDSolve[
         Join[initialcond, syshelp], {x2, y2, p2} , {t, 0, tmax},
         AccuracyGoal -> 10, PrecisionGoal -> 15], 
      {a, 0.01, 1, 0.15}, {b, 0.01, 1 - a, 0.15}];

functions =
    Map[{y2[t] + p2[t]/2, p2[t]*Sqrt[3]/2} /. # &, Flatten[solhelp, 2]];

ParametricPlot[Evaluate[functions], {t, 0, tmax},
    PlotRange -> {{0, 1}, {0, 1}}, AspectRatio -> Automatic]

Third day with Mathematica...

2

2 Answers

3
votes

You could find a map from the triangle in the 2D plot to the one in 3D using FindGeometricTransformation and use that in ParametricPlot3D to plot your function, e.g.

corners2D = {{0, 0}, {1, 0}, {1/2, 1}};
corners3D = {{1, 0, 0}, {0, 1, 0}, {0, 0, 1}};

fun[pts1_, pts2_] := FindGeometricTransform[Append[pts2, Mean[pts2]], 
   PadRight[#, 3] & /@ Append[pts1, Mean[pts1]], 
  "Transformation" -> "Affine"][[2]]

ParametricPlot3D[Evaluate[fun[corners2D, corners3D][{##, 0}] & @@@ functions], 
  {t, 0, tmax}, PlotRange -> {{0, 1}, {0, 1}, {0, 1}}]

Mathematica graphics

2
votes

Since your solution has the property that x2[t]+y2[t]+p2[t]==1 it should be enough to plot something like:

functions3D = Map[{x2[t], y2[t], p2[t]} /. # &, Flatten[solhelp, 2]];

ParametricPlot3D[Evaluate[functions3D], {t, 0, tmax}, 
 PlotRange -> {{0, 1}, {0, 1}, {0, 1}}]

enter image description here