function valf = F(xr,xp,xs,zeta,lc)
%
% Return the value of the function F defined by
%              F:= (xs-zeta)+ beta *(3*(z.r)^2/r^5-1./r^3);
%
global beta
global r0
global Z
global g
lgeta = beta/(g*r0^4);

r    = sqrt(xr^2+xp^2+(xs+lc)^2);
zr   = Z(1)*xr+Z(2)*xp+Z(3)*(xs+lc);
valf = (xs-zeta)+ lgeta*(3*(zr^2)/(r^5) - 1/r^3);
   
