calculation in a residue system
This page tells on Decimal BASIC Ver.7 and Ver.8.
a×b mod m is expressed as MOD(a*b, m) in BASIC.
Care that a*b does not exceed the upper limit on which the number can be accurately expressed, which is 16 digits in the decimal mode.
It is okay if a<1E8 and b<1E8 in the decimal mode. When you exceed it, calculate in rational mode.
Calculating an mod m with MOD (a^n,m ) increases the likelihood that a^n will exceed the limit at which the number can be accurately expressed.
For example, in decimal mode
10 PRINT MOD(9^88,89) 20 END
gives 41, but it is 1 correctly. You will get the correct answer in the rational number mode.
When the number of modulo is 1E8 or less, ，you can ask for the correct value in decimal mode by calculating the power multiplication in the residue system by repeating the multiplication.
The following program calculates the correct answer even in decimal mode.
100 DECLARE FUNCTION PowerMod 110 PRINT PowerMod(9,88,89) 120 END 130 EXTERNAL FUNCTION PowerMod(a,n,m) 140 DECLARE NUMERIC x,i 150 LET x=1 160 FOR i=1 TO n 170 LET x=MOD(x*a,m) 180 NEXT i 190 LET PowerMod=x 200 END FUNCTION
However, the more repetitions you have, the longer it will take. Calculating PowerMod (2, p-1, p) with p=99999989 took about 40 seconds in decimal mode.
a2, a4, a8, a16, .. .. can be quickly calculated by iterating over squaring.
Taking advantage of this, for example, Power calculation can be sped up by transforming it into as a37=a1+4+32=a1a4a32
This calculation method is particularly useful for calculating powers in residue systems. The next program will give you the answer in almost an instant.
100 DECLARE FUNCTION PowerMod 110 LET p=99999989 120 PRINT PowerMod(2,p-1,p) 130 END 140 EXTERNAL FUNCTION PowerMod(a,n,m) 150 DECLARE NUMERIC x,i,y 160 LET x=1 170 LET y=a 180 DO WHILE n>0 190 IF MOD(n,2)=1 THEN LET x=MOD(x*y,m) 200 LET y=MOD(y^2,m) 210 LET n=INT(n/2) 220 LOOP 230 LET PowerMod=x 240 END FUNCTION
If you want to modulo a number with an even larger number of digits, calculate in the rational number mode.
100 OPTION ARITHMETIC RATIONAL 110 DECLARE FUNCTION PowerMod 120 LET a=2 130 LET p=2^32+1 140 PRINT PowerMod(a,p-1,p) 150 END 800 EXTERNAL FUNCTION PowerMod(a,n,m) 810 OPTION ARITHMETIC RATIONAL 820 DECLARE NUMERIC x,i,y 830 LET x=1 840 LET y=a 850 DO WHILE n>0 860 IF MOD(n,2)=1 THEN LET x=MOD(x*y,m) 870 LET y=MOD(y^2,m) 880 LET n=INT(n/2) 890 LOOP 900 LET PowerMod=x 910 END FUNCTION
When p=232+1, if a=2 then ap-1≡1 (mod p), but if a=3 then ap-1≡1 (mod p) ，thus 232+1 is not a prime number (by Fermat's Little Theorem).
The reciprocal of a in the residue system of modulo m is x such that ax≡1 (mod m).
This means the existence of an integer k that satisfies the equation ax=1+km, it can be found from integers x and k that satisfy the eqation ax-mk=1.
Such x and k can be found as letting k=-y in the solution x and y for the linear indeterminate equation ax+by=1, but it is not necessary to find k, so you can get the reciprocal by solving the linear indeterminate equation ax+by=1 and find x.
The function InvMod (a,m ) defined in line 800 of the following program seeks the reciprocal of a in the residue system modulo m. It generates an exception for EXTYPE 999 when it does not exist.
The external sub-program Euclid (a,b,x,y,g) asks for an integer x,y that satisfy g=GCD(a,b) and ax+by=g.
100 OPTION ARITHMETIC RATIONAL 110 DECLARE EXTERNAL FUNCTION InvMod 120 DECLARE NUMERIC a,p 130 INPUT a,p 140 PRINT InvMod(a,p) 150 END 1000 EXTERNAL FUNCTION InvMod(a,m) 1010 OPTION ARITHMETIC RATIONAL 1020 DECLARE EXTERNAL FUNCTION Euclid 1030 DECLARE NUMERIC x,y,g 1040 CALL Euclid(a,m,x,y,g) 1050 ! ax+my=g, g=GCD(a,m) 1060 IF g<>1 THEN CAUSE EXCEPTION 999 1070 LET InvMod=MOD(x,m) 1080 END FUNCTION 2000 EXTERNAL SUB Euclid(a,b,x,y,g) 2010 OPTION ARITHMETIC RATIONAL 2020 DECLARE NUMERIC q,r,u,v 2030 LET q=INT(a/b) 2040 LET r=MOD(a,b) 2050 IF r=0 THEN 2060 LET x=0 2070 LET y=1 2080 LET g=b 2090 ELSE 2100 CALL Euclid(b,r,u,v,g) 2110 LET x=v 2120 LET y=u-q*v 2130 END IF 2140 END SUB