Computation of the Lagrangian points
Consider an object moving in the plane and subject to a set of forces. Our goal in this application consists in determining the equilibrium positions of this object.
- Question :
- Write the code to represent the following kinds of forces :
- an elastic force (with zero natural length)
- a centrifugal force
$$f_c : \left[ \begin{array}{c} x \\\ y \end{array} \right] \rightarrow \left[ \begin{array}{c} k (x-x_0) \\\ k (y-y_0) \end{array} \right] $$
- a gravitational force
$$f_g : \left[ \begin{array}{c} x \\\ y \end{array} \right] \rightarrow \left[ \begin{array}{c} -k \cdot \frac{x-x_0}{((x-x_0)^2+(y-y_0)^2)^{3/2}} \\\ -k \cdot \frac{y-y_0}{((x-x_0)^2+(y-y_0)^2)^{3/2}} \end{array} \right] $$
Remarks:
- It is strongly recommended to use functional programming techniques to represent the forces and their Jacobian matrices.
- For each of these cases, it should be possible to parameterize the force by a constant $k$ symbolizing its intensity, as well as a central point $(x_0,y_0)$ from which the force is issued.
- Question :
- Use the Newton-Raphson method to obtain the equilibrium points in the following case :
- Two gravitational forces with respective coefficients $1$ (resp. $0.01$) originating from $[0,0]$ (resp. $[1,0]$);
- And a centrifugal force centered on the barycenter of the two masses, with coefficient $1$.
As a means of verification, and with the conditions above, at $U = [1.5, 0]$, one gets the following values for the total force and its Jacobian :
$f(U) = \left[ \begin{array}{c} 1.00565457 \\\ 0 \end{array} \right]$ and $df(U) = \left[ \begin{array}{cc} 1.75259259 & 0. \\\ 0. & 0.6237037 \end{array}\right]$.