%%This M-file will generate the phase portrait for the Inverse Harmonic
%%Oscillator Potential V=-1/2*Kx^2.

%%Total Energy E is constant:

clf;
hold on;
axis([-3 3 -3 3]);

for E=-4:.25:4

    if E >= 0
        x=[-5:0.01:-0.01,0.01:0.01:5];
        dx=gradient(x);
        v=sqrt(2*E-x.^2-1./x.^2);       
        dv=gradient(v);
        plot(x,v,'black');
        plot(x,-v,'black');
        
        for i=1:length(v)
            if imag(v(i))==0
                h=quiver(x(i),v(i),dx(i),dv(i)); 
                i=quiver(-x(i),v(i),dx(i),-dv(i));
                set(i,'LineWidth',3);
                set(i,'Color','black');
                set(h,'LineWidth',3);
                set(h,'MarkerSize',10000);
                set(h,'Color','black');
            break    
            end
        end 
    end
    if E <0
        xmin=sqrt(-2*E);
        x=(xmin+.00001:0.01:4*xmin);
        v=sqrt(2*E-x.^2-1./x.^2);
        plot(x,v);
        plot(-x,v);
        plot(x,-v);
        plot(-x,-v);
    end    

end

title('Phase Portrait for Potential V(x)=1