Energy Surfaces

* Türkçe’si için

In CR3BP, the third body (m3) is actually moving in a field which is also called effective potential, and the potential function (-U) gives this effective potential surface, which is shown in below;

As it is seen, L1, L2 and L3 are on the saddle in this surface, which means the motion are unstable in there. But L4 and L5 are on the peaks, which is actually dip in possitive U plot (+U), which means the motion are stable around there. The Matlab code of effective potantial is in below;

function pot_fun

mue=0.2; %---mass ratio

[X,Y] = meshgrid(-2:.01:2);
hold on
%potential function
 surf(X,Y,Z,'FaceColor',[1 0.99 1],'EdgeColor','none') %[1 0.99 1]
 camlight left; lighting phong
zlim([-3.8 -2.83])
alpha(0.5)% --- for transpancy
view([60 80])
grid on

A horizontal cut of the plot of the effective potential (2U-C) gives Hill’s region, where the third body can only move in there at a given energy level C (Jacobi constant). Hill’s region is also called as zero velocity curves (or surfece in 3D), which is also obtined from Jacobi integral function (V²=2U-C), making velocity (V) zero.

In this figure, the three realms; namely the inner region, the neck region and the outher region that are connected by the neck region. For other values of the energy (C), one or more of these realms may be prohibited due to conservation of energy; such as, the necks may close off. The Matlab code of Hill’s region (or can be said as zero velocity curves in 2D or or zero velocity surfaces in 3Dis below. Copy-paste the code in your matlab work folder, then write in the commant line as follow;

  • for 2 dimensional plots; zrvs(‘2d’, μ, C)
  • for 3 dimensional plots; zrvs(‘3d’, μ, C, transparancy)

Where μ is mass ratio, C is jacobi constant. Transparancy is set from 0 to 1; 0 for clear, 1 for opaque (so set as 0.5). Such as zrvs(‘3d’,0.01215,3.16, 0.5).

function zrvs(dim,mue,C,tra)
% plot the zero velocity surface
%dim : '2d' or '3d'
%mue : mass ratio
%C : jacobi cost.
%tra : transparancy
%exp: zrvs('3d',0.01215,3.16,0.6)
hold on
switch dim 
 case '3d'
[X,Y,Z] = meshgrid(d,d,d);

r1=((X + mue).^2+Y.^2+Z.^2).^(1/2);
r2=((X + mue-1).^2+Y.^2+Z.^2).^(1/2);
 V=2*[(1/2)*(X.^2 + Y.^2)+(1-mue)./r1+mue./r2];
% you may change the color in below 
h = patch(isosurface(X,Y,Z,V,C),'FaceColor','white','EdgeColor','none');
camlight; lighting phong
axis off
axis equal
 case '2d'
syms X Y
axis equal
grid on;

Zero Velocity Surfaces

******** ho ho ho *****


Leave a Reply

Fill in your details below or click an icon to log in: Logo

You are commenting using your account. Log Out /  Change )

Google+ photo

You are commenting using your Google+ account. Log Out /  Change )

Twitter picture

You are commenting using your Twitter account. Log Out /  Change )

Facebook photo

You are commenting using your Facebook account. Log Out /  Change )

Connecting to %s