restart: prodFirstPrimes := proc(N) # This function returns the product of the primes less than N local x,i,myPrime; i := 1: x := 1: myPrime := 2: while(myPrime < N) do i := 1: while (myPrime^i <= N) do i := i + 1: od: i := i - 1: x := x * myPrime^i: myPrime := nextprime(myPrime): od: x end proc;pollardPMinusOneMethod := proc(n, B1, B2) # This function uses the Pollard P-1 method of factoring to try and factor n # This method works only when n has a prime factor p such that p-1 has small prime # factors local a,k,b,d,i,myPrime,testFactor,currentB,t,c,index,gap,randomNumber: i := 0: d := 1: k := prodFirstPrimes(B1): randomNumber := rand(2..n-2): while ((d = 1) and (i < 10)) do # these first 4 lines are what is called the p-1 first stage in the literature # relatively simple a := randomNumber(): b := a&^k mod n: d := gcd(b - 1, n): i := i + 1: # this if block is the second stage p-1 factoring method. A little more # complex, yet is a more sophisticated method. if (d = 1) then myPrime := nextPrime(B1): t := 1: currentB := (b&^(myPrime)) mod N: c := 10: while ((t <= floor(B2/c)*c + 1) and (d = 1)) do testFactor := 1: for index from 0 to c do testFactor := testFactor * (currentB - 1): gap := nextPrime(myPrime) - myPrime: currentB := (b&^(gap)*currentB) mod N od: d = igcd(testFactor, N): t := t + c: od: end if: od: d end proc;lucasVFunction := proc(P,N,r) # this function computes V_r(P) mod N where V_n is the Lucas Function # defined by a^n + b^n where a,b are roots of the polynomial # x^2 - px + 1 local i, V_Zero, V_bigtemp, V_smalltemp, V_big, V_small,numBits,binaryR: # this is initialized to v_0 V_small := 2: # this is initialized to v_1 V_big := P: # convert r to binary and then to a string so we may access the bits individually binaryR := convert(convert(r, binary), string); numBits := evalf(log(r)/log(2)): if (ceil(numBits) > numBits) then numBits := numBits + 1: end if: # this loop calculates V_r(P) using an algorithm that is polynomial # in the number of bits of r, using the following recurrences: # V_(2f-1) = V_f*V_(f-1) - P mod N # V_(2f) = V_f*V_f - 2 mod N # V_(2f-1) = P*V_f^2 - V_f*V_(f-1) - P mod N for i from 2 to numBits do V_smalltemp := V_small: V_bigtemp := V_big: if (binaryR[i] = "1") then V_big := (P*V_bigtemp^2 - V_bigtemp*V_smalltemp - P) mod N: V_small := (V_bigtemp^2 - 2) mod N: else V_big := (V_bigtemp^2 - 2) mod N: V_small := (V_bigtemp*V_smalltemp - P) mod N: end if: od: V_big end proc;williamsPPlusOneMethod := proc(n,B) # This function uses the William's p+1 method of factoring an integer n. # This method works when n has a prime divisor p such that p+1 has small # prime divisors. This function returns a factor. local k,pnot,loopCount, result: loopCount := 1: # this finds the pnot such that (pnot^2 - 4) = 1 pnot := 3; while (igcd(pnot^2 - 4, n) <> 1) do pnot := pnot + 1: od: k := prodFirstPrimes(B): result := igcd(lucasVFunction(pnot,n,k) - 2, n): # loop until we have found a factor. Don't perform this # loop more than 10 times though. while ((result = 1) and (loopCount < 10)) do loopCount := loopCount + 1: # get the next pnot pnot := pnot + 1: while (igcd(pnot^2 - 4, n) <> 1) do pnot := pnot + 1: od: result := igcd(lucasVFunction(pnot,n,k) - 2, n): od: result end proc;getKBitPrime := proc(k) # This function returns a prime number of at least k bits local a,i: i := 0: a := rand(2^(k-1)..2^(k-1) + 2^(k-2) - 1)(): if (type(a, even)) then a := a + 1; end if; while (not(isprime(a))) do i := i + 1: a := a + 2: od: # note the return value is a pair. i is the number of times isprime() is called # and a is the prime of at least k bits [i,a] end proc;getKBitStrongPrime := proc(k) # This function returns a prime number of at least k bits that is deemed # a strong prime. This means that p is such that t = p-1 has a large prime # factor, p+1 has a large prime factor, and t's large prime factor has # another large prime factor. local r,s,t,p,pnot,L,i,retList,j: # get a large prime s of at least floor(k/2) bits retList := getKBitPrime(floor(k/2)): s := retList[2]: i := retList[1]: # get a large prime t of at least floor(k/2) bits retList := getKBitPrime(floor(k/2)): t := retList[2]: i := i + retList[1]: # this algorithm creates a prime r that is equivalent to 1 mod t(i.e. t divides r-1) L := floor(k/2): r := 2*L*t + 1: while (not(isprime(r))) do i := i + 1: L := L + 1; r := 2*L*t + 1: od: # The paper states that a prime that such that p+1 has factor s and p-1 has factor r # must be of the form p = pnot + 2*j*r*s where pnot is defined below. pnot := (s&^(r-1) - r&^(s-1)) mod r*s: if (type(pnot, even)) then pnot := pnot + r*s end if: j := 1; p := pnot + 2*j*r*s; while (not(isprime(p))) do i := i + 1: j := j + 1; p := pnot + 2*j*r*s; od: # Note that the return value is again a list, with the last value is the desired prime p [i,r,s,t,p] end proc;getRSA_e_d := proc(phi,start) local x,i,e,d: for i from start to phi-1 do x := i: if (igcdex(i,phi,'d','blah') = 1) then e := i: break: end if: od: [e mod phi,d mod phi] end proc;encoder := (P,e,n) -> (P&^e) mod n: decoder := (C,d,n) -> (C&^d) mod n:encodeRSA := proc(e,n,message) local blockSize, encMessage, messageLength, numBlocks, encMessage2, theBlocks, theEncodedBlocks, encodedMessage, i, j, morbledMessage; blockSize := floor(evalf(log(sqrt(n))/log(256)))-1; encMessage := convert(message, bytes); messageLength := nops(encMessage); numBlocks := ceil(messageLength/blockSize); encMessage2 := array(1..blockSize*numBlocks); for i from 1 to messageLength do encMessage2[i] := encMessage[i]: od: for i from messageLength+1 to blockSize*numBlocks do encMessage2[i] := 32: od: theBlocks := array(1..numBlocks); for i from 1 to numBlocks do theBlocks[i] := encMessage2[(i-1)*blockSize+1]: for j from 1 to blockSize-1 do theBlocks[i] := theBlocks[i] + (256^j)*encMessage2[(i-1)*blockSize+1+j]: od: od: theEncodedBlocks := convert(map(encoder, theBlocks,e,n),array); encodedMessage := map(convert,map(convert,theEncodedBlocks,base,256),bytes); morbledMessage :=``; for i from 1 to numBlocks do morbledMessage := cat(morbledMessage, encodedMessage[i]): od: [morbledMessage, theEncodedBlocks] end proc;decodeRSA := proc(theEncodedBlocks, d, n) local blockSize,decodedMessage, demorbledMessage,i; blockSize := max($op(2,theEncodedBlocks())); decodedMessage := map(convert,map(convert,map(decoder,theEncodedBlocks,d,n),base,256),bytes); demorbledMessage :=``; for i from 1 to blockSize do demorbledMessage := cat(demorbledMessage, decodedMessage[i]): od: demorbledMessage end proc;

If we want to use RSA, the first thing we should do is create our keys. So we'll need to choose our n = p*q. This is accomplished by the getKBitStrongPrime(k) procedure. Here k is a lower bound on the # of digits in the prime. Since we want p and q to be far apart, we'll make p 512 digits and q 550 digits.p := getKBitStrongPrime(512)[5]; q := getKBitStrongPrime(550)[5]; n := p*q; phi := (p-1)*(q-1); One suggestion that the original authors make is that gcd(p-1,q-1) should be fairly small. At this point we're not sure why, but we can check that this is the case and choose different primes if necessary. gcd(p-1,q-1); Now that we have strong primes, we can set up our public key (e,n) and private key (d,n) where ed=1 (mod phi). Our algorithm works by searching for a number e such that gcd(e,phi)=1. We must first pick a random place to start. start:= rand((phi/2)+1..(phi-1))(): if (type(start, even)) then start:=start+1; end if: start:=start: key:= getRSA_e_d(phi,start): e:= key[1]; d:=key[2]; Now we're ready to encrypt! Let's say someone wants to send us a secret message. They just need to look up our public key (e,n) and run the encryption function. message:="Hi guys. RSA is so very awesome...isn't it? I just wanted to say good job on your project and keep up the excellent work. Your's Truly, Ronald Linn Rivest."; encodeResult := encodeRSA(e,n,message): codeBlocks:= encodeResult[2]: encodeResult[1];We just recieved this encrypted message from our friend R.L. Rivest at M.I.T.. To decode his message, we simply break it into blocks of the appropriate size and apply the decoding algorithm. decodeResult := decodeRSA(codeBlocks,d,n);