Why acos fails?

10 views (last 30 days)
OldCar
OldCar on 25 Oct 2015
Answered: Star Strider on 25 Oct 2015
I have a Space Mechanics problem: I have a position vector and a velocity vector of a circular orbit with i=0. I have a problem with the function acos. In fact acos(something) is not 0 when something is 1. the code is
mi=398600.4;
x=7000;
y=0;
z=0;
u=0;
v=sqrt(mi/7000) ;
w=0;
r=sqrt(x^2+y^2+z^2);
V=sqrt((u^2+v^2+w^2));
a=mi*r/(2*mi-r*V^2);
e=sqrt((((V^2-mi/r)*x-(x*u+y*v+z*w)*u)/mi)^2+...
(((V^2-mi/r)*y-(x*u+y*v+z*w)*v)/mi)^2+...
(((V^2-mi/r)*z-(x*u+y*v+z*w)*w)/mi)^2);
if e>=1
error('l''orbita non è una ellisse')
end
E=atan2(e*sqrt(1/mi*a)*(x*u+y*v+z*w),e*(1-r/a));
%tau=-(a^3/mi)^(1/3)*(E-e*sin(E)); %tempo dal passaggio al perigeo
ni=2*atan(sqrt((1+e)/(1-e))*tan(E/2));
H=sqrt(mi*a*(1-e^2));
i=acos((x*v-y*u)/H);
%if i<10^-5
% error('l''orbita è poco inclinata e bisogna corregere con un trick')
%end
raan=atan2((y*w-z*v)/(H*sin(i)),(x*w-z*u)/(H*sin(i)));
ot=atan2(z/(r*sin(i)),y/r*sin(raan)+x/r*cos(raan)); %omega+teta
omegap=ot-ni;
(x*v-y*u)/H is equal to 1, but acos((x*v-y*u)/H)=0,00000000000000 + 2,10734242554470e-08i Why????

Accepted Answer

Star Strider
Star Strider on 25 Oct 2015
(x*v-y*u)/H is equal to 1
It is probably a bit greater than 1, because of floating-point approximation error. Thus the complex result.

More Answers (0)

Categories

Find more on Programming in Help Center and File Exchange

Community Treasure Hunt

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!