1
$\begingroup$

So I went through the physics of a baseball and I got the exact same equations as shown here on page 8. I input all of the same parameters into my code, and I don't get the same thing. Any idea where I am going wrong?

ClearAll[t, x, y, z];
parms = {Cd -> .3, Cm -> 1, ωx -> -1500, ωy -> 0, ωz -> 0, m -> .142, ρ -> 1.225, A -> .608, R -> .22};
term = Sqrt[x'[t]^2 + y'[t]^2 + z'[t]^2];
eq1 = m x''[t] == -(1/2 ρ*A*Cd*x'[t]*term) + (4 π*ρ*R^3*(ωy*z'[t] - ωz*y'[t]));
eq2 = m y''[t] == -(1/2 ρ*A*Cd*y'[t]*term) + (4 π*ρ*R^3*(ωz*x'[t] - ωx*z'[t]))-9.81*m;
eq3 = m z''[t] == -(1/2 ρ*A*Cd*z'[t]*term) + (4 π*ρ*R^3*(ωx*y'[t] - ωy*x'[t]));
ic1 = {x'[0] == 0, x[0] == 0};
ic2 = {y'[0] == 0, y[0] == 1.6};
ic3 = {z'[0] == 90, z[0] == 0};
sol = NDSolve[{eq1, eq2, eq3, ic1, ic2, ic3} /. parms, {x[t], y[t], z[t]}, {t, 0, 1}];
ParametricPlot[Evaluate[{z[t], y[t]} /. sol], {t, 0, .01}, PlotRange -> 1.8]

This should be the same as the bottom image on page 13.

$\endgroup$
  • $\begingroup$ Subscript[[Omega], y] is missing the assignment rule $\endgroup$ – Dr. belisarius Dec 9 '13 at 2:37
  • $\begingroup$ it should be 0 there $\endgroup$ – yankeefan11 Dec 9 '13 at 2:38
  • $\begingroup$ can't read your code with all this subscript this and subscript that. Sorry. Can't read code, can't help. $\endgroup$ – Nasser Dec 9 '13 at 3:07
  • $\begingroup$ When I copied it from mathematica that is how it copied over :( $\endgroup$ – yankeefan11 Dec 9 '13 at 3:11
  • 1
    $\begingroup$ I know that offcourse. But the problem is that when copying it back from here to my notebook, the subscripts stay in the long form. Hence impossible to read the code. That is why it is not good idea to use it subscripts in code. Become hard to move around as plain text code. $\endgroup$ – Nasser Dec 9 '13 at 3:16
5
$\begingroup$

As I said in a comment, you have some errors on your constants:

ClearAll[t, x, y, z];
parms = {Cd -> .3, Cm -> 1, ωx -> 0, ωy -> 0, ωz -> -1500/60, m -> .142, ρ -> 1.225, 
         A -> Pi .03^2, R -> .03};
term = Sqrt[x'[t]^2 + y'[t]^2 + z'[t]^2];
eq1 = m x''[t] == -(1/2 ρ*A*Cd*x'[t]*term) + (4 π*ρ*R^3*(ωy*z'[t] - ωz*y'[t]));
eq2 = m y''[t] == -(1/2 ρ*A*Cd*y'[t]*term) + (4 π*ρ*R^3*(ωz*x'[t] - ωx*z'[t])) - 9.81*m;
eq3 = m z''[t] == -(1/2 ρ*A*Cd*z'[t]*term) + (4 π*ρ*R^3*(ωx*y'[t] - ωy*x'[t]));
ic1 = {x'[0] == 90, x[0] == 0};
ic2 = {y'[0] == 0,  y[0] == 1.6};
ic3 = {z'[0] == 0,  z[0] == 0};
sol = NDSolve[{eq1, eq2, eq3, ic1, ic2, ic3} /. parms, {x[t], y[t], z[t]}, {t, 0, 3}];
ParametricPlot[Evaluate[{x[t], y[t]} /. sol], {t, 0, .2}, AspectRatio -> 1,
               PlotRange -> {{0, 18}, {0, 1.7}}, GridLines -> Automatic]

Mathematica graphics

| improve this answer | |
$\endgroup$
  • 1
    $\begingroup$ area = pi * r^3? :-) $\endgroup$ – Peltio Dec 9 '13 at 5:40
  • 3
    $\begingroup$ @Peltio You know, it's the New Kind of Geometry :D $\endgroup$ – Dr. belisarius Dec 9 '13 at 5:44
  • 4
    $\begingroup$ @belisarius You all gonna end in Wolfram Hell for those jokes :P $\endgroup$ – Kuba Dec 9 '13 at 5:49
  • 5
    $\begingroup$ @Kuba You don't know how many times a day I feel like being already there :) $\endgroup$ – Dr. belisarius Dec 9 '13 at 6:10
  • 1
    $\begingroup$ @belisarius alternating between limbo, purgatory and hell, tho. $\endgroup$ – Yves Klett Dec 9 '13 at 10:21

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.