why does factor(1) return 1? shouldn't it return null?

I'm using factor() to find the prime factors of a range of numbers.

As far as I understand it (which may be very limited as I am new to this!), factor(n) should accept a value of n=1 as 1 is a "real, nonnegative integer scalar", but should only return prime factors.

I believe that 1 is not regarded as a prime.

Given this, shouldn't the function return an empty/null row vector?

Thanks.

Answer 1

[dpb]

help factor

 factor - Prime factors

    Syntax
      f = factor(n)

    Input Arguments
      n - Input value
        real, nonnegative integer scalar

    Examples
      openExample('matlab/PrimeFactorsofDoubleIntegerValueExample')
      openExample('matlab/PrimeFactorsofUnsignedIntegerValueExample')

    See also isprime, primes

    Introduced in MATLAB before R2006a
    Documentation for factor
       doc factor

    Other uses of factor

       sym/factor

From factor, _"F = factor(x) returns all irreducible factors of x in vector F. If x is an integer, factor returns the prime factorization of x. "_

The two sentences are conflicting; the first description is what is actually the implementation while the second is true only for integers >=2. I think the doc should be clarified.

For your purposes, use

v=0:10;         % a range of integers
pfactors=arrayfun(@primes,v,'uni',0)

pfactors = 1×11 cell array
    {1×0 double}    {1×0 double}    {[2]}    {[2 3]}    {[2 3]}    {[2 3 5]}    {[2 3 5]}    {[2 3 5 7]}    {[2 3 5 7]}    {[2 3 5 7]}    {[2 3 5 7]}

Oooops! The above doesn't do what you want/need, either, since primes returns all primes less than the number, not the factors in the number.

You'll have to make up your own version to fix up the shortcut factor takes to maintain multiplication rule.

pfactors=arrayfun(@factors,v,'uni',0)

pfactors = 1×11 cell array
    {0×0 double}    {0×0 double}    {[2]}    {[3]}    {[2 2]}    {[5]}    {[2 3]}    {[7]}    {[2 2 2]}    {[3 3]}    {[2 5]}

function f = factors(n)
%FACTOR Prime factors.
%   FACTOR(N) returns a vector containing the prime factors of N
%   Unlike MATHWORKS version 0, 1 return the empty result, not the value

  f=factor(n);      % let the builtin do its thing first
  if n<2            % fixup the nonmath definition ones
    f=[];
    return
  end
end

Answer 2

[Steven Lord]

Technically, yes. If you want prod(factor(x)) to return x and the output of factor to only include prime numbers, having factor(1) return the empty matrix would be consistent.

But there would still be an inconsistency -- what should factor(0) return? 0 isn't prime. Having it error feels like it could be technically correct but practically annoying. [The fundamental theorem of arithmetic covers integer values greater than 1, so technically it doesn't apply to either of those cases.]

Doing a little archaeology and going back to the earliest version of factor.m that I have access to (over 30 years old) I believe factor(1) would actually error (with an incorrect error message.)

Almost 27 years ago, factor.m was modified to just return the input if that input is non-negative and is less than 4. That code is still present in factor.m today. This optimization works for n = 3 and n = 2. It's this that causes factor to return 1 for n = 1 and 0 for n = 0. [Negative values are checked for and cause factor to error out on an earlier line of code.]

dbtype 29:34 factor.m

29    if n < 4
30       f = floor(n); 
31       return
32    else
33       f = [];
34    end

We could change how factor.m behaves in the n = 1 case (and in the n = 0 case as well.) But given how long this behavior has been in the product I'd be hesitant to change it without a good use case and an analysis of the potential impact of the incompatibility.

[John D'Errico]

@Steven Lord - I agree, that factor as applied to 0 or 1 should remain as they are. But I do think the help should reflect that.

[Steven Lord]

That's fair. They are edge cases and as such it's reasonable to ask that they be documented.

Answer 3

[John D'Errico]

Yes, I agree, in that

f1 = factor(1)

f1 = 1

could arguably return empty. 1 is not indeed prime, and so the result is inconsistent with the documentation.

I'd also note that

f0 = factor(0)

f0 = 0

fails to return a set of prime factors. Again, 0 is not prime.

However, you could also argue that 1 is indeed a "factor" of the number 1. In fact, it is the only such number that divides 1, and is not itself reducible. It is just that 1 is technically not a prime factor of 1, which lacks any truly prime factors at all. It gets even more difficult to decide what should be returned for the factor(s) of 0. But the help for factor allows any non-negative arguement.

It is also true that the product of the vector returned by factor reproduces the input. And so we see that

prod(f0)

ans = 0

prod(f1)

ans = 1

do indeed reproduce the original number as factored. And this could be a useful property of factor, relied upon in legacy code, going back as much as 40 years now. Yes, factor has been doing this for a long time. And any change to factor would potentially cause legacy code to fail were they to arbitrarily change the code. That would mean any change to factor would need to happen with a long lead time. They would probably keep reminding us for multiple years to change any code that relies on factor(1).

Personally, I think the argument about what should be returned as the factors of either 0 or 1 as being more a pedantic thing, up there with the number of angels that fit on the head of a pin. I don't think it merits a need to change the code, because of the issues that might potentially cause.

Better I think is to just introduce a documentation change, indicating that factor(1) will return 1 as a special case (the same applies to factor(0)). That can happen quite quickly, probably seeing such a change reflected in R2027, or even sooner. A change to the documentation IS DEFINITELY a good idea here. I would recommend this be introduced to support. I'll post it as a request today.

[dpb]

I don't suppose the archeology uncovered any other rationale notes besides the shortcut, but I presumed that keeping the product of the results of factor() reproducing the inputs was the motivation. I agree also on only a documentation update as noted before. I hadn't ever thought of it specifically before as I don't know that had ever had reason to inspect factor(1) before...

[John D'Errico]

Funny thing is, I don't think it ever bothered me either. Returning a factorization that contains 1 (or 0) just seems to make more sense. Do I really want factor to generate an error message for a boundary case like this? Not really. And even returning an empty vector for factor(1) or factor(0), that now would fail the requirement that prod(factor(n))==N. So I sort of like the return of a non-prime for those boundary cases.

And that means, in my eyes, really the only right thing to do is to update the documentation. Just have it say that well, yeah, factor does not return primes for these special cases, but what can we do there?