Animation of Differential Equations from NDSolve with ParametricPlot3D and Evaluate

Question

I have a system of differential equations (referred to as "s") and use NDSolve to obtain the solution. I substitute the interpolated functions for the original functions in order to represent the trajectories of two objects in 3 dimensions. Here is the system that I am working with:

(I added this as an edit after I got some feedback that some more information about what I was attempting to animate would be useful )

sysstars = {((sx[t]^2 + sy[t]^2 + sz[t]^2)^(3/2))*(1 + 
      2*(mstar1/mstar2) + (mstar1/mstar2)^2)*sx''[t] == 
   g*mstar2*(-sx[t]), ((sx[t]^2 + sy[t]^2 + sz[t]^2)^(3/2))*(1 + 
      2*(mstar1/mstar2) + (mstar1/mstar2)^2)*sy''[t] == 
   g*mstar2*(-sy[t]), ((sx[t]^2 + sy[t]^2 + sz[t]^2)^(3/2))*(1 + 
      2*(mstar1/mstar2) + (mstar1/mstar2)^2)*sz''[t] == 
   g*mstar2*(-sz[t]), sx[0] == tpostar1[[1]], sy[0] == tpostar1[[2]], 
  sz[0] == tpostar1[[3]], sx'[0] == v0com1[[1]], 
  sy'[0] == v0com1[[2]], 
  sz'[0] == 
   v0com1[[3]], ((qx[t]^2 + qy[t]^2 + qz[t]^2)^(3/2))*(1 + 
      2*(mstar2/mstar1) + (mstar2/mstar1)^2)*qx''[t] == 
   g*mstar1*(-qx[t]), ((qx[t]^2 + qy[t]^2 + qz[t]^2)^(3/2))*(1 + 
      2*(mstar2/mstar1) + (mstar2/mstar1)^2)*qy''[t] == 
   g*mstar2*(-qy[t]), ((qx[t]^2 + qy[t]^2 + qz[t]^2)^(3/2))*(1 + 
      2*(mstar2/mstar1) + (mstar2/mstar1)^2)*qz''[t] == 
   g*mstar2*(-qz[t]), qx[0] == tpostar2[[1]], qy[0] == tpostar2[[2]], 
  qz[0] == tpostar2[[3]], qx'[0] == v0com2[[1]], 
  qy'[0] == v0com2[[2]], qz'[0] == v0com2[[3]]}

s = NDSolve[
  sysstars, {sx[t], sy[t], sz[t], qx[t], qy[t], qz[t]}, {t, 0, 
   8*10^7}, MaxSteps -> 10^8]  

And here is how I am plotting them:

ParametricPlot3D[
 Evaluate[{sx[t], sy[t], sz[t], qx[t], qy[t], qz[t]} /. s, {t, 0, 
   10^7}], PlotRange -> All]
ParametricPlot3D[
 Evaluate[{qx[t] + comx[t], qy[t] + comy[t], qz[t] + comz[t], 
    sx[t] + comx[t], sy[t] + comy[t], sz[t] + comz[t]} /. s, {t, 0, 
   3*10^7}], PlotRange -> All,]

(The second Parametric Plot above is where the center of mass vector is added back to show the system from outside the center of mass inertial frame. So the first one plots two objects on initial trajectories that then go into orbit around the center of mass, with the center of mass treated as static)

The results of this are good - static lines in 3D, with the sensible paths. Now, I would like to animate this so that, say two spheres, move as a function of t. I have tried a few on-line suggestions using Manipulate[Show[..., but without success so far. I believe that the difficulty may relate to the use of Evaluate inside ParametricPlot3D - in conjunction with the animation techniques. It does not seem possible to get rid of Evaluate with ParamatricPlot3D, and still get ParametricPlot3D to work with the interpolated functions which get substituted to generate the representation. There must be a way to do this in Mathematica - the program is way to good not to have a good way to do such an obvious thing with solutions, but I am stumped.

Suggestions greatly appreciated (after too many hours) !

Delphi

Answer 1

Edit

Please note that your eq. could be written in a compact form:

eq[s_List, m1_, m2_, tPos_List, v0Com_List, g_] :=
  With[{st = Through[s[t]]},
   ((Norm[st]^3 (1 + 2 m1/m2 + (m1/m2)^2) #''[t] == -g m2 #[t]) & /@ s)
    ~Join~
    Thread[Through[s[0]] == tPos]
    ~Join~
    (Thread[D[st, t] == v0Com] /. t -> 0)
   ];

And solve them as

s = NDSolve[Join[
    eq[{sx, sy, sz}, 1, 1.7, {1, -1, 1}, {1, 1, 1}, 10],
    eq[{qx, qy, qz}, 1.7, 1, -{1, -1, 1}, -{1, 1, 1}, 10]],
   {sx, sy, sz, qx, qy, qz}, {t, 0, 8}, MaxSteps -> 10^3];

---

For Example

s = Table[k[j] -> Interpolation[Table[{i, Sin[i j]}, {i, 6}]], {j, 1, 6}]; 
Manipulate[
 ParametricPlot3D[{{k[4][t], k[5][t], k[6][t]}, {k[1][t], k[2][t], k[3][t]}} /. s, {t, 1, a}, 
  PlotStyle -> {Red, Green}, 
  Evaluated -> True, 
  PlotRange -> {{-1.5, 1.5}, {-1.5, 1.5}, {-1.5, 1.5}}], {a, 1.01, 6}]

Edit

Solving your equations:

(*Some constants*)
g = 10;
mstar1 = 1;
mstar2 = 1.5;
tpostar1 = {1, -1, 1};
tpostar2 = -{1, -1, 1};
v0com2 = {1, 1, 1};
v0com1 = {-1, -1, -1};

Your def for sysstars = ... comes here ....

s = NDSolve[sysstars, {sx[t], sy[t], sz[t], qx[t], qy[t], qz[t]},{t, 0, 8},MaxSteps -> 10^3]

Manipulate[
 ParametricPlot3D[{{sx[t], sy[t], sz[t]}, {qx[t], qy[t], qz[t]}} /. s[[1]], {t, 1, a},
  PlotStyle -> {Red, Green},
  PlotRange -> 5 {{-1, 1}, {-1, 1}, {-1, 1}}, 
  Evaluated -> True],
 {a, 0, 10}]

Edit

Adding my orbiting balls:

tt = Manipulate[Show@{
    ParametricPlot3D[{{sx[t], sy[t], sz[t]}, {qx[t], qy[t], qz[t]}} /. s[[1]], {t, 1, a},
     PlotStyle -> {Red, Green},
     PlotRange -> 8.5 {{-1, 1}, {-1, 1}, {-1, 1}}, Evaluated -> True],
    Graphics3D[{Sphere[{sx[t], sy[t], sz[t]}], 
                Sphere[{qx[t], qy[t], qz[t]}]} /. s[[1]]] /. t -> a},
  {a, 0, 10}]

Answer 2

Without knowing anything about your functions, I can see that the code you showed isn't quite correct: The Evaluate should wrap only the first plot argument and not the second one containing the plot range. Also, Evaluate can't be followed by { as you typed - it has to be [. Finally, the grouping of the vector components has to be in triplets as shown in the answer by belisarius.

Moreover, when making any kind of trajectory animation you have to fix the PlotRange and AspectRatio of your plot frames to a constant value, otherwise there will be no stationary reference frame to which the movement of your objects can be compared. This means you can't use PlotRange -> All here - except in the first plot over the full time range, to get an idea of the final plot range you're going to need.

After that, to make the animation, you should replace this by a fixed PlotRange as belisarius also shows.