Equation of Motion For a Specific Gravitation Question: Is harmonic motion involved? If yes/no, how and why? This is a question that I have thought about myself, but I'm struggling a bit regarding how to answer it.
This is my question:
Imagine we have a universe with 3 bodies.  Out of these 3, 2 are identical and massive, and 1 is small and insignificant.  The 2 identical celestial bodies are literally identical in every way possible -- density, mass, volume, etc.  They are also fixed in space and unable to move.  If their centers are perfectly aligned vertically (although this doesn't really matter since we're in space) and separated by distance $d$, then their Lagrangian point will be at $\frac{d}{2}$ since the bodies are identical.
Now, let's make another axis from this Lagrangian point that is orthogonal to the vertical axis (which $d$ falls on), and call it the horizontal axis.  In other words, let's just say this produces some $xy$-plane.
If our third body, which we may assume to be a perfect sphere and insignificant mass, is placed on this horizontal axis $x$ units of distance away from the origin of this $xy$-plane (which is the Lagrangian point), how will it's motion behave?  Will it be harmonic and form some sort of space-pendulum or will it eventually fall into the Lagrangian point?
And in either case, what will be its equation of motion?  If it falls into the Lagrangian point, depending on $x$ and other parameters, how long will it take to fall in?
I have included a diagram below.

EDITED:
So basically, you can start out with vector summation of the forces and get to here:
$$2\mathbf{F_{G}}\cos(\theta)=m \cdot \mathbf{a}$$
$\theta$ is basically the angle above the horizontal axis.
$$ \frac{2GMm}{\frac{d^2}{4}+\Big(x(t)\Big)^2} \cdot \frac{x(t)}{\sqrt{\frac{d^2}{4}+\Big(x(t)\Big)^2}} = m \cdot a(t)$$
$$\frac{2GM\cdot x(t)}{\bigg (\frac{d^2}{4}+\Big(x(t)\Big)^2 \bigg) \sqrt{\frac{d^2}{4}+\Big(x(t)\Big)^2}} = a(t)$$
But now I get kind of stuck.
 A: You can tell immediately that it’s not harmonic motion: that requires a restoring force proportional to displacement, hence rising with distance from the origin. Here, the restoring force gets smaller as it moves farther away: larger oscillations will feel less force, be restored more slowly, and take longer. 
That’s not the behavior of a harmonic oscillator. 
A: Okay, so I think I've figured out the answer myself, haha.  Just reiterating the work already done above.
$$2\mathbf{F_{G}}\cos(\theta)=m \cdot \mathbf{a}$$
$$ \frac{2GMm}{\frac{d^2}{4}+\Big(x(t)\Big)^2} \cdot \frac{x(t)}{\sqrt{\frac{d^2}{4}+\Big(x(t)\Big)^2}} = m \cdot a(t)$$
$$\frac{2GM\cdot x(t)}{\bigg (\frac{d^2}{4}+\Big(x(t)\Big)^2 \bigg) \sqrt{\frac{d^2}{4}+\Big(x(t)\Big)^2}} = a(t)$$
$$\left(\frac{2GM\cdot x}{\bigg (\frac{d^2}{4}+x^2 \bigg) \sqrt{\frac{d^2}{4}+x^2}}  \right) \cdot v\,dt = \left ( \frac{dv}{dt} \right) v\,dt$$
$$\left(\frac{2GM\cdot x}{\bigg (\frac{d^2}{4}+x^2 \bigg) \sqrt{\frac{d^2}{4}+x^2}}  \right) \cdot \left ( \frac{dx}{dt} \right)\,dt = v \, dv$$
$$\frac{2GM\cdot x}{\bigg (\frac{d^2}{4}+x^2 \bigg) \sqrt{\frac{d^2}{4}+x^2}} \,dx = v\,dv$$
$$\int_{x_0}^{0}\frac{2GM\cdot x}{\bigg (\frac{d^2}{4}+x^2 \bigg) \sqrt{\frac{d^2}{4}+x^2}} \,dx = \int_{v_0}^{v_f,L}v\,dv$$
$$-\frac{2GM}{\sqrt{x^2+\frac{d^2}{4}}} \Bigg |_{x_0}^0=\frac{v^2}{2}  \Bigg |_{v_0}^{v_{f,L}}$$
$$-\frac{4GM}{d}+\frac{2GM}{\sqrt{x_{0}^{2}+\frac{d^2}{4}}}=\frac{v_{f,L}^{2}}{2}$$
Assuming that $v_0$ is zero.
$$v_{f,L}=\sqrt{\frac{4GM}{\sqrt{x_{0}^{2}+\frac{d^2}{4}}}-\frac{8GM}{d}}$$
Where $v_{f,L}$ is the final velocity when it reaches the Lagrangian point.
Now, it will pass the Lagrangian point and go some distance $u$.  As a result, we notice the same situation but it is now reflected on the vertical axis.
$$-2\mathbf{F_{G}}\cos(\alpha)=m \cdot \mathbf{a}$$
$$ \frac{2GMm}{\frac{d^2}{4}+u^2} \cdot \frac{u}{\sqrt{\frac{d^2}{4}+u^2}} = m \cdot a$$
$$\frac{2GM\cdot u}{\bigg (\frac{d^2}{4}+u^2 \bigg) \sqrt{\frac{d^2}{4}+u^2}} = a$$
Since we are getting the same thing, we can conclude that by going through the same integration method, we will get:
$$-\left(-\frac{2GM}{\sqrt{u^2+\frac{d^2}{4}}}\right) \Bigg |_{0}^{u_f}=\frac{v^2}{2}  \Bigg |_{v_{f,L}}^{v_f}$$
$$\frac{2GM}{\sqrt{u_{f}^{2}+\frac{d^2}{4}}}-\frac{4GM}{d}=\frac{v_{f}^{2}}{2}-\frac{v_{f,L}^{2}}{2}$$
$$\frac{2GM}{\sqrt{u_{f}^{2}+\frac{d^2}{4}}}-\frac{4GM}{d}=\frac{v_{f}^{2}}{2}+\frac{4GM}{d}-\frac{2GM}{\sqrt{x_{0}^{2}+\frac{d^2}{4}}}$$
$$v_f=\sqrt{\frac{4GM}{\sqrt{u_{f}^{2}+\frac{d^2}{4}}}+\frac{4GM}{\sqrt{x_{0}^{2}+\frac{d^2}{4}}}-\frac{16GM}{d}}$$
Since we want to calculate it's ultimate final speed, $u_f \to \infty$ such that it's deceleration will approach zero.
$$v_f=\sqrt{\frac{4GM}{\sqrt{x_{0}^{2}+\frac{d^2}{4}}}-\frac{16GM}{d}}$$
And with some simplification:
$$v_f = 2\sqrt{\frac{2GM\, \left (d-2\sqrt{4x_{0}^{2}+d^2} \right)}{d\,\sqrt{4x_{0}^{2}+d^2}}}$$
Is this correct?
A: Well, first let us have an intuitive idea as to what is going on in the system. 
Let us fix some co-ordinates :
The Lagrange point might be the Origin $(0,0)$. The massive masses are then located at $(+\frac d2,0)$ and $(-\frac d2,0)$ respectively. The initial position of the small mass be $(-A,0)$.
Now, the body is gravitationally attracted by both the masses, but only their cosine components act along the line joining the body and Lagrange point. The sine components cancel out. Net effect is an attractive force towards the Lagrange point. When the body is moving towards the Lagrange point, the distance between the bodies is decreasing and so the force is increasing, and thus the acceleration. At exactly the Lagrange point, the forces balance out. But due to inertia, the body moves out of the point, and is now attracted in exactly the opposite direction to the former case. The whole situation is now reversed, and the body decelerates till it reaches $(-A,0)$. At that point, the velocity is $0$ and the above steps are repeated again in exactly the same fashion except the direction reversed. 
Thus, we have an oscillatory motion with mean position $(0,0)$ (or the Lagrange point if you like) and Amplitude A. 
Now, we want the maths to reflect this.
We will use the following variables in our work :
$$M=Mass\ of\ massive\ bodies$$
$$m=Mass\ of\ small\ boy$$
$$A=Amplitude\ of\ motion$$
$$d=Distance\ between\ massive\ bodies$$
$$\theta=Angle\ between\ x-axis\ and\ the\ line\ joining\ small\ body\ to\ massive$$ $$body\ when\ former\ is\ at\ a\ distance\ x\ from\ the\ origin$$
$$\theta_0=Angle\ between\ x-axis\ and\ the\ line\ joining\ small\ body\ to\ massive$$ $$body\ when\ former\ is\ at\ a\ distance\ A\ from\ the\ origin$$
We start with Newton's second law of motion : 
$$m\frac{d^2x}{d\ t^2} = \frac{-2GMm}{x^2+\frac d2^2}cos\theta$$
This is the required equation of motion
Eliminating $\theta$, doing some algebra, writing $\frac{d^2x}{d\ t^2}$ as $v\frac{dv}{dx}$ and solving the equation gives :
$$\frac{dx}{dt}=-\sqrt{4GM(\frac 1{\sqrt{x^2+(\frac d2)^2}}-\frac 1{\sqrt{A^2+(\frac d2)^2}})}$$
We now bring back the angles. We have $x=\frac 12d\ cot \theta$, $\sqrt{x^2+(\frac d2)^2}=\frac d2 cosec\ \theta$ and $\sqrt{A^2+(\frac d2)^2}=\frac d2 cosec\ \theta_0$.
Putting everything in, and using substitutions $\alpha_0=\frac{\pi}{2}-\theta_0$ and a similar one for $\theta$. We also use $sin \frac{\alpha}{2}=k\ sin\phi$ where $k=sin \frac{\alpha_0}{2}$.
Finally we get the time period of oscillation $T$ as : 
$$T=\sqrt{\frac{2d^3}{GM}}\int^\frac{\pi}{2}_0 \frac{d\phi}{1-4k^2sin^\phi(1-k^2sin^2\phi)\sqrt{1-k^2sin^2\phi}}$$
This is an elliptic integral and cannot be expressed as an expression in elementary functions.
If the oscillations are small, then $k=0$ and thus we have after solving:
$$T=\pi \sqrt{\frac{d^3}{2GM}}$$
This is the required time period of small oscillations
If any point is unclear, feel free to ask. I will be happy to help.
Cheers!
