CODEWITHSUNDEEP
matlabfor-loopprimes
Vikram Kaushik
Vikram Kaushik

Reputation: 51

Program to find prime numbers in MATLAB

I wrote some code to display the prime numbers between 2 and a user-chosen number, following the pseudocode on wikipedia. I am not sure why this does not work, as I have the increments correct following Erastothenes' Sieve. Please help me.

I have tried changing the bounds but this did not work.

There are no errors, but it returns the wrong output. If I enter 10, it returns 2, 3, 4, 5, 6, 7, 8, 9, 10.

n=input("Enter an upper limit: ");
nums= 2:n;
p=2;
for i = p:sqrt(n)
    for j = (i^2):i:sqrt(n)
        nums(j) = 0;
    end
end
for k = 1:n-1
    if nums(k) ~= 0
        disp(nums(k))
    end
end

Upvotes: 2

Views: 5219

Answers (3)

Salma Bendary
Salma Bendary

Reputation: 1

n=input('Enter a number:')
F=factor(n);
if F==n
    disp('Prime number.')
else
    disp('Not a prime number.')
end

Upvotes: -2

Wolfie
Wolfie

Reputation: 30101

You can use the primes function in MATLAB for this

N = 10;        % upper limit
p = primes(N); % List of all primes up to (and including) N

With one step less automation, you could use another in-built isprime

p = 1:N;                  % List of all numbers up to N
p( ~isprime( p ) ) = [];  % Remove non-primes

Finally without using built-ins, we can address your code! I assume you're referring to this pseudocode for the Sieve of Eratosthenes on Wikipedia.

 Input: an integer n > 1.

 Let A be an array of Boolean values, indexed by integers 2 to n,
 initially all set to true.

 for i = 2, 3, 4, ..., not exceeding √n:
   if A[i] is true:
     for j = i2, i2+i, i2+2i, i2+3i, ..., not exceeding n:
       A[j] := false.

 Output: all i such that A[i] is true.

I'll follow it step by step, pointing out differences to your code:

n = 10;
A = [false; true(n-1,1)];   % array of true Booleans, first element (1) is not prime
% But I've included a first element to make indexing easier.
% In your code, you were using index 'i' which was incorrect, as your array started at 2.
% Two options: (1) take my approach and pad the array
%              (2) take your approach and using indices i-1 and j-1
for ii = 2:sqrt(n)
    if A(ii) == true        % YOU WERE MISSING THIS STEP!
        for jj = ii^2:ii:n   % YOU ONLY LOOPED UNTIL SQRT(n)!
            A(jj) = false;
        end
    end
end

p = find(A);
disp(p)

This outputs the expected values.

Note that, at the end of the manual looping method, A is equivalent to isprime(1:n), mirroring my earlier suggestions.

Upvotes: 5

obchardon
obchardon

Reputation: 10792

There is two mistakes in your code:

  1. The multiple should be check until n and not sqrt(n)

  2. Since your nums vector start with 2 and not 1, if you want to access the right value you need to use nums(j-1) = 0

So:

n=100
nums= 2:n;
p=2;
for i = p:sqrt(n)
    for j = (i^2):i:n
        nums(j-1) = 0;
    end
end
for k = 1:n-1
    if nums(k) ~= 0
        disp(nums(k))
    end
end

Noticed that you can skip one for loop using a modulo, it's probably not faster than the previous solution since this code create a logical index that include each prime that already been found.

n   = 100
nums= 2:n;

for i = 2:sqrt(n) 
    nums(mod(nums,i)==0 & nums != i) = [];
end

nums.'

I simply delete the value in nums that can be divided by x but not x.

Upvotes: 2

Related Questions