function [V] = potentia(l, r)
% Schwarzschild space-like geodesic
% Potential V for parametr l = L/M (L denotes angular momentum
% about an axis normal to the invariant plane, M is mass)
% r is in mass units.
V = (1-2./r) .*(l^2./r.^2-1);