I am trying to find the 5 Lagrange points of the Three-Body system by plotting the given potential function in Matlab. The only thing is that I'm not very good at programming. Any assistance would be greatly appreciated. What I want to know is how to make this code give me a decent contour plot:
function Lagrange(a)
x = ( -10000: 1 : 10000);
y = ( -10000: 1 : 10000);
Potential = zeros(length(x));
for i = 1: length(x)
for j = 1 : length(y)
Potential(i,j) = ( 1 - a ) / sqrt( ( x(i) - a )^2 + y(j)^2) + a / sqrt( ( x(i) + 1 - a )^2 + y(j)^2 ) + ( x(i)^2 + y(j)^2 ) / 2 ;
end
j = 1;
end
contour(Potential);
xlabel('X axis');
ylabel('Y axis');
zlabel('Z axis');
The way the three-body problem is set up, the distance coordinates are normalized to a. Thus, you should pick x and y to be more like:
x = linspace(-1.5, 1.5, 1000);
y = linspace(-1.5, 1.5, 1000);
For the contour plot, you can use meshgrid, which allows you to avoid that for loop and plot a little easier:
[X, Y] = meshgrid(x, y);
For the potential, try plotting 2U - this is called the Jacobi constant and is a bit more informative.
U = (1-a)./sqrt(Y.^2 + (X + a).^2) + ...
a./sqrt(Y.^2 + (X + a - 1).^2) + ...
0.5*(Y.^2 + X.^2);
Z = 2*U;
Finally, you'll need contours. You'll want to tweak these for your plot, but I used something like
c = [2.988:0.05:3.1, 3.2:0.2:5];
for the Earth-Moon system. Now, to plot, simply use contourf as follows:
figure
contourf(X, Y, Z, c)
colorbar
Also note that you can solve for the Lagrange points themselves analytically using the equations of motion - you may consider plotting these too, since the contours will only converge on the points but will never hit them.