Looking for a sieving program

[m_f_h]

Originally Posted by jasong

Since there doesn't seem to be a program that fits my needs, I guess I might as well reveal what I wanted to do. I wanted to find the lowest 114-digit number with two brilliant factors. A brilliant number is a number where all the prime factors have the same number of digits. So far, digits 1 through 113 are spoken for.

What you say is equivalent to :
1) size(N) = 114
2) N = B1*B2*F3
3) B1 = F11 * F12 *...* F1r with all F1i of the same size
4) B2 = F21 * F22 *...* F2s with all F2j of the same size ?

I think it should read:
2') N = B1*B2
3') size(B1) = size(B2)
4') B1 and B2 prime

Jason, I suggest you to use PARI rather than BASH to tackle the problem.
PARI is free and quite fast. e.g. to get a list of candidates, you can do

n=10^113-1;
while( n+=2, if( matsize(factor( n, 2^62))[1]==1, print( n )))

This will print out about 10 values/second (on a PC having 2 mprime in background...) for numbers having no factor <2^62.
In one day you have a million candidates...
To trial factor those, you can use mprime or so.
(But I suppose some big crunchers are already somewhat ahead on that way...)

[Jens K Andersen]

Originally Posted by jasong

I might as well reveal what I wanted to do. I wanted to find the lowest 114-digit number with two brilliant factors.

Saying what you actually want is usually a good idea. You may get much better help and others don't spend time on irrelevant answers. I have already sieved this to 3*10^11 last year. Candidates for this and other brilliant numbers are in http://mersenneforum.org/showpost.ph...4&postcount=40. The used sieve is not public.

[m_f_h]

Originally Posted by Jens K Andersen

I have already sieved this to 3*10^11 last year. Candidates for this and other brilliant numbers are in http://mersenneforum.org/showpost.ph...4&postcount=40.

With my primitive PARI 1-liner you get about 10 candidates/second which have no factor < 2^62 = 4.6e18 (on a 64 bit machine).
The density of them is about 4% i.e. 10 candidates = 250 k's eliminated so you get all k<5e5 in about 2000 sec, i.e. less than one hour.
So this is really the trivial part of the problem; I suppose it's somewhat different with stage 2 and 3 (ECM, NFS)...

[m_f_h]

However, I wonder if there could not be a more intelligent method to find minimal brilliant numbers.
AFAICS, the problem consists, for given N (e.g. N=10^113), in finding a prime p in [ sqrt(N/10), sqrt(N) ] such that nextprime( N/p ) - N/p is as small as possible; in particular less than one, i.e. otherwise stated, such that N/p is very close to ceiling( N/p ) which must be a prime. (Very close means of the order of 1/sqrt( N )).
Couldn't this reasoning allow for quick elimination of prime factor candidates ?

Also, can one accelerate the trial factoring process by using the fact that we're only interested in prime factors in a "small" (in terms of orders of magnitude.. say better "specific") interval (sqrt(N/10) ... sqrt(N)) ?
AFAIK, usual factoring methods don't use that.

Other "brainstorming" ideas : write the factor p starting from its highest bit (which we know (almost?) since our interval is of ("relative") size sqrt(10) so there's less than a factor 2 to each side when we're in the (geometric) middle), and then do a "binary path search" from high bits downward to low bits, by eliminating a branch as soon as possible (maybe using methods like minimax, alpha-beta pruning etc.)
(Well, there's quite a couple of bits to go downto, but...)

[Jens K Andersen]

Originally Posted by m_f_h

n=10^113-1;
while( n+=2, if( matsize(factor( n, 2^62))[1]==1, print( n )))

This will print out about 10 values/second (on a PC having 2 mprime in background...) for numbers having no factor <2^62.

Factoring to 2^62 in less than a second? Not likely (except algorithms for special large forms with known factor properties).

My PARI/GP 2.3.1 says:

? ?factor
factor(x,{lim}): factorization of x. lim is optional and can be set whenever x is of (possibly recursive) rational type. If lim is set return partial factorization, using primes up to lim (up to primelimit if lim=0)

On my 32-bit Windows, factor(n,2^31-1) is the highest allowed:

? factor(10^113+33963,2^62)
*** integer too big: factor(10^113+33963,2^62)
^-----
? factor(10^113+33963,2^31)
*** integer too big: factor(10^113+33963,2^31)
^-----
? factor(10^113+33963,2^31-1)
%48 =
[10000000000000000000000000000000000000000000000000000000000
0000000000000000000000000000000000000000000000000033963 1]

But it doesn't actually factor to 2^31-1. It misses this factor from my sieve:

? (10^113+33963)%1983963043
%49 = 0

My default primelimit is 500000:

? default(primelimit)
%50 = 500000

I guess that if x > primelimit, then factor(n,x) is only guaranteed to find factors below primelimit. In my limited experiments, it found some factors a little above primelimit, but not a lot above, except for perfect powers:

? factor(500009*500029,2^31-1)
%70 =
[500009 1]

[500029 1]

? factor(600011*600043,2^31-1)
%71 =
[360032400473 1]

? factor((600011*600043)^2,2^31-1)
%72 =
[360032400473 2]

Anyway, it's right that trial factoring is a minor part of a brilliant search. I also said that in http://mersenneforum.org/showpost.ph...4&postcount=40

[Citrix]

I have posted a program here (long time back) called plus.zip, it can sieve to 500M. Hope that helps?
http://mersenneforum.org/showthread....8754#post78754

[m_f_h]

Originally Posted by Jens K Andersen

Factoring to 2^62 in less than a second? Not likely
My PARI/GP 2.3.1 says:
? ?factor
factor(x,{lim}): factorization of x. lim is optional and can be set whenever x is of (possibly recursive) rational type. If lim is set return partial factorization, using primes up to lim (up to primelimit if lim=0)
On my 32-bit Windows, factor(n,2^31-1) is the highest allowed:
(...)
I guess that if x > primelimit, then factor(n,x) is only guaranteed to find factors below primelimit. In my limited experiments, it found some factors a little above primelimit, but not a lot above, except for perfect powers:

Indeed on 32 bit machines it's 2^31, on 64 bit machines it's 2^62.
But you may be right on the primelimit (I just checked and it's 5e5 for me also).
So the PARI documentation is wrong (at least unprecise and misleading - just like Maple's doc on factor(), see the other thread). They should add "..lim = 0 or lim > primelimit."

So I apologize for what I said. (However, 5e5 is still better than 1000 as somebody suggested earlier...) I'll see what will be the speed dropdown if I increase this to a somewhat higher value. (It will be less than 2^62 though....)

PS: just set primelimit to 10^7, then
(15:04) gp > factor((nextprime(9*10^7)*nextprime(10^7+200))^2,2^31)
%384 =
[10000223 2]
[90000049 2]
(15:05) gp > factor((nextprime(3*10^7)*nextprime(10^7+200)),2^31)
%389 =
[10000223 1]
[30000001 1]
(15:05) gp > factor((nextprime(3*10^7)*nextprime(2*10^7+200)),2^31)
%390 =
[600006410000213 1]

So you seem right with "...guaranteed to ... primelimit..." :-( !