41 digit exponentiation modulo P*Q

 a 98[97]-step program for an hp21s pocket calculator

Using double-precision arithmetic, [m].[n] designates m@10+n@0 --
also written as concatenated decimal, "m,n" (zeros un-suppressed) --
and similarly [[21-digits]] means [upper-11-digits].[lower-10-digits] --
j:k denotes separate sequential process parts chained (all carries kept) --

 STO 1,2 'modulator' [5@10].[@10] maximum STO 5,6 'root' [5@10].[@10] maximum; [0].[5@10] works STO 7,8 'result' [0].[1] initially 'X' 'exponent' [2@11] maximum: chains (binary) low-order first XEQ D 'exponentiate' XEQ C 'multiply' XEQ A 'adjust' also 9,4 'multiplier' [10@10].[@10] maximum; [0].[9@10] works also root = 'multiplicand' also 3,0 'product'

 ```LBL D IP X=0? RTN ÷ 2 INPUT FRAC X=0? GTO 6 RCL 8 STO 4 RCL 7 STO 9 XEQ C RCL 0 STO 8 RCL 3 STO 7 -------- LBL 6 RCL 6 STO 4 RCL 5 STO 9 XEQ C ``` ```RCL 0 STO 6 RCL 3 STO 5 SWAP GTO D -------- LBL B RCL 9 X=0? RTN XEQ A 1 0 STO × 0 STO × 3 STO × 4 STO × 9 RCL 9 IP STO - 9 × ( × RCL 6 ) ``` ```STO + 0 RCL 5 = STO + 3 GTO B -------- LBL C 0 STO 0 STO 3 E 1 1 +/- STO + 4 STO × 9 STO × 4 XEQ B RCL 4 STO 9 XEQ B RCL 6 STO - 0 RCL 5 STO - 3 1 ``` ```0 STO ÷ 3 STO ÷ 0 -------- LBL A RCL 3 ÷ RCL 1 = IP × ( × RCL 2 ) STO - 0 RCL 1 = STO - 3 RCL 0 IP STO - 0 STO + 3 RTN -------- SHOW d96C ```
M = (mp-mq/Q mod P) * Q + mq : mp,mq ≡ M mod P,Q primes; Q(/Q) ≡ 1 mod P ;
M^97 mod P*Q : P ~ 20.65 digits ; Q ~ 20.35 digits ; P ~ 2*Q ; P*Q ~ 41 digits
/Q ≡ Q^(P-2) mod P ≡ Q^(ØP-1) mod P

Example: Find a prime (qualifying a prime-suspect: not a final proof of primality)

 STO 1,2 'prime-suspect' [5@10].[@10] maximum STO 7,8 'result' [0].[1] initially STO 5,6 'chain-link' =2^33~35 [[2^37]] maximum; [0].[2^35] works (estimating '1+2'/'5+6' for major remainder-fraction>0.5) STO 9,4 'integer-part' =[IP['1+2'/'5+6']]; [0].[9@10] works XEQ C 'multiply' STO 5,6 'result' [0].[2] trially, or other - or m- '-3-0' or '1-3'+'2-0' positive 'product' modulo 'modulator' ÷ RCL 8 (or ×10^10) integer form XEQ D 'exponentiate' (about 6 minutes) RCL 1+2 estimate of 'prime-suspect' ÷ .[2^same 33~35] IP 'integer-part' of 'prime-suspect'/2^same 33~35 XEQ D 'exponentiate' (continue, about 6 minutes) RCL 7,8 Check 'result' for original 5,6 trial (congruent)

Sample values:
P = 4,9567823169,6710863373 =(33bit)= 5,7704540865:0,4969761293 (chain pair)
P = 4,9658723516,8246629391 =(34bit)= 2,8905181398:1,0816390159 (chain pair)
P = 4,8653265932,1434373841 =(35bit)= 1,4159964028:1,9062947537 (chain pair)
P*Q = 0,9999999999,9999999988,6922049508,2506023793
P = 1,3210304190,5569697827 =(33bit)= 1,5378818137:0,6476402723
Q = 0,7569848396,9417076059 =(34bit)= 0,4406231686:0,9541311835
KEY97 = 0,7134209118:4928465501 = 0,6128238969,4198475357
/97 mod ØP ≡ 0,7081812555,7624992649
/97 mod ØQ ≡ 0,2029031529,0771587397
/Q ≡ 1,0938094288,6891802521
/P ≡ 0,1302035459,2185788965

Ø is the totient function, Ø(P*Q^n) = (P-1)*(Q-1)*Q^(n-1)
@ delimits the radix-10 characteristic (base ten exponent), 2@3 = 2*10^3 = 2000

[A similar program should be writable for an RPN calculator; the hp21s was at hand]

A premise discovery under the title,