function feq=orbY(yi,C,xp) global GEO wt mu [GU U]=Gpot(0,yi,0,GEO,mu); feq=1/2*wt^2*yi^2+U-C-1/2*xp^2; end