www.gusucode.com > MATLAB2008应用程序接口编程技术源码程序 > MATLAB2008应用程序接口编程技术源码程序/code/第6章/6.6/myprimes.m

    function p = myprimes(n)
if length(n)~=1
error(‘N must be a scalar’);
end
if n<2
    p = zeros(1,0);
    return;
end
p = 1:2:n;
q = length(p);
p(1) = 2;
for k = 3:2:sqrt(n)
    if p((k+1)/2)
        p(((k*k+1)/2:k:q = 0;
    end
end
p = (p(p>0));