Lucas-Lehmer Test

The Lucas-Lehmer Test consists of a recursive algorithm performed p-1 times, where p is a prime exponent of the Mersenne Number to be tested, 2^p-1
The algorithm is defined as: S₁=4, S₂ = (S₁)² - 2, ..., S_P-1 = (S_p-2)² - 2
When the last result in the sequence is congruent to zero Mod(2^p-1), 2^p-1 is a prime number.
Here is an example:
p=7, 2^p-1 = 127

S₁=4, S₂=14, S₃=194, S₄=37634, S₅=1416317954, S₆=2005956546822746114

As one can see the numbers become quite large quickly. However, using modulo 127 on each result, the algorithm is simplified:

S₁=4, S₂=14, S₃=67, S₄=42, S₅=-16, S₆=0

Another interesting result is that S_p-2 = ±2^(p+1)/2 Mod(2^p-1)

There are several ways to generate the sequence of numbers for the Lucas-Lehmer Test. They do not seem to be efficient for large values of p, however they are interesting and may lead to a faster method for computing the last result.
A Binet formula for each result exists:
S_k = (2 + sqrt(3))^2(k-1) + (2 - sqrt(3))^2(k-1)
Example: k=1, S₁ = (2 + sqrt(3))^2(1-1) + (2 - sqrt(3))^2(1-1)
S₁ = (2 + sqrt(3)) + (2 - sqrt(3)) = 4
Example: k=3, S₃ = (2 + sqrt(3))^2(3-1) + (2 - sqrt(3))^2(3-1)
S₃ = (2 + sqrt(3))⁴ + (2 - sqrt(3))⁴
S₃ = (4 + 4sqrt(3)+3)² + (4 - 4sqrt(3)+3)²
S₃ = (7 + 4sqrt(3))² + (7 - 4sqrt(3))²
S₃ = (49 + 56sqrt(3)+48) + (49 - 56sqrt(3)+48)
S₃ = 194

Substituting 2=a, sqrt(3)=b, S_k = (a + b)^2(k-1) + (a - b)^2(k-1)
S₁ = 2a = 4
S₂ = 2(a² + b²) = 14
S₃ = 2(a⁴ + 6a²b² + b⁴) = 194
S₄ = 2(a⁸ + 28a⁶b²+ 70a⁴b⁴+ 28a²b⁶ + b⁸) = 37634

Now to simplify this further we employ Pascal's Triangle. Pascal's Triangle is an easy way to compute the coefficients of the binomial expansion of (a + b)ⁿ. The first and last coefficient are always 1. The intermidiate coefficients are generated by adding the two coefficients in the previous row.

							1		1									(a + b)1
						1		2		1								(a + b)2
					1		3		3		1							(a + b)3
				1		4		6		4		1						(a + b)4
			1		5		10		10		5		1					(a + b)5
		1		6		15		20		15		6		1				(a + b)6
	1		7		21		35		35		21		7		1			(a + b)7
1		8		28		56		70		56		28		8		1		(a + b)8	= a8b0 + 8a7b1 + 28a6b2 + 56a5b3 + 70a4b4 + 56a3b5 + 28a2b6 + 8a1b7 + a0b8

Since the Binet formula is raised to a power of 2, these are the only rows of Pascal's Trinagle we are interested in, 1,2,4,8,16,32, etc. Plus, we are only interested in every-other coefficient. Now we have another way of computing the last result in the algorithm for the Lucas-Lehmer Test. S₄, above, is equal to, 2 times the binomial expansion of every-other term.
The coefficients can also be computed using the Combination of n things taken r at a time: C(n,r) = n!/((n-r)!(r!)). For (a + b)⁸ the coefficients are C(8,8), C(8,7), C(8,6), C(8,5), C(8,4), C(8,3), C(8,2), C(8,1), C(8,0).

Here's another interesting way to use Pascal's Trinagle. Each row represents the number eleven in any base taken to a power. Lets do base 10 as an example, row 2 is 11² or 121, next row is 11³, 1331, etc. However, we have to do a little manipulating when the number in a cell is greater than the base. For example: row 8, 1 8 28 56 70 56 28 8 1. Here we have to use a way to get the correct digits. We divide-up and leave the moduloas the digit.

Step 1:								8	1
Step 2: 28 mod(10)							8	8	1
Step 3: 20/10 added to next cell						58	8	8	1
Step 4: 58 mod(10)						8	8	8	1
Step 5: 50/10 added to next cell					75	8	8	8	1
Step 6: 75 mod(10)					5	8	8	8	1
Step 7: 70/10 added to next cell				63	5	8	8	8	1
Step 8: 63 mod(10)				3	5	8	8	8	1
Step 9: 60/10 added to next cell			34	3	5	8	8	8	1
Step 10: 34 mod(10)			4	3	5	8	8	8	1
Step 11: 30/10 added to next cell		11	4	3	5	8	8	8	1
Step 12: 11 mod(10)		1	4	3	5	8	8	8	1
Step 13: 10/10 added to next cell	2	1	4	3	5	8	8	8	1

We're finished 11⁸ = 214358881

This works for any base where the representation of the desired number is 11. In binary that would be a 3 in base 10. The same process is applied. However, we would Mod(2) and divide by 2. So here it is for base 2:

Step 1:													1
Step 2: 8 mod(2)												0	1
Step 3: 8/2 added to next cell											32	0	1
Step 4: 32 mod(2)											0	0	1
Step 5: 32/2 added to next cell										72	0	0	1
Step 6: 72 mod(2)										0	0	0	1
Step 7: 72/2 added to next cell									106	0	0	0	1
Step 8: 106 mod(2)									0	0	0	0	1
Step 9: 106/2 added to next cell								109	0	0	0	0	1
Step 10: 109 mod(2)								1	0	0	0	0	1
Step 11: 108/2 added to next cell							82	1	0	0	0	0	1
Step 12: 82 mod(2)							0	1	0	0	0	0	1
Step 13: 82/2 added to next cell						49	0	1	0	0	0	0	1
Step 14: 49 mod(2)						1	0	1	0	0	0	0	1
Step 15: 48/2 added to next cell					25	1	0	1	0	0	0	0	1
Step 16: 25 mod(2)					1	1	0	1	0	0	0	0	1
Step 17: 24/2 added to next cell				12	1	1	0	1	0	0	0	0	1
Step 18: 12 mod(2)				0	1	1	0	1	0	0	0	0	1
Step 19: 12/2 added to next cell			6	0	1	1	0	1	0	0	0	0	1
Step 20: 6 mod(2)			0	0	1	1	0	1	0	0	0	0	1
Step 21: 6/2 added to next cell		3	0	0	1	1	0	1	0	0	0	0	1
Step 22: 3 mod(2)		1	0	0	1	1	0	1	0	0	0	0	1
Step 23: 2/2 added to next cell	1	1	0	0	1	1	0	1	0	0	0	0	1