[QUOTE=Jim White;113449]
It runs fastest of all so far, but I think that's mainly because it has a faster Pell solver - it uses various 'tricks' to speed up the Pell equation solving process (eg going directly from mid-period to fundamental solution).
[/QUOTE]

Thanks.
Yesterday I've figured out some tricks. I'll try to rewrite my program.

Do you know that period-length can be odd? All of the posted period lengths in the maximal solutions are even but there are examples when the period length of the continued fraction expansion of sqrt(n) is odd ( where n is even ).

Yes. Ignoring the case period = 1, odd periods are particularly rare when you have even [I]D's [/I](as we do currently, under the Lehmer scheme).

Here I've reported over various runs, the maximum odd period length that turned up, and the maximum [I]D [/I]which had an odd period:
[code]
p= 83 <none>

p= 89
Odd periods = 37
Max length pl = 31, D = 22952210
Max D pl = 27, D = 416392730

p= 97
Odd periods = 60
Max length pl = 31, D = 1715930
MAx D pl = 17, D = 3774000501602

p=101
Odd periods = 72
Max length pl = 31, D = 20402
Max D pl = 29, D = 23010564730570

p=103 <none>

p=107 <none>
[/code]

An odd period can generally be treated as just the mid-point of an even period, I think.

Pell CF Tricks
 

Regarding Pell equation (CF development, etc) efficiency, the restriction of the CF development to the mid-point is definitely worth implementing, even in our relatively short periods - it's still effective for period lengths of say, 8 or more.

This is generally known as "the mid-point trick", it seems

Another trick is one that seems peculiar to this problem - if we are going to reject period lengths that exceed some max [B]L[/B], and we know that the great majority of cases will be thus rejected, then it makes very good sense to separate the CF expansion into two stages:
[INDENT][B]1. Periodicity + CF Denominators[/B] - compute and store just the CF denominators, aborting any cases where periodicity is not detected in less than [B]L [/B]steps[/INDENT][INDENT][B]2. Convergents [/B]- for acceptable cases only, generate the required convergents (via the recurrence relation multipliers given by the stored CF denominators)[/INDENT]Combine with mid-point shortcut in stage [B]2, [/B]and you have program [tex]\sp[/tex] in a nutshell. That's why it's quite fast, even though some aspects are not fully optimised.

We can also use size bounding, of course, during stage 2. However, it looks like what happened when I applying period length bounds to [I]your [/I]code, they didn't make any real difference, as the size bound usually triggered first.

With the period-bounding method described here, it's just the opposite - by rejecting periods > [B]L[/B][I],[/I]we have already filtered out most of those that would run over size (and we have not wasted much computational effort on them at all - [I]this[/I] is the main efficiency gain)

So we only develop the convergents for those guys we know that have short periods, and we can still keep an eye on the size, although we expect to find such oversize cases less often.

Finger trouble
 

[quote=Jim White;113567]We can also use size bounding, of course, during stage 2.

[I]However, [U]it looks like what happened when I applying[/U] period length bounds to your code, they didn't make any real difference, as the size bound usually triggered first.[/I]

With the period-bounding method described here, it's just the opposite - ... ...
[/quote]

Egad! Must have fallen asleep at the keyboard (again!)

I'm sorry! That should read:
[quote=Jim White;113567]We can also use size bounding, of course, during stage 2.

If I add period length bound checks to Robert's original program, however, they don't make any real difference - his existing size bound check is almost always more effective (catches exceptions earlier) than the priod length bound.

With the period-bounding method described here, it's just the opposite - [/quote]

Prime Combinations - Ordering Choices
 

The set of [I]D [/I]values that we have to generate can be done in various combinatorial ways. Two examples:
[INDENT][B]1. Recursive Generation[/B] - eg: 3 -> 3,5 -> 3,5,7 -> 3,5,7,11 -> 3,5,11 -> 3,7 - > 3,7,11 -> 3,11 ==> 5 -> ... (a pre-order traversal of the combination tree)[/INDENT][INDENT][B]2. Generating Same-size Selections [/B]- this involves generating all combinations of size [B][I]k [/I][/B]for k = 1, 2, ... . So pass 1 = {3 -> 5 -> 7 -> 11}, pass 2 = {3,5 -> 3,7 -> 3,11 -> 5,7 -> 5, 11 -> 7, 11} and so on[/INDENT]The first method is very much simpler to code, and neither one makes any significant difference to running time in this problem.

But the second method has at least one interesting feature - it can locate the majority of the targets in a flash!

For example, with [B]n = 30 (p = 113)[/B], we are combining the 28 primes {3, 5, ... 109} to form [I]D' [/I](which we then combine with 2 and 113 in the various appropriate ways), of which there are [tex]2^{28}[/tex] possible combinations (well, give or take 1).

Here's a log from a sample run showing the results for various [B]k :[/B]

[code]
[U][B][FONT=Courier New] k C(n,k) Solns Total run time[/FONT][/B][/U]
[FONT=Courier New] 1 28, 19, 23, [/FONT]
[FONT=Courier New] 2 378, 72, 95, [/FONT]
[FONT=Courier New] 3 ** 3276, 191, 286, [/FONT]
[FONT=Courier New] 4 ***** 20475, 415, 701, .1s [/FONT]
[FONT=Courier New] 5 ****** 98280, 676, 1377, .4s[/FONT]
[FONT=Courier New] 6 ********* 376740, 878, 2255, 1.7s[/FONT]
[FONT=Courier New] 7 ********* 1184040, 934, 3189, 5.6s[/FONT]
[FONT=Courier New] 8 ********* 3108105, 819, 4008, 15.9s[/FONT]
[FONT=Courier New] 9 ***** 6906900, 575, 4583, 39.7s[/FONT]
[FONT=Courier New]10 **** 13123110, 357, 4940, 92.7s[/FONT]
[FONT=Courier New]11 ** 21474180, 193, 5133, 181.8s[/FONT]
[FONT=Courier New]12 * 30421755, 91, 5224, 311.4s[/FONT]
[FONT=Courier New]13 37442160, 36, 5260, 471.8s[/FONT]
[FONT=Courier New]14 40116600, 18, 5278, 643.8s[/FONT]
[FONT=Courier New]15 37442160, 3, 5281, 804.5s[/FONT]
[FONT=Courier New]16 30421755, 2, 5283, 934.9s[/FONT]
[FONT=Courier New]17 21474180, 1, 5284, 1026.8s[/FONT]
[FONT=Courier New]18 13123110, 0, 5284, 1083.0s[/FONT]
[FONT=Courier New]19 6906900, 0, 5284, 1112.7s[/FONT]
[FONT=Courier New]..... [/FONT]
[/code]
So it took just 16 seconds to find 4000 of the 5300 odd solutions, then nearly 20 minutes to find the the rest ...

The number of possible combinations that have to be considered is listed, and this gives you yet another view of the "fundamental" problem ...

The "peak" solution pass seems to be around [B]n/4 [/B]for all cases I've looked at so far.

If you write a generic [B]k [/B]from [B]n [/B]selection generator, then you can do the selection passes in any order you like, of course ....

Some known tuning settings
 

These settings might be useful if you want to save time on a particular run:

[code]
[U][B]n p D'2 S'2 v(D') v(S') L(D)[/B][/U]
25 97 64 57 13 17 26
26 101 66 63 14 17 28
27 103 72 76 15 18 28
28 107 74 60 15 19 28
29 109 77 68 16 20 32
30 113 87 69 17 20 36
31 127 81 67 16 20 32
32 131 88 73 16 21 30
33 137 88 78 17 21 30
[/code]

[B]n = [/B]index([B]p[/B])

[B]D'2 [/B]= max [I]reduced [/I]producer size (bits). The reduced determinants [I]D' [/I]are the odd square-free values. ie. subsets of {3, 5, ... [I]p'[/I]}, where [I]p'[/I] is the [B]n-1[/B]'th prime. This bound can save time by helping to eliminate many of the possible divisor combinations that are supposed to be tested.

[B]S'2 [/B]= max solution size (bits). If [B]S' [/B]is the max [I]S [/I]for which [I]S(S-1) [/I]is a hit, then [B]S'2 = [/B][tex]log_2(2S-1)[/tex]. Like knowing max period length, this allows Pell equation solution development to be abandoned early when known to be futile.

[B][I]v[/I](D') [/B]= [I]v(N) [/I]is the number of prime divisors of N. [B][I]v[/I](D') [/B]is the highest value for any (reduced) producer [I][B]D'[/B][/I]. In other words, it's pointless combining any more than this when making any [I]D'[/I].

[B][I]v[/I](S') [/B]= the highest number of divisors of any [I]S(S-1) [/I]pair. Includes all divisors, so includes 2 and [B]p. [/B]

[B]L(D) [/B]= the largest period length of any producer.

[B]D'2 [/B]and [B]L(D) [/B]are very effective bounds (when known!)


As an example, a current program I have takes over 1 hour to run (70 mins) for the case[B] n=[/B]32, [B]p[/B]=131, using a reasonably good period bound (these seem very predictable, eg [B]n + 10 [/B]would usually do the job).

Once I have a sharper value for the period bound, and also know the size of the largest producer I need to generate, re-running the job is quicker, completing in 40 minutes.

[SIZE=3][FONT=Times New Roman]I am also interested in this topic.[/FONT][/SIZE]
[FONT=Times New Roman][SIZE=3]I searched all p34-smooth number pair by my c language program and get 63428 solutions. My arithmetic is base on a conjecture, let (n,n+1) is a p34 –smooth number pair, the max of n should at most with 128 bits. The following is the largest 16 solutions (n) that satisfy with both n and n+1 is [/SIZE][/FONT][FONT=Arial CYR]139-smooth number[/FONT]

[SIZE=3][FONT=Times New Roman]107109384837019296875[/FONT][/SIZE]
[SIZE=3][FONT=Times New Roman]118991993751522079210[/FONT][/SIZE]
[SIZE=3][FONT=Times New Roman]152295745769656587384[/FONT][/SIZE]
[SIZE=3][FONT=Times New Roman]205142063213188103639[/FONT][/SIZE]
[SIZE=3][FONT=Times New Roman]507882950033870979071[/FONT][/SIZE]
[SIZE=3][FONT=Times New Roman]1043073004436787852800[/FONT][/SIZE]
[SIZE=3][FONT=Times New Roman]1338614896790283750000[/FONT][/SIZE]
[SIZE=3][FONT=Times New Roman]1492860563840742187104[/FONT][/SIZE]
[SIZE=3][FONT=Times New Roman]1805692829069444659200[/FONT][/SIZE]
[SIZE=3][FONT=Times New Roman]1811103583756769555714[/FONT][/SIZE]
[SIZE=3][FONT=Times New Roman]2187239579855392547199[/FONT][/SIZE]
[SIZE=3][FONT=Times New Roman]2888494687736348753790[/FONT][/SIZE]
[SIZE=3][FONT=Times New Roman]3482435534325338530939[/FONT][/SIZE]
[SIZE=3][FONT=Times New Roman]4114304445616636016031[/FONT][/SIZE]
[SIZE=3][FONT=Times New Roman]19316158377073923834000[/FONT][/SIZE]
[FONT=Times New Roman][SIZE=3]124225935845233319439173[/SIZE][/FONT]

I google number [FONT=Arial CYR]124225935845233319439173 in internet and found a table as bellow, of course, it is spectacular result. The table give the max p_n smooth nubmer pair and n up to 40. I don't understand how to search such result in very very big range. [/FONT]
[FONT=Arial CYR][FONT=Times New Roman][SIZE=3]I believe your state is reasonable. “[/SIZE][/FONT][FONT=Arial]There is a bigger 157-limit comma”[/FONT]
[FONT=Arial]This table can be found from [FONT=Times New Roman][URL="http://www.primefan.ru/stuff/math/maxs.xls"][COLOR=#800080]www.primefan.ru/stuff/math/maxs.xls[/COLOR][/URL][/FONT][/FONT][/FONT]

[QUOTE]
Index p k=1 (A002072)
1 2 1
2 3 8
3 5 80
4 7 4374
5 11 9800
6 13 123200
7 17 336140
8 19 11859210
9 23 11859210
10 29 177182720
11 31 1611308699
12 37 3463199999
13 41 63927525375
14 43 421138799639
15 47 1109496723125
16 53 1453579866024
17 59 20628591204480
18 61 31887350832896
19 67 31887350832896
20 71 119089041053696
21 73 2286831727304144
22 79 9591468737351909375
23 83 9591468737351909375
24 89 9591468737351909375
25 97 9591468737351909375
26 101 9591468737351909375
27 103 19316158377073923834000
28 107 19316158377073923834000
29 109 19316158377073923834000
30 113 19316158377073923834000
31 127 19316158377073923834000
32 131 19316158377073923834000
33 137 124225935845233319439173
34 139 124225935845233319439173
35 149 124225935845233319439173
36 151 124225935845233319439173
37 157 6318268857746831540296874
38 163 6318268857746831540296874
39 167 6318268857746831540296874
40 173 6318268857746831540296874

[/QUOTE]

You (or the guy behind of primefan.ru) should report that [URL="http://oeis.org/A002072/list"]here[/URL].
My feeling is that the list would not be very difficult to extend, I will give it a try in pari if I find some free time over teh weekend.

Basically, total 4 OEIS is related with this topic. They are A117581, A002072, A002071 and A145606.
let both m and m+1 is P(n)-smooth number, then A002072= max(m) for each of P(n) smooth number pair, and A117581=max(m+1), so A117581 [i] = A002072[i]+1.
The A002071 is number count for all of m which satisfy m and m+1 is P(n) smooth number.
The A145606 denote that largest number x such that x and x+1 are prime(n)-smooth but not prime(n-1)-smooth.
It is said in both A14506 and A002072 "Jim White, Results to P = 173". If the OEIS owner thinks Jim White's result and computation method is exactly correct, then all of these 4 serial should extend to 40 term. But I think that because Jim White's computation method base on a conjecture, the upperbound should less than a const and this upperbound can not been prove by mathematical theorem, so only a part of serial be given.

Jim White, If I am wrong, please correct me.

[QUOTE=Jim White;113248]It makes sense to start a specific thread for this topic, I think ..

The discussion developed originally in this thread
[URL]http://www.mersenneforum.org/showthread.php?t=5647[/URL]

I'll start this thread by reproducing my first posting there, which has some useful data, and also give a simple definition of the problem for newcomers:


[B][U]Størmer's Problem (1897)[/U][/B]

[SIZE=3][COLOR=black]Write [I][B]gpd[/B][/I]([I][B]n[/B][/I]) for the greatest prime divisor of [I][B]n[/B][/I][/COLOR][/SIZE]


[COLOR=black][SIZE=3]Given some prime [B]p[/B], [/SIZE][SIZE=3]for which [I][B]n[/B][/I] do we have [I][B]gpd[/B][/I]([B]n[/B]([B]n[/B]-[B]1[/B])) <= [B]p[/B]? [/SIZE][/COLOR]


We can refer to the set of all such values of [B][I]n [/I][/B]for given [B]p [/B]as [FONT=Courier New][SIZE=3][COLOR=#000080][FONT=Verdana][SIZE=2][COLOR=black][B]S(p)[/B], and [/COLOR][/SIZE][/FONT][/COLOR][/SIZE][/FONT][SIZE=3][COLOR=#000080][SIZE=2][COLOR=black][B]S'[/B]([B]p[/B])[/COLOR][/SIZE] [/COLOR][/SIZE]for the [I]maximum[/I] value [B][I]n [/I][/B]for which [B][I]gpd[/I][/B]([B]n[/B]([B]n[/B]-[B]1[/B])) [B]=[/B] [B]p [/B](note that these sets are demonstrably [I]finite, [/I]and so [B]S'[/B]([B]p[/B]) always exists).

The pairs ([B]n[/B],[B] n-1[/B]) are obviously representable as (twice the value of) triangular numbers, via their product, and are also known in musical interval theory as [I]Superparticular Ratios, [/I]when represented as rational quotients [B]n / n-1 [/B](these ratios are highly composite, and are as close to possible to 1)

=======================================================

OEIS sequences [B]A002072[/B] and [B]A117581[/B] both list the sequence [B]S'(p)[/B], only the first one lists the lower value of each pair, ie [B]S'(p) - 1[/B]

Using the higher value is preferable, I think, since the only published tables, the 1964 results of Dick Lehmer ("On a Problem of Størmer", [I]Illinois J. Math[/I], [B]8[/B], 1964, pp 57-79), for [B]p [/B]up to 41,use that convention.

Also his proofs are given in the same terms, and it corresponds to the convention used by the muso's.

Both sequences as recorded in OEIS have an inconvenient feature, though - there are duplicate values in the sequence, in order (presumably) to keep the sequence increasing.

But this omits useful information. If we ask what is [B][I]S'[/I](23)[/B], the [I]greatest[/I] value of S for which [I]gpd [/I](S(S-1)) = 23?

The answer is 5142501, but this is [I]less[/I] than [B][I]S'[/I](19) = [/B]11859211, so it does not appear in the sequence, the latter value is just repeated.

Here is a list which applies the definition of [B][I]S'[/I](p) [/B]strictly, and which extends the sequence to the 35th prime (149):

[code]
[SIZE=3][FONT=Courier New][COLOR=navy]N p S'(p) log2(S')[/COLOR][/FONT][/SIZE]
[SIZE=3][FONT=Courier New][COLOR=navy]=============================================[/COLOR][/FONT][/SIZE]
[SIZE=3][FONT=Courier New][COLOR=navy] 1. 2 2 1[/COLOR][/FONT][/SIZE]
[SIZE=3][FONT=Courier New][COLOR=navy] 2. 3 9 3.1699[/COLOR][/FONT][/SIZE]
[SIZE=3][FONT=Courier New][COLOR=navy] 3. 5 81[/COLOR][/FONT][/SIZE]
[SIZE=3][FONT=Courier New][COLOR=navy] 4. 7 4375[/COLOR][/FONT][/SIZE]
[SIZE=3][FONT=Courier New][COLOR=navy] 5. 11 9801[/COLOR][/FONT][/SIZE]
[SIZE=3][FONT=Courier New][COLOR=navy] 6. 13 123201[/COLOR][/FONT][/SIZE]
[SIZE=3][FONT=Courier New][COLOR=navy] 7. 17 336141[/COLOR][/FONT][/SIZE]
[SIZE=3][FONT=Courier New][COLOR=navy] 8. 19 11859211 23.4995[/COLOR][/FONT][/SIZE]
[SIZE=3][FONT=Courier New][COLOR=navy] 9. 23 5142501 22.2940[/COLOR][/FONT][/SIZE]
[SIZE=3][FONT=Courier New][COLOR=navy]10. 29 177182721 27.4077[/COLOR][/FONT][/SIZE]
[SIZE=3][FONT=Courier New][COLOR=navy]11. 31 1611308700 30.5856 [/COLOR][/FONT][/SIZE]
[SIZE=3][FONT=Courier New][COLOR=navy]12. 37 3463200000 31.6895[/COLOR][/FONT][/SIZE]
[SIZE=3][FONT=Courier New][COLOR=navy]13. 41 63927525376 35.8957[/COLOR][/FONT][/SIZE]
[SIZE=3][FONT=Courier New][COLOR=navy]14. 43 421138799640 38.6155[/COLOR][/FONT][/SIZE]
[SIZE=3][FONT=Courier New][COLOR=navy]15. 47 1109496723126 40.0103[/COLOR][/FONT][/SIZE]
[SIZE=3][FONT=Courier New][COLOR=navy]16. 53 1453579866025 40.4027[/COLOR][/FONT][/SIZE]
[SIZE=3][FONT=Courier New][COLOR=navy]17. 59 20628591204481 44.2297[/COLOR][/FONT][/SIZE]
[SIZE=3][FONT=Courier New][COLOR=navy]18. 61 31887350832897 44.8580[/COLOR][/FONT][/SIZE]
[SIZE=3][FONT=Courier New][COLOR=navy]19. 67 12820120234376 43.5435[/COLOR][/FONT][/SIZE]
[SIZE=3][FONT=Courier New][COLOR=navy]20. 71 119089041053697 46.7090[/COLOR][/FONT][/SIZE]
[SIZE=3][FONT=Courier New][COLOR=navy]21. 73 2286831727304145 51.0223[/COLOR][/FONT][/SIZE]
[SIZE=3][FONT=Courier New][COLOR=navy]22. 79 9591468737351909376 63.0565 [/COLOR][/FONT][/SIZE]
[SIZE=3][FONT=Courier New][COLOR=navy]23. 83 17451620110781857 53.9542 [/COLOR][/FONT][/SIZE]
[SIZE=3][FONT=Courier New][COLOR=navy]24. 89 166055401586083681 57.2044[/COLOR][/FONT][/SIZE]
[SIZE=3][FONT=Courier New][COLOR=navy]25. 97 49956990469100001 55.4715[/COLOR][/FONT][/SIZE]
[SIZE=3][FONT=Courier New][COLOR=navy]26. 101 4108258965739505500 61.8332[/COLOR][/FONT][/SIZE]
[SIZE=3][FONT=Courier New][COLOR=navy]27. 103 19316158377073923834001 74.0322 [/COLOR][/FONT][/SIZE]
[SIZE=3][FONT=Courier New][COLOR=navy]28. 107 386539843111191225 58.4234[/COLOR][/FONT][/SIZE]
[SIZE=3][FONT=Courier New][COLOR=navy]29. 109 90550606380841216611 66.2954 [/COLOR][/FONT][/SIZE]
[SIZE=3][FONT=Courier New][COLOR=navy]30. 113 205142063213188103640 67.4752[/COLOR][/FONT][/SIZE]
[SIZE=3][FONT=Courier New][COLOR=navy]31. 127 53234795127882729825 65.5290[/COLOR][/FONT][/SIZE]
[SIZE=3][FONT=Courier New][COLOR=navy]32. 131 4114304445616636016032 71.8011[/COLOR][/FONT][/SIZE]
[SIZE=3][FONT=Courier New][COLOR=navy]33. 137 124225935845233319439174 76.7173[/COLOR][/FONT][/SIZE]
[SIZE=3][FONT=Courier New][COLOR=navy]34. 139 3482435534325338530940 71.5606[/COLOR][/FONT][/SIZE]
[SIZE=3][FONT=Courier New][COLOR=navy]35. 149 6418869735252139369210 72.4428[/COLOR][/FONT][/SIZE]
[/code]


Jim White
Math. Sciences Institute
Australian National University[/QUOTE]

brute force is basically:

[CODE]gpd(y) = factor(y)[#factor(y)[,1],1][/CODE] and a infinite loop of values to check, a less brute force is knowing that for all all primes q<p , all values q*x for x<p are in the set also all powers q^y for q<p are in the set which lead to an infinite set.