# Thread: Differential equation problem (originally kinematics)

1. ## Differential equation problem (originally kinematics)

A point mass ‘A’ is kept at the origin. Another B is kept at the x axis at x = H. Another C is kept at distance H from origin and distance h from B. A, B and C thus form an isosceles triangle with vertex A. Given that A's mass >> B's mass >> C's mass; Newton’s gravitational law governs the bodies. C will have motion because of both A and B while B will have motion because of A only.
Find:
The equation of path of C in any coordinate system given that B and C have zero initial velocities.

This problem after a bit of physics will reduce to a problem of differential equations. I applied r, $\vartheta$ coordinate system and got two differential equations in r, $\vartheta$ and both their single and double time derivatives. A very complex eqn indeed!

Solve or suggest a way to solve this problem, please!

2. Numerically. It's easy to code Mathematica to solve coupled sets. If you want, post the equations along with the appropriate initial conditions, parameter values and I'll show you how to set up the NDSolve function to solve it.

3. The initial velocity of both B and C is 0. A is fixed to ground.

$\displaystyle r_{b/c}=\sqrt{r_b^2\,+\,r_c^2\,-2bc\cos \theta}\quad(1)$

$\displaystyle \frac{r_{b/c}}{\sin \theta}\,=\,\frac{r_b}{\sin \phi}\quad(2)$
These two equations define the angle $\displaystyle \phi$ and $\displaystyle r_b/c$ for us.

Proceeding...

$\displaystyle -\frac{k_b}{r_{b/c}^2}\,\left(1-\,\frac{r_b^2}{\sin^2 \theta r_{b/c}^2}\right)\,=\,\ddot{r_c}\,-\,r_c\dot{\theta^2}\quad(3)$
$\displaystyle -\left(\frac{k_b}{r_{b/c}^2}\,\frac{r_b}{\sin \theta r_{b/c}}\,+\,\frac{k_a}{r_c^2}\right) \,= \,r_c\ddot{\theta}\,+\,2\dot{r_c}\dot{\theta}\quad (4)$

$\displaystyle r_{b/c}$ is defined by first two equations. Eqations 3 and 4 then need be solved for $\displaystyle r_c$ and $\displaystyle \theta$.

When we get $\displaystyle r_c$ and $\displaystyle \theta$, we shall have got the eqn of path of C.

Method applied to reach till the equations 3 and 4:

I used polar coordinate systems with the point A as origin. I found out radial and angular components of the accleration of C and compard them with the standard expressions with the same as given here:
Polar coordinate system - Wikipedia, the free encyclopedia

I don't have any intermediate values so I dont know how yu gonna go for a numerical method. Anyways, if you can solve it by whatever means, I will be very happy and grateful.

Thanks.

4. Here's some observations:

1. Need to better define what all the parameters and variables are.

2. You have an extraneous symbol: $\displaystyle \phi$. Either it's a constant or some external function of t.

3. I don't see the significance of equation (2): $\displaystyle r_{b/c}$ is defined in (1) and that can be substituted into (3) and (4).

4. Make up values for the parameters and initial conditions:

$\displaystyle r_c(0)=1;\quad r'_c(0)=1$

$\displaystyle \theta(0)=\pi/2;\quad \theta'(0)=2$

and any other values that need to be initialized.

5. Once you have done this, then it's a breeze to solve it numerically in Mathematica and plot the results

5. ## Parametric values!

$\displaystyle r_c(t=0)=132.44 \times 10^{24} metres\quad \dot {r_c}(t=0)=0$

$\displaystyle \theta(t=0)=1.297176448\times10^{-13}\deg \quad \dot {\theta}(t=0)=0$

In equations 3 and 4 below in the page, $\displaystyle k_b = G\times m_a$
Similarly for $\displaystyle k_a$.

masses: m_a=1.6 * 10^55; m_b=2 * 10^30; m_c=6 * 10^24.

This way we quantify $\displaystyle k_a$ and $\displaystyle k_b$ in equations, 3 and 4.

Further, the expression for $\displaystyle r_b$ in the two equations will actually bring along a $\displaystyle t$ (time) term, too complicating matters for us. that is, $\displaystyle r_b=((132.44 \times 10^{24})^{1.5}\, -\, k_a t)^{2/3}$
Do you think it can be solved still? *dejected*

The expression for r_{b/c} (equation 1) also contains a r_b and we will have to put above equation in this equation.

After doing all this we have two equations (3 and 4 ie) which will thus have one $\displaystyle t$ term too along with those visible.
Well?

6. Taureau, I'm optimistic it can be solved but one has to be Zen about it: you don't just wip it out. You know that. Rather you approach the solution asymptotically, gradually building upon it's construct.

Also, during the initial phase, don't use all those complicated initial values. Just use simple values like 1, 2, 100, -5, etc. Yea, I know, that may result in quite different solutions, but then you gradually ramp up those values towards your actual values.

The important thing is to "rough it in" even if it's wrong. You can clean it up later. I've done that below although I suspect you won't be impressed. It's a start though and does yield a result (no syntax, logic errors). That's a big start. Now you can start "cleaning it up".

Keep in mind that as you ramp up the parameters, you may reach a point in which the software has problems doing the numerical calculations. That's a hazard you have to accept. You can always attempt to ramp up the precision and numerical accuracy.

I dropped all the subscripts in the code below as they interfered with the coding: I used "rbc" for r_{bc}, "rb" for r_b, and th(t) for $\displaystyle \theta(t)$ and so on.

Note when I do get (pseudo) results for $\displaystyle r_c(t)$ and $\displaystyle \theta(t)$, I assume these are in polar coordinates. I then plotted the trajectory as $\displaystyle r_c(t)\cos(\theta(t)),r_c(t)\sin(\theta(t))$. Not sure that's what you want.

I plotted the results for this run below. Don't laugh. I realize it's not even close. Remember: Zen.

Code:
g = 0.0001;
ma = 1000;
mb = 1000;
b = 1;
c = 1;
ka = g ma;
kb = g mb;
rb[t_] := (100^(1/5) - ka t)^(2/3);

rbc[t_] := Sqrt[rb[t]^2 + rc[t]^2 - 2 b c Cos[th[t]]]

eqn1 = rc''[t] -
rc[t] th'[t]^2 == -(kb/
rbc[t]^2) (1 - rb[t]^2/(Sin[th[t]]^2 rbc[t]^2));

eqn2 = rc[t] th''[t] +
2 rc'[t] th'[t] ==
-(kb/rbc[t]^2 rb[t]/(Sin[th[t]] rbc[t]) + ka/
rc[t]^2);

sol = NDSolve[{eqn1, eqn2, rc[0] == 100,
rc'[0] == 0, th[0] == 0.1,
th'[0] == 0}, {rc, th}, {t, 0, 10}]
ParametricPlot[
Evaluate[{rc[t] Cos[th[t]], rc[t] Sin[th[t]]} /.
Flatten[sol]], {t,
0, 10}]

7. Dear shawsend,

You correctly anticipated that those are polar coordinates I was working with and the graph I want is properly drawn in $\displaystyle (x,y)$ where x takes the value $\displaystyle r cos\theta$, and y takes the value, $\displaystyle r sin \theta$.

I much like the asymptotic approach you suggest (and will make it a point to follow it in future problems!).
I see that you have provided me with the code that you have employed in Mathematika to solve the equations. Much appreciated.

Albeit I do not have Mathematika on my computer, I will try and obtain the software and then work on it to get my graphic solution. Thanks for your cooperation, Shawsend, much appreciate it

8. I will try and stay in touch and notify you if my 'asymptotic' approach gives me any result. If not, I will ask for your help again!

9. Originally Posted by taureau20
Albeit I do not have Mathematika on my computer, I will try and obtain the software and then work on it to get my graphic solution. Thanks for your cooperation, Shawsend, much appreciate it
Hey, just go to your local university or just pick one at random and show up in the math department and start asking people there (again randomly) that you have a differential equation you can't solve and would like to run it on Mathematica. They love doing math there and would welcome the challenge.

Also, Mathematica (or other) in my opinion is indispensable in doing mathematics. You've got to find a way to get it and start learning how to use it.

10. Hi!

I have installed Mathematica on my system. Now, when I use the NDSolve function or anything whose output is supposed to be an image displayed, it doesn't display that image (could be a plot, etc.) rather it just says, -Graphics- and leaves it at that. Could you tell me how to get round it?
I write code on Mathematical 6 Kernel window; my OS is Vista.
Thanks!