in reply to Re^5: OT: Finding Factor Closest To Square Root
in thread OT: Finding Factor Closest To Square Root

The original question started "given N and N's prime factorization", so it seems inappropriate to take into account the time taken to factorize.

In this case though, note that 988041964007 = 7 x 11087 x 12731023 - you are losing accuracy by throwing out the BigInts.

In general if the factorization of N is taking longer than it takes to do sqrt(N) trial divisions then the factorization algorithm is pretty bad, since "trial division up to the square root" is the classical example of a slow algorithm.

Many of the high precision maths libraries out there do factorization - for example the pari library (Math::Pari) has a factorint() function:

  use Math::Pari qw/ factorint /;
  my $n = Math::Pari->new("988041964007");
  my $f = factorint($n);
  my $count = $#{ $f->[0] };
  print join(" * ", map {
    my($prime, $power) = ($f->[0][$_], $f->[1][$_]);
    $power > 1 ? "$prime ^ $power" : $prime
  } 0 .. $count), "\n";
7 * 11087 * 12731023

[download]

Note that this algorithm does not guarantee primality: for large factors you need to verify with an isprime() check.

Hugo

Replies are listed 'Best First'.
Re^7: OT: Finding Factor Closest To Square Root
by BrowserUk (Patriarch) on Feb 21, 2005 at 17:56 UTC
    so it seems inappropriate to take into account the time taken to factorize

    Inappropriate maybe, but it's bloody boring waiting for it :)

    The other part of the argument against doing the factorising, is I haven't thought of, nor seen, a way to use them, that arrives at the result more quickly.

    As the factor you are looking for can be any combination (small c), of the prime factors of N, I haven't seen any algorithm that doesn't require an exhaustive search of those combinations plus trial division to determine if you've found the answer. The result is that the descending, linear search is easier to code and runs more quickly.

    I have produced the results for all the integers 1-2 million, and 3-4 million, 10-11 million and 100-101 million and then atttempted to find some sort of pattern to the results. There may be one there, but my meagre memory of numerical analysis is not enough to devine it.

    28.x% of of the low ranges are primes. This falls of very slowly as the scale of the range increases, but it's only dropped to 27.76% by the time you get the the 100-101 million range. This is worst case for the linear search. If the factorisation ran more quickly, then it would allow you to skip the search in those cases.

    Whether that would result in an overall speed up, given that you would have to factor all the numbers and (so far) fall back to a linear search for those Ns that are not primes, I think is dubious.

    As the range gets higher, the frequency of the primes gets less, so the benefit of performing the factorisation reduces as the cost of the linear search increases, but also the cost of the factorisation increase and much more rapidly that the cost of the linear search.

    I will install Math::Pari. I think it would be interesting to see if there is any pattern to the prime factors of the factor.

    This has probably already been analysed to death sometime, but if it has, I'm not typing in the right keywords to locate it.

    Examine what is said, not who speaks.
    Silence betokens consent.
    Love the truth but pardon error.
[reply]
      If you install Math::Pari, try the function sigma($N,0). That computes the total number of divisors of $N. It requires factorization to compute, but Math::Pari has very fast factorization algorithms.

      If N = (p1^a1)*(p2^a2)*...(pR^aR) then sigma(N,0) = (a1+1)*(a2+1)*...(aR + 1)

      This number is usually very small, even for huge N.

      If sigma($N,0) is reasonable, you can call divisors($N), which returns a list of all the divisors in order. A simple binary search will find the desired number.

      This should be incredibly faster than brute-force countdown searches for large N.

[reply]
        but Math::Pari has very fast factorization algorithms.

        Indeed it does--blindly fast. And that completely changes the perspective on using the factors in some sort of solution. I'm still finding my way around the huge set of functions in Math::Pari. If POD allowed tables, the looong list of functions could be made more compact and easier to navigate. Even so, I doubt I would have picked out sigma as being useful in this context from the description: 

        sigma(x,{k = 1})

sum of the k^{th} powers of the positive divisors 
of |x|. x must be of type integer.

The library syntax is sumdiv(x) ( = sigma(x)) 
or gsumdivk(x,k) ( = sigma(x,k)), 
where k is a C long integer.

        [download]

        Mmm'kay! :)

        Do you have the time/inclination to post code for your solution?

        Examine what is said, not who speaks.
        Silence betokens consent.
        Love the truth but pardon error.
[reply]
[d/l]

In Section Seekers of Perl Wisdom