for p
    for q
        B=union(B,p/q);
end

function Y=F(n)
    Y=%e.^(-0.1*n).*sin((2*n+1)*%pi/2);
endfunction
X=1:100;
Z=F(X);
plot(X,Z,'b*')