Modular inverse: Difference between revisions
Realize in F# |
m →{{header|Wren}}: Minor tidy |
||
(125 intermediate revisions by 46 users not shown) | |||
Line 17: | Line 17: | ||
Either by implementing the algorithm, by using a dedicated library or by using a built-in function in |
Either by implementing the algorithm, by using a dedicated library or by using a built-in function in |
||
your language, compute the modular inverse of 42 modulo 2017. |
your language, compute the modular inverse of 42 modulo 2017. |
||
<br><br> |
|||
=={{header|11l}}== |
|||
{{trans|C}} |
|||
<syntaxhighlight lang="11l">F mul_inv(=a, =b) |
|||
V b0 = b |
|||
V x0 = 0 |
|||
V x1 = 1 |
|||
I b == 1 {R 1} |
|||
L a > 1 |
|||
V q = a I/ b |
|||
(a, b) = (b, a % b) |
|||
(x0, x1) = (x1 - q * x0, x0) |
|||
I x1 < 0 {x1 += b0} |
|||
R x1 |
|||
print(mul_inv(42, 2017))</syntaxhighlight> |
|||
{{out}} |
|||
<pre> |
|||
1969 |
|||
</pre> |
|||
=={{header|8th}}== |
=={{header|8th}}== |
||
<syntaxhighlight lang="forth"> |
|||
<lang Forth> |
|||
\ return "extended gcd" of a and b; The result satisfies the equation: |
\ return "extended gcd" of a and b; The result satisfies the equation: |
||
\ a*x + b*y = gcd(a,b) |
\ a*x + b*y = gcd(a,b) |
||
Line 46: | Line 69: | ||
42 2017 n:invmod . cr bye |
42 2017 n:invmod . cr bye |
||
</syntaxhighlight> |
|||
</lang> |
|||
{{out}} |
{{out}} |
||
<pre> |
<pre> |
||
1969 |
1969 |
||
</pre> |
|||
=={{header|Action!}}== |
|||
<syntaxhighlight lang="action!">INT FUNC ModInverse(INT a,b) |
|||
INT t,nt,r,nr,q,tmp |
|||
IF b<0 THEN b=-b FI |
|||
IF a<0 THEN a=b-(-a MOD b) FI |
|||
t=0 nt=1 |
|||
r=b nr=a MOD b |
|||
WHILE nr#0 |
|||
DO |
|||
q=r/nr |
|||
tmp=nt nt=t-q*nt t=tmp |
|||
tmp=nr nr=r-q*nr r=tmp |
|||
OD |
|||
IF r>1 THEN |
|||
RETURN (-1) |
|||
FI |
|||
IF t<0 THEN |
|||
t==+b |
|||
FI |
|||
RETURN (t) |
|||
PROC Test(INT a,b) |
|||
INT res |
|||
res=ModInverse(a,b) |
|||
IF res>=0 THEN |
|||
PrintF("%I MODINV %I=%I%E",a,b,res) |
|||
ELSE |
|||
PrintF("%I MODINV %I has no result%E",a,b) |
|||
FI |
|||
RETURN |
|||
PROC Main() |
|||
Test(42,2017) |
|||
Test(40,1) |
|||
Test(52,-217) |
|||
Test(-486,217) |
|||
Test(40,2018) |
|||
RETURN</syntaxhighlight> |
|||
{{out}} |
|||
[https://gitlab.com/amarok8bit/action-rosetta-code/-/raw/master/images/Modular_inverse.png Screenshot from Atari 8-bit computer] |
|||
<pre> |
|||
42 MODINV 2017=1969 |
|||
40 MODINV 1=0 |
|||
52 MODINV -217=96 |
|||
-486 MODINV 217=121 |
|||
40 MODINV 2018 has no result |
|||
</pre> |
</pre> |
||
=={{header|Ada}}== |
=={{header|Ada}}== |
||
<syntaxhighlight lang="ada"> |
|||
<lang Ada> |
|||
with Ada.Text_IO;use Ada.Text_IO; |
with Ada.Text_IO;use Ada.Text_IO; |
||
procedure modular_inverse is |
procedure modular_inverse is |
||
Line 62: | Line 135: | ||
-- To calculate the inverse we do as if we would calculate the GCD with the Euclid extended algorithm |
-- To calculate the inverse we do as if we would calculate the GCD with the Euclid extended algorithm |
||
-- (but we just keep the coefficient on a) |
-- (but we just keep the coefficient on a) |
||
function inverse (a, b, u, v : Integer) return Integer is |
function inverse (a, b, u, v : Integer) return Integer is |
||
(if b=0 then u else inverse (b, a mod b, v, u-(v*a)/b)); |
|||
begin |
begin |
||
return inverse (a, n, 1, 0); |
return inverse (a, n, 1, 0); |
||
Line 72: | Line 146: | ||
Put_Line (inv_mod (42,77)'img); |
Put_Line (inv_mod (42,77)'img); |
||
exception when others => Put_Line ("The inverse doesn't exist."); |
exception when others => Put_Line ("The inverse doesn't exist."); |
||
end |
end modular_inverse; |
||
</syntaxhighlight> |
|||
</lang> |
|||
=={{header|ALGOL 68}}== |
=={{header|ALGOL 68}}== |
||
< |
<syntaxhighlight lang="algol68"> |
||
BEGIN |
BEGIN |
||
PROC modular inverse = (INT a, m) INT : |
PROC modular inverse = (INT a, m) INT : |
||
Line 112: | Line 186: | ||
CO |
CO |
||
END |
END |
||
</syntaxhighlight> |
|||
</lang> |
|||
{{out}} |
{{out}} |
||
<pre>42 ^ -1 (mod 2017) = 1969 |
<pre>42 ^ -1 (mod 2017) = 1969 |
||
</pre> |
</pre> |
||
=={{header|Arturo}}== |
|||
{{trans|D}} |
|||
<syntaxhighlight lang="rebol">modInverse: function [a,b][ |
|||
if b = 1 -> return 1 |
|||
b0: b x0: 0 x1: 1 |
|||
z: a |
|||
while [z > 1][ |
|||
q: z / b t: b |
|||
b: z % b z: t |
|||
t: x0 x0: x1 - q * x0 |
|||
x1: t |
|||
] |
|||
(x1 < 0) ? -> x1 + b0 |
|||
-> x1 |
|||
] |
|||
print modInverse 42 2017</syntaxhighlight> |
|||
{{out}} |
|||
<pre>1969</pre> |
|||
=={{header|ATS}}== |
|||
===Using allocated memory=== |
|||
In addition to solving the task, I demonstrate some aspects of call-by-reference in ATS. In particular, ATS can distinguish ''at compile time'' between uninitialized and initialized variables. |
|||
The code is written as templates that will expand to code for any of the (non-dependent) signed integer types. In the main program, I use <code>llint</code> (typically exactly the same as <code>long long int</code> in C). |
|||
The return value of <code>inverse</code> is a linear optional value. It is allocated in the heap, used once, and freed. |
|||
<syntaxhighlight lang="ats"> |
|||
(* |
|||
Using the algorithm described at |
|||
https://en.wikipedia.org/w/index.php?title=Extended_Euclidean_algorithm&oldid=1135569411#Modular_integers |
|||
*) |
|||
#include "share/atspre_staload.hats" |
|||
fn {tk : tkind} |
|||
division_with_nonnegative_remainder |
|||
(n : g0int tk, d : g0int tk, |
|||
(* q and r are called by reference, and start out |
|||
uninitialized. *) |
|||
q : &g0int tk? >> g0int tk, |
|||
r : &g0int tk? >> g0int tk) |
|||
: void = |
|||
let |
|||
(* The C optimizer most likely will reduce these these two |
|||
divisions to just one. They are simply synonyms for C '/' and |
|||
'%', and perform division that rounds the quotient towards |
|||
zero. *) |
|||
val q0 = g0int_div (n, d) |
|||
val r0 = g0int_mod (n, d) |
|||
in |
|||
(* The following calculation results in 'floor division', if the |
|||
divisor is positive, or 'ceiling division', if the divisor is |
|||
negative. This choice of method results in the remainder never |
|||
being negative. *) |
|||
if isgtez n || iseqz r0 then |
|||
(q := q0; r := r0) |
|||
else if isltz d then |
|||
(q := succ q0; r := r0 - d) |
|||
else |
|||
(q := pred q0; r := r0 + d) |
|||
end |
|||
fn {tk : tkind} |
|||
inverse (a : g0int tk, n : g0int tk) : Option_vt (g0int tk) = |
|||
let |
|||
typedef integer = g0int tk |
|||
fun |
|||
loop (t : integer, newt : integer, |
|||
r : integer, newr : integer) : Option_vt integer = |
|||
if iseqz newr then |
|||
begin |
|||
if r > g0i2i 1 then |
|||
None_vt () |
|||
else if t < g0i2i 0 then |
|||
Some_vt (t + n) |
|||
else |
|||
Some_vt t |
|||
end |
|||
else |
|||
let |
|||
(* These become C variables. *) |
|||
var quotient : g0int tk? |
|||
var remainder : g0int tk? |
|||
(* Show the type AT COMPILE TIME. *) |
|||
prval _ = $showtype quotient |
|||
prval _ = $showtype remainder |
|||
val () = |
|||
division_with_nonnegative_remainder |
|||
(r, newr, quotient, remainder) |
|||
(* THE TYPES WILL HAVE CHANGED, because the storage is |
|||
initialized by the call to |
|||
division_with_nonnegative_remainder. *) |
|||
prval _ = $showtype quotient |
|||
prval _ = $showtype remainder |
|||
val t = newt |
|||
and newt = t - (quotient * newt) |
|||
and r = newr |
|||
and newr = remainder |
|||
in |
|||
loop (t, newt, r, newr) |
|||
end |
|||
in |
|||
loop (g0i2i 0, g0i2i 1, n, a) |
|||
end |
|||
implement |
|||
main0 () = |
|||
case+ inverse (42LL, 2017LL) of |
|||
| ~ None_vt () => println! "There is no inverse." |
|||
| ~ Some_vt value => println! value |
|||
</syntaxhighlight> |
|||
{{out}} |
|||
First compile the program: |
|||
<pre>$ patscc -DATS_MEMALLOC_LIBC -g -O2 modular-inverse.dats |
|||
**SHOWTYPE[UP]**(/home/trashman/src/chemoelectric/rosettacode-contributions/modular-inverse.dats: 1786(line=62, offs=31) -- 1794(line=62, offs=39)): S2Etop(knd=0; S2Eapp(S2Ecst(g0int_t0ype); S2Evar(tk(8481)))): S2RTbas(S2RTBASimp(1; t@ype)) |
|||
**SHOWTYPE[UP]**(/home/trashman/src/chemoelectric/rosettacode-contributions/modular-inverse.dats: 1825(line=63, offs=31) -- 1834(line=63, offs=40)): S2Etop(knd=0; S2Eapp(S2Ecst(g0int_t0ype); S2Evar(tk(8481)))): S2RTbas(S2RTBASimp(1; t@ype)) |
|||
**SHOWTYPE[UP]**(/home/trashman/src/chemoelectric/rosettacode-contributions/modular-inverse.dats: 2137(line=72, offs=31) -- 2145(line=72, offs=39)): S2Eapp(S2Ecst(g0int_t0ype); S2Evar(tk(8481))): S2RTbas(S2RTBASimp(1; t@ype)) |
|||
**SHOWTYPE[UP]**(/home/trashman/src/chemoelectric/rosettacode-contributions/modular-inverse.dats: 2176(line=73, offs=31) -- 2185(line=73, offs=40)): S2Eapp(S2Ecst(g0int_t0ype); S2Evar(tk(8481))): S2RTbas(S2RTBASimp(1; t@ype)) |
|||
</pre> |
|||
You may notice there is a subtle change in the type of <code>quotient</code> and <code>remainder</code>, once they have been initialized. ATS is making it safe ''not'' to initialize variables (with phony values) when you declare them. |
|||
Now run the program: |
|||
<pre>$ ./a.out |
|||
1969</pre> |
|||
===Safely avoiding the need for an allocator=== |
|||
Here I demonstrate an optional value that requires no runtime overhead at all, but which is safe. If there is no inverse, then the compiler knows that <code>inverse_value</code> is still uninitialized, and will not let you use its value. (Try it and see.) |
|||
<syntaxhighlight lang="ats"> |
|||
(* |
|||
Using the algorithm described at |
|||
https://en.wikipedia.org/w/index.php?title=Extended_Euclidean_algorithm&oldid=1135569411#Modular_integers |
|||
*) |
|||
#include "share/atspre_staload.hats" |
|||
fn {tk : tkind} |
|||
division_with_nonnegative_remainder |
|||
(n : g0int tk, d : g0int tk, |
|||
(* q and r are called by reference, and start out |
|||
uninitialized. *) |
|||
q : &g0int tk? >> g0int tk, |
|||
r : &g0int tk? >> g0int tk) |
|||
: void = |
|||
let |
|||
(* The C optimizer most likely will reduce these these two |
|||
divisions to just one. They are simply synonyms for C '/' and |
|||
'%', and perform division that rounds the quotient towards |
|||
zero. *) |
|||
val q0 = g0int_div (n, d) |
|||
val r0 = g0int_mod (n, d) |
|||
in |
|||
(* The following calculation results in 'floor division', if the |
|||
divisor is positive, or 'ceiling division', if the divisor is |
|||
negative. This choice of method results in the remainder never |
|||
being negative. *) |
|||
if isgtez n || iseqz r0 then |
|||
(q := q0; r := r0) |
|||
else if isltz d then |
|||
(q := succ q0; r := r0 - d) |
|||
else |
|||
(q := pred q0; r := r0 + d) |
|||
end |
|||
fn {tk : tkind} |
|||
inverse (a : g0int tk, n : g0int tk, |
|||
inverse_exists : &bool? >> bool exists, |
|||
inverse_value : &g0int tk? >> opt (g0int tk, exists)) |
|||
: #[exists: bool] void = |
|||
let |
|||
typedef integer = g0int tk |
|||
fun |
|||
loop (t : integer, newt : integer, |
|||
r : integer, newr : integer, |
|||
inverse_exists : &bool? >> bool exists, |
|||
inverse_value : &g0int tk? >> opt (g0int tk, exists)) |
|||
: #[exists: bool] void = |
|||
if iseqz newr then |
|||
begin |
|||
if r > g0i2i 1 then |
|||
let |
|||
val () = inverse_exists := false |
|||
prval () = opt_none inverse_value |
|||
in |
|||
end |
|||
else if t < g0i2i 0 then |
|||
let |
|||
val () = inverse_exists := true |
|||
val () = inverse_value := t + n |
|||
prval () = opt_some inverse_value |
|||
in |
|||
end |
|||
else |
|||
let |
|||
val () = inverse_exists := true |
|||
val () = inverse_value := t |
|||
prval () = opt_some inverse_value |
|||
in |
|||
end |
|||
end |
|||
else |
|||
let |
|||
(* These become C variables. *) |
|||
var quotient : g0int tk? |
|||
var remainder : g0int tk? |
|||
val () = |
|||
division_with_nonnegative_remainder |
|||
(r, newr, quotient, remainder) |
|||
val t = newt |
|||
and newt = t - (quotient * newt) |
|||
and r = newr |
|||
and newr = remainder |
|||
in |
|||
loop (t, newt, r, newr, inverse_exists, inverse_value) |
|||
end |
|||
in |
|||
loop (g0i2i 0, g0i2i 1, n, a, inverse_exists, inverse_value) |
|||
end |
|||
implement |
|||
main0 () = |
|||
let |
|||
var inverse_exists : bool? |
|||
var inverse_value : llint? |
|||
in |
|||
inverse (42LL, 2017LL, inverse_exists, inverse_value); |
|||
if inverse_exists then |
|||
let |
|||
prval () = opt_unsome inverse_value |
|||
in |
|||
println! inverse_value |
|||
end |
|||
else |
|||
let |
|||
prval () = opt_unnone inverse_value |
|||
in |
|||
println! "There is no inverse." |
|||
end |
|||
end |
|||
</syntaxhighlight> |
|||
{{out}} |
|||
There is no need to tell patscc what allocator to use, because none is used. |
|||
<pre>$ patscc -g -O2 modular-inverse-noheap.dats && ./a.out |
|||
1969</pre> |
|||
=={{header|AutoHotkey}}== |
=={{header|AutoHotkey}}== |
||
Translation of [http://rosettacode.org/wiki/Modular_inverse#C C]. |
Translation of [http://rosettacode.org/wiki/Modular_inverse#C C]. |
||
< |
<syntaxhighlight lang="autohotkey">MsgBox, % ModInv(42, 2017) |
||
ModInv(a, b) { |
ModInv(a, b) { |
||
Line 137: | Line 481: | ||
x1 += b0 |
x1 += b0 |
||
return x1 |
return x1 |
||
}</ |
}</syntaxhighlight> |
||
{{out}} |
{{out}} |
||
<pre>1969</pre> |
<pre>1969</pre> |
||
=={{header|AWK}}== |
=={{header|AWK}}== |
||
<syntaxhighlight lang="awk"> |
|||
<lang AWK> |
|||
# syntax: GAWK -f MODULAR_INVERSE.AWK |
# syntax: GAWK -f MODULAR_INVERSE.AWK |
||
# converted from C |
# converted from C |
||
Line 170: | Line 514: | ||
return(x1) |
return(x1) |
||
} |
} |
||
</syntaxhighlight> |
|||
</lang> |
|||
{{out}} |
{{out}} |
||
<pre> |
<pre> |
||
1969 |
1969 |
||
</pre> |
</pre> |
||
=={{header|BASIC}}== |
|||
==={{header|ASIC}}=== |
|||
{{trans|Pascal}} |
|||
<syntaxhighlight lang="basic"> |
|||
REM Modular inverse |
|||
E = 42 |
|||
T = 2017 |
|||
GOSUB CalcModInv: |
|||
PRINT ModInv |
|||
END |
|||
CalcModInv: |
|||
REM Increments E Step times until Bal is greater than T |
|||
REM Repeats until Bal = 1 (MOD = 1) and returns Count |
|||
REM Bal will not be greater than T + E |
|||
D = 0 |
|||
IF E < T THEN |
|||
Bal = E |
|||
Count = 1 |
|||
Loop: |
|||
Step = T - Bal |
|||
Step = Step / E |
|||
Step = Step + 1 |
|||
REM So ... Step = (T - Bal) / E + 1 |
|||
StepTimesE = Step * E |
|||
Bal = Bal + StepTimesE |
|||
Count = Count + Step |
|||
Bal = Bal - T |
|||
IF Bal <> 1 THEN Loop: |
|||
D = Count |
|||
ENDIF |
|||
ModInv = D |
|||
RETURN |
|||
</syntaxhighlight> |
|||
{{out}} |
|||
<pre> |
|||
1969 |
|||
</pre> |
|||
==={{header|BASIC256}}=== |
|||
<syntaxhighlight lang="basic">print multInv(42, 2017) |
|||
end |
|||
function multInv(a,b) |
|||
x0 = 0 |
|||
b0 = b |
|||
multInv = 1 |
|||
if b = 1 then return |
|||
while a > 1 |
|||
q = a / b |
|||
t = b |
|||
b = a mod b |
|||
a = t |
|||
t = x0 |
|||
x0 = multInv - q * x0 |
|||
multInv = int(t) |
|||
end while |
|||
if multInv < 0 then return multInv + b0 |
|||
end function</syntaxhighlight> |
|||
{{out}} |
|||
<pre>1969</pre> |
|||
==={{header|Chipmunk Basic}}=== |
|||
{{works with|Chipmunk Basic|3.6.4}} |
|||
{{works with|QBasic}} |
|||
{{trans|BASIC256}} |
|||
<syntaxhighlight lang="qbasic">10 CLS |
|||
20 CALL modularinverse(42, 2017) |
|||
30 CALL modularinverse(40, 1) |
|||
40 END |
|||
50 SUB modularinverse(e,t) |
|||
60 d = 0 |
|||
70 IF e < t THEN |
|||
80 b = e |
|||
90 c = 1 |
|||
100 WHILE b > 1 |
|||
110 s = INT(((t-b)/e)+1) |
|||
120 b = b+s*e |
|||
130 c = c+s |
|||
140 b = b-t |
|||
150 WEND |
|||
160 d = c |
|||
170 ENDIF |
|||
180 m = d |
|||
190 PRINT m |
|||
200 END SUB</syntaxhighlight> |
|||
==={{header|Minimal BASIC}}=== |
|||
{{trans|Pascal}} |
|||
{{works with|Applesoft BASIC}} |
|||
{{works with|Commodore BASIC|3.5}} |
|||
{{works with|Nascom ROM BASIC|4.7}} |
|||
<syntaxhighlight lang="gwbasic"> |
|||
10 REM Modular inverse |
|||
20 LET E = 42 |
|||
30 LET T = 2017 |
|||
40 GOSUB 500 |
|||
50 PRINT M |
|||
60 END |
|||
490 REM Calculate modular inverse |
|||
500 LET D = 0 |
|||
510 IF E >= T THEN 600 |
|||
520 LET B = E |
|||
530 LET C = 1 |
|||
540 LET S1 = INT((T-B)/E)+1 |
|||
550 LET B = B+S1*E |
|||
560 LET C = C+S1 |
|||
570 LET B = B-T |
|||
580 IF B <> 1 THEN 540 |
|||
590 LET D = C |
|||
600 LET M = D |
|||
610 RETURN |
|||
</syntaxhighlight> |
|||
==={{header|QBasic}}=== |
|||
The [[#Chipmunk_Basic|Chipmunk Basic]] solution works without any changes. |
|||
<syntaxhighlight lang="qbasic"></syntaxhighlight> |
|||
==={{header|True BASIC}}=== |
|||
<syntaxhighlight lang="qbasic">SUB modularinverse(e,t) |
|||
LET d = 0 |
|||
IF e < t then |
|||
LET b = e |
|||
LET c = 1 |
|||
DO WHILE b > 1 |
|||
LET s = int(((t-b)/e)+1) |
|||
LET b = b+s*e |
|||
LET c = c+s |
|||
LET b = b-t |
|||
LOOP |
|||
LET d = c |
|||
END IF |
|||
LET m = d |
|||
PRINT m |
|||
END SUB |
|||
CALL modularinverse(42,2017) |
|||
CALL modularinverse(40,1) |
|||
END</syntaxhighlight> |
|||
=={{header|Batch File}}== |
=={{header|Batch File}}== |
||
Line 180: | Line 664: | ||
{{trans|Perl}} |
{{trans|Perl}} |
||
< |
<syntaxhighlight lang="dos">@echo off |
||
setlocal enabledelayedexpansion |
setlocal enabledelayedexpansion |
||
%== Calls the "function" ==% |
%== Calls the "function" ==% |
||
Line 222: | Line 706: | ||
if !t! lss 0 set /a t+=b |
if !t! lss 0 set /a t+=b |
||
set %3=!t! |
set %3=!t! |
||
goto :EOF</ |
goto :EOF</syntaxhighlight> |
||
{{Out}} |
{{Out}} |
||
<pre>1969 |
<pre>1969 |
||
Line 229: | Line 713: | ||
121 |
121 |
||
-1</pre> |
-1</pre> |
||
=={{header|BCPL}}== |
|||
<syntaxhighlight lang="bcpl">get "libhdr" |
|||
let mulinv(a, b) = |
|||
b<0 -> mulinv(a, -b), |
|||
a<0 -> mulinv(b - (-a rem b), b), |
|||
valof |
|||
$( let t, nt, r, nr = 0, 1, b, a rem b |
|||
until nr = 0 |
|||
$( let tmp, q = ?, r / nr |
|||
tmp := nt ; nt := t - q*nt ; t := tmp |
|||
tmp := nr ; nr := r - q*nr ; r := tmp |
|||
$) |
|||
resultis r>1 -> -1, |
|||
t<0 -> t + b, |
|||
t |
|||
$) |
|||
let show(a, b) be |
|||
$( let mi = mulinv(a, b) |
|||
test mi>=0 |
|||
do writef("%N, %N -> %N*N", a, b, mi) |
|||
or writef("%N, %N -> no inverse*N", a, b) |
|||
$) |
|||
let start() be |
|||
$( show(42, 2017) |
|||
show(40, 1) |
|||
show(52, -217) |
|||
show(-486, 217) |
|||
show(40, 2018) |
|||
$)</syntaxhighlight> |
|||
{{out}} |
|||
<pre>42, 2017 -> 1969 |
|||
40, 1 -> 0 |
|||
52, -217 -> 96 |
|||
-486, 217 -> 121 |
|||
40, 2018 -> no inverse</pre> |
|||
=={{header|Bracmat}}== |
=={{header|Bracmat}}== |
||
{{trans|Julia}} |
{{trans|Julia}} |
||
< |
<syntaxhighlight lang="bracmat">( ( mod-inv |
||
= a b b0 x0 x1 q |
= a b b0 x0 x1 q |
||
. !arg:(?a.?b) |
. !arg:(?a.?b) |
||
Line 247: | Line 770: | ||
) |
) |
||
& out$(mod-inv$(42.2017)) |
& out$(mod-inv$(42.2017)) |
||
};</ |
};</syntaxhighlight> |
||
Output |
Output |
||
<pre>1969</pre> |
<pre>1969</pre> |
||
=={{header|C}}== |
=={{header|C}}== |
||
< |
<syntaxhighlight lang="c">#include <stdio.h> |
||
int mul_inv(int a, int b) |
int mul_inv(int a, int b) |
||
Line 271: | Line 794: | ||
printf("%d\n", mul_inv(42, 2017)); |
printf("%d\n", mul_inv(42, 2017)); |
||
return 0; |
return 0; |
||
}</ |
}</syntaxhighlight> |
||
The above method has some problems. Most importantly, when given a pair (a,b) with no solution, it generates an FP exception. When given <tt>b=1</tt>, it returns 1 which is not a valid result mod 1. When given negative a or b the results are incorrect. The following generates results that should match Pari/GP for numbers in the int range. |
The above method has some problems. Most importantly, when given a pair (a,b) with no solution, it generates an FP exception. When given <tt>b=1</tt>, it returns 1 which is not a valid result mod 1. When given negative a or b the results are incorrect. The following generates results that should match Pari/GP for numbers in the int range. |
||
{{trans|Perl}} |
{{trans|Perl}} |
||
< |
<syntaxhighlight lang="c">#include <stdio.h> |
||
int mul_inv(int a, int b) |
int mul_inv(int a, int b) |
||
Line 299: | Line 822: | ||
printf("%d\n", mul_inv(40, 2018)); |
printf("%d\n", mul_inv(40, 2018)); |
||
return 0; |
return 0; |
||
}</ |
}</syntaxhighlight> |
||
{{out}} |
{{out}} |
||
<pre> |
<pre> |
||
Line 308: | Line 831: | ||
-1 |
-1 |
||
</pre> |
</pre> |
||
=={{header|C sharp}}== |
|||
<syntaxhighlight lang="csharp">public class Program |
|||
{ |
|||
static void Main() |
|||
{ |
|||
System.Console.WriteLine(42.ModInverse(2017)); |
|||
} |
|||
} |
|||
public static class IntExtensions |
|||
{ |
|||
public static int ModInverse(this int a, int m) |
|||
{ |
|||
if (m == 1) return 0; |
|||
int m0 = m; |
|||
(int x, int y) = (1, 0); |
|||
while (a > 1) { |
|||
int q = a / m; |
|||
(a, m) = (m, a % m); |
|||
(x, y) = (y, x - q * y); |
|||
} |
|||
return x < 0 ? x + m0 : x; |
|||
} |
|||
}</syntaxhighlight> |
|||
=={{header|C++}}== |
=={{header|C++}}== |
||
===Iterative implementation=== |
|||
{{trans|C}} |
{{trans|C}} |
||
< |
<syntaxhighlight lang="cpp">#include <iostream> |
||
using namespace std; |
|||
int mul_inv(int a, int b) |
int mul_inv(int a, int b) |
||
Line 329: | Line 878: | ||
int main(void) { |
int main(void) { |
||
cout<<mul_inv(42, 2017)<<endl; |
std::cout << mul_inv(42, 2017) << std::endl; |
||
return 0; |
return 0; |
||
}</syntaxhighlight> |
|||
} |
|||
</lang> |
|||
Recursive implementation |
===Recursive implementation=== |
||
< |
<syntaxhighlight lang="cpp">#include <iostream> |
||
short ObtainMultiplicativeInverse(int a, int b, int s0 = 1, int s1 = 0) |
short ObtainMultiplicativeInverse(int a, int b, int s0 = 1, int s1 = 0) |
||
Line 346: | Line 894: | ||
std::cout << ObtainMultiplicativeInverse(42, 2017) << std::endl; |
std::cout << ObtainMultiplicativeInverse(42, 2017) << std::endl; |
||
return 0; |
return 0; |
||
}</syntaxhighlight> |
|||
} |
|||
</lang> |
|||
=={{header|Clojure}}== |
=={{header|Clojure}}== |
||
< |
<syntaxhighlight lang="lisp">(ns test-p.core |
||
(:require [clojure.math.numeric-tower :as math])) |
(:require [clojure.math.numeric-tower :as math])) |
||
Line 394: | Line 941: | ||
(println (mul_inv 40 2018)) |
(println (mul_inv 40 2018)) |
||
</syntaxhighlight> |
|||
</lang> |
|||
'''Output:''' |
'''Output:''' |
||
Line 404: | Line 951: | ||
No inverse since gcd is: 2 |
No inverse since gcd is: 2 |
||
</pre> |
</pre> |
||
=={{header|CLU}}== |
|||
{{trans|Perl}} |
|||
<syntaxhighlight lang="clu">mul_inv = proc (a, b: int) returns (int) signals (no_inverse) |
|||
if b<0 then b := -b end |
|||
if a<0 then a := b - (-a // b) end |
|||
t: int := 0 |
|||
nt: int := 1 |
|||
r: int := b |
|||
nr: int := a // b |
|||
while nr ~= 0 do |
|||
q: int := r / nr |
|||
t, nt := nt, t - q*nt |
|||
r, nr := nr, r - q*nr |
|||
end |
|||
if r>1 then signal no_inverse end |
|||
if t<0 then t := t+b end |
|||
return(t) |
|||
end mul_inv |
|||
start_up = proc () |
|||
pair = struct[a, b: int] |
|||
tests: sequence[pair] := sequence[pair]$ |
|||
[pair${a: 42, b: 2017}, |
|||
pair${a: 40, b: 1}, |
|||
pair${a: 52, b: -217}, |
|||
pair${a: -486, b: 217}, |
|||
pair${a: 40, b: 2018}] |
|||
po: stream := stream$primary_output() |
|||
for test: pair in sequence[pair]$elements(tests) do |
|||
stream$puts(po, int$unparse(test.a) || ", " |
|||
|| int$unparse(test.b) || " -> ") |
|||
stream$putl(po, int$unparse(mul_inv(test.a, test.b))) |
|||
except when no_inverse: |
|||
stream$putl(po, "no modular inverse") |
|||
end |
|||
end |
|||
end start_up</syntaxhighlight> |
|||
{{out}} |
|||
<pre>42, 2017 -> 1969 |
|||
40, 1 -> 0 |
|||
52, -217 -> 96 |
|||
-486, 217 -> 121 |
|||
40, 2018 -> no modular inverse</pre> |
|||
=={{header|Comal}}== |
|||
<syntaxhighlight lang="comal">0010 FUNC mulinv#(a#,b#) CLOSED |
|||
0020 IF b#<0 THEN b#:=-b# |
|||
0030 IF a#<0 THEN a#:=b#-(-a# MOD b#) |
|||
0040 t#:=0;nt#:=1;r#:=b#;nr#:=a# MOD b# |
|||
0050 WHILE nr#<>0 DO |
|||
0060 q#:=r# DIV nr# |
|||
0070 tmp#:=nt#;nt#:=t#-q#*nt#;t#:=tmp# |
|||
0080 tmp#:=nr#;nr#:=r#-q#*nr#;r#:=tmp# |
|||
0090 ENDWHILE |
|||
0100 IF r#>1 THEN RETURN -1 |
|||
0110 IF t#<0 THEN t#:+b# |
|||
0120 RETURN t# |
|||
0130 ENDFUNC mulinv# |
|||
0140 // |
|||
0150 WHILE NOT EOD DO |
|||
0160 READ a#,b# |
|||
0170 PRINT a#,", ",b#," -> ",mulinv#(a#,b#) |
|||
0180 ENDWHILE |
|||
0190 END |
|||
0200 // |
|||
0210 DATA 42,2017,40,1,52,-217,-486,217,40,2018</syntaxhighlight> |
|||
{{out}} |
|||
<pre>42, 2017 -> 1969 |
|||
40, 1 -> 0 |
|||
52, -217 -> 96 |
|||
-486, 217 -> 121 |
|||
40, 2018 -> -1</pre> |
|||
=={{header|Common Lisp}}== |
=={{header|Common Lisp}}== |
||
< |
<syntaxhighlight lang="lisp"> |
||
;; |
;; |
||
;; Calculates the GCD of a and b based on the Extended Euclidean Algorithm. The function also returns |
;; Calculates the GCD of a and b based on the Extended Euclidean Algorithm. The function also returns |
||
Line 430: | Line 1,051: | ||
(unless (= 1 r) (error "invmod: Values ~a and ~a are not coprimes." a m)) |
(unless (= 1 r) (error "invmod: Values ~a and ~a are not coprimes." a m)) |
||
s)) |
s)) |
||
</syntaxhighlight> |
|||
</lang> |
|||
{{out}} |
{{out}} |
||
Line 441: | Line 1,062: | ||
1969 |
1969 |
||
</pre> |
</pre> |
||
=={{header|Cowgol}}== |
|||
<syntaxhighlight lang="cowgol">include "cowgol.coh"; |
|||
sub mulinv(a: int32, b: int32): (t: int32) is |
|||
if b<0 then b := -b; end if; |
|||
if a<0 then a := b - (-a % b); end if; |
|||
t := 0; |
|||
var nt: int32 := 1; |
|||
var r := b; |
|||
var nr := a % b; |
|||
while nr != 0 loop |
|||
var q := r / nr; |
|||
var tmp := nt; nt := t - q*nt; t := tmp; |
|||
tmp := nr; nr := r - q*nr; r := tmp; |
|||
end loop; |
|||
if r>1 then t := -1; |
|||
elseif t<0 then t := t + b; |
|||
end if; |
|||
end sub; |
|||
record Pair is |
|||
a: int32; |
|||
b: int32; |
|||
end record; |
|||
var data: Pair[] := { |
|||
{42, 2017}, |
|||
{40, 1}, |
|||
{52, -217}, |
|||
{-486, 217}, |
|||
{40, 2018} |
|||
}; |
|||
var i: @indexof data := 0; |
|||
while i < @sizeof data loop |
|||
print_i32(data[i].a as uint32); |
|||
print(", "); |
|||
print_i32(data[i].b as uint32); |
|||
print(" -> "); |
|||
var mi := mulinv(data[i].a, data[i].b); |
|||
if mi<0 |
|||
then print("no inverse"); |
|||
else print_i32(mi as uint32); |
|||
end if; |
|||
print_nl(); |
|||
i := i + 1; |
|||
end loop;</syntaxhighlight> |
|||
{{out}} |
|||
<pre>42, 2017 -> 1969 |
|||
40, 1 -> 0 |
|||
52, 4294967079 -> 96 |
|||
4294966810, 217 -> 121 |
|||
40, 2018 -> no inverse</pre> |
|||
=={{header|Craft Basic}}== |
|||
<syntaxhighlight lang="basic">let e = 42 |
|||
let t = 2017 |
|||
gosub modularinverse |
|||
end |
|||
sub modularinverse |
|||
let d = 0 |
|||
if e < t then |
|||
let b = e |
|||
let c = 1 |
|||
do |
|||
let s = int(((t - b) / e) + 1) |
|||
let b = b + s * e |
|||
let c = c + s |
|||
let b = b - t |
|||
loop b <> 1 |
|||
let d = c |
|||
endif |
|||
let m = d |
|||
print m |
|||
return</syntaxhighlight> |
|||
{{out| Output}}<pre>1969</pre> |
|||
=={{header|Crystal}}== |
|||
{{trans|Ruby}} |
|||
<syntaxhighlight lang="ruby">def modinv(a0, m0) |
|||
return 1 if m0 == 1 |
|||
a, m = a0, m0 |
|||
x0, inv = 0, 1 |
|||
while a > 1 |
|||
inv -= (a // m) * x0 |
|||
a, m = m, a % m |
|||
x0, inv = inv, x0 |
|||
end |
|||
inv += m0 if inv < 0 |
|||
inv |
|||
end</syntaxhighlight> |
|||
{{out}} |
|||
<pre> |
|||
> modinv(42,2017) |
|||
=> 1969</pre> |
|||
=={{header|D}}== |
=={{header|D}}== |
||
{{trans|C}} |
{{trans|C}} |
||
< |
<syntaxhighlight lang="d">T modInverse(T)(T a, T b) pure nothrow { |
||
if (b == 1) |
if (b == 1) |
||
return 1; |
return 1; |
||
Line 466: | Line 1,198: | ||
import std.stdio; |
import std.stdio; |
||
writeln(modInverse(42, 2017)); |
writeln(modInverse(42, 2017)); |
||
}</ |
}</syntaxhighlight> |
||
{{out}} |
{{out}} |
||
<pre>1969</pre> |
<pre>1969</pre> |
||
=={{header|dc}}== |
|||
{{trans|C}} |
|||
This solution prints the inverse <code>u</code> only if it exists (<code>a*u = 1 mod m</code>). |
|||
<syntaxhighlight lang="dc"> |
|||
dc -e "[m=]P?dsm[a=]P?dsa1sv[dsb~rsqlbrldlqlv*-lvsdsvd0<x]dsxxldd[dlmr+]sx0>xdla*lm%[p]sx1=x" |
|||
</syntaxhighlight> |
|||
If <code>~</code> is not implemented, it can be replaced by <code>SdSnlnld/LnLd%</code>. |
|||
Replace <code>[p]sx1=x</code> at the end by <code>[pq]sx1=x16i6E6F7420636F7072696D65P</code> if an error message "not coprime" is desired. |
|||
{{out}} |
|||
<pre> |
|||
m=2 800^1+ |
|||
a=37 |
|||
342411551958695219479776173037037562556082184118925013641969995739234\ |
|||
344644689214483533004909620355470582887300743869103978073598454778206\ |
|||
829469635119691272637318902731800747596752668736012071540136041369140\ |
|||
1228044652005748974399041408477572 |
|||
</pre> |
|||
<pre> |
|||
m=2017 |
|||
a=42 |
|||
1969 |
|||
</pre> |
|||
<pre> |
|||
m=42 |
|||
a=7 |
|||
</pre> |
|||
=={{header|Delphi}}== |
|||
See [[#Pascal]]. |
|||
=={{header|Draco}}== |
|||
<syntaxhighlight lang="draco">proc mulinv(int a, b) int: |
|||
int t, nt, r, nr, q, tmp; |
|||
if b<0 then b := -b fi; |
|||
if a<0 then a := b - (-a % b) fi; |
|||
t := 0; nt := 1; r := b; nr := a % b; |
|||
while nr /= 0 do |
|||
q := r / nr; |
|||
tmp := nt; nt := t - q*nt; t := tmp; |
|||
tmp := nr; nr := r - q*nr; r := tmp |
|||
od; |
|||
if r>1 then -1 |
|||
elif t<0 then t+b |
|||
else t |
|||
fi |
|||
corp |
|||
proc show(int a, b) void: |
|||
int mi; |
|||
mi := mulinv(a, b); |
|||
if mi>=0 |
|||
then writeln(a:5, ", ", b:5, " -> ", mi:5) |
|||
else writeln(a:5, ", ", b:5, " -> no inverse") |
|||
fi |
|||
corp |
|||
proc main() void: |
|||
show(42, 2017); |
|||
show(40, 1); |
|||
show(52, -217); |
|||
show(-486, 217); |
|||
show(40, 2018) |
|||
corp</syntaxhighlight> |
|||
{{out}} |
|||
<pre> 42, 2017 -> 1969 |
|||
40, 1 -> 0 |
|||
52, -217 -> 96 |
|||
-486, 217 -> 121 |
|||
40, 2018 -> no inverse</pre> |
|||
=={{header|EasyLang}}== |
|||
{{trans|AWK}} |
|||
<syntaxhighlight lang="easylang"> |
|||
func mod_inv a b . |
|||
b0 = b |
|||
x1 = 1 |
|||
if b = 1 |
|||
return 1 |
|||
. |
|||
while a > 1 |
|||
q = a div b |
|||
t = b |
|||
b = a mod b |
|||
a = t |
|||
t = x0 |
|||
x0 = x1 - q * x0 |
|||
x1 = t |
|||
. |
|||
if x1 < 0 |
|||
x1 += b0 |
|||
. |
|||
return x1 |
|||
. |
|||
print mod_inv 42 2017 |
|||
</syntaxhighlight> |
|||
=={{header|EchoLisp}}== |
=={{header|EchoLisp}}== |
||
< |
<syntaxhighlight lang="scheme"> |
||
(lib 'math) ;; for egcd = extended gcd |
(lib 'math) ;; for egcd = extended gcd |
||
Line 482: | Line 1,311: | ||
(mod-inv 42 666) |
(mod-inv 42 666) |
||
🔴 error: not-coprimes (42 666) |
🔴 error: not-coprimes (42 666) |
||
</syntaxhighlight> |
|||
</lang> |
|||
=={{header|Elixir}}== |
=={{header|Elixir}}== |
||
{{trans|Ruby}} |
{{trans|Ruby}} |
||
< |
<syntaxhighlight lang="elixir">defmodule Modular do |
||
def extended_gcd(a, b) do |
def extended_gcd(a, b) do |
||
{last_remainder, last_x} = extended_gcd(abs(a), abs(b), 1, 0, 0, 1) |
{last_remainder, last_x} = extended_gcd(abs(a), abs(b), 1, 0, 0, 1) |
||
Line 506: | Line 1,335: | ||
end |
end |
||
IO.puts Modular.inverse(42,2017)</ |
IO.puts Modular.inverse(42,2017)</syntaxhighlight> |
||
{{out}} |
{{out}} |
||
Line 514: | Line 1,343: | ||
=={{header|ERRE}}== |
=={{header|ERRE}}== |
||
< |
<syntaxhighlight lang="erre">PROGRAM MOD_INV |
||
!$INTEGER |
!$INTEGER |
||
Line 540: | Line 1,369: | ||
MUL_INV(40,2018->T) PRINT(T) |
MUL_INV(40,2018->T) PRINT(T) |
||
END PROGRAM |
END PROGRAM |
||
</syntaxhighlight> |
|||
</lang> |
|||
{{out}} |
{{out}} |
||
<pre> |
<pre> |
||
Line 551: | Line 1,380: | ||
=={{header|F_Sharp|F#}}== |
=={{header|F_Sharp|F#}}== |
||
< |
<syntaxhighlight lang="fsharp"> |
||
//Calculate the |
// Calculate the inverse of a (mod m) |
||
// See here for eea specs: |
|||
let MI n g = |
|||
// https://en.wikipedia.org/wiki/Extended_Euclidean_algorithm |
|||
let rec fN n i g e l a = |
|||
let modInv m a = |
|||
match e with |
|||
let rec eea t t' r r' = |
|||
match r' with |
|||
| 0 -> t |
|||
fN e l a (n-o*e) (i-o*l) (g-o*a) |
|||
| _ -> |
|||
let div = r/r' |
|||
</lang> |
|||
eea t' (t - div * t') r' (r - div * r') |
|||
(m + eea 0 1 m a) % m |
|||
</syntaxhighlight> |
|||
{{in}} |
|||
<pre> |
|||
// Inverse of 347 (mod 29) |
|||
modInv 29 347 |
|||
</pre> |
|||
{{out}} |
{{out}} |
||
<pre> |
<pre> |
||
28 |
|||
MI 2017 42 -> 1969 |
|||
</pre> |
</pre> |
||
=={{header|Factor}}== |
=={{header|Factor}}== |
||
<lang>USE: math.functions |
<syntaxhighlight lang="text">USE: math.functions |
||
42 2017 mod-inv</ |
42 2017 mod-inv</syntaxhighlight> |
||
{{out}} |
|||
<pre> |
|||
1969 |
|||
</pre> |
|||
=={{header|Forth}}== |
|||
ANS Forth with double-number word set |
|||
<syntaxhighlight lang="forth"> |
|||
: invmod { a m | v b c -- inv } |
|||
m to v |
|||
1 to c |
|||
0 to b |
|||
begin a |
|||
while v a / >r |
|||
c b s>d c s>d r@ 1 m*/ d- d>s to c to b |
|||
a v s>d a s>d r> 1 m*/ d- d>s to a to v |
|||
repeat b m mod dup to b 0< |
|||
if m b + else b then ; |
|||
</syntaxhighlight> |
|||
ANS Forth version without locals |
|||
<syntaxhighlight lang="forth"> |
|||
: modinv ( a m - inv) |
|||
dup 1- \ a m (m != 1)? |
|||
if \ a m |
|||
tuck 1 0 \ m0 a m 1 0 |
|||
begin \ m0 a m inv x0 |
|||
2>r over 1 > \ m0 a m (a > 1)? R: inv x0 |
|||
while \ m0 a m R: inv x0 |
|||
tuck /mod \ m0 m (a mod m) (a/m) R: inv x0 |
|||
r> tuck * \ m0 a' m' x0 (a/m)*x0 R: inv |
|||
r> swap - \ m0 a' m' x0 (inv-q) R: |
|||
repeat \ m0 a' m' inv' x0' |
|||
2drop \ m0 R: inv x0 |
|||
2r> drop \ m0 inv R: |
|||
dup 0< \ m0 inv (inv < 0)? |
|||
if over + then \ m0 (inv + m0) |
|||
then \ x inv' |
|||
nip \ inv |
|||
; |
|||
</syntaxhighlight> |
|||
<pre> |
|||
42 2017 invmod . 1969 |
|||
42 2017 modinv . 1969 |
|||
</pre> |
|||
=={{header|Fortran}}== |
|||
<syntaxhighlight lang="fortran"> |
|||
program modular_inverse_task |
|||
implicit none |
|||
write (*,*) inverse (42, 2017) |
|||
contains |
|||
! Returns -1 if there is no inverse. I assume n > 0. The algorithm |
|||
! is described at |
|||
! https://en.wikipedia.org/w/index.php?title=Extended_Euclidean_algorithm&oldid=1135569411#Modular_integers |
|||
function inverse (a, n) result (inverse_value) |
|||
integer, intent(in) :: a, n |
|||
integer :: inverse_value |
|||
integer :: t, newt |
|||
integer :: r, newr |
|||
integer :: quotient, remainder, tmp |
|||
if (n <= 0) error stop |
|||
t = 0; newt = 1 |
|||
r = n; newr = a |
|||
do while (newr /= 0) |
|||
remainder = modulo (r, newr) ! Floor division. |
|||
quotient = (r - remainder) / newr |
|||
tmp = newt; newt = t - (quotient * newt); t = tmp |
|||
r = newr; newr = remainder |
|||
end do |
|||
if (r > 1) then |
|||
inverse_value = -1 |
|||
else if (t < 0) then |
|||
inverse_value = t + n |
|||
else |
|||
inverse_value = t |
|||
end if |
|||
end function inverse |
|||
end program modular_inverse_task |
|||
</syntaxhighlight> |
|||
{{out}} |
|||
<pre>$ gfortran -Wall -Wextra modular_inverse_task.f90 && ./a.out |
|||
1969 |
|||
</pre> |
|||
=={{header|FreeBASIC}}== |
|||
<syntaxhighlight lang="freebasic">' version 10-07-2018 |
|||
' compile with: fbc -s console |
|||
Type ext_euclid |
|||
Dim As Integer a, b |
|||
End Type |
|||
' "Table method" aka "The Magic Box" |
|||
Function magic_box(x As Integer, y As Integer) As ext_euclid |
|||
Dim As Integer a(1 To 128), b(1 To 128), d(1 To 128), k(1 To 128) |
|||
a(1) = 1 : b(1) = 0 : d(1) = x |
|||
a(2) = 0 : b(2) = 1 : d(2) = y : k(2) = x \ y |
|||
Dim As Integer i = 2 |
|||
While Abs(d(i)) <> 1 |
|||
i += 1 |
|||
a(i) = a(i -2) - k(i -1) * a(i -1) |
|||
b(i) = b(i -2) - k(i -1) * b(i -1) |
|||
d(i) = d(i -2) Mod d(i -1) |
|||
k(i) = d(i -1) \ d(i) |
|||
'Print a(i),b(i),d(i),k(i) |
|||
If d(i -1) Mod d(i) = 0 Then Exit While |
|||
Wend |
|||
If d(i) = -1 Then ' -1 * (ab + by) = -1 * -1 ==> -ab -by = 1 |
|||
a(i) = -a(i) |
|||
b(i) = -b(i) |
|||
End If |
|||
Function = Type( a(i), b(i) ) |
|||
End Function |
|||
' ------=< MAIN >=------ |
|||
Dim As Integer x, y, gcd |
|||
Dim As ext_euclid result |
|||
Do |
|||
Read x, y |
|||
If x = 0 AndAlso y = 0 Then Exit Do |
|||
result = magic_box(x, y) |
|||
With result |
|||
gcd = .a * x + .b * y |
|||
Print "a * "; Str(x); " + b * "; Str(y); |
|||
Print " = GCD("; Str(x); ", "; Str(y); ") ="; gcd |
|||
If gcd > 1 Then |
|||
Print "No solution, numbers are not coprime" |
|||
Else |
|||
Print "a = "; .a; ", b = ";.b |
|||
Print "The Modular inverse of "; x; " modulo "; y; " = "; |
|||
While .a < 0 : .a += IIf(y > 0, y, -y) : Wend |
|||
Print .a |
|||
'Print "The Modular inverse of "; y; " modulo "; x; " = "; |
|||
'While .b < 0 : .b += IIf(x > 0, x, -x) : Wend |
|||
'Print .b |
|||
End if |
|||
End With |
|||
Print |
|||
Loop |
|||
Data 42, 2017 |
|||
Data 40, 1 |
|||
Data 52, -217 |
|||
Data -486, 217 |
|||
Data 40, 2018 |
|||
Data 0, 0 |
|||
' empty keyboard buffer |
|||
While Inkey <> "" : Wend |
|||
Print : Print "hit any key to end program" |
|||
Sleep |
|||
End</syntaxhighlight> |
|||
{{out}} |
|||
<pre>a * 42 + b * 2017 = GCD(42, 2017) = 1 |
|||
a = -48, b = 1 |
|||
The Modular inverse of 42 modulo 2017 = 1969 |
|||
a * 40 + b * 1 = GCD(40, 1) = 1 |
|||
a = 0, b = 1 |
|||
The Modular inverse of 40 modulo 1 = 0 |
|||
a * 52 + b * -217 = GCD(52, -217) = 1 |
|||
a = 96, b = 23 |
|||
The Modular inverse of 52 modulo -217 = 96 |
|||
a * -486 + b * 217 = GCD(-486, 217) = 1 |
|||
a = -96, b = -215 |
|||
The Modular inverse of -486 modulo 217 = 121 |
|||
a * 40 + b * 2018 = GCD(40, 2018) = 2 |
|||
No solution, numbers are not coprime</pre> |
|||
=={{header|Frink}}== |
|||
<syntaxhighlight lang="frink">println[modInverse[42, 2017]]</syntaxhighlight> |
|||
{{out}} |
{{out}} |
||
<pre> |
<pre> |
||
Line 574: | Line 1,602: | ||
=={{header|FunL}}== |
=={{header|FunL}}== |
||
< |
<syntaxhighlight lang="funl">import integers.egcd |
||
def modinv( a, m ) = |
def modinv( a, m ) = |
||
Line 585: | Line 1,613: | ||
if res < 0 then res + m else res |
if res < 0 then res + m else res |
||
println( modinv(42, 2017) )</ |
println( modinv(42, 2017) )</syntaxhighlight> |
||
{{out}} |
{{out}} |
||
Line 591: | Line 1,619: | ||
<pre> |
<pre> |
||
1969 |
1969 |
||
</pre> |
|||
=={{header|Forth}}== |
|||
ANS Forth with double-number word set |
|||
<lang forth> |
|||
: invmod { a m | v b c -- inv } |
|||
m to v |
|||
1 to c |
|||
0 to b |
|||
begin a |
|||
while v a / >r |
|||
c b s>d c s>d r@ 1 m*/ d- d>s to c to b |
|||
a v s>d a s>d r> 1 m*/ d- d>s to a to v |
|||
repeat b m mod dup to b 0< |
|||
if m b + else b then ; |
|||
</lang> |
|||
<pre> |
|||
42 2017 invmod . 1969 |
|||
</pre> |
</pre> |
||
=={{header|Go}}== |
=={{header|Go}}== |
||
The standard library function uses the extended Euclidean algorithm internally. |
The standard library function uses the extended Euclidean algorithm internally. |
||
< |
<syntaxhighlight lang="go">package main |
||
import ( |
import ( |
||
Line 625: | Line 1,635: | ||
k := new(big.Int).ModInverse(a, m) |
k := new(big.Int).ModInverse(a, m) |
||
fmt.Println(k) |
fmt.Println(k) |
||
}</ |
}</syntaxhighlight> |
||
{{out}} |
{{out}} |
||
<pre> |
<pre> |
||
1969 |
1969 |
||
</pre> |
|||
=={{header|GW-BASIC}}== |
|||
{{trans|Pascal}} |
|||
{{works with|PC-BASIC|any}} |
|||
<syntaxhighlight lang="qbasic"> |
|||
10 ' Modular inverse |
|||
20 LET E% = 42 |
|||
30 LET T% = 2017 |
|||
40 GOSUB 1000 |
|||
50 PRINT MODINV% |
|||
60 END |
|||
990 ' increments e stp (step) times until bal is greater than t |
|||
992 ' repeats until bal = 1 (mod = 1) and returns count |
|||
994 ' bal will not be greater than t + e |
|||
1000 LET D% = 0 |
|||
1010 IF E% >= T% THEN GOTO 1140 |
|||
1020 LET BAL% = E% |
|||
1025 ' At least one iteration is necessary |
|||
1030 LET STP% = ((T% - BAL%) \ E%) + 1 |
|||
1040 LET BAL% = BAL% + STP% * E% |
|||
1050 LET COUNT% = 1 + STP% |
|||
1060 LET BAL% = BAL% - T% |
|||
1070 WHILE BAL% <> 1 |
|||
1080 LET STP% = ((T% - BAL%) \ E%) + 1 |
|||
1090 LET BAL% = BAL% + STP% * E% |
|||
1100 LET COUNT% = COUNT% + STP% |
|||
1110 LET BAL% = BAL% - T% |
|||
1120 WEND |
|||
1130 LET D% = COUNT% |
|||
1140 LET MODINV% = D% |
|||
1150 RETURN |
|||
</syntaxhighlight> |
|||
{{out}} |
|||
<pre> |
|||
1969 |
|||
</pre> |
</pre> |
||
=={{header|Haskell}}== |
=={{header|Haskell}}== |
||
< |
<syntaxhighlight lang="haskell">-- Given a and m, return Just x such that ax = 1 mod m. |
||
-- If there is no such x return Nothing. |
|||
-- such that ax + by = g, where g = gcd(a,b). Note that x or y may be negative. |
|||
modInv :: Int -> Int -> Maybe Int |
|||
gcdExt a 0 = (1, 0, a) |
|||
modInv a m |
|||
gcdExt a b = let (q, r) = a `quotRem` b |
|||
| 1 == g = Just (mkPos i) |
|||
| otherwise = Nothing |
|||
in (t, s - q * t, g) |
|||
where |
|||
(i, _, g) = gcdExt a m |
|||
mkPos x |
|||
| x < 0 = x + m |
|||
| otherwise = x |
|||
-- Extended Euclidean algorithm. |
|||
-- Given a and m, return Just x such that ax = 1 mod m. If there is no such x |
|||
-- return |
-- Given non-negative a and b, return x, y and g |
||
-- such that ax + by = g, where g = gcd(a,b). |
|||
-- Note that x or y may be negative. |
|||
in if g == 1 then Just (mkPos i) else Nothing |
|||
gcdExt :: Int -> Int -> (Int, Int, Int) |
|||
where mkPos x = if x < 0 then x + m else x |
|||
gcdExt a 0 = (1, 0, a) |
|||
gcdExt a b = |
|||
let (q, r) = a `quotRem` b |
|||
(s, t, g) = gcdExt b r |
|||
in (t, s - q * t, g) |
|||
main |
main :: IO () |
||
print |
main = mapM_ print [2 `modInv` 4, 42 `modInv` 2017]</syntaxhighlight> |
||
print $ 42 `modInv` 2017</lang> |
|||
{{out}} |
{{out}} |
||
<pre>Nothing |
<pre>Nothing |
||
Line 654: | Line 1,710: | ||
=={{header|Icon}} and {{header|Unicon}}== |
=={{header|Icon}} and {{header|Unicon}}== |
||
{{trans|C}} |
{{trans|C}} |
||
< |
<syntaxhighlight lang="unicon">procedure main(args) |
||
a := integer(args[1]) | 42 |
a := integer(args[1]) | 42 |
||
b := integer(args[2]) | 2017 |
b := integer(args[2]) | 2017 |
||
Line 669: | Line 1,725: | ||
} |
} |
||
return if (x1 > 0) then x1 else x1+b0 |
return if (x1 > 0) then x1 else x1+b0 |
||
end</ |
end</syntaxhighlight> |
||
{{out}} |
{{out}} |
||
Line 680: | Line 1,736: | ||
Adding a coprime test: |
Adding a coprime test: |
||
< |
<syntaxhighlight lang="unicon">link numbers |
||
procedure main(args) |
procedure main(args) |
||
Line 698: | Line 1,754: | ||
} |
} |
||
return if (x1 > 0) then x1 else x1+b0 |
return if (x1 > 0) then x1 else x1+b0 |
||
end</ |
end</syntaxhighlight> |
||
=={{header|IS-BASIC}}== |
|||
<syntaxhighlight lang="is-basic">100 PRINT MODINV(42,2017) |
|||
120 DEF MODINV(A,B) |
|||
130 LET B=ABS(B) |
|||
140 IF A<0 THEN LET A=B-MOD(-A,B) |
|||
150 LET T=0:LET NT=1:LET R=B:LET NR=MOD(A,B) |
|||
160 DO WHILE NR<>0 |
|||
170 LET Q=INT(R/NR) |
|||
180 LET TMP=NT:LET NT=T-Q*NT:LET T=TMP |
|||
190 LET TMP=NR:LET NR=R-Q*NR:LET R=TMP |
|||
200 LOOP |
|||
210 IF R>1 THEN |
|||
220 LET MODINV=-1 |
|||
230 ELSE IF T<0 THEN |
|||
240 LET MODINV=T+B |
|||
250 ELSE |
|||
260 LET MODINV=T |
|||
270 END IF |
|||
280 END DEF</syntaxhighlight> |
|||
=={{header|J}}== |
=={{header|J}}== |
||
'''Solution''':< |
'''Solution''':<syntaxhighlight lang="j"> modInv =: dyad def 'x y&|@^ <: 5 p: y'"0</syntaxhighlight> |
||
'''Example''':< |
'''Example''':<syntaxhighlight lang="j"> 42 modInv 2017 |
||
1969</ |
1969</syntaxhighlight> |
||
'''Notes''': |
'''Notes''': |
||
Line 714: | Line 1,790: | ||
=={{header|Java}}== |
=={{header|Java}}== |
||
The <code>BigInteger</code> library has a method for this: |
The <code>BigInteger</code> library has a method for this: |
||
< |
<syntaxhighlight lang="java">System.out.println(BigInteger.valueOf(42).modInverse(BigInteger.valueOf(2017)));</syntaxhighlight> |
||
{{out}} |
{{out}} |
||
<pre>1969</pre> |
<pre>1969</pre> |
||
Alternatively, working from first principles. |
|||
<syntaxhighlight lang="java"> |
|||
public final class ModularInverse { |
|||
public static void main(String[] aArgs) { |
|||
System.out.println(inverseModulus(42, 2017)); |
|||
} |
|||
private static Egcd extendedGCD(int aOne, int aTwo) { |
|||
if ( aOne == 0 ) { |
|||
return new Egcd(aTwo, 0, 1); |
|||
} |
|||
Egcd value = extendedGCD(aTwo % aOne, aOne); |
|||
return new Egcd(value.aG, value.aY - ( aTwo / aOne ) * value.aX, value.aX); |
|||
} |
|||
private static int inverseModulus(int aNumber, int aModulus) { |
|||
Egcd value = extendedGCD(aNumber, aModulus); |
|||
return ( value.aG == 1 ) ? ( value.aX + aModulus ) % aModulus : 0; |
|||
} |
|||
private static record Egcd(int aG, int aX, int aY) {} |
|||
} |
|||
</syntaxhighlight> |
|||
{{ out }} |
|||
<pre> |
|||
1969 |
|||
</pre> |
|||
=={{header|JavaScript}}== |
=={{header|JavaScript}}== |
||
Using brute force. |
Using brute force. |
||
< |
<syntaxhighlight lang="javascript">var modInverse = function(a, b) { |
||
a %= b; |
a %= b; |
||
for (var x = 1; x < b; x++) { |
for (var x = 1; x < b; x++) { |
||
Line 727: | Line 1,832: | ||
} |
} |
||
} |
} |
||
}</ |
}</syntaxhighlight> |
||
=={{header|jq}}== |
|||
{{works with|jq}} |
|||
'''Works with gojq, the Go implementation of jq''' |
|||
<syntaxhighlight lang="jq"># Integer division: |
|||
# If $j is 0, then an error condition is raised; |
|||
# otherwise, assuming infinite-precision integer arithmetic, |
|||
# if the input and $j are integers, then the result will be an integer. |
|||
def idivide($j): |
|||
. as $i |
|||
| ($i % $j) as $mod |
|||
| ($i - $mod) / $j ; |
|||
# the multiplicative inverse of . modulo $n |
|||
def modInv($n): |
|||
if $n == 1 then 1 |
|||
else . as $this |
|||
| { r : $n, |
|||
t : 0, |
|||
newR: length, # abs |
|||
newT: 1} |
|||
| until(.newR == 0; |
|||
.newR as $newR |
|||
| (.r | idivide($newR)) as $q |
|||
| {r : $newR, |
|||
t : .newT, |
|||
newT: (.t - $q * .newT), |
|||
newR: (.r - $q * $newR) } ) |
|||
| if (.r|length) != 1 then "\($this) and \($n) are not co-prime." | error |
|||
else .t |
|||
| if . < 0 then . + $n |
|||
elif $this < 0 then - . |
|||
else . |
|||
end |
|||
end |
|||
end ; |
|||
# Example: |
|||
42 | modInv(2017) |
|||
</syntaxhighlight> |
|||
{{out}} |
|||
<pre> |
|||
1969 |
|||
</pre> |
|||
=={{header|Julia}}== |
=={{header|Julia}}== |
||
{{works with|Julia|1.2}} |
|||
===Built-in=== |
===Built-in=== |
||
Julia includes a built-in function for this: |
Julia includes a built-in function for this: |
||
<lang |
<syntaxhighlight lang="julia">invmod(a, b)</syntaxhighlight> |
||
===C translation=== |
===C translation=== |
||
Line 738: | Line 1,891: | ||
The following code works on any integer type. |
The following code works on any integer type. |
||
To maximize performance, we ensure (via a promotion rule) that the operands are the same type (and use built-ins <code>zero(T)</code> and <code>one(T)</code> for initialization of temporary variables to ensure that they remain of the same type throughout execution). |
To maximize performance, we ensure (via a promotion rule) that the operands are the same type (and use built-ins <code>zero(T)</code> and <code>one(T)</code> for initialization of temporary variables to ensure that they remain of the same type throughout execution). |
||
< |
<syntaxhighlight lang="julia">function modinv(a::T, b::T) where T <: Integer |
||
b0 = b |
b0 = b |
||
x0, x1 = zero(T), one(T) |
x0, x1 = zero(T), one(T) |
||
Line 748: | Line 1,901: | ||
x1 < 0 ? x1 + b0 : x1 |
x1 < 0 ? x1 + b0 : x1 |
||
end |
end |
||
modinv(a::Integer, b::Integer) = modinv(promote(a,b)...)</ |
modinv(a::Integer, b::Integer) = modinv(promote(a,b)...)</syntaxhighlight> |
||
{{out}} |
{{out}} |
||
Line 760: | Line 1,913: | ||
=={{header|Kotlin}}== |
=={{header|Kotlin}}== |
||
< |
<syntaxhighlight lang="scala">// version 1.0.6 |
||
import java.math.BigInteger |
import java.math.BigInteger |
||
Line 768: | Line 1,921: | ||
val m = BigInteger.valueOf(2017) |
val m = BigInteger.valueOf(2017) |
||
println(a.modInverse(m)) |
println(a.modInverse(m)) |
||
}</ |
}</syntaxhighlight> |
||
{{out}} |
{{out}} |
||
Line 774: | Line 1,927: | ||
1969 |
1969 |
||
</pre> |
</pre> |
||
=={{header|Lambdatalk}}== |
|||
{{trans|Phix}} |
|||
<syntaxhighlight lang="scheme"> |
|||
{def mulinv |
|||
{def mulinv.loop |
|||
{lambda {:t :nt :r :nr} |
|||
{if {not {= :nr 0}} |
|||
then {mulinv.loop :nt |
|||
{- :t {* {floor {/ :r :nr}} :nt}} |
|||
:nr |
|||
{- :r {* {floor {/ :r :nr}} :nr}} } |
|||
else {cons :t :r} }}} |
|||
{lambda {:a :n} |
|||
{let { {:a :a} {:n :n} |
|||
{:cons {mulinv.loop 0 |
|||
1 |
|||
{if {< :n 0} then {- :n} else :n} |
|||
{if {< :a 0} then {- :n {% {- :a} :n}} else :a}}} |
|||
} {if {> {cdr :cons} 1} |
|||
then not invertible |
|||
else {if {< {car :cons} 0} |
|||
then {+ {car :cons} :n} |
|||
else {car :cons} }}}}} |
|||
-> mulinv |
|||
{mulinv 42 2017} |
|||
-> 1969 |
|||
{mulinv 40 1} |
|||
-> 0 |
|||
{mulinv 52 -217} |
|||
-> 96 |
|||
{mulinv -486 217} |
|||
-> 121 |
|||
{mulinv 40 218} |
|||
-> not invertible |
|||
</syntaxhighlight> |
|||
=={{header|m4}}== |
|||
{{trans|Mercury}} |
|||
Note that <code>$0</code> is the name of the macro being evaluated. Therefore, in the following, <code>_$0</code> is the name of another macro, the same as the name of the first macro, except for an underscore prepended. This is a common idiom. |
|||
The core of the program is <code>__inverse</code> recursively calling itself. |
|||
m4 is a macro-preprocessor, and so some of the following is there simply to keep things from being echoed to the output. :) |
|||
<syntaxhighlight lang="m4"> |
|||
divert(-1) |
|||
# I assume non-negative arguments. The algorithm is described at |
|||
# https://en.wikipedia.org/w/index.php?title=Extended_Euclidean_algorithm&oldid=1135569411#Modular_integers |
|||
define(`inverse',`_$0(eval(`$1'), eval(`$2'))') |
|||
define(`_inverse',`_$0($2, 0, 1, $2, $1)') |
|||
define(`__inverse', |
|||
`dnl n = $1, t = $2, newt = $3, r = $4, newr = $5 |
|||
ifelse(eval($5 != 0), 1, `$0($1, $3, |
|||
eval($2 - (($4 / $5) * $3)), |
|||
$5,eval($4 % $5))', |
|||
eval($4 > 1), 1, `no inverse', |
|||
eval($2 < 0), 1, eval($2 + $1), |
|||
$2)') |
|||
divert`'dnl |
|||
inverse(42, 2017) |
|||
</syntaxhighlight> |
|||
{{out}} |
|||
<pre>$ m4 modular-inverse-task.m4 |
|||
1969</pre> |
|||
=={{header|MAD}}== |
|||
<syntaxhighlight lang="mad"> NORMAL MODE IS INTEGER |
|||
INTERNAL FUNCTION(AA, BB) |
|||
ENTRY TO MULINV. |
|||
A = AA |
|||
B = BB |
|||
WHENEVER B.L.0, B = -B |
|||
WHENEVER A.L.0, A = B - (-(A-A/B*B)) |
|||
T = 0 |
|||
NT = 1 |
|||
R = B |
|||
NR = A-A/B*B |
|||
LOOP WHENEVER NR.NE.0 |
|||
Q = R/NR |
|||
TMP = NT |
|||
NT = T - Q*NT |
|||
T = TMP |
|||
TMP = NR |
|||
NR = R - Q*NR |
|||
R = TMP |
|||
TRANSFER TO LOOP |
|||
END OF CONDITIONAL |
|||
WHENEVER R.G.1, FUNCTION RETURN -1 |
|||
WHENEVER T.L.0, T = T+B |
|||
FUNCTION RETURN T |
|||
END OF FUNCTION |
|||
INTERNAL FUNCTION(AA, BB) |
|||
VECTOR VALUES FMT = $I5,2H, ,I5,2H: ,I5*$ |
|||
ENTRY TO SHOW. |
|||
PRINT FORMAT FMT, AA, BB, MULINV.(AA, BB) |
|||
END OF FUNCTION |
|||
SHOW.(42,2017) |
|||
SHOW.(40,1) |
|||
SHOW.(52,-217) |
|||
SHOW.(-486,217) |
|||
SHOW.(40,2018) |
|||
END OF PROGRAM</syntaxhighlight> |
|||
{{out}} |
|||
<pre> 42, 2017: 1969 |
|||
40, 1: 0 |
|||
52, -217: 96 |
|||
-486, 217: 121 |
|||
40, 2018: -1</pre> |
|||
=={{header|Maple}}== |
=={{header|Maple}}== |
||
<syntaxhighlight lang="maple"> |
|||
<lang Maple> |
|||
1/42 mod 2017; |
1/42 mod 2017; |
||
</syntaxhighlight> |
|||
</lang> |
|||
{{out}} |
{{out}} |
||
<pre> |
<pre> |
||
Line 785: | Line 2,057: | ||
</pre> |
</pre> |
||
=={{header|Mathematica}}== |
=={{header|Mathematica}}/{{header|Wolfram Language}}== |
||
Use the built-in function <code>ModularInverse</code>: |
|||
<lang Mathematica>modInv[a_, m_] := |
|||
Block[{x,k}, x /. FindInstance[a x == 1 + k m, {x, k}, Integers]]</lang> |
|||
<syntaxhighlight lang="mathematica">ModularInverse[a, m]</syntaxhighlight> |
|||
Another way by using the built-in function <code>PowerMod</code> : |
|||
For example: |
|||
<lang Mathematica>PowerMod[a,-1,m]</lang> |
|||
<pre>ModularInverse[42, 2017] |
|||
For example : |
|||
1969</pre> |
|||
<pre>modInv[42, 2017] |
|||
=={{header|Maxima}}== |
|||
{1969} |
|||
Using built-in function inv_mod |
|||
<syntaxhighlight lang="maxima"> |
|||
inv_mod(42,2017); |
|||
</syntaxhighlight> |
|||
{{out}} |
|||
<pre> |
|||
1969 |
|||
</pre> |
|||
=={{header|Mercury}}== |
|||
PowerMod[42, -1, 2017] |
|||
{{works with|Mercury|22.01.1}} |
|||
<syntaxhighlight lang="mercury"> |
|||
1969</pre> |
|||
%%% -*- mode: mercury; prolog-indent-width: 2; -*- |
|||
%%% |
|||
%%% Compile with: |
|||
%%% mmc --make --use-subdirs modular_inverse_task |
|||
%%% |
|||
:- module modular_inverse_task. |
|||
:- interface. |
|||
:- import_module io. |
|||
:- pred main(io::di, io::uo) is det. |
|||
:- implementation. |
|||
:- import_module exception. |
|||
:- import_module int. |
|||
%% inverse(A, N, Inverse). I assume N > 0, and throw an exception if |
|||
%% it is not. The predicate fails if there is no inverse (and thus is |
|||
%% "semidet"). The algorithm is described at |
|||
%% https://en.wikipedia.org/w/index.php?title=Extended_Euclidean_algorithm&oldid=1135569411#Modular_integers |
|||
:- pred inverse(int::in, int::in, int::out) is semidet. |
|||
inverse(A, N, Inverse) :- |
|||
if (N =< 0) then throw(domain_error("inverse")) |
|||
else inverse_(N, 0, 1, N, A, Inverse). |
|||
:- pred inverse_(int::in, int::in, int::in, int::in, int::in, |
|||
int::out) is semidet. |
|||
inverse_(N, T, NewT, R, NewR, Inverse) :- |
|||
if (NewR \= 0) |
|||
then (Quotient = div(R, NewR), % Floor division. |
|||
inverse_(N, |
|||
NewT, T - (Quotient * NewT), |
|||
NewR, R - (Quotient * NewR), |
|||
Inverse)) % Tail recursion. |
|||
else (R =< 1, % R =< 1 FAILS if R > 1. |
|||
(if (T < 0) then Inverse = T + N else Inverse = T)). |
|||
main(!IO) :- |
|||
if inverse(42, 2017, Inverse) |
|||
then (print(Inverse, !IO), nl(!IO)) |
|||
else (print("There is no inverse.", !IO), nl(!IO)). |
|||
:- end_module modular_inverse_task. |
|||
</syntaxhighlight> |
|||
{{out}} |
|||
<pre>$ mmc --make --use-subdirs modular_inverse_task && ./modular_inverse_task |
|||
Making Mercury/int3s/modular_inverse_task.int3 |
|||
Making Mercury/ints/modular_inverse_task.int |
|||
Making Mercury/cs/modular_inverse_task.c |
|||
Making Mercury/os/modular_inverse_task.o |
|||
Making modular_inverse_task |
|||
1969 |
|||
</pre> |
|||
=={{header|МК-61/52}}== |
=={{header|МК-61/52}}== |
||
<lang>П1 П2 <-> П0 0 П5 1 П6 ИП1 1 |
<syntaxhighlight lang="text">П1 П2 <-> П0 0 П5 1 П6 ИП1 1 |
||
- x=0 14 С/П ИП0 1 - /-/ x<0 50 |
- x=0 14 С/П ИП0 1 - /-/ x<0 50 |
||
ИП0 ИП1 / [x] П4 ИП1 П3 ИП0 ^ ИП1 |
ИП0 ИП1 / [x] П4 ИП1 П3 ИП0 ^ ИП1 |
||
/ [x] ИП1 * - П1 ИП3 П0 ИП5 П3 |
/ [x] ИП1 * - П1 ИП3 П0 ИП5 П3 |
||
ИП6 ИП4 ИП5 * - П5 ИП3 П6 БП 14 |
ИП6 ИП4 ИП5 * - П5 ИП3 П6 БП 14 |
||
ИП6 x<0 55 ИП2 + С/П</ |
ИП6 x<0 55 ИП2 + С/П</syntaxhighlight> |
||
=={{header|Modula-2}}== |
|||
{{trans|C}} |
|||
<syntaxhighlight lang="modula-2">MODULE ModularInverse; |
|||
FROM InOut IMPORT WriteString, WriteInt, WriteLn; |
|||
TYPE Data = RECORD x : INTEGER; |
|||
y : INTEGER |
|||
END; |
|||
VAR c : INTEGER; |
|||
ab : ARRAY [1..5] OF Data; |
|||
PROCEDURE mi(VAR a, b : INTEGER): INTEGER; |
|||
VAR t, nt, r, nr, q, tmp : INTEGER; |
|||
BEGIN |
|||
b := ABS(b); |
|||
IF a < 0 THEN a := b - (-a MOD b) END; |
|||
t := 0; nt := 1; r := b; nr := a MOD b; |
|||
WHILE (nr # 0) DO |
|||
q := r / nr; |
|||
tmp := nt; nt := t - q * nt; t := tmp; |
|||
tmp := nr; nr := r - q * nr; r := tmp; |
|||
END; |
|||
IF (r > 1) THEN RETURN -1 END; |
|||
IF (t < 0) THEN RETURN t + b END; |
|||
RETURN t; |
|||
END mi; |
|||
BEGIN |
|||
ab[1].x := 42; ab[1].y := 2017; |
|||
ab[2].x := 40; ab[2].y := 1; |
|||
ab[3].x := 52; ab[3].y := -217; |
|||
ab[4].x := -486; ab[4].y := 217; |
|||
ab[5].x := 40; ab[5].y := 2018; |
|||
WriteLn; |
|||
WriteString("Modular inverse"); |
|||
WriteLn; |
|||
FOR c := 1 TO 5 DO |
|||
WriteInt(ab[c].x, 6); WriteString(", "); |
|||
WriteInt(ab[c].y, 6); WriteString(" = "); |
|||
WriteInt(mi(ab[c].x, ab[c].y),6); |
|||
WriteLn; |
|||
END; |
|||
END ModularInverse.</syntaxhighlight> |
|||
{{out}} |
|||
<pre>Modular inverse |
|||
42, 2017 = 1969 |
|||
40, 1 = 0 |
|||
52, -217 = 96 |
|||
-486, 217 = 121 |
|||
40, 2018 = -1</pre> |
|||
=={{header|newLISP}}== |
|||
<syntaxhighlight lang="newlisp"> |
|||
(define (modular-multiplicative-inverse a n) |
|||
(if (< n 0) |
|||
(setf n (abs n))) |
|||
(if (< a 0) |
|||
(setf a (- n (% (- 0 a) n)))) |
|||
(setf t 0) |
|||
(setf nt 1) |
|||
(setf r n) |
|||
(setf nr (mod a n)) |
|||
(while (not (zero? nr)) |
|||
(setf q (int (div r nr))) |
|||
(setf tmp nt) |
|||
(setf nt (sub t (mul q nt))) |
|||
(setf t tmp) |
|||
(setf tmp nr) |
|||
(setf nr (sub r (mul q nr))) |
|||
(setf r tmp)) |
|||
(if (> r 1) |
|||
(setf retvalue nil)) |
|||
(if (< t 0) |
|||
(setf retvalue (add t n)) |
|||
(setf retvalue t)) |
|||
retvalue) |
|||
(println (modular-multiplicative-inverse 42 2017)) |
|||
</syntaxhighlight> |
|||
Output: |
|||
<pre> |
|||
1969 |
|||
</pre> |
|||
=={{header|Nim}}== |
=={{header|Nim}}== |
||
{{trans|C}} |
{{trans|C}} |
||
< |
<syntaxhighlight lang="nim">proc modInv(a0, b0: int): int = |
||
var (a, b, x0) = (a0, b0, 0) |
var (a, b, x0) = (a0, b0, 0) |
||
result = 1 |
result = 1 |
||
if b == 1: return |
if b == 1: return |
||
while a > 1: |
while a > 1: |
||
result = result - (a div b) * x0 |
|||
a = a mod b |
a = a mod b |
||
swap a, b |
swap a, b |
||
result = result - q * x0 |
|||
swap x0, result |
swap x0, result |
||
if result < 0: result += b0 |
if result < 0: result += b0 |
||
echo |
echo modInv(42, 2017)</syntaxhighlight> |
||
{{out}} |
|||
<pre>1969</pre> |
|||
=={{header|Oberon-2}}== |
|||
{{works with|Oxford Oberon-2 Compiler|3.3.0}} |
|||
<syntaxhighlight lang="oberon"> |
|||
(*-*- mode: indented-text; tab-width: 2; -*-*) |
|||
MODULE modularInverseInOberon2; |
|||
IMPORT Out; |
|||
(* Division with a non-negative remainder. This will work no matter |
|||
how your compiler handles DIV (and mine seems not to do what the |
|||
Oberon-2 specification says). *) |
|||
PROCEDURE euclidDiv (x, y : INTEGER) : INTEGER; |
|||
VAR q : INTEGER; |
|||
BEGIN |
|||
IF 0 <= y THEN (* Do floor division. *) |
|||
IF 0 <= x THEN |
|||
q := x DIV y |
|||
ELSE |
|||
q := -((-x) DIV y); |
|||
IF (-x) MOD y # 0 THEN q := q - 1 END |
|||
END; |
|||
ELSE (* Do ceiling division. *) |
|||
IF 0 <= x THEN |
|||
q := -(x DIV (-y)) |
|||
ELSE |
|||
q := ((-x) DIV (-y)); |
|||
IF (-x) MOD (-y) # 0 THEN q := q + 1 END |
|||
END |
|||
END; |
|||
RETURN q |
|||
END euclidDiv; |
|||
(* I have added this unit test because, earlier, I posted a buggy |
|||
version of euclidDiv. *) |
|||
PROCEDURE testEuclidDiv; |
|||
VAR x, y, q, r : INTEGER; |
|||
BEGIN |
|||
FOR x := -100 TO 100 DO |
|||
FOR y := -100 TO 100 DO |
|||
IF y # 0 THEN |
|||
q := euclidDiv (x, y); |
|||
r := x - (q * y); |
|||
IF (r < 0) OR (ABS (y) <= r) THEN |
|||
(* A remainder was outside the expected range. *) |
|||
Out.String ("euclidDiv fails its test") |
|||
END |
|||
END |
|||
END |
|||
END |
|||
END testEuclidDiv; |
|||
PROCEDURE inverse (a, n : INTEGER) : INTEGER; |
|||
VAR t, newt : INTEGER; |
|||
VAR r, newr : INTEGER; |
|||
VAR quotient : INTEGER; |
|||
VAR tmp : INTEGER; |
|||
BEGIN |
|||
t := 0; newt := 1; |
|||
r := n; newr := a; |
|||
WHILE newr # 0 DO |
|||
quotient := euclidDiv (r, newr); |
|||
tmp := newt; newt := t - (quotient * newt); t := tmp; |
|||
tmp := newr; newr := r - (quotient * newr); r := tmp |
|||
END; |
|||
IF r > 1 THEN |
|||
t := -1 |
|||
ELSIF t < 0 THEN |
|||
t := t + n |
|||
END; |
|||
RETURN t |
|||
END inverse; |
|||
BEGIN |
|||
testEuclidDiv; |
|||
Out.Int (inverse (42, 2017), 0); |
|||
Out.Ln |
|||
END modularInverseInOberon2. |
|||
</syntaxhighlight> |
|||
{{out}} |
|||
<pre>$ obc modularInverseInOberon2.Mod && ./a.out |
|||
1969 |
|||
</pre> |
|||
=={{header|ObjectIcon}}== |
|||
{{trans|Oberon-2}} |
|||
{{trans|Mercury}} |
|||
<syntaxhighlight lang="objecticon"> |
|||
# -*- ObjectIcon -*- |
|||
import exception |
|||
import io |
|||
procedure main () |
|||
test_euclid_div () |
|||
io.write (inverse (42, 2017)) |
|||
end |
|||
procedure inverse (a, n) # FAILS if there is no inverse. |
|||
local t, newt, r, newr, quotient, tmp |
|||
if n <= 0 then throw ("non-positive modulus") |
|||
t := 0; newt := 1 |
|||
r := n; newr := a |
|||
while newr ~= 0 do |
|||
{ |
|||
quotient := euclid_div (r, newr) |
|||
tmp := newt; newt := t - (quotient * newt); t := tmp |
|||
tmp := newr; newr := r - (quotient * newr); r := tmp |
|||
} |
|||
r <= 1 | fail |
|||
return (if t < 0 then t + n else t) |
|||
end |
|||
procedure euclid_div (x, y) |
|||
# This kind of integer division always gives a remainder between 0 |
|||
# and abs(y)-1, inclusive. Thus the remainder is always a LEAST |
|||
# RESIDUE modulo abs(y). (If y is a positive modulus, then only the |
|||
# floor division branch is used.) |
|||
return \ |
|||
if 0 <= y then # Do floor division. |
|||
(if 0 <= x then x / y |
|||
else if (-x) % y = 0 then -((-x) / y) |
|||
else -((-x) / y) - 1) |
|||
else # Do ceiling division. |
|||
(if 0 <= x then -(x / (-y)) |
|||
else if (-x) % (-y) = 0 then ((-x) / (-y)) |
|||
else ((-x) / (-y)) + 1) |
|||
end |
|||
procedure test_euclid_div () |
|||
local x, y, q, r |
|||
every x := -100 to 100 do |
|||
every y := -100 to 100 & y ~= 0 do |
|||
{ |
|||
q := euclid_div (x, y) |
|||
r := x - (q * y) |
|||
if r < 0 | abs (y) <= r then |
|||
# A remainder was outside the expected range. |
|||
throw ("Test of euclid_div failed.") |
|||
} |
|||
end |
|||
</syntaxhighlight> |
|||
{{out}} |
|||
<pre>$ oiscript modular-inverse-task-OI.icn |
|||
1969</pre> |
|||
=={{header|OCaml}}== |
=={{header|OCaml}}== |
||
==={{trans|C}}=== |
==={{trans|C}}=== |
||
< |
<syntaxhighlight lang="ocaml">let mul_inv a = function 1 -> 1 | b -> |
||
let rec aux a b x0 x1 = |
let rec aux a b x0 x1 = |
||
if a <= 1 then x1 else |
if a <= 1 then x1 else |
||
Line 835: | Line 2,413: | ||
in |
in |
||
let x = aux a b 0 1 in |
let x = aux a b 0 1 in |
||
if x < 0 then x + b else x</ |
if x < 0 then x + b else x</syntaxhighlight> |
||
Testing: |
Testing: |
||
Line 845: | Line 2,423: | ||
==={{trans|Haskell}}=== |
==={{trans|Haskell}}=== |
||
< |
<syntaxhighlight lang="ocaml">let rec gcd_ext a = function |
||
| 0 -> (1, 0, a) |
| 0 -> (1, 0, a) |
||
| b -> |
| b -> |
||
Line 855: | Line 2,433: | ||
match gcd_ext a m with |
match gcd_ext a m with |
||
| i, _, 1 -> mk_pos i |
| i, _, 1 -> mk_pos i |
||
| _ -> failwith "mod_inv"</ |
| _ -> failwith "mod_inv"</syntaxhighlight> |
||
Testing: |
Testing: |
||
Line 867: | Line 2,445: | ||
Usage : a modulus invmod |
Usage : a modulus invmod |
||
< |
<syntaxhighlight lang="oforth">// euclid ( a b -- u v r ) |
||
// Return r = gcd(a, b) and (u, v) / r = au + bv |
// Return r = gcd(a, b) and (u, v) / r = au + bv |
||
Line 888: | Line 2,466: | ||
: invmod(a, modulus) |
: invmod(a, modulus) |
||
a modulus euclid 1 == ifFalse: [ drop drop null return ] |
a modulus euclid 1 == ifFalse: [ drop drop null return ] |
||
drop dup 0 < ifTrue: [ modulus + ] ;</ |
drop dup 0 < ifTrue: [ modulus + ] ;</syntaxhighlight> |
||
{{out}} |
{{out}} |
||
Line 895: | Line 2,473: | ||
1969 |
1969 |
||
</pre> |
</pre> |
||
=={{header|Owl Lisp}}== |
|||
{{works with|Owl Lisp|0.2.1}} |
|||
<syntaxhighlight lang="scheme"> |
|||
(import (owl math) |
|||
(owl math extra)) |
|||
(define (euclid-quotient x y) |
|||
(if (<= 0 y) |
|||
(cond ((<= 0 x) (quotient x y)) |
|||
((zero? (remainder (negate x) y)) |
|||
(negate (quotient (negate x) y))) |
|||
(else (- (negate (quotient (negate x) y)) 1))) |
|||
(cond ((<= 0 x) (negate (quotient x (negate y)))) |
|||
((zero? (remainder (negate x) (negate y))) |
|||
(quotient (negate x) (negate y))) |
|||
(else (+ (quotient (negate x) (negate y)) 1))))) |
|||
;; A unit test of euclid-quotient. |
|||
(let repeat ((x -100) |
|||
(y -100)) |
|||
(cond ((= x 101) #t) |
|||
((= y 0) (repeat x (+ y 1))) |
|||
((= y 101) (repeat (+ x 1) -100)) |
|||
(else (let* ((q (euclid-quotient x y)) |
|||
(r (- x (* q y)))) |
|||
(cond ((< r 0) (display "negative remainder\n")) |
|||
((<= (abs y) r) (display "remainder too large\n")) |
|||
(else (repeat x (+ y 1)))))))) |
|||
(define (inverse a n) |
|||
(let repeat ((t 0) (newt 1) |
|||
(r n) (newr a)) |
|||
(cond ((not (zero? newr)) |
|||
(let ((quotient (euclid-quotient r newr))) |
|||
(repeat newt (- t (* quotient newt)) |
|||
newr (- r (* quotient newr))))) |
|||
((< 1 r) #f) ; The inverse does not exist. |
|||
((negative? t) (+ t n)) |
|||
(else t)))) |
|||
(display (inverse 42 2017)) |
|||
(newline) |
|||
</syntaxhighlight> |
|||
{{out}} |
|||
<pre>$ ol modular-inverse-task-Owl.scm |
|||
1969</pre> |
|||
=={{header|PARI/GP}}== |
=={{header|PARI/GP}}== |
||
<lang |
<syntaxhighlight lang="parigp">Mod(1/42,2017)</syntaxhighlight> |
||
=={{header|Pascal}}== |
=={{header|Pascal}}== |
||
<syntaxhighlight lang="pascal"> |
|||
<lang Pascal> |
|||
// increments e step times until bal is greater than t |
// increments e step times until bal is greater than t |
||
// repeats until bal = 1 (mod = 1) and returns count |
// repeats until bal = 1 (mod = 1) and returns count |
||
Line 924: | Line 2,551: | ||
end; |
end; |
||
modInv := d; |
modInv := d; |
||
end;</ |
end;</syntaxhighlight> |
||
Testing: |
Testing: |
||
Line 935: | Line 2,562: | ||
=={{header|Perl}}== |
=={{header|Perl}}== |
||
Various CPAN modules can do this, such as: |
Various CPAN modules can do this, such as: |
||
< |
<syntaxhighlight lang="perl">use bigint; say 42->bmodinv(2017); |
||
# or |
# or |
||
use Math::ModInt qw/mod/; say mod(42, 2017)->inverse->residue; |
use Math::ModInt qw/mod/; say mod(42, 2017)->inverse->residue; |
||
Line 943: | Line 2,570: | ||
use Math::GMP qw/:constant/; say 42->bmodinv(2017); |
use Math::GMP qw/:constant/; say 42->bmodinv(2017); |
||
# or |
# or |
||
use ntheory qw/invmod/; say invmod(42, 2017);</ |
use ntheory qw/invmod/; say invmod(42, 2017);</syntaxhighlight> |
||
or we can write our own: |
or we can write our own: |
||
< |
<syntaxhighlight lang="perl">sub invmod { |
||
my($a,$n) = @_; |
my($a,$n) = @_; |
||
my($t,$nt,$r,$nr) = (0, 1, $n, $a % $n); |
my($t,$nt,$r,$nr) = (0, 1, $n, $a % $n); |
||
Line 959: | Line 2,586: | ||
} |
} |
||
say invmod(42,2017);</ |
say invmod(42,2017);</syntaxhighlight> |
||
'''Notes''': Special cases to watch out for include (1) where the inverse doesn't exist, such as invmod(14,28474), which should return undef or raise an exception, not return a wrong value. (2) the high bit of a or n is set, e.g. invmod(11,2**63), (3) negative first arguments, e.g. invmod(-11,23). |
'''Notes''': Special cases to watch out for include (1) where the inverse doesn't exist, such as invmod(14,28474), which should return undef or raise an exception, not return a wrong value. (2) the high bit of a or n is set, e.g. invmod(11,2**63), (3) negative first arguments, e.g. invmod(-11,23). |
||
The modules and code above handle these cases, |
The modules and code above handle these cases, |
||
but some other language implementations for this task do not. |
but some other language implementations for this task do not. |
||
=={{header| |
=={{header|Phix}}== |
||
{{trans|C}} |
|||
<lang perl6>sub inverse($n, :$modulo) { |
|||
<!--<syntaxhighlight lang="phix">(phixonline)--> |
|||
my ($c, $d, $uc, $vc, $ud, $vd) = ($n % $modulo, $modulo, 1, 0, 0, 1); |
|||
<span style="color: #008080;">function</span> <span style="color: #000000;">mul_inv</span><span style="color: #0000FF;">(</span><span style="color: #004080;">integer</span> <span style="color: #000000;">a</span><span style="color: #0000FF;">,</span> <span style="color: #000000;">n</span><span style="color: #0000FF;">)</span> |
|||
my $q; |
|||
<span style="color: #008080;">if</span> <span style="color: #000000;">n</span><span style="color: #0000FF;"><</span><span style="color: #000000;">0</span> <span style="color: #008080;">then</span> <span style="color: #000000;">n</span> <span style="color: #0000FF;">=</span> <span style="color: #0000FF;">-</span><span style="color: #000000;">n</span> <span style="color: #008080;">end</span> <span style="color: #008080;">if</span> |
|||
while $c != 0 { |
|||
<span style="color: #008080;">if</span> <span style="color: #000000;">a</span><span style="color: #0000FF;"><</span><span style="color: #000000;">0</span> <span style="color: #008080;">then</span> <span style="color: #000000;">a</span> <span style="color: #0000FF;">=</span> <span style="color: #000000;">n</span> <span style="color: #0000FF;">-</span> <span style="color: #7060A8;">mod</span><span style="color: #0000FF;">(-</span><span style="color: #000000;">a</span><span style="color: #0000FF;">,</span><span style="color: #000000;">n</span><span style="color: #0000FF;">)</span> <span style="color: #008080;">end</span> <span style="color: #008080;">if</span> |
|||
($q, $c, $d) = ($d div $c, $d % $c, $c); |
|||
<span style="color: #004080;">integer</span> <span style="color: #000000;">t</span> <span style="color: #0000FF;">=</span> <span style="color: #000000;">0</span><span style="color: #0000FF;">,</span> <span style="color: #000000;">nt</span> <span style="color: #0000FF;">=</span> <span style="color: #000000;">1</span><span style="color: #0000FF;">,</span> |
|||
($uc, $vc, $ud, $vd) = ($ud - $q*$uc, $vd - $q*$vc, $uc, $vc); |
|||
<span style="color: #000000;">r</span> <span style="color: #0000FF;">=</span> <span style="color: #000000;">n</span><span style="color: #0000FF;">,</span> <span style="color: #000000;">nr</span> <span style="color: #0000FF;">=</span> <span style="color: #000000;">a</span><span style="color: #0000FF;">;</span> |
|||
} |
|||
<span style="color: #008080;">while</span> <span style="color: #000000;">nr</span><span style="color: #0000FF;">!=</span><span style="color: #000000;">0</span> <span style="color: #008080;">do</span> |
|||
return $ud % $modulo; |
|||
<span style="color: #004080;">integer</span> <span style="color: #000000;">q</span> <span style="color: #0000FF;">=</span> <span style="color: #7060A8;">floor</span><span style="color: #0000FF;">(</span><span style="color: #000000;">r</span><span style="color: #0000FF;">/</span><span style="color: #000000;">nr</span><span style="color: #0000FF;">)</span> |
|||
} |
|||
<span style="color: #0000FF;">{</span><span style="color: #000000;">t</span><span style="color: #0000FF;">,</span> <span style="color: #000000;">nt</span><span style="color: #0000FF;">}</span> <span style="color: #0000FF;">=</span> <span style="color: #0000FF;">{</span><span style="color: #000000;">nt</span><span style="color: #0000FF;">,</span> <span style="color: #000000;">t</span><span style="color: #0000FF;">-</span><span style="color: #000000;">q</span><span style="color: #0000FF;">*</span><span style="color: #000000;">nt</span><span style="color: #0000FF;">}</span> |
|||
<span style="color: #0000FF;">{</span><span style="color: #000000;">r</span><span style="color: #0000FF;">,</span> <span style="color: #000000;">nr</span><span style="color: #0000FF;">}</span> <span style="color: #0000FF;">=</span> <span style="color: #0000FF;">{</span><span style="color: #000000;">nr</span><span style="color: #0000FF;">,</span> <span style="color: #000000;">r</span><span style="color: #0000FF;">-</span><span style="color: #000000;">q</span><span style="color: #0000FF;">*</span><span style="color: #000000;">nr</span><span style="color: #0000FF;">}</span> |
|||
say inverse 42, :modulo(2017)</lang> |
|||
<span style="color: #008080;">end</span> <span style="color: #008080;">while</span> |
|||
<span style="color: #008080;">if</span> <span style="color: #000000;">r</span><span style="color: #0000FF;">></span><span style="color: #000000;">1</span> <span style="color: #008080;">then</span> <span style="color: #008080;">return</span> <span style="color: #008000;">"a is not invertible"</span> <span style="color: #008080;">end</span> <span style="color: #008080;">if</span> |
|||
<span style="color: #008080;">if</span> <span style="color: #000000;">t</span><span style="color: #0000FF;"><</span><span style="color: #000000;">0</span> <span style="color: #008080;">then</span> <span style="color: #000000;">t</span> <span style="color: #0000FF;">+=</span> <span style="color: #000000;">n</span> <span style="color: #008080;">end</span> <span style="color: #008080;">if</span> |
|||
<span style="color: #008080;">return</span> <span style="color: #000000;">t</span> |
|||
<span style="color: #008080;">end</span> <span style="color: #008080;">function</span> |
|||
<span style="color: #0000FF;">?</span><span style="color: #000000;">mul_inv</span><span style="color: #0000FF;">(</span><span style="color: #000000;">42</span><span style="color: #0000FF;">,</span><span style="color: #000000;">2017</span><span style="color: #0000FF;">)</span> |
|||
<span style="color: #0000FF;">?</span><span style="color: #000000;">mul_inv</span><span style="color: #0000FF;">(</span><span style="color: #000000;">40</span><span style="color: #0000FF;">,</span> <span style="color: #000000;">1</span><span style="color: #0000FF;">)</span> |
|||
<span style="color: #0000FF;">?</span><span style="color: #000000;">mul_inv</span><span style="color: #0000FF;">(</span><span style="color: #000000;">52</span><span style="color: #0000FF;">,</span> <span style="color: #0000FF;">-</span><span style="color: #000000;">217</span><span style="color: #0000FF;">)</span> <span style="color: #000080;font-style:italic;">/* Pari semantics for negative modulus */</span> |
|||
<span style="color: #0000FF;">?</span><span style="color: #000000;">mul_inv</span><span style="color: #0000FF;">(-</span><span style="color: #000000;">486</span><span style="color: #0000FF;">,</span> <span style="color: #000000;">217</span><span style="color: #0000FF;">)</span> |
|||
<span style="color: #0000FF;">?</span><span style="color: #000000;">mul_inv</span><span style="color: #0000FF;">(</span><span style="color: #000000;">40</span><span style="color: #0000FF;">,</span> <span style="color: #000000;">2018</span><span style="color: #0000FF;">)</span> |
|||
<!--</syntaxhighlight>--> |
|||
{{out}} |
|||
<pre> |
|||
1969 |
|||
0 |
|||
96 |
|||
121 |
|||
"a is not invertible" |
|||
</pre> |
|||
=={{header|PHP}}== |
=={{header|PHP}}== |
||
'''Algorithm Implementation''' |
'''Algorithm Implementation''' |
||
< |
<syntaxhighlight lang="php"><?php |
||
function invmod($a,$n){ |
function invmod($a,$n){ |
||
if ($n < 0) $n = -$n; |
if ($n < 0) $n = -$n; |
||
Line 994: | Line 2,641: | ||
} |
} |
||
printf("%d\n", invmod(42, 2017)); |
printf("%d\n", invmod(42, 2017)); |
||
?></ |
?></syntaxhighlight> |
||
{{Out}} |
{{Out}} |
||
<pre>1969</pre> |
<pre>1969</pre> |
||
Line 1,000: | Line 2,647: | ||
=={{header|PicoLisp}}== |
=={{header|PicoLisp}}== |
||
{{trans|C}} |
{{trans|C}} |
||
< |
<syntaxhighlight lang="picolisp">(de modinv (A B) |
||
(let (B0 B X0 0 X1 1 Q 0 T1 0) |
(let (B0 B X0 0 X1 1 Q 0 T1 0) |
||
(while (< 1 A) |
(while (< 1 A) |
||
Line 1,016: | Line 2,663: | ||
(modinv 42 2017) ) |
(modinv 42 2017) ) |
||
(bye)</ |
(bye)</syntaxhighlight> |
||
=={{header|PL/I}}== |
=={{header|PL/I}}== |
||
{{trans|REXX}} |
{{trans|REXX}} |
||
< |
<syntaxhighlight lang="pli">*process source attributes xref or(!); |
||
/*-------------------------------------------------------------------- |
/*-------------------------------------------------------------------- |
||
* 13.07.2015 Walter Pachl |
* 13.07.2015 Walter Pachl |
||
Line 1,052: | Line 2,699: | ||
Return(d); |
Return(d); |
||
End; |
End; |
||
End;</ |
End;</syntaxhighlight> |
||
{{out}} |
{{out}} |
||
<pre>modular inverse of 42 by 2017 ---> 1969</pre> |
<pre>modular inverse of 42 by 2017 ---> 1969</pre> |
||
=={{header|PowerShell}}== |
=={{header|PowerShell}}== |
||
< |
<syntaxhighlight lang="powershell">function invmod($a,$n){ |
||
if ([int]$n -lt 0) {$n = -$n} |
if ([int]$n -lt 0) {$n = -$n} |
||
if ([int]$a -lt 0) {$a = $n - ((-$a) % $n)} |
if ([int]$a -lt 0) {$a = $n - ((-$a) % $n)} |
||
Line 1,079: | Line 2,726: | ||
} |
} |
||
invmod 42 2017</ |
invmod 42 2017</syntaxhighlight> |
||
{{Out}} |
{{Out}} |
||
<pre>PS> .\INVMOD.PS1 |
<pre>PS> .\INVMOD.PS1 |
||
1969 |
1969 |
||
PS> </pre> |
PS> </pre> |
||
=={{header|Prolog}}== |
|||
<syntaxhighlight lang="prolog"> |
|||
egcd(_, 0, 1, 0) :- !. |
|||
egcd(A, B, X, Y) :- |
|||
divmod(A, B, Q, R), |
|||
egcd(B, R, S, X), |
|||
Y is S - Q*X. |
|||
modinv(A, B, N) :- |
|||
egcd(A, B, X, Y), |
|||
A*X + B*Y =:= 1, |
|||
N is X mod B. |
|||
</syntaxhighlight> |
|||
{{Out}} |
|||
<pre> |
|||
?- modinv(42, 2017, N). |
|||
N = 1969. |
|||
?- modinv(42, 64, X). |
|||
false. |
|||
</pre> |
|||
=={{header|PureBasic}}== |
|||
Using brute force. |
|||
<syntaxhighlight lang="purebasic">EnableExplicit |
|||
Declare main() |
|||
Declare.i mi(a.i, b.i) |
|||
If OpenConsole("MODULAR-INVERSE") |
|||
main() : Input() : End |
|||
EndIf |
|||
Macro ModularInverse(a, b) |
|||
PrintN(~"\tMODULAR-INVERSE(" + RSet(Str(a),5) + "," + |
|||
RSet(Str(b),5)+") = " + |
|||
RSet(Str(mi(a, b)),5)) |
|||
EndMacro |
|||
Procedure main() |
|||
ModularInverse(42, 2017) ; = 1969 |
|||
ModularInverse(40, 1) ; = 0 |
|||
ModularInverse(52, -217) ; = 96 |
|||
ModularInverse(-486, 217) ; = 121 |
|||
ModularInverse(40, 2018) ; = -1 |
|||
EndProcedure |
|||
Procedure.i mi(a.i, b.i) |
|||
Define x.i = 1, |
|||
y.i = Int(Abs(b)), |
|||
r.i = 0 |
|||
If y = 1 : ProcedureReturn 0 : EndIf |
|||
While x < y |
|||
r = (a * x) % b |
|||
If r = 1 Or (y + r) = 1 |
|||
Break |
|||
EndIf |
|||
x + 1 |
|||
Wend |
|||
If x > y - 1 : x = -1 : EndIf |
|||
ProcedureReturn x |
|||
EndProcedure</syntaxhighlight> |
|||
{{out}} |
|||
<pre> MODULAR-INVERSE( 42, 2017) = 1969 |
|||
MODULAR-INVERSE( 40, 1) = 0 |
|||
MODULAR-INVERSE( 52, -217) = 96 |
|||
MODULAR-INVERSE( -486, 217) = 121 |
|||
MODULAR-INVERSE( 40, 2018) = -1</pre> |
|||
=={{header|Python}}== |
=={{header|Python}}== |
||
===Builtin function=== |
|||
Since python3.8, builtin function "pow" can be used directly to compute modular inverses by specifying an exponent of -1: |
|||
<syntaxhighlight lang="python">>>> pow(42, -1, 2017) |
|||
1969 |
|||
</syntaxhighlight> |
|||
===Iteration and error-handling=== |
|||
Implementation of this [http://en.wikipedia.org/wiki/Extended_Euclidean_algorithm#Iterative_method_2 pseudocode] with [https://en.wikibooks.org/wiki/Algorithm_Implementation/Mathematics/Extended_Euclidean_algorithm#Modular_inverse this]. |
Implementation of this [http://en.wikipedia.org/wiki/Extended_Euclidean_algorithm#Iterative_method_2 pseudocode] with [https://en.wikibooks.org/wiki/Algorithm_Implementation/Mathematics/Extended_Euclidean_algorithm#Modular_inverse this]. |
||
< |
<syntaxhighlight lang="python">>>> def extended_gcd(aa, bb): |
||
lastremainder, remainder = abs(aa), abs(bb) |
lastremainder, remainder = abs(aa), abs(bb) |
||
x, lastx, y, lasty = 0, 1, 1, 0 |
x, lastx, y, lasty = 0, 1, 1, 0 |
||
Line 1,104: | Line 2,827: | ||
>>> modinv(42, 2017) |
>>> modinv(42, 2017) |
||
1969 |
1969 |
||
>>> </ |
>>> </syntaxhighlight> |
||
===Recursion and an option type=== |
|||
Or, using functional composition as an alternative to iterative mutation, and wrapping the resulting value in an option type, to allow for the expression of computations which establish the '''absence''' of a modular inverse: |
|||
<syntaxhighlight lang="python">from functools import (reduce) |
|||
from itertools import (chain) |
|||
# modInv :: Int -> Int -> Maybe Int |
|||
def modInv(a): |
|||
return lambda m: ( |
|||
lambda ig=gcdExt(a)(m): ( |
|||
lambda i=ig[0]: ( |
|||
Just(i + m if 0 > i else i) if 1 == ig[2] else ( |
|||
Nothing() |
|||
) |
|||
) |
|||
)() |
|||
)() |
|||
# gcdExt :: Int -> Int -> (Int, Int, Int) |
|||
def gcdExt(x): |
|||
def go(a, b): |
|||
if 0 == b: |
|||
return (1, 0, a) |
|||
else: |
|||
(q, r) = divmod(a, b) |
|||
(s, t, g) = go(b, r) |
|||
return (t, s - q * t, g) |
|||
return lambda y: go(x, y) |
|||
# TEST --------------------------------------------------- |
|||
# Numbers between 2010 and 2015 which do yield modular inverses for 42: |
|||
# main :: IO () |
|||
def main(): |
|||
print ( |
|||
mapMaybe( |
|||
lambda y: bindMay(modInv(42)(y))( |
|||
lambda mInv: Just((y, mInv)) |
|||
) |
|||
)( |
|||
enumFromTo(2010)(2025) |
|||
) |
|||
) |
|||
# -> [(2011, 814), (2015, 48), (2017, 1969), (2021, 1203)] |
|||
# GENERIC ABSTRACTIONS ------------------------------------ |
|||
# enumFromTo :: Int -> Int -> [Int] |
|||
def enumFromTo(m): |
|||
return lambda n: list(range(m, 1 + n)) |
|||
# bindMay (>>=) :: Maybe a -> (a -> Maybe b) -> Maybe b |
|||
def bindMay(m): |
|||
return lambda mf: ( |
|||
m if m.get('Nothing') else mf(m.get('Just')) |
|||
) |
|||
# Just :: a -> Maybe a |
|||
def Just(x): |
|||
return {'type': 'Maybe', 'Nothing': False, 'Just': x} |
|||
# mapMaybe :: (a -> Maybe b) -> [a] -> [b] |
|||
def mapMaybe(mf): |
|||
return lambda xs: reduce( |
|||
lambda a, x: maybe(a)(lambda j: a + [j])(mf(x)), |
|||
xs, |
|||
[] |
|||
) |
|||
# maybe :: b -> (a -> b) -> Maybe a -> b |
|||
def maybe(v): |
|||
return lambda f: lambda m: v if m.get('Nothing') else ( |
|||
f(m.get('Just')) |
|||
) |
|||
# Nothing :: Maybe a |
|||
def Nothing(): |
|||
return {'type': 'Maybe', 'Nothing': True} |
|||
# MAIN --- |
|||
main()</syntaxhighlight> |
|||
{{Out}} |
|||
<pre> |
|||
[(2011, 814), (2015, 48), (2017, 1969), (2021, 1203)]</pre> |
|||
=={{header|Quackery}}== |
|||
{{trans|Forth}} |
|||
<syntaxhighlight lang="quackery"> [ dup 1 != if |
|||
[ tuck 1 0 |
|||
[ swap temp put |
|||
temp put |
|||
over 1 > while |
|||
tuck /mod swap |
|||
temp take tuck * |
|||
temp take swap - |
|||
again ] |
|||
2drop |
|||
temp release |
|||
temp take |
|||
dup 0 < if |
|||
[ over + ] ] |
|||
nip ] is modinv ( n n --> n ) |
|||
42 2017 modinv echo</syntaxhighlight> |
|||
{{out}} |
|||
<pre>1969</pre> |
|||
===Using Extended Euclidean Algorithm=== |
|||
Handles negative args. Returns -1 for non-coprime args. |
|||
<syntaxhighlight lang="quackery"> [ dup 0 = iff |
|||
[ 2drop 1 0 ] |
|||
done |
|||
dup unrot /mod |
|||
dip swap recurse |
|||
tuck 2swap * |
|||
dip swap - ] is egcd ( n n --> n n ) |
|||
[ dup 0 < if negate |
|||
over 0 < if |
|||
[ swap negate |
|||
over tuck mod |
|||
- swap ] |
|||
dup rot 2dup egcd |
|||
2swap 2over rot * |
|||
unrot * + 1 != iff |
|||
[ drop 2drop -1 ] |
|||
done |
|||
nip swap mod ] is modinv ( n n --> n ) |
|||
say " 42 2017 modinv --> " 42 2017 modinv echo cr ( --> 1969 ) |
|||
say " 40 1 modinv --> " 40 1 modinv echo cr ( --> 0 ) |
|||
say " 52 -217 modinv --> " 52 -217 modinv echo cr ( --> 96 ) |
|||
say "-486 217 modinv --> " -486 217 modinv echo cr ( --> 121 ) |
|||
say " 40 2018 modinv --> " 40 2018 modinv echo cr ( --> -1 )</syntaxhighlight> |
|||
{{out}} |
|||
<pre> 42 2017 modinv --> 1969 |
|||
40 1 modinv --> 0 |
|||
52 -217 modinv --> 96 |
|||
-486 217 modinv --> 121 |
|||
40 2018 modinv --> -1 |
|||
</pre> |
|||
=={{header|Racket}}== |
=={{header|Racket}}== |
||
< |
<syntaxhighlight lang="racket"> |
||
(require math) |
(require math) |
||
(modular-inverse 42 2017) |
(modular-inverse 42 2017) |
||
</syntaxhighlight> |
|||
</lang> |
|||
{{out}} |
{{out}} |
||
< |
<syntaxhighlight lang="racket"> |
||
1969 |
1969 |
||
</syntaxhighlight> |
|||
</lang> |
|||
=={{header|Raku}}== |
|||
(formerly Perl 6) |
|||
<syntaxhighlight lang="raku" line>sub inverse($n, :$modulo) { |
|||
my ($c, $d, $uc, $vc, $ud, $vd) = ($n % $modulo, $modulo, 1, 0, 0, 1); |
|||
my $q; |
|||
while $c != 0 { |
|||
($q, $c, $d) = ($d div $c, $d % $c, $c); |
|||
($uc, $vc, $ud, $vd) = ($ud - $q*$uc, $vd - $q*$vc, $uc, $vc); |
|||
} |
|||
return $ud % $modulo; |
|||
} |
|||
say inverse 42, :modulo(2017)</syntaxhighlight> |
|||
=={{header|REXX}}== |
=={{header|REXX}}== |
||
< |
<syntaxhighlight lang="rexx">/*REXX program calculates and displays the modular inverse of an integer X modulo Y.*/ |
||
parse arg x y . /*obtain two integers from the C.L. */ |
parse arg x y . /*obtain two integers from the C.L. */ |
||
if x=='' | x=="," then x= 42 /*Not specified? Then use the default.*/ |
if x=='' | x=="," then x= 42 /*Not specified? Then use the default.*/ |
||
Line 1,124: | Line 3,024: | ||
exit /*stick a fork in it, we're all done. */ |
exit /*stick a fork in it, we're all done. */ |
||
/*──────────────────────────────────────────────────────────────────────────────────────*/ |
/*──────────────────────────────────────────────────────────────────────────────────────*/ |
||
modInv: parse arg a,b 1 ob; z= |
modInv: parse arg a,b 1 ob; z= 0 /*B & OB are obtained from the 2nd arg.*/ |
||
$= 1 /*initialize modular inverse to unity. */ |
|||
$=1 |
|||
if b\=1 |
if b\=1 then do while a>1 |
||
parse value a/b a//b b z with q b a t |
|||
z= $ - q * z |
|||
$= trunc(t) |
|||
end /*while*/ |
|||
return $</lang> |
|||
if $<0 then $= $ + ob /*Negative? Then add the original B. */ |
|||
'''output''' when using the default inputs of: <tt> 42 2017 </tt> |
|||
return $</syntaxhighlight> |
|||
{{out|output|text= when using the default inputs of: <tt> 42 2017 </tt>}} |
|||
<pre> |
<pre> |
||
modular inverse of 42 by 2017 ───► 1969 |
modular inverse of 42 by 2017 ───► 1969 |
||
Line 1,138: | Line 3,040: | ||
=={{header|Ring}}== |
=={{header|Ring}}== |
||
< |
<syntaxhighlight lang="ring"> |
||
see "42 %! 2017 = " + multInv(42, 2017) + nl |
see "42 %! 2017 = " + multInv(42, 2017) + nl |
||
Line 1,157: | Line 3,059: | ||
if multInv < 0 multInv = multInv + b0 ok |
if multInv < 0 multInv = multInv + b0 ok |
||
return multInv |
return multInv |
||
</syntaxhighlight> |
|||
</lang> |
|||
Output: |
Output: |
||
<pre> |
<pre> |
||
42 %! 2017 = 1969 |
42 %! 2017 = 1969 |
||
</pre> |
|||
=={{header|RPL}}== |
|||
Using complex numbers allows to ‘parallelize’ calculations and keeps the stack depth low: never more than 4 levels despite the simultaneous use of 6 variables: r, r’, u, u’, q - and b for the final touch. |
|||
{{works with|Halcyon Calc|4.2.7}} |
|||
{| class="wikitable" |
|||
! RPL code |
|||
! Comment |
|||
|- |
|||
| |
|||
≪ |
|||
DUP ROT 1 R→C ROT 0 R→C |
|||
'''WHILE''' DUP RE '''REPEAT''' |
|||
OVER RE OVER RE / FLOOR |
|||
OVER * NEG ROT + |
|||
'''END''' |
|||
DROP C→R ROT MOD |
|||
SWAP 1 == SWAP 0 IFTE |
|||
≫ ‘'''MODINV'''’ STO |
|||
| |
|||
'''MODINV''' ''( a b -- x )'' with ax = 1 mod b |
|||
3: b 2: (r,u)←(a,1) 1:(r',u')←(b,0) |
|||
While r' ≠ 0 |
|||
q ← r // r' |
|||
(r - q*r', u - q*u') |
|||
Forget (r',u') and calculate u mod b |
|||
Test r and return zero if a and b are not co-prime |
|||
|} |
|||
{{in}} |
|||
<pre> |
|||
123 456 MODINV |
|||
42 2017 MODINV |
|||
</pre> |
|||
{{out}} |
|||
<pre> |
|||
2: 1969 |
|||
1: 0 |
|||
</pre> |
</pre> |
||
=={{header|Ruby}}== |
=={{header|Ruby}}== |
||
< |
<syntaxhighlight lang="ruby">#based on pseudo code from http://en.wikipedia.org/wiki/Extended_Euclidean_algorithm#Iterative_method_2 and from translating the python implementation. |
||
def extended_gcd(a, b) |
def extended_gcd(a, b) |
||
last_remainder, remainder = a.abs, b.abs |
last_remainder, remainder = a.abs, b.abs |
||
Line 1,183: | Line 3,124: | ||
end |
end |
||
x % et |
x % et |
||
end</ |
end</syntaxhighlight> |
||
<pre |
<pre> |
||
> invmod(42,2017) |
|||
=> 1969</pre> |
=> 1969</pre> |
||
Simplified equivalent implementation |
|||
<syntaxhighlight lang="ruby"> |
|||
def modinv(a, m) # compute a^-1 mod m if possible |
|||
raise "NO INVERSE - #{a} and #{m} not coprime" unless a.gcd(m) == 1 |
|||
return m if m == 1 |
|||
m0, inv, x0 = m, 1, 0 |
|||
while a > 1 |
|||
inv -= (a / m) * x0 |
|||
a, m = m, a % m |
|||
inv, x0 = x0, inv |
|||
end |
|||
inv += m0 if inv < 0 |
|||
inv |
|||
end |
|||
</syntaxhighlight> |
|||
<pre> |
|||
> modinv(42,2017) |
|||
=> 1969 |
|||
</pre> |
|||
The OpenSSL module has modular inverse built in: |
|||
<syntaxhighlight lang="ruby">require 'openssl' |
|||
p OpenSSL::BN.new(42).mod_inverse(2017).to_i</syntaxhighlight> |
|||
=={{header|Run BASIC}}== |
=={{header|Run BASIC}}== |
||
< |
<syntaxhighlight lang="runbasic">print multInv(42, 2017) |
||
end |
end |
||
Line 1,206: | Line 3,171: | ||
if multInv < 0 then multInv = multInv + b0 |
if multInv < 0 then multInv = multInv + b0 |
||
[endFun] |
[endFun] |
||
end function</ |
end function</syntaxhighlight> |
||
{{out}} |
{{out}} |
||
Line 1,214: | Line 3,179: | ||
=={{header|Rust}}== |
=={{header|Rust}}== |
||
< |
<syntaxhighlight lang="rust">fn mod_inv(a: isize, module: isize) -> isize { |
||
let mut mn = (module, a); |
let mut mn = (module, a); |
||
let mut xy = (0, 1); |
let mut xy = (0, 1); |
||
Line 1,231: | Line 3,196: | ||
fn main() { |
fn main() { |
||
println!("{}", mod_inv(42, 2017)) |
println!("{}", mod_inv(42, 2017)) |
||
}</ |
}</syntaxhighlight> |
||
{{out}} |
|||
<pre> |
|||
1969 |
|||
</pre> |
|||
Alternative implementation |
|||
<syntaxhighlight lang="rust"> |
|||
fn modinv(a0: isize, m0: isize) -> isize { |
|||
if m0 == 1 { return 1 } |
|||
let (mut a, mut m, mut x0, mut inv) = (a0, m0, 0, 1); |
|||
while a > 1 { |
|||
inv -= (a / m) * x0; |
|||
a = a % m; |
|||
std::mem::swap(&mut a, &mut m); |
|||
std::mem::swap(&mut x0, &mut inv); |
|||
} |
|||
if inv < 0 { inv += m0 } |
|||
inv |
|||
} |
|||
fn main() { |
|||
println!("{}", modinv(42, 2017)) |
|||
}</syntaxhighlight> |
|||
{{out}} |
{{out}} |
||
Line 1,240: | Line 3,232: | ||
=={{header|Scala}}== |
=={{header|Scala}}== |
||
Based on the ''Handbook of Applied Cryptography'', Chapter 2. See http://cacr.uwaterloo.ca/hac/ . |
Based on the ''Handbook of Applied Cryptography'', Chapter 2. See http://cacr.uwaterloo.ca/hac/ . |
||
< |
<syntaxhighlight lang="scala"> |
||
def gcdExt(u: Int, v: Int): (Int, Int, Int) = { |
def gcdExt(u: Int, v: Int): (Int, Int, Int) = { |
||
@tailrec |
@tailrec |
||
Line 1,255: | Line 3,247: | ||
val (i, j, g) = gcdExt(a, m) |
val (i, j, g) = gcdExt(a, m) |
||
if (g == 1) Option(if (i < 0) i + m else i) else Option.empty |
if (g == 1) Option(if (i < 0) i + m else i) else Option.empty |
||
}</ |
}</syntaxhighlight> |
||
Translated from C++ (on this page) |
Translated from C++ (on this page) |
||
< |
<syntaxhighlight lang="scala"> |
||
def modInv(a: Int, m: Int, x:Int = 1, y:Int = 0) : Int = if (m == 0) x else modInv(m, a%m, y, x - y*(a/m)) |
def modInv(a: Int, m: Int, x:Int = 1, y:Int = 0) : Int = if (m == 0) x else modInv(m, a%m, y, x - y*(a/m)) |
||
</syntaxhighlight> |
|||
</lang> |
|||
{{out}} |
{{out}} |
||
Line 1,277: | Line 3,269: | ||
If a and b are not coprime (gcd(a, b) <> 1) the exception RANGE_ERROR is raised. |
If a and b are not coprime (gcd(a, b) <> 1) the exception RANGE_ERROR is raised. |
||
< |
<syntaxhighlight lang="seed7">const func bigInteger: modInverse (in var bigInteger: a, |
||
in var bigInteger: b) is func |
in var bigInteger: b) is func |
||
result |
result |
||
Line 1,319: | Line 3,311: | ||
raise RANGE_ERROR; |
raise RANGE_ERROR; |
||
end if; |
end if; |
||
end func;</ |
end func;</syntaxhighlight> |
||
Original source: [http://seed7.sourceforge.net/algorith/math.htm#modInverse] |
Original source: [http://seed7.sourceforge.net/algorith/math.htm#modInverse] |
||
Line 1,325: | Line 3,317: | ||
=={{header|Sidef}}== |
=={{header|Sidef}}== |
||
Built-in: |
Built-in: |
||
<lang |
<syntaxhighlight lang="ruby">say 42.modinv(2017)</syntaxhighlight> |
||
Algorithm implementation: |
Algorithm implementation: |
||
< |
<syntaxhighlight lang="ruby">func invmod(a, n) { |
||
var (t, nt, r, nr) = (0, 1, n, a % n) |
var (t, nt, r, nr) = (0, 1, n, a % n) |
||
while (nr != 0) { |
while (nr != 0) { |
||
Line 1,340: | Line 3,332: | ||
} |
} |
||
say invmod(42, 2017)</ |
say invmod(42, 2017)</syntaxhighlight> |
||
{{out}} |
{{out}} |
||
<pre> |
<pre> |
||
1969 |
1969 |
||
</pre> |
</pre> |
||
=={{header|Swift}}== |
|||
<syntaxhighlight lang="swift">extension BinaryInteger { |
|||
@inlinable |
|||
public func modInv(_ mod: Self) -> Self { |
|||
var (m, n) = (mod, self) |
|||
var (x, y) = (Self(0), Self(1)) |
|||
while n != 0 { |
|||
(x, y) = (y, x - (m / n) * y) |
|||
(m, n) = (n, m % n) |
|||
} |
|||
while x < 0 { |
|||
x += mod |
|||
} |
|||
return x |
|||
} |
|||
} |
|||
print(42.modInv(2017))</syntaxhighlight> |
|||
{{out}} |
|||
<pre>1969</pre> |
|||
=={{header|Tcl}}== |
=={{header|Tcl}}== |
||
{{trans|Haskell}} |
{{trans|Haskell}} |
||
< |
<syntaxhighlight lang="tcl">proc gcdExt {a b} { |
||
if {$b == 0} { |
if {$b == 0} { |
||
return [list 1 0 $a] |
return [list 1 0 $a] |
||
Line 1,364: | Line 3,383: | ||
while {$i < 0} {incr i $m} |
while {$i < 0} {incr i $m} |
||
return $i |
return $i |
||
}</ |
}</syntaxhighlight> |
||
Demonstrating |
Demonstrating |
||
< |
<syntaxhighlight lang="tcl">puts "42 %! 2017 = [modInv 42 2017]" |
||
catch { |
catch { |
||
puts "2 %! 4 = [modInv 2 4]" |
puts "2 %! 4 = [modInv 2 4]" |
||
} msg; puts $msg</ |
} msg; puts $msg</syntaxhighlight> |
||
{{out}} |
{{out}} |
||
<pre> |
<pre> |
||
42 %! 2017 = 1969 |
42 %! 2017 = 1969 |
||
no inverse exists of 2 %! 4 |
no inverse exists of 2 %! 4 |
||
</pre> |
|||
=={{header|Tiny BASIC}}== |
|||
2017 causes integer overflow, so I'll do the inverse of 42 modulo 331 instead. |
|||
<syntaxhighlight lang="tinybasic"> |
|||
PRINT "Modular inverse." |
|||
PRINT "What is the modulus?" |
|||
INPUT M |
|||
PRINT "What number is to be inverted?" |
|||
INPUT X |
|||
PRINT "Solution is:" |
|||
10 LET A = A + 1 |
|||
GOTO 20 |
|||
15 IF B = 1 THEN PRINT A |
|||
IF B = 1 THEN END |
|||
IF A = M-1 THEN PRINT "nonexistent" |
|||
IF A = M-1 THEN END |
|||
GOTO 10 |
|||
20 LET B = A*X |
|||
30 IF B < M THEN GOTO 15 |
|||
LET B = B - M |
|||
GOTO 30 |
|||
</syntaxhighlight> |
|||
{{out}} |
|||
<pre>Modular inverse. |
|||
What is the modulus? |
|||
331 |
|||
What number is to be inverted? |
|||
42 |
|||
Solution is: |
|||
134</pre> |
|||
Another version: |
|||
{{trans|GW-BASIC}} |
|||
<syntaxhighlight lang="tinybasic"> |
|||
REM Modular inverse |
|||
LET E=42 |
|||
LET T=2017 |
|||
GOSUB 100 |
|||
PRINT M |
|||
END |
|||
REM Increments E S (step) times until B is greater than T |
|||
REM Repeats until B = 1 (MOD = 1) and C (count) |
|||
REM B will not be greater than T + E |
|||
100 LET D=0 |
|||
IF E>=T THEN GOTO 130 |
|||
LET B=E |
|||
REM At least one iteration is necessary |
|||
LET S=((T-B)/E)+1 |
|||
LET B=B+S*E |
|||
LET C=1+S |
|||
LET B=B-T |
|||
110 IF B=1 THEN GOTO 120 |
|||
LET S=((T-B)/E)+1 |
|||
LET B=B+S*E |
|||
LET C=C+S |
|||
LET B=B-T |
|||
GOTO 110 |
|||
120 LET D=C |
|||
130 LET M=D |
|||
RETURN |
|||
</syntaxhighlight> |
|||
{{out}} |
|||
<pre> |
|||
1969 |
|||
</pre> |
</pre> |
||
=={{header|tsql}}== |
=={{header|tsql}}== |
||
< |
<syntaxhighlight lang="tsql">;WITH Iterate(N,A,B,X0,X1) |
||
AS ( |
AS ( |
||
SELECT |
SELECT |
||
Line 1,410: | Line 3,495: | ||
) |
) |
||
SELECT * |
SELECT * |
||
FROM ModularInverse</ |
FROM ModularInverse</syntaxhighlight> |
||
== {{header|TypeScript}} == |
|||
{{trans|Pascal}} |
|||
<syntaxhighlight lang="javascript"> |
|||
// Modular inverse |
|||
function modInv(e: number, t: number): number { |
|||
var d = 0; |
|||
if (e < t) { |
|||
var count = 1; |
|||
var bal = e; |
|||
do { |
|||
var step = Math.floor((t - bal) / e) + 1; |
|||
bal += step * e; |
|||
count += step; |
|||
bal -= t; |
|||
} while (bal != 1); |
|||
d = count; |
|||
} |
|||
return d; |
|||
} |
|||
console.log(`${modInv(42, 2017)}`); // 1969 |
|||
</syntaxhighlight> |
|||
{{out}} |
|||
<pre> |
|||
1969 |
|||
</pre> |
|||
=={{header|uBasic/4tH}}== |
=={{header|uBasic/4tH}}== |
||
{{trans|C}} |
{{trans|C}} |
||
<lang>Print FUNC(_MulInv(42, 2017)) |
<syntaxhighlight lang="text">Print FUNC(_MulInv(42, 2017)) |
||
End |
End |
||
Line 1,438: | Line 3,551: | ||
If g@ < 0 Then g@ = g@ + c@ |
If g@ < 0 Then g@ = g@ + c@ |
||
Return (g@)</ |
Return (g@)</syntaxhighlight> |
||
{{trans|Perl}} |
{{trans|Perl}} |
||
<lang>Print FUNC(_mul_inv(42, 2017)) |
<syntaxhighlight lang="text">Print FUNC(_mul_inv(42, 2017)) |
||
Print FUNC(_mul_inv(40, 1)) |
Print FUNC(_mul_inv(40, 1)) |
||
Print FUNC(_mul_inv(52, -217)) |
Print FUNC(_mul_inv(52, -217)) |
||
Line 1,463: | Line 3,576: | ||
If (e@ > 1) Return (-1) ' No inverse' |
If (e@ > 1) Return (-1) ' No inverse' |
||
If (c@ < 0) c@ = c@ + b@ |
If (c@ < 0) c@ = c@ + b@ |
||
Return (c@)</ |
Return (c@)</syntaxhighlight> |
||
{{out}} |
{{out}} |
||
<pre>1969 |
<pre>1969 |
||
Line 1,472: | Line 3,585: | ||
0 OK, 0:156</pre> |
0 OK, 0:156</pre> |
||
=={{header|UNIX Shell}}== |
|||
{{works with|Bourne Again SHell}} |
|||
{{works with|Korn Shell}} |
|||
{{works with|Zsh}} |
|||
{{trans|PowerShell}}<syntaxhighlight lang="sh">function invmod { |
|||
typeset -i a=$1 n=$2 |
|||
if (( n < 0 )); then (( n = -n )); fi |
|||
if (( a < 0 )); then (( a = n - (-a) % n )); fi |
|||
typeset -i t=0 nt=1 r=n nr q tmp |
|||
(( nr = a % n )) |
|||
while (( nr )); do |
|||
(( q = r/nr )) |
|||
(( tmp = nt )) |
|||
(( nt = t - q*nt )) |
|||
(( t = tmp )) |
|||
(( tmp = nr )) |
|||
(( nr = r - q*nr )) |
|||
(( r = tmp )) |
|||
done |
|||
if (( r > 1 )); then |
|||
return 1 |
|||
fi |
|||
while (( t < 0 )); do (( t += n )); done |
|||
printf '%s\n' "$t" |
|||
} |
|||
invmod 42 2017</syntaxhighlight> |
|||
{{Out}} |
|||
<pre>1969</pre> |
|||
=={{header|VBA}}== |
|||
{{trans|Phix}}<syntaxhighlight lang="vb"> |
|||
Private Function mul_inv(a As Long, n As Long) As Variant |
|||
If n < 0 Then n = -n |
|||
If a < 0 Then a = n - ((-a) Mod n) |
|||
Dim t As Long: t = 0 |
|||
Dim nt As Long: nt = 1 |
|||
Dim r As Long: r = n |
|||
Dim nr As Long: nr = a |
|||
Dim q As Long |
|||
Do While nr <> 0 |
|||
q = r \ nr |
|||
tmp = t |
|||
t = nt |
|||
nt = tmp - q * nt |
|||
tmp = r |
|||
r = nr |
|||
nr = tmp - q * nr |
|||
Loop |
|||
If r > 1 Then |
|||
mul_inv = "a is not invertible" |
|||
Else |
|||
If t < 0 Then t = t + n |
|||
mul_inv = t |
|||
End If |
|||
End Function |
|||
Public Sub mi() |
|||
Debug.Print mul_inv(42, 2017) |
|||
Debug.Print mul_inv(40, 1) |
|||
Debug.Print mul_inv(52, -217) '/* Pari semantics for negative modulus */ |
|||
Debug.Print mul_inv(-486, 217) |
|||
Debug.Print mul_inv(40, 2018) |
|||
End Sub</syntaxhighlight>{{out}} |
|||
<pre> 1969 |
|||
0 |
|||
96 |
|||
121 |
|||
a is not invertible</pre> |
|||
=={{header|V (Vlang)}}== |
|||
<syntaxhighlight lang="Zig"> |
|||
fn main() { |
|||
println("42 %! 2017 = ${mult_inv(42, 2017)}") |
|||
} |
|||
fn mult_inv(aa int, bb int) int { |
|||
mut a, mut b := aa, bb |
|||
mut x0, mut t := 0, 0 |
|||
mut b0 := b |
|||
mut x1 := 1 |
|||
if b == 1 {return 1} |
|||
for a > 1 { |
|||
q := a / b |
|||
t = b |
|||
b = a % b |
|||
a = t |
|||
t = x0 |
|||
x0 = x1 - q * x0 |
|||
x1 = t |
|||
} |
|||
if x1 < 0 {x1 += b0} |
|||
return x1 |
|||
} |
|||
</syntaxhighlight> |
|||
{{out}} |
|||
<pre> |
|||
42 %! 2017 = 1969 |
|||
</pre> |
|||
=={{header|Wren}}== |
|||
{{libheader|Wren-big}} |
|||
<syntaxhighlight lang="wren">import "./big" for BigInt |
|||
var a = BigInt.new(42) |
|||
var b = BigInt.new(2017) |
|||
System.print(a.modInv(b))</syntaxhighlight> |
|||
{{out}} |
|||
<pre> |
|||
1969 |
|||
</pre> |
|||
=={{header|XPL0}}== |
=={{header|XPL0}}== |
||
< |
<syntaxhighlight lang="xpl0">code IntOut=11, Text=12; |
||
int X; |
int X; |
||
def A=42, M=2017; |
def A=42, M=2017; |
||
Line 1,480: | Line 3,708: | ||
if rem(A*X/M) = 1 then [IntOut(0, X); exit]; |
if rem(A*X/M) = 1 then [IntOut(0, X); exit]; |
||
Text(0, "Does not exist"); |
Text(0, "Does not exist"); |
||
]</ |
]</syntaxhighlight> |
||
{{out}} |
{{out}} |
||
<pre> |
<pre> |
||
Line 1,488: | Line 3,716: | ||
=={{header|zkl}}== |
=={{header|zkl}}== |
||
{{trans|Haskell}} |
{{trans|Haskell}} |
||
< |
<syntaxhighlight lang="zkl">fcn gcdExt(a,b){ |
||
if(b==0) return(1,0,a); |
if(b==0) return(1,0,a); |
||
q,r:=a.divr(b); s,t,g:=gcdExt(b,r); return(t,s-q*t,g); |
q,r:=a.divr(b); s,t,g:=gcdExt(b,r); return(t,s-q*t,g); |
||
} |
} |
||
fcn modInv(a,m){i,_,g:=gcdExt(a,m); if(g==1) {if(i<0)i+m} else Void}</ |
fcn modInv(a,m){i,_,g:=gcdExt(a,m); if(g==1) {if(i<0)i+m} else Void}</syntaxhighlight> |
||
divr(a,b) is [integer] (a/b,remainder) |
divr(a,b) is [integer] (a/b,remainder) |
||
{{out}} |
{{out}} |
Latest revision as of 10:25, 4 January 2024
![Task](http://static.miraheze.org/rosettacodewiki/thumb/b/ba/Rcode-button-task-crushed.png/64px-Rcode-button-task-crushed.png)
You are encouraged to solve this task according to the task description, using any language you may know.
From Wikipedia:
In modular arithmetic, the modular multiplicative inverse of an integer a modulo m is an integer x such that
Or in other words, such that:
It can be shown that such an inverse exists if and only if a and m are coprime, but we will ignore this for this task.
- Task
Either by implementing the algorithm, by using a dedicated library or by using a built-in function in
your language, compute the modular inverse of 42 modulo 2017.
11l
F mul_inv(=a, =b)
V b0 = b
V x0 = 0
V x1 = 1
I b == 1 {R 1}
L a > 1
V q = a I/ b
(a, b) = (b, a % b)
(x0, x1) = (x1 - q * x0, x0)
I x1 < 0 {x1 += b0}
R x1
print(mul_inv(42, 2017))
- Output:
1969
8th
\ return "extended gcd" of a and b; The result satisfies the equation:
\ a*x + b*y = gcd(a,b)
: n:xgcd \ a b -- gcd x y
dup 0 n:= if
1 swap \ -- a 1 0
else
tuck n:/mod
-rot recurse
tuck 4 roll
n:* n:neg n:+
then ;
\ Return modular inverse of n modulo mod, or null if it doesn't exist (n and mod
\ not coprime):
: n:invmod \ n mod -- invmod
dup >r
n:xgcd rot 1 n:= not if
2drop null
else
drop dup 0 n:< if r@ n:+ then
then
rdrop ;
42 2017 n:invmod . cr bye
- Output:
1969
Action!
INT FUNC ModInverse(INT a,b)
INT t,nt,r,nr,q,tmp
IF b<0 THEN b=-b FI
IF a<0 THEN a=b-(-a MOD b) FI
t=0 nt=1
r=b nr=a MOD b
WHILE nr#0
DO
q=r/nr
tmp=nt nt=t-q*nt t=tmp
tmp=nr nr=r-q*nr r=tmp
OD
IF r>1 THEN
RETURN (-1)
FI
IF t<0 THEN
t==+b
FI
RETURN (t)
PROC Test(INT a,b)
INT res
res=ModInverse(a,b)
IF res>=0 THEN
PrintF("%I MODINV %I=%I%E",a,b,res)
ELSE
PrintF("%I MODINV %I has no result%E",a,b)
FI
RETURN
PROC Main()
Test(42,2017)
Test(40,1)
Test(52,-217)
Test(-486,217)
Test(40,2018)
RETURN
- Output:
Screenshot from Atari 8-bit computer
42 MODINV 2017=1969 40 MODINV 1=0 52 MODINV -217=96 -486 MODINV 217=121 40 MODINV 2018 has no result
Ada
with Ada.Text_IO;use Ada.Text_IO;
procedure modular_inverse is
-- inv_mod calculates the inverse of a mod n. We should have n>0 and, at the end, the contract is a*Result=1 mod n
-- If this is false then we raise an exception (don't forget the -gnata option when you compile
function inv_mod (a : Integer; n : Positive) return Integer with post=> (a * inv_mod'Result) mod n = 1 is
-- To calculate the inverse we do as if we would calculate the GCD with the Euclid extended algorithm
-- (but we just keep the coefficient on a)
function inverse (a, b, u, v : Integer) return Integer is
(if b=0 then u else inverse (b, a mod b, v, u-(v*a)/b));
begin
return inverse (a, n, 1, 0);
end inv_mod;
begin
-- This will output -48 (which is correct)
Put_Line (inv_mod (42,2017)'img);
-- The further line will raise an exception since the GCD will not be 1
Put_Line (inv_mod (42,77)'img);
exception when others => Put_Line ("The inverse doesn't exist.");
end modular_inverse;
ALGOL 68
BEGIN
PROC modular inverse = (INT a, m) INT :
BEGIN
PROC extended gcd = (INT x, y) []INT :
CO
Algol 68 allows us to return three INTs in several ways. A [3]INT
is used here but it could just as well be a STRUCT.
CO
BEGIN
INT v := 1, a := 1, u := 0, b := 0, g := x, w := y;
WHILE w>0
DO
INT q := g % w, t := a - q * u;
a := u; u := t;
t := b - q * v;
b := v; v := t;
t := g - q * w;
g := w; w := t
OD;
a PLUSAB (a < 0 | u | 0);
(a, b, g)
END;
[] INT egcd = extended gcd (a, m);
(egcd[3] > 1 | 0 | egcd[1] MOD m)
END;
printf (($"42 ^ -1 (mod 2017) = ", g(0)$, modular inverse (42, 2017)))
CO
Note that if ϕ(m) is known, then a^-1 = a^(ϕ(m)-1) mod m which
allows an alternative implementation in terms of modular
exponentiation but, in general, this requires the factorization of
m. If m is prime the factorization is trivial and ϕ(m) = m-1.
2017 is prime which may, or may not, be ironic within the context
of the Rosetta Code conditions.
CO
END
- Output:
42 ^ -1 (mod 2017) = 1969
Arturo
modInverse: function [a,b][
if b = 1 -> return 1
b0: b x0: 0 x1: 1
z: a
while [z > 1][
q: z / b t: b
b: z % b z: t
t: x0 x0: x1 - q * x0
x1: t
]
(x1 < 0) ? -> x1 + b0
-> x1
]
print modInverse 42 2017
- Output:
1969
ATS
Using allocated memory
In addition to solving the task, I demonstrate some aspects of call-by-reference in ATS. In particular, ATS can distinguish at compile time between uninitialized and initialized variables.
The code is written as templates that will expand to code for any of the (non-dependent) signed integer types. In the main program, I use llint
(typically exactly the same as long long int
in C).
The return value of inverse
is a linear optional value. It is allocated in the heap, used once, and freed.
(*
Using the algorithm described at
https://en.wikipedia.org/w/index.php?title=Extended_Euclidean_algorithm&oldid=1135569411#Modular_integers
*)
#include "share/atspre_staload.hats"
fn {tk : tkind}
division_with_nonnegative_remainder
(n : g0int tk, d : g0int tk,
(* q and r are called by reference, and start out
uninitialized. *)
q : &g0int tk? >> g0int tk,
r : &g0int tk? >> g0int tk)
: void =
let
(* The C optimizer most likely will reduce these these two
divisions to just one. They are simply synonyms for C '/' and
'%', and perform division that rounds the quotient towards
zero. *)
val q0 = g0int_div (n, d)
val r0 = g0int_mod (n, d)
in
(* The following calculation results in 'floor division', if the
divisor is positive, or 'ceiling division', if the divisor is
negative. This choice of method results in the remainder never
being negative. *)
if isgtez n || iseqz r0 then
(q := q0; r := r0)
else if isltz d then
(q := succ q0; r := r0 - d)
else
(q := pred q0; r := r0 + d)
end
fn {tk : tkind}
inverse (a : g0int tk, n : g0int tk) : Option_vt (g0int tk) =
let
typedef integer = g0int tk
fun
loop (t : integer, newt : integer,
r : integer, newr : integer) : Option_vt integer =
if iseqz newr then
begin
if r > g0i2i 1 then
None_vt ()
else if t < g0i2i 0 then
Some_vt (t + n)
else
Some_vt t
end
else
let
(* These become C variables. *)
var quotient : g0int tk?
var remainder : g0int tk?
(* Show the type AT COMPILE TIME. *)
prval _ = $showtype quotient
prval _ = $showtype remainder
val () =
division_with_nonnegative_remainder
(r, newr, quotient, remainder)
(* THE TYPES WILL HAVE CHANGED, because the storage is
initialized by the call to
division_with_nonnegative_remainder. *)
prval _ = $showtype quotient
prval _ = $showtype remainder
val t = newt
and newt = t - (quotient * newt)
and r = newr
and newr = remainder
in
loop (t, newt, r, newr)
end
in
loop (g0i2i 0, g0i2i 1, n, a)
end
implement
main0 () =
case+ inverse (42LL, 2017LL) of
| ~ None_vt () => println! "There is no inverse."
| ~ Some_vt value => println! value
- Output:
First compile the program:
$ patscc -DATS_MEMALLOC_LIBC -g -O2 modular-inverse.dats **SHOWTYPE[UP]**(/home/trashman/src/chemoelectric/rosettacode-contributions/modular-inverse.dats: 1786(line=62, offs=31) -- 1794(line=62, offs=39)): S2Etop(knd=0; S2Eapp(S2Ecst(g0int_t0ype); S2Evar(tk(8481)))): S2RTbas(S2RTBASimp(1; t@ype)) **SHOWTYPE[UP]**(/home/trashman/src/chemoelectric/rosettacode-contributions/modular-inverse.dats: 1825(line=63, offs=31) -- 1834(line=63, offs=40)): S2Etop(knd=0; S2Eapp(S2Ecst(g0int_t0ype); S2Evar(tk(8481)))): S2RTbas(S2RTBASimp(1; t@ype)) **SHOWTYPE[UP]**(/home/trashman/src/chemoelectric/rosettacode-contributions/modular-inverse.dats: 2137(line=72, offs=31) -- 2145(line=72, offs=39)): S2Eapp(S2Ecst(g0int_t0ype); S2Evar(tk(8481))): S2RTbas(S2RTBASimp(1; t@ype)) **SHOWTYPE[UP]**(/home/trashman/src/chemoelectric/rosettacode-contributions/modular-inverse.dats: 2176(line=73, offs=31) -- 2185(line=73, offs=40)): S2Eapp(S2Ecst(g0int_t0ype); S2Evar(tk(8481))): S2RTbas(S2RTBASimp(1; t@ype))
You may notice there is a subtle change in the type of quotient
and remainder
, once they have been initialized. ATS is making it safe not to initialize variables (with phony values) when you declare them.
Now run the program:
$ ./a.out 1969
Safely avoiding the need for an allocator
Here I demonstrate an optional value that requires no runtime overhead at all, but which is safe. If there is no inverse, then the compiler knows that inverse_value
is still uninitialized, and will not let you use its value. (Try it and see.)
(*
Using the algorithm described at
https://en.wikipedia.org/w/index.php?title=Extended_Euclidean_algorithm&oldid=1135569411#Modular_integers
*)
#include "share/atspre_staload.hats"
fn {tk : tkind}
division_with_nonnegative_remainder
(n : g0int tk, d : g0int tk,
(* q and r are called by reference, and start out
uninitialized. *)
q : &g0int tk? >> g0int tk,
r : &g0int tk? >> g0int tk)
: void =
let
(* The C optimizer most likely will reduce these these two
divisions to just one. They are simply synonyms for C '/' and
'%', and perform division that rounds the quotient towards
zero. *)
val q0 = g0int_div (n, d)
val r0 = g0int_mod (n, d)
in
(* The following calculation results in 'floor division', if the
divisor is positive, or 'ceiling division', if the divisor is
negative. This choice of method results in the remainder never
being negative. *)
if isgtez n || iseqz r0 then
(q := q0; r := r0)
else if isltz d then
(q := succ q0; r := r0 - d)
else
(q := pred q0; r := r0 + d)
end
fn {tk : tkind}
inverse (a : g0int tk, n : g0int tk,
inverse_exists : &bool? >> bool exists,
inverse_value : &g0int tk? >> opt (g0int tk, exists))
: #[exists: bool] void =
let
typedef integer = g0int tk
fun
loop (t : integer, newt : integer,
r : integer, newr : integer,
inverse_exists : &bool? >> bool exists,
inverse_value : &g0int tk? >> opt (g0int tk, exists))
: #[exists: bool] void =
if iseqz newr then
begin
if r > g0i2i 1 then
let
val () = inverse_exists := false
prval () = opt_none inverse_value
in
end
else if t < g0i2i 0 then
let
val () = inverse_exists := true
val () = inverse_value := t + n
prval () = opt_some inverse_value
in
end
else
let
val () = inverse_exists := true
val () = inverse_value := t
prval () = opt_some inverse_value
in
end
end
else
let
(* These become C variables. *)
var quotient : g0int tk?
var remainder : g0int tk?
val () =
division_with_nonnegative_remainder
(r, newr, quotient, remainder)
val t = newt
and newt = t - (quotient * newt)
and r = newr
and newr = remainder
in
loop (t, newt, r, newr, inverse_exists, inverse_value)
end
in
loop (g0i2i 0, g0i2i 1, n, a, inverse_exists, inverse_value)
end
implement
main0 () =
let
var inverse_exists : bool?
var inverse_value : llint?
in
inverse (42LL, 2017LL, inverse_exists, inverse_value);
if inverse_exists then
let
prval () = opt_unsome inverse_value
in
println! inverse_value
end
else
let
prval () = opt_unnone inverse_value
in
println! "There is no inverse."
end
end
- Output:
There is no need to tell patscc what allocator to use, because none is used.
$ patscc -g -O2 modular-inverse-noheap.dats && ./a.out 1969
AutoHotkey
Translation of C.
MsgBox, % ModInv(42, 2017)
ModInv(a, b) {
if (b = 1)
return 1
b0 := b, x0 := 0, x1 :=1
while (a > 1) {
q := a // b
, t := b
, b := Mod(a, b)
, a := t
, t := x0
, x0 := x1 - q * x0
, x1 := t
}
if (x1 < 0)
x1 += b0
return x1
}
- Output:
1969
AWK
# syntax: GAWK -f MODULAR_INVERSE.AWK
# converted from C
BEGIN {
printf("%s\n",mod_inv(42,2017))
exit(0)
}
function mod_inv(a,b, b0,t,q,x0,x1) {
b0 = b
x0 = 0
x1 = 1
if (b == 1) {
return(1)
}
while (a > 1) {
q = int(a / b)
t = b
b = int(a % b)
a = t
t = x0
x0 = x1 - q * x0
x1 = t
}
if (x1 < 0) {
x1 += b0
}
return(x1)
}
- Output:
1969
BASIC
ASIC
REM Modular inverse
E = 42
T = 2017
GOSUB CalcModInv:
PRINT ModInv
END
CalcModInv:
REM Increments E Step times until Bal is greater than T
REM Repeats until Bal = 1 (MOD = 1) and returns Count
REM Bal will not be greater than T + E
D = 0
IF E < T THEN
Bal = E
Count = 1
Loop:
Step = T - Bal
Step = Step / E
Step = Step + 1
REM So ... Step = (T - Bal) / E + 1
StepTimesE = Step * E
Bal = Bal + StepTimesE
Count = Count + Step
Bal = Bal - T
IF Bal <> 1 THEN Loop:
D = Count
ENDIF
ModInv = D
RETURN
- Output:
1969
BASIC256
print multInv(42, 2017)
end
function multInv(a,b)
x0 = 0
b0 = b
multInv = 1
if b = 1 then return
while a > 1
q = a / b
t = b
b = a mod b
a = t
t = x0
x0 = multInv - q * x0
multInv = int(t)
end while
if multInv < 0 then return multInv + b0
end function
- Output:
1969
Chipmunk Basic
10 CLS
20 CALL modularinverse(42, 2017)
30 CALL modularinverse(40, 1)
40 END
50 SUB modularinverse(e,t)
60 d = 0
70 IF e < t THEN
80 b = e
90 c = 1
100 WHILE b > 1
110 s = INT(((t-b)/e)+1)
120 b = b+s*e
130 c = c+s
140 b = b-t
150 WEND
160 d = c
170 ENDIF
180 m = d
190 PRINT m
200 END SUB
Minimal BASIC
10 REM Modular inverse
20 LET E = 42
30 LET T = 2017
40 GOSUB 500
50 PRINT M
60 END
490 REM Calculate modular inverse
500 LET D = 0
510 IF E >= T THEN 600
520 LET B = E
530 LET C = 1
540 LET S1 = INT((T-B)/E)+1
550 LET B = B+S1*E
560 LET C = C+S1
570 LET B = B-T
580 IF B <> 1 THEN 540
590 LET D = C
600 LET M = D
610 RETURN
QBasic
The Chipmunk Basic solution works without any changes.
True BASIC
SUB modularinverse(e,t)
LET d = 0
IF e < t then
LET b = e
LET c = 1
DO WHILE b > 1
LET s = int(((t-b)/e)+1)
LET b = b+s*e
LET c = c+s
LET b = b-t
LOOP
LET d = c
END IF
LET m = d
PRINT m
END SUB
CALL modularinverse(42,2017)
CALL modularinverse(40,1)
END
Batch File
Based from C's second implementation
@echo off
setlocal enabledelayedexpansion
%== Calls the "function" ==%
call :ModInv 42 2017 result
echo !result!
call :ModInv 40 1 result
echo !result!
call :ModInv 52 -217 result
echo !result!
call :ModInv -486 217 result
echo !result!
call :ModInv 40 2018 result
echo !result!
pause>nul
exit /b 0
%== The "function" ==%
:ModInv
set a=%1
set b=%2
if !b! lss 0 (set /a b=-b)
if !a! lss 0 (set /a a=b - ^(-a %% b^))
set t=0&set nt=1&set r=!b!&set /a nr=a%%b
:while_loop
if !nr! neq 0 (
set /a q=r/nr
set /a tmp=nt
set /a nt=t - ^(q*nt^)
set /a t=tmp
set /a tmp=nr
set /a nr=r - ^(q*nr^)
set /a r=tmp
goto while_loop
)
if !r! gtr 1 (set %3=-1&goto :EOF)
if !t! lss 0 set /a t+=b
set %3=!t!
goto :EOF
- Output:
1969 0 96 121 -1
BCPL
get "libhdr"
let mulinv(a, b) =
b<0 -> mulinv(a, -b),
a<0 -> mulinv(b - (-a rem b), b),
valof
$( let t, nt, r, nr = 0, 1, b, a rem b
until nr = 0
$( let tmp, q = ?, r / nr
tmp := nt ; nt := t - q*nt ; t := tmp
tmp := nr ; nr := r - q*nr ; r := tmp
$)
resultis r>1 -> -1,
t<0 -> t + b,
t
$)
let show(a, b) be
$( let mi = mulinv(a, b)
test mi>=0
do writef("%N, %N -> %N*N", a, b, mi)
or writef("%N, %N -> no inverse*N", a, b)
$)
let start() be
$( show(42, 2017)
show(40, 1)
show(52, -217)
show(-486, 217)
show(40, 2018)
$)
- Output:
42, 2017 -> 1969 40, 1 -> 0 52, -217 -> 96 -486, 217 -> 121 40, 2018 -> no inverse
Bracmat
( ( mod-inv
= a b b0 x0 x1 q
. !arg:(?a.?b)
& ( !b:1
| (!b.0.1):(?b0.?x0.?x1)
& whl
' ( !a:>1
& div$(!a.!b):?q
& (!b.mod$(!a.!b)):(?a.?b)
& (!x1+-1*!q*!x0.!x0):(?x0.?x1)
)
& (!x:>0|!x1+!b0)
)
)
& out$(mod-inv$(42.2017))
};
Output
1969
C
#include <stdio.h>
int mul_inv(int a, int b)
{
int b0 = b, t, q;
int x0 = 0, x1 = 1;
if (b == 1) return 1;
while (a > 1) {
q = a / b;
t = b, b = a % b, a = t;
t = x0, x0 = x1 - q * x0, x1 = t;
}
if (x1 < 0) x1 += b0;
return x1;
}
int main(void) {
printf("%d\n", mul_inv(42, 2017));
return 0;
}
The above method has some problems. Most importantly, when given a pair (a,b) with no solution, it generates an FP exception. When given b=1, it returns 1 which is not a valid result mod 1. When given negative a or b the results are incorrect. The following generates results that should match Pari/GP for numbers in the int range.
#include <stdio.h>
int mul_inv(int a, int b)
{
int t, nt, r, nr, q, tmp;
if (b < 0) b = -b;
if (a < 0) a = b - (-a % b);
t = 0; nt = 1; r = b; nr = a % b;
while (nr != 0) {
q = r/nr;
tmp = nt; nt = t - q*nt; t = tmp;
tmp = nr; nr = r - q*nr; r = tmp;
}
if (r > 1) return -1; /* No inverse */
if (t < 0) t += b;
return t;
}
int main(void) {
printf("%d\n", mul_inv(42, 2017));
printf("%d\n", mul_inv(40, 1));
printf("%d\n", mul_inv(52, -217)); /* Pari semantics for negative modulus */
printf("%d\n", mul_inv(-486, 217));
printf("%d\n", mul_inv(40, 2018));
return 0;
}
- Output:
1969 0 96 121 -1
C#
public class Program
{
static void Main()
{
System.Console.WriteLine(42.ModInverse(2017));
}
}
public static class IntExtensions
{
public static int ModInverse(this int a, int m)
{
if (m == 1) return 0;
int m0 = m;
(int x, int y) = (1, 0);
while (a > 1) {
int q = a / m;
(a, m) = (m, a % m);
(x, y) = (y, x - q * y);
}
return x < 0 ? x + m0 : x;
}
}
C++
Iterative implementation
#include <iostream>
int mul_inv(int a, int b)
{
int b0 = b, t, q;
int x0 = 0, x1 = 1;
if (b == 1) return 1;
while (a > 1) {
q = a / b;
t = b, b = a % b, a = t;
t = x0, x0 = x1 - q * x0, x1 = t;
}
if (x1 < 0) x1 += b0;
return x1;
}
int main(void) {
std::cout << mul_inv(42, 2017) << std::endl;
return 0;
}
Recursive implementation
#include <iostream>
short ObtainMultiplicativeInverse(int a, int b, int s0 = 1, int s1 = 0)
{
return b==0? s0: ObtainMultiplicativeInverse(b, a%b, s1, s0 - s1*(a/b));
}
int main(int argc, char* argv[])
{
std::cout << ObtainMultiplicativeInverse(42, 2017) << std::endl;
return 0;
}
Clojure
(ns test-p.core
(:require [clojure.math.numeric-tower :as math]))
(defn extended-gcd
"The extended Euclidean algorithm--using Clojure code from RosettaCode for Extended Eucliean
(see http://en.wikipedia.orwiki/Extended_Euclidean_algorithm)
Returns a list containing the GCD and the Bézout coefficients
corresponding to the inputs with the result: gcd followed by bezout coefficients "
[a b]
(cond (zero? a) [(math/abs b) 0 1]
(zero? b) [(math/abs a) 1 0]
:else (loop [s 0
s0 1
t 1
t0 0
r (math/abs b)
r0 (math/abs a)]
(if (zero? r)
[r0 s0 t0]
(let [q (quot r0 r)]
(recur (- s0 (* q s)) s
(- t0 (* q t)) t
(- r0 (* q r)) r))))))
(defn mul_inv
" Get inverse using extended gcd. Extended GCD returns
gcd followed by bezout coefficients. We want the 1st coefficients
(i.e. second of extend-gcd result). We compute mod base so result
is between 0..(base-1) "
[a b]
(let [b (if (neg? b) (- b) b)
a (if (neg? a) (- b (mod (- a) b)) a)
egcd (extended-gcd a b)]
(if (= (first egcd) 1)
(mod (second egcd) b)
(str "No inverse since gcd is: " (first egcd)))))
(println (mul_inv 42 2017))
(println (mul_inv 40 1))
(println (mul_inv 52 -217))
(println (mul_inv -486 217))
(println (mul_inv 40 2018))
Output:
1969 0 96 121 No inverse since gcd is: 2
CLU
mul_inv = proc (a, b: int) returns (int) signals (no_inverse)
if b<0 then b := -b end
if a<0 then a := b - (-a // b) end
t: int := 0
nt: int := 1
r: int := b
nr: int := a // b
while nr ~= 0 do
q: int := r / nr
t, nt := nt, t - q*nt
r, nr := nr, r - q*nr
end
if r>1 then signal no_inverse end
if t<0 then t := t+b end
return(t)
end mul_inv
start_up = proc ()
pair = struct[a, b: int]
tests: sequence[pair] := sequence[pair]$
[pair${a: 42, b: 2017},
pair${a: 40, b: 1},
pair${a: 52, b: -217},
pair${a: -486, b: 217},
pair${a: 40, b: 2018}]
po: stream := stream$primary_output()
for test: pair in sequence[pair]$elements(tests) do
stream$puts(po, int$unparse(test.a) || ", "
|| int$unparse(test.b) || " -> ")
stream$putl(po, int$unparse(mul_inv(test.a, test.b)))
except when no_inverse:
stream$putl(po, "no modular inverse")
end
end
end start_up
- Output:
42, 2017 -> 1969 40, 1 -> 0 52, -217 -> 96 -486, 217 -> 121 40, 2018 -> no modular inverse
Comal
0010 FUNC mulinv#(a#,b#) CLOSED
0020 IF b#<0 THEN b#:=-b#
0030 IF a#<0 THEN a#:=b#-(-a# MOD b#)
0040 t#:=0;nt#:=1;r#:=b#;nr#:=a# MOD b#
0050 WHILE nr#<>0 DO
0060 q#:=r# DIV nr#
0070 tmp#:=nt#;nt#:=t#-q#*nt#;t#:=tmp#
0080 tmp#:=nr#;nr#:=r#-q#*nr#;r#:=tmp#
0090 ENDWHILE
0100 IF r#>1 THEN RETURN -1
0110 IF t#<0 THEN t#:+b#
0120 RETURN t#
0130 ENDFUNC mulinv#
0140 //
0150 WHILE NOT EOD DO
0160 READ a#,b#
0170 PRINT a#,", ",b#," -> ",mulinv#(a#,b#)
0180 ENDWHILE
0190 END
0200 //
0210 DATA 42,2017,40,1,52,-217,-486,217,40,2018
- Output:
42, 2017 -> 1969 40, 1 -> 0 52, -217 -> 96 -486, 217 -> 121 40, 2018 -> -1
Common Lisp
;;
;; Calculates the GCD of a and b based on the Extended Euclidean Algorithm. The function also returns
;; the Bézout coefficients s and t, such that gcd(a, b) = as + bt.
;;
;; The algorithm is described on page http://en.wikipedia.org/wiki/Extended_Euclidean_algorithm#Iterative_method_2
;;
(defun egcd (a b)
(do ((r (cons b a) (cons (- (cdr r) (* (car r) q)) (car r))) ; (r+1 r) i.e. the latest is first.
(s (cons 0 1) (cons (- (cdr s) (* (car s) q)) (car s))) ; (s+1 s)
(u (cons 1 0) (cons (- (cdr u) (* (car u) q)) (car u))) ; (t+1 t)
(q nil))
((zerop (car r)) (values (cdr r) (cdr s) (cdr u))) ; exit when r+1 = 0 and return r s t
(setq q (floor (/ (cdr r) (car r)))))) ; inside loop; calculate the q
;;
;; Calculates the inverse module for a = 1 (mod m).
;;
;; Note: The inverse is only defined when a and m are coprimes, i.e. gcd(a, m) = 1.”
;;
(defun invmod (a m)
(multiple-value-bind (r s k) (egcd a m)
(unless (= 1 r) (error "invmod: Values ~a and ~a are not coprimes." a m))
s))
- Output:
* (invmod 42 2017) -48 * (mod -48 2017) 1969
Cowgol
include "cowgol.coh";
sub mulinv(a: int32, b: int32): (t: int32) is
if b<0 then b := -b; end if;
if a<0 then a := b - (-a % b); end if;
t := 0;
var nt: int32 := 1;
var r := b;
var nr := a % b;
while nr != 0 loop
var q := r / nr;
var tmp := nt; nt := t - q*nt; t := tmp;
tmp := nr; nr := r - q*nr; r := tmp;
end loop;
if r>1 then t := -1;
elseif t<0 then t := t + b;
end if;
end sub;
record Pair is
a: int32;
b: int32;
end record;
var data: Pair[] := {
{42, 2017},
{40, 1},
{52, -217},
{-486, 217},
{40, 2018}
};
var i: @indexof data := 0;
while i < @sizeof data loop
print_i32(data[i].a as uint32);
print(", ");
print_i32(data[i].b as uint32);
print(" -> ");
var mi := mulinv(data[i].a, data[i].b);
if mi<0
then print("no inverse");
else print_i32(mi as uint32);
end if;
print_nl();
i := i + 1;
end loop;
- Output:
42, 2017 -> 1969 40, 1 -> 0 52, 4294967079 -> 96 4294966810, 217 -> 121 40, 2018 -> no inverse
Craft Basic
let e = 42
let t = 2017
gosub modularinverse
end
sub modularinverse
let d = 0
if e < t then
let b = e
let c = 1
do
let s = int(((t - b) / e) + 1)
let b = b + s * e
let c = c + s
let b = b - t
loop b <> 1
let d = c
endif
let m = d
print m
return
- Output:
1969
Crystal
def modinv(a0, m0)
return 1 if m0 == 1
a, m = a0, m0
x0, inv = 0, 1
while a > 1
inv -= (a // m) * x0
a, m = m, a % m
x0, inv = inv, x0
end
inv += m0 if inv < 0
inv
end
- Output:
> modinv(42,2017) => 1969
D
T modInverse(T)(T a, T b) pure nothrow {
if (b == 1)
return 1;
T b0 = b,
x0 = 0,
x1 = 1;
while (a > 1) {
immutable q = a / b;
auto t = b;
b = a % b;
a = t;
t = x0;
x0 = x1 - q * x0;
x1 = t;
}
return (x1 < 0) ? (x1 + b0) : x1;
}
void main() {
import std.stdio;
writeln(modInverse(42, 2017));
}
- Output:
1969
dc
This solution prints the inverse u
only if it exists (a*u = 1 mod m
).
dc -e "[m=]P?dsm[a=]P?dsa1sv[dsb~rsqlbrldlqlv*-lvsdsvd0<x]dsxxldd[dlmr+]sx0>xdla*lm%[p]sx1=x"
If ~
is not implemented, it can be replaced by SdSnlnld/LnLd%
.
Replace [p]sx1=x
at the end by [pq]sx1=x16i6E6F7420636F7072696D65P
if an error message "not coprime" is desired.
- Output:
m=2 800^1+ a=37 342411551958695219479776173037037562556082184118925013641969995739234\ 344644689214483533004909620355470582887300743869103978073598454778206\ 829469635119691272637318902731800747596752668736012071540136041369140\ 1228044652005748974399041408477572
m=2017 a=42 1969
m=42 a=7
Delphi
See #Pascal.
Draco
proc mulinv(int a, b) int:
int t, nt, r, nr, q, tmp;
if b<0 then b := -b fi;
if a<0 then a := b - (-a % b) fi;
t := 0; nt := 1; r := b; nr := a % b;
while nr /= 0 do
q := r / nr;
tmp := nt; nt := t - q*nt; t := tmp;
tmp := nr; nr := r - q*nr; r := tmp
od;
if r>1 then -1
elif t<0 then t+b
else t
fi
corp
proc show(int a, b) void:
int mi;
mi := mulinv(a, b);
if mi>=0
then writeln(a:5, ", ", b:5, " -> ", mi:5)
else writeln(a:5, ", ", b:5, " -> no inverse")
fi
corp
proc main() void:
show(42, 2017);
show(40, 1);
show(52, -217);
show(-486, 217);
show(40, 2018)
corp
- Output:
42, 2017 -> 1969 40, 1 -> 0 52, -217 -> 96 -486, 217 -> 121 40, 2018 -> no inverse
EasyLang
func mod_inv a b .
b0 = b
x1 = 1
if b = 1
return 1
.
while a > 1
q = a div b
t = b
b = a mod b
a = t
t = x0
x0 = x1 - q * x0
x1 = t
.
if x1 < 0
x1 += b0
.
return x1
.
print mod_inv 42 2017
EchoLisp
(lib 'math) ;; for egcd = extended gcd
(define (mod-inv x m)
(define-values (g inv q) (egcd x m))
(unless (= 1 g) (error 'not-coprimes (list x m) ))
(if (< inv 0) (+ m inv) inv))
(mod-inv 42 2017) → 1969
(mod-inv 42 666)
🔴 error: not-coprimes (42 666)
Elixir
defmodule Modular do
def extended_gcd(a, b) do
{last_remainder, last_x} = extended_gcd(abs(a), abs(b), 1, 0, 0, 1)
{last_remainder, last_x * (if a < 0, do: -1, else: 1)}
end
defp extended_gcd(last_remainder, 0, last_x, _, _, _), do: {last_remainder, last_x}
defp extended_gcd(last_remainder, remainder, last_x, x, last_y, y) do
quotient = div(last_remainder, remainder)
remainder2 = rem(last_remainder, remainder)
extended_gcd(remainder, remainder2, x, last_x - quotient*x, y, last_y - quotient*y)
end
def inverse(e, et) do
{g, x} = extended_gcd(e, et)
if g != 1, do: raise "The maths are broken!"
rem(x+et, et)
end
end
IO.puts Modular.inverse(42,2017)
- Output:
1969
ERRE
PROGRAM MOD_INV
!$INTEGER
PROCEDURE MUL_INV(A,B->T)
LOCAL NT,R,NR,Q,TMP
IF B<0 THEN B=-B
IF A<0 THEN A=B-(-A MOD B)
T=0 NT=1 R=B NR=A MOD B
WHILE NR<>0 DO
Q=R DIV NR
TMP=NT NT=T-Q*NT T=TMP
TMP=NR NR=R-Q*NR R=TMP
END WHILE
IF (R>1) THEN T=-1 EXIT PROCEDURE ! NO INVERSE
IF (T<0) THEN T+=B
END PROCEDURE
BEGIN
MUL_INV(42,2017->T) PRINT(T)
MUL_INV(40,1->T) PRINT(T)
MUL_INV(52,-217->T) PRINT(T) ! pari semantics for negative modulus
MUL_INV(-486,217->T) PRINT(T)
MUL_INV(40,2018->T) PRINT(T)
END PROGRAM
- Output:
1969 0 96 121 -1
F#
// Calculate the inverse of a (mod m)
// See here for eea specs:
// https://en.wikipedia.org/wiki/Extended_Euclidean_algorithm
let modInv m a =
let rec eea t t' r r' =
match r' with
| 0 -> t
| _ ->
let div = r/r'
eea t' (t - div * t') r' (r - div * r')
(m + eea 0 1 m a) % m
- Input:
// Inverse of 347 (mod 29) modInv 29 347
- Output:
28
Factor
USE: math.functions
42 2017 mod-inv
- Output:
1969
Forth
ANS Forth with double-number word set
: invmod { a m | v b c -- inv }
m to v
1 to c
0 to b
begin a
while v a / >r
c b s>d c s>d r@ 1 m*/ d- d>s to c to b
a v s>d a s>d r> 1 m*/ d- d>s to a to v
repeat b m mod dup to b 0<
if m b + else b then ;
ANS Forth version without locals
: modinv ( a m - inv)
dup 1- \ a m (m != 1)?
if \ a m
tuck 1 0 \ m0 a m 1 0
begin \ m0 a m inv x0
2>r over 1 > \ m0 a m (a > 1)? R: inv x0
while \ m0 a m R: inv x0
tuck /mod \ m0 m (a mod m) (a/m) R: inv x0
r> tuck * \ m0 a' m' x0 (a/m)*x0 R: inv
r> swap - \ m0 a' m' x0 (inv-q) R:
repeat \ m0 a' m' inv' x0'
2drop \ m0 R: inv x0
2r> drop \ m0 inv R:
dup 0< \ m0 inv (inv < 0)?
if over + then \ m0 (inv + m0)
then \ x inv'
nip \ inv
;
42 2017 invmod . 1969 42 2017 modinv . 1969
Fortran
program modular_inverse_task
implicit none
write (*,*) inverse (42, 2017)
contains
! Returns -1 if there is no inverse. I assume n > 0. The algorithm
! is described at
! https://en.wikipedia.org/w/index.php?title=Extended_Euclidean_algorithm&oldid=1135569411#Modular_integers
function inverse (a, n) result (inverse_value)
integer, intent(in) :: a, n
integer :: inverse_value
integer :: t, newt
integer :: r, newr
integer :: quotient, remainder, tmp
if (n <= 0) error stop
t = 0; newt = 1
r = n; newr = a
do while (newr /= 0)
remainder = modulo (r, newr) ! Floor division.
quotient = (r - remainder) / newr
tmp = newt; newt = t - (quotient * newt); t = tmp
r = newr; newr = remainder
end do
if (r > 1) then
inverse_value = -1
else if (t < 0) then
inverse_value = t + n
else
inverse_value = t
end if
end function inverse
end program modular_inverse_task
- Output:
$ gfortran -Wall -Wextra modular_inverse_task.f90 && ./a.out 1969
FreeBASIC
' version 10-07-2018
' compile with: fbc -s console
Type ext_euclid
Dim As Integer a, b
End Type
' "Table method" aka "The Magic Box"
Function magic_box(x As Integer, y As Integer) As ext_euclid
Dim As Integer a(1 To 128), b(1 To 128), d(1 To 128), k(1 To 128)
a(1) = 1 : b(1) = 0 : d(1) = x
a(2) = 0 : b(2) = 1 : d(2) = y : k(2) = x \ y
Dim As Integer i = 2
While Abs(d(i)) <> 1
i += 1
a(i) = a(i -2) - k(i -1) * a(i -1)
b(i) = b(i -2) - k(i -1) * b(i -1)
d(i) = d(i -2) Mod d(i -1)
k(i) = d(i -1) \ d(i)
'Print a(i),b(i),d(i),k(i)
If d(i -1) Mod d(i) = 0 Then Exit While
Wend
If d(i) = -1 Then ' -1 * (ab + by) = -1 * -1 ==> -ab -by = 1
a(i) = -a(i)
b(i) = -b(i)
End If
Function = Type( a(i), b(i) )
End Function
' ------=< MAIN >=------
Dim As Integer x, y, gcd
Dim As ext_euclid result
Do
Read x, y
If x = 0 AndAlso y = 0 Then Exit Do
result = magic_box(x, y)
With result
gcd = .a * x + .b * y
Print "a * "; Str(x); " + b * "; Str(y);
Print " = GCD("; Str(x); ", "; Str(y); ") ="; gcd
If gcd > 1 Then
Print "No solution, numbers are not coprime"
Else
Print "a = "; .a; ", b = ";.b
Print "The Modular inverse of "; x; " modulo "; y; " = ";
While .a < 0 : .a += IIf(y > 0, y, -y) : Wend
Print .a
'Print "The Modular inverse of "; y; " modulo "; x; " = ";
'While .b < 0 : .b += IIf(x > 0, x, -x) : Wend
'Print .b
End if
End With
Print
Loop
Data 42, 2017
Data 40, 1
Data 52, -217
Data -486, 217
Data 40, 2018
Data 0, 0
' empty keyboard buffer
While Inkey <> "" : Wend
Print : Print "hit any key to end program"
Sleep
End
- Output:
a * 42 + b * 2017 = GCD(42, 2017) = 1 a = -48, b = 1 The Modular inverse of 42 modulo 2017 = 1969 a * 40 + b * 1 = GCD(40, 1) = 1 a = 0, b = 1 The Modular inverse of 40 modulo 1 = 0 a * 52 + b * -217 = GCD(52, -217) = 1 a = 96, b = 23 The Modular inverse of 52 modulo -217 = 96 a * -486 + b * 217 = GCD(-486, 217) = 1 a = -96, b = -215 The Modular inverse of -486 modulo 217 = 121 a * 40 + b * 2018 = GCD(40, 2018) = 2 No solution, numbers are not coprime
Frink
println[modInverse[42, 2017]]
- Output:
1969
FunL
import integers.egcd
def modinv( a, m ) =
val (g, x, _) = egcd( a, m )
if g != 1 then error( a + ' and ' + m + ' not coprime' )
val res = x % m
if res < 0 then res + m else res
println( modinv(42, 2017) )
- Output:
1969
Go
The standard library function uses the extended Euclidean algorithm internally.
package main
import (
"fmt"
"math/big"
)
func main() {
a := big.NewInt(42)
m := big.NewInt(2017)
k := new(big.Int).ModInverse(a, m)
fmt.Println(k)
}
- Output:
1969
GW-BASIC
10 ' Modular inverse
20 LET E% = 42
30 LET T% = 2017
40 GOSUB 1000
50 PRINT MODINV%
60 END
990 ' increments e stp (step) times until bal is greater than t
992 ' repeats until bal = 1 (mod = 1) and returns count
994 ' bal will not be greater than t + e
1000 LET D% = 0
1010 IF E% >= T% THEN GOTO 1140
1020 LET BAL% = E%
1025 ' At least one iteration is necessary
1030 LET STP% = ((T% - BAL%) \ E%) + 1
1040 LET BAL% = BAL% + STP% * E%
1050 LET COUNT% = 1 + STP%
1060 LET BAL% = BAL% - T%
1070 WHILE BAL% <> 1
1080 LET STP% = ((T% - BAL%) \ E%) + 1
1090 LET BAL% = BAL% + STP% * E%
1100 LET COUNT% = COUNT% + STP%
1110 LET BAL% = BAL% - T%
1120 WEND
1130 LET D% = COUNT%
1140 LET MODINV% = D%
1150 RETURN
- Output:
1969
Haskell
-- Given a and m, return Just x such that ax = 1 mod m.
-- If there is no such x return Nothing.
modInv :: Int -> Int -> Maybe Int
modInv a m
| 1 == g = Just (mkPos i)
| otherwise = Nothing
where
(i, _, g) = gcdExt a m
mkPos x
| x < 0 = x + m
| otherwise = x
-- Extended Euclidean algorithm.
-- Given non-negative a and b, return x, y and g
-- such that ax + by = g, where g = gcd(a,b).
-- Note that x or y may be negative.
gcdExt :: Int -> Int -> (Int, Int, Int)
gcdExt a 0 = (1, 0, a)
gcdExt a b =
let (q, r) = a `quotRem` b
(s, t, g) = gcdExt b r
in (t, s - q * t, g)
main :: IO ()
main = mapM_ print [2 `modInv` 4, 42 `modInv` 2017]
- Output:
Nothing Just 1969
Icon and Unicon
procedure main(args)
a := integer(args[1]) | 42
b := integer(args[2]) | 2017
write(mul_inv(a,b))
end
procedure mul_inv(a,b)
if b == 1 then return 1
(b0 := b, x0 := 0, x1 := 1)
while a > 1 do {
q := a/b
(t := b, b := a%b, a := t)
(t := x0, x0 := x1-q*x0, x1 := t)
}
return if (x1 > 0) then x1 else x1+b0
end
- Output:
->mi 1969 ->
Adding a coprime test:
link numbers
procedure main(args)
a := integer(args[1]) | 42
b := integer(args[2]) | 2017
write(mul_inv(a,b))
end
procedure mul_inv(a,b)
if b == 1 then return 1
if gcd(a,b) ~= 1 then return "not coprime"
(b0 := b, x0 := 0, x1 := 1)
while a > 1 do {
q := a/b
(t := b, b := a%b, a := t)
(t := x0, x0 := x1-q*x0, x1 := t)
}
return if (x1 > 0) then x1 else x1+b0
end
IS-BASIC
100 PRINT MODINV(42,2017)
120 DEF MODINV(A,B)
130 LET B=ABS(B)
140 IF A<0 THEN LET A=B-MOD(-A,B)
150 LET T=0:LET NT=1:LET R=B:LET NR=MOD(A,B)
160 DO WHILE NR<>0
170 LET Q=INT(R/NR)
180 LET TMP=NT:LET NT=T-Q*NT:LET T=TMP
190 LET TMP=NR:LET NR=R-Q*NR:LET R=TMP
200 LOOP
210 IF R>1 THEN
220 LET MODINV=-1
230 ELSE IF T<0 THEN
240 LET MODINV=T+B
250 ELSE
260 LET MODINV=T
270 END IF
280 END DEF
J
Solution:
modInv =: dyad def 'x y&|@^ <: 5 p: y'"0
Example:
42 modInv 2017
1969
Notes:
- Calculates the modular inverse as a^( totient(m) - 1 ) mod m.
- 5 p: y is Euler's totient function of y.
- J has a fast implementation of modular exponentiation (which avoids the exponentiation altogether), invoked with the form m&|@^ (hence, we use explicitly-named arguments for this entry, as opposed to the "variable free" tacit style: the m&| construct must freeze the value before it can be used but we want to use different values in that expression at different times...).
Java
The BigInteger
library has a method for this:
System.out.println(BigInteger.valueOf(42).modInverse(BigInteger.valueOf(2017)));
- Output:
1969
Alternatively, working from first principles.
public final class ModularInverse {
public static void main(String[] aArgs) {
System.out.println(inverseModulus(42, 2017));
}
private static Egcd extendedGCD(int aOne, int aTwo) {
if ( aOne == 0 ) {
return new Egcd(aTwo, 0, 1);
}
Egcd value = extendedGCD(aTwo % aOne, aOne);
return new Egcd(value.aG, value.aY - ( aTwo / aOne ) * value.aX, value.aX);
}
private static int inverseModulus(int aNumber, int aModulus) {
Egcd value = extendedGCD(aNumber, aModulus);
return ( value.aG == 1 ) ? ( value.aX + aModulus ) % aModulus : 0;
}
private static record Egcd(int aG, int aX, int aY) {}
}
- Output:
1969
JavaScript
Using brute force.
var modInverse = function(a, b) {
a %= b;
for (var x = 1; x < b; x++) {
if ((a*x)%b == 1) {
return x;
}
}
}
jq
Works with gojq, the Go implementation of jq
# Integer division:
# If $j is 0, then an error condition is raised;
# otherwise, assuming infinite-precision integer arithmetic,
# if the input and $j are integers, then the result will be an integer.
def idivide($j):
. as $i
| ($i % $j) as $mod
| ($i - $mod) / $j ;
# the multiplicative inverse of . modulo $n
def modInv($n):
if $n == 1 then 1
else . as $this
| { r : $n,
t : 0,
newR: length, # abs
newT: 1}
| until(.newR == 0;
.newR as $newR
| (.r | idivide($newR)) as $q
| {r : $newR,
t : .newT,
newT: (.t - $q * .newT),
newR: (.r - $q * $newR) } )
| if (.r|length) != 1 then "\($this) and \($n) are not co-prime." | error
else .t
| if . < 0 then . + $n
elif $this < 0 then - .
else .
end
end
end ;
# Example:
42 | modInv(2017)
- Output:
1969
Julia
Built-in
Julia includes a built-in function for this:
invmod(a, b)
C translation
The following code works on any integer type.
To maximize performance, we ensure (via a promotion rule) that the operands are the same type (and use built-ins zero(T)
and one(T)
for initialization of temporary variables to ensure that they remain of the same type throughout execution).
function modinv(a::T, b::T) where T <: Integer
b0 = b
x0, x1 = zero(T), one(T)
while a > 1
q = div(a, b)
a, b = b, a % b
x0, x1 = x1 - q * x0, x0
end
x1 < 0 ? x1 + b0 : x1
end
modinv(a::Integer, b::Integer) = modinv(promote(a,b)...)
- Output:
julia> invmod(42, 2017) 1969 julia> modinv(42, 2017) 1969
Kotlin
// version 1.0.6
import java.math.BigInteger
fun main(args: Array<String>) {
val a = BigInteger.valueOf(42)
val m = BigInteger.valueOf(2017)
println(a.modInverse(m))
}
- Output:
1969
Lambdatalk
{def mulinv
{def mulinv.loop
{lambda {:t :nt :r :nr}
{if {not {= :nr 0}}
then {mulinv.loop :nt
{- :t {* {floor {/ :r :nr}} :nt}}
:nr
{- :r {* {floor {/ :r :nr}} :nr}} }
else {cons :t :r} }}}
{lambda {:a :n}
{let { {:a :a} {:n :n}
{:cons {mulinv.loop 0
1
{if {< :n 0} then {- :n} else :n}
{if {< :a 0} then {- :n {% {- :a} :n}} else :a}}}
} {if {> {cdr :cons} 1}
then not invertible
else {if {< {car :cons} 0}
then {+ {car :cons} :n}
else {car :cons} }}}}}
-> mulinv
{mulinv 42 2017}
-> 1969
{mulinv 40 1}
-> 0
{mulinv 52 -217}
-> 96
{mulinv -486 217}
-> 121
{mulinv 40 218}
-> not invertible
m4
Note that $0
is the name of the macro being evaluated. Therefore, in the following, _$0
is the name of another macro, the same as the name of the first macro, except for an underscore prepended. This is a common idiom.
The core of the program is __inverse
recursively calling itself.
m4 is a macro-preprocessor, and so some of the following is there simply to keep things from being echoed to the output. :)
divert(-1)
# I assume non-negative arguments. The algorithm is described at
# https://en.wikipedia.org/w/index.php?title=Extended_Euclidean_algorithm&oldid=1135569411#Modular_integers
define(`inverse',`_$0(eval(`$1'), eval(`$2'))')
define(`_inverse',`_$0($2, 0, 1, $2, $1)')
define(`__inverse',
`dnl n = $1, t = $2, newt = $3, r = $4, newr = $5
ifelse(eval($5 != 0), 1, `$0($1, $3,
eval($2 - (($4 / $5) * $3)),
$5,eval($4 % $5))',
eval($4 > 1), 1, `no inverse',
eval($2 < 0), 1, eval($2 + $1),
$2)')
divert`'dnl
inverse(42, 2017)
- Output:
$ m4 modular-inverse-task.m4 1969
MAD
NORMAL MODE IS INTEGER
INTERNAL FUNCTION(AA, BB)
ENTRY TO MULINV.
A = AA
B = BB
WHENEVER B.L.0, B = -B
WHENEVER A.L.0, A = B - (-(A-A/B*B))
T = 0
NT = 1
R = B
NR = A-A/B*B
LOOP WHENEVER NR.NE.0
Q = R/NR
TMP = NT
NT = T - Q*NT
T = TMP
TMP = NR
NR = R - Q*NR
R = TMP
TRANSFER TO LOOP
END OF CONDITIONAL
WHENEVER R.G.1, FUNCTION RETURN -1
WHENEVER T.L.0, T = T+B
FUNCTION RETURN T
END OF FUNCTION
INTERNAL FUNCTION(AA, BB)
VECTOR VALUES FMT = $I5,2H, ,I5,2H: ,I5*$
ENTRY TO SHOW.
PRINT FORMAT FMT, AA, BB, MULINV.(AA, BB)
END OF FUNCTION
SHOW.(42,2017)
SHOW.(40,1)
SHOW.(52,-217)
SHOW.(-486,217)
SHOW.(40,2018)
END OF PROGRAM
- Output:
42, 2017: 1969 40, 1: 0 52, -217: 96 -486, 217: 121 40, 2018: -1
Maple
1/42 mod 2017;
- Output:
1969
Mathematica /Wolfram Language
Use the built-in function ModularInverse
:
ModularInverse[a, m]
For example:
ModularInverse[42, 2017] 1969
Maxima
Using built-in function inv_mod
inv_mod(42,2017);
- Output:
1969
Mercury
%%% -*- mode: mercury; prolog-indent-width: 2; -*-
%%%
%%% Compile with:
%%% mmc --make --use-subdirs modular_inverse_task
%%%
:- module modular_inverse_task.
:- interface.
:- import_module io.
:- pred main(io::di, io::uo) is det.
:- implementation.
:- import_module exception.
:- import_module int.
%% inverse(A, N, Inverse). I assume N > 0, and throw an exception if
%% it is not. The predicate fails if there is no inverse (and thus is
%% "semidet"). The algorithm is described at
%% https://en.wikipedia.org/w/index.php?title=Extended_Euclidean_algorithm&oldid=1135569411#Modular_integers
:- pred inverse(int::in, int::in, int::out) is semidet.
inverse(A, N, Inverse) :-
if (N =< 0) then throw(domain_error("inverse"))
else inverse_(N, 0, 1, N, A, Inverse).
:- pred inverse_(int::in, int::in, int::in, int::in, int::in,
int::out) is semidet.
inverse_(N, T, NewT, R, NewR, Inverse) :-
if (NewR \= 0)
then (Quotient = div(R, NewR), % Floor division.
inverse_(N,
NewT, T - (Quotient * NewT),
NewR, R - (Quotient * NewR),
Inverse)) % Tail recursion.
else (R =< 1, % R =< 1 FAILS if R > 1.
(if (T < 0) then Inverse = T + N else Inverse = T)).
main(!IO) :-
if inverse(42, 2017, Inverse)
then (print(Inverse, !IO), nl(!IO))
else (print("There is no inverse.", !IO), nl(!IO)).
:- end_module modular_inverse_task.
- Output:
$ mmc --make --use-subdirs modular_inverse_task && ./modular_inverse_task Making Mercury/int3s/modular_inverse_task.int3 Making Mercury/ints/modular_inverse_task.int Making Mercury/cs/modular_inverse_task.c Making Mercury/os/modular_inverse_task.o Making modular_inverse_task 1969
МК-61/52
П1 П2 <-> П0 0 П5 1 П6 ИП1 1
- x=0 14 С/П ИП0 1 - /-/ x<0 50
ИП0 ИП1 / [x] П4 ИП1 П3 ИП0 ^ ИП1
/ [x] ИП1 * - П1 ИП3 П0 ИП5 П3
ИП6 ИП4 ИП5 * - П5 ИП3 П6 БП 14
ИП6 x<0 55 ИП2 + С/П
Modula-2
MODULE ModularInverse;
FROM InOut IMPORT WriteString, WriteInt, WriteLn;
TYPE Data = RECORD x : INTEGER;
y : INTEGER
END;
VAR c : INTEGER;
ab : ARRAY [1..5] OF Data;
PROCEDURE mi(VAR a, b : INTEGER): INTEGER;
VAR t, nt, r, nr, q, tmp : INTEGER;
BEGIN
b := ABS(b);
IF a < 0 THEN a := b - (-a MOD b) END;
t := 0; nt := 1; r := b; nr := a MOD b;
WHILE (nr # 0) DO
q := r / nr;
tmp := nt; nt := t - q * nt; t := tmp;
tmp := nr; nr := r - q * nr; r := tmp;
END;
IF (r > 1) THEN RETURN -1 END;
IF (t < 0) THEN RETURN t + b END;
RETURN t;
END mi;
BEGIN
ab[1].x := 42; ab[1].y := 2017;
ab[2].x := 40; ab[2].y := 1;
ab[3].x := 52; ab[3].y := -217;
ab[4].x := -486; ab[4].y := 217;
ab[5].x := 40; ab[5].y := 2018;
WriteLn;
WriteString("Modular inverse");
WriteLn;
FOR c := 1 TO 5 DO
WriteInt(ab[c].x, 6); WriteString(", ");
WriteInt(ab[c].y, 6); WriteString(" = ");
WriteInt(mi(ab[c].x, ab[c].y),6);
WriteLn;
END;
END ModularInverse.
- Output:
Modular inverse 42, 2017 = 1969 40, 1 = 0 52, -217 = 96 -486, 217 = 121 40, 2018 = -1
newLISP
(define (modular-multiplicative-inverse a n)
(if (< n 0)
(setf n (abs n)))
(if (< a 0)
(setf a (- n (% (- 0 a) n))))
(setf t 0)
(setf nt 1)
(setf r n)
(setf nr (mod a n))
(while (not (zero? nr))
(setf q (int (div r nr)))
(setf tmp nt)
(setf nt (sub t (mul q nt)))
(setf t tmp)
(setf tmp nr)
(setf nr (sub r (mul q nr)))
(setf r tmp))
(if (> r 1)
(setf retvalue nil))
(if (< t 0)
(setf retvalue (add t n))
(setf retvalue t))
retvalue)
(println (modular-multiplicative-inverse 42 2017))
Output:
1969
Nim
proc modInv(a0, b0: int): int =
var (a, b, x0) = (a0, b0, 0)
result = 1
if b == 1: return
while a > 1:
result = result - (a div b) * x0
a = a mod b
swap a, b
swap x0, result
if result < 0: result += b0
echo modInv(42, 2017)
- Output:
1969
Oberon-2
(*-*- mode: indented-text; tab-width: 2; -*-*)
MODULE modularInverseInOberon2;
IMPORT Out;
(* Division with a non-negative remainder. This will work no matter
how your compiler handles DIV (and mine seems not to do what the
Oberon-2 specification says). *)
PROCEDURE euclidDiv (x, y : INTEGER) : INTEGER;
VAR q : INTEGER;
BEGIN
IF 0 <= y THEN (* Do floor division. *)
IF 0 <= x THEN
q := x DIV y
ELSE
q := -((-x) DIV y);
IF (-x) MOD y # 0 THEN q := q - 1 END
END;
ELSE (* Do ceiling division. *)
IF 0 <= x THEN
q := -(x DIV (-y))
ELSE
q := ((-x) DIV (-y));
IF (-x) MOD (-y) # 0 THEN q := q + 1 END
END
END;
RETURN q
END euclidDiv;
(* I have added this unit test because, earlier, I posted a buggy
version of euclidDiv. *)
PROCEDURE testEuclidDiv;
VAR x, y, q, r : INTEGER;
BEGIN
FOR x := -100 TO 100 DO
FOR y := -100 TO 100 DO
IF y # 0 THEN
q := euclidDiv (x, y);
r := x - (q * y);
IF (r < 0) OR (ABS (y) <= r) THEN
(* A remainder was outside the expected range. *)
Out.String ("euclidDiv fails its test")
END
END
END
END
END testEuclidDiv;
PROCEDURE inverse (a, n : INTEGER) : INTEGER;
VAR t, newt : INTEGER;
VAR r, newr : INTEGER;
VAR quotient : INTEGER;
VAR tmp : INTEGER;
BEGIN
t := 0; newt := 1;
r := n; newr := a;
WHILE newr # 0 DO
quotient := euclidDiv (r, newr);
tmp := newt; newt := t - (quotient * newt); t := tmp;
tmp := newr; newr := r - (quotient * newr); r := tmp
END;
IF r > 1 THEN
t := -1
ELSIF t < 0 THEN
t := t + n
END;
RETURN t
END inverse;
BEGIN
testEuclidDiv;
Out.Int (inverse (42, 2017), 0);
Out.Ln
END modularInverseInOberon2.
- Output:
$ obc modularInverseInOberon2.Mod && ./a.out 1969
ObjectIcon
# -*- ObjectIcon -*-
import exception
import io
procedure main ()
test_euclid_div ()
io.write (inverse (42, 2017))
end
procedure inverse (a, n) # FAILS if there is no inverse.
local t, newt, r, newr, quotient, tmp
if n <= 0 then throw ("non-positive modulus")
t := 0; newt := 1
r := n; newr := a
while newr ~= 0 do
{
quotient := euclid_div (r, newr)
tmp := newt; newt := t - (quotient * newt); t := tmp
tmp := newr; newr := r - (quotient * newr); r := tmp
}
r <= 1 | fail
return (if t < 0 then t + n else t)
end
procedure euclid_div (x, y)
# This kind of integer division always gives a remainder between 0
# and abs(y)-1, inclusive. Thus the remainder is always a LEAST
# RESIDUE modulo abs(y). (If y is a positive modulus, then only the
# floor division branch is used.)
return \
if 0 <= y then # Do floor division.
(if 0 <= x then x / y
else if (-x) % y = 0 then -((-x) / y)
else -((-x) / y) - 1)
else # Do ceiling division.
(if 0 <= x then -(x / (-y))
else if (-x) % (-y) = 0 then ((-x) / (-y))
else ((-x) / (-y)) + 1)
end
procedure test_euclid_div ()
local x, y, q, r
every x := -100 to 100 do
every y := -100 to 100 & y ~= 0 do
{
q := euclid_div (x, y)
r := x - (q * y)
if r < 0 | abs (y) <= r then
# A remainder was outside the expected range.
throw ("Test of euclid_div failed.")
}
end
- Output:
$ oiscript modular-inverse-task-OI.icn 1969
OCaml
let mul_inv a = function 1 -> 1 | b ->
let rec aux a b x0 x1 =
if a <= 1 then x1 else
if b = 0 then failwith "mul_inv" else
aux b (a mod b) (x1 - (a / b) * x0) x0
in
let x = aux a b 0 1 in
if x < 0 then x + b else x
Testing:
# mul_inv 42 2017 ;; - : int = 1969
let rec gcd_ext a = function
| 0 -> (1, 0, a)
| b ->
let s, t, g = gcd_ext b (a mod b) in
(t, s - (a / b) * t, g)
let mod_inv a m =
let mk_pos x = if x < 0 then x + m else x in
match gcd_ext a m with
| i, _, 1 -> mk_pos i
| _ -> failwith "mod_inv"
Testing:
# mod_inv 42 2017 ;; - : int = 1969
Oforth
Usage : a modulus invmod
// euclid ( a b -- u v r )
// Return r = gcd(a, b) and (u, v) / r = au + bv
: euclid(a, b)
| q u u1 v v1 |
b 0 < ifTrue: [ b neg ->b ]
a 0 < ifTrue: [ b a neg b mod - ->a ]
1 dup ->u ->v1
0 dup ->v ->u1
while(b) [
b a b /mod ->q ->b ->a
u1 u u1 q * - ->u1 ->u
v1 v v1 q * - ->v1 ->v
]
u v a ;
: invmod(a, modulus)
a modulus euclid 1 == ifFalse: [ drop drop null return ]
drop dup 0 < ifTrue: [ modulus + ] ;
- Output:
42 2017 invmod println 1969
Owl Lisp
(import (owl math)
(owl math extra))
(define (euclid-quotient x y)
(if (<= 0 y)
(cond ((<= 0 x) (quotient x y))
((zero? (remainder (negate x) y))
(negate (quotient (negate x) y)))
(else (- (negate (quotient (negate x) y)) 1)))
(cond ((<= 0 x) (negate (quotient x (negate y))))
((zero? (remainder (negate x) (negate y)))
(quotient (negate x) (negate y)))
(else (+ (quotient (negate x) (negate y)) 1)))))
;; A unit test of euclid-quotient.
(let repeat ((x -100)
(y -100))
(cond ((= x 101) #t)
((= y 0) (repeat x (+ y 1)))
((= y 101) (repeat (+ x 1) -100))
(else (let* ((q (euclid-quotient x y))
(r (- x (* q y))))
(cond ((< r 0) (display "negative remainder\n"))
((<= (abs y) r) (display "remainder too large\n"))
(else (repeat x (+ y 1))))))))
(define (inverse a n)
(let repeat ((t 0) (newt 1)
(r n) (newr a))
(cond ((not (zero? newr))
(let ((quotient (euclid-quotient r newr)))
(repeat newt (- t (* quotient newt))
newr (- r (* quotient newr)))))
((< 1 r) #f) ; The inverse does not exist.
((negative? t) (+ t n))
(else t))))
(display (inverse 42 2017))
(newline)
- Output:
$ ol modular-inverse-task-Owl.scm 1969
PARI/GP
Mod(1/42,2017)
Pascal
// increments e step times until bal is greater than t
// repeats until bal = 1 (mod = 1) and returns count
// bal will not be greater than t + e
function modInv(e, t : integer) : integer;
var
d : integer;
bal, count, step : integer;
begin
d := 0;
if e < t then
begin
count := 1;
bal := e;
repeat
step := ((t-bal) DIV e)+1;
bal := bal + step * e;
count := count + step;
bal := bal - t;
until bal = 1;
d := count;
end;
modInv := d;
end;
Testing:
Writeln(modInv(42,2017));
- Output:
1969
Perl
Various CPAN modules can do this, such as:
use bigint; say 42->bmodinv(2017);
# or
use Math::ModInt qw/mod/; say mod(42, 2017)->inverse->residue;
# or
use Math::Pari qw/PARI lift/; say lift PARI "Mod(1/42,2017)";
# or
use Math::GMP qw/:constant/; say 42->bmodinv(2017);
# or
use ntheory qw/invmod/; say invmod(42, 2017);
or we can write our own:
sub invmod {
my($a,$n) = @_;
my($t,$nt,$r,$nr) = (0, 1, $n, $a % $n);
while ($nr != 0) {
# Use this instead of int($r/$nr) to get exact unsigned integer answers
my $quot = int( ($r - ($r % $nr)) / $nr );
($nt,$t) = ($t-$quot*$nt,$nt);
($nr,$r) = ($r-$quot*$nr,$nr);
}
return if $r > 1;
$t += $n if $t < 0;
$t;
}
say invmod(42,2017);
Notes: Special cases to watch out for include (1) where the inverse doesn't exist, such as invmod(14,28474), which should return undef or raise an exception, not return a wrong value. (2) the high bit of a or n is set, e.g. invmod(11,2**63), (3) negative first arguments, e.g. invmod(-11,23). The modules and code above handle these cases, but some other language implementations for this task do not.
Phix
function mul_inv(integer a, n) if n<0 then n = -n end if if a<0 then a = n - mod(-a,n) end if integer t = 0, nt = 1, r = n, nr = a; while nr!=0 do integer q = floor(r/nr) {t, nt} = {nt, t-q*nt} {r, nr} = {nr, r-q*nr} end while if r>1 then return "a is not invertible" end if if t<0 then t += n end if return t end function ?mul_inv(42,2017) ?mul_inv(40, 1) ?mul_inv(52, -217) /* Pari semantics for negative modulus */ ?mul_inv(-486, 217) ?mul_inv(40, 2018)
- Output:
1969 0 96 121 "a is not invertible"
PHP
Algorithm Implementation
<?php
function invmod($a,$n){
if ($n < 0) $n = -$n;
if ($a < 0) $a = $n - (-$a % $n);
$t = 0; $nt = 1; $r = $n; $nr = $a % $n;
while ($nr != 0) {
$quot= intval($r/$nr);
$tmp = $nt; $nt = $t - $quot*$nt; $t = $tmp;
$tmp = $nr; $nr = $r - $quot*$nr; $r = $tmp;
}
if ($r > 1) return -1;
if ($t < 0) $t += $n;
return $t;
}
printf("%d\n", invmod(42, 2017));
?>
- Output:
1969
PicoLisp
(de modinv (A B)
(let (B0 B X0 0 X1 1 Q 0 T1 0)
(while (< 1 A)
(setq
Q (/ A B)
T1 B
B (% A B)
A T1
T1 X0
X0 (- X1 (* Q X0))
X1 T1 ) )
(if (lt0 X1) (+ X1 B0) X1) ) )
(println
(modinv 42 2017) )
(bye)
PL/I
*process source attributes xref or(!);
/*--------------------------------------------------------------------
* 13.07.2015 Walter Pachl
*-------------------------------------------------------------------*/
minv: Proc Options(main);
Dcl (x,y) Bin Fixed(31);
x=42;
y=2017;
Put Edit('modular inverse of',x,' by ',y,' ---> ',modinv(x,y))
(Skip,3(a,f(4)));
modinv: Proc(a,b) Returns(Bin Fixed(31));
Dcl (a,b,ob,ox,d,t) Bin Fixed(31);
ob=b;
ox=0;
d=1;
If b=1 Then;
Else Do;
Do While(a>1);
q=a/b;
r=mod(a,b);
a=b;
b=r;
t=ox;
ox=d-q*ox;
d=t;
End;
End;
If d<0 Then
d=d+ob;
Return(d);
End;
End;
- Output:
modular inverse of 42 by 2017 ---> 1969
PowerShell
function invmod($a,$n){
if ([int]$n -lt 0) {$n = -$n}
if ([int]$a -lt 0) {$a = $n - ((-$a) % $n)}
$t = 0
$nt = 1
$r = $n
$nr = $a % $n
while ($nr -ne 0) {
$q = [Math]::truncate($r/$nr)
$tmp = $nt
$nt = $t - $q*$nt
$t = $tmp
$tmp = $nr
$nr = $r - $q*$nr
$r = $tmp
}
if ($r -gt 1) {return -1}
if ($t -lt 0) {$t += $n}
return $t
}
invmod 42 2017
- Output:
PS> .\INVMOD.PS1 1969 PS>
Prolog
egcd(_, 0, 1, 0) :- !.
egcd(A, B, X, Y) :-
divmod(A, B, Q, R),
egcd(B, R, S, X),
Y is S - Q*X.
modinv(A, B, N) :-
egcd(A, B, X, Y),
A*X + B*Y =:= 1,
N is X mod B.
- Output:
?- modinv(42, 2017, N). N = 1969. ?- modinv(42, 64, X). false.
PureBasic
Using brute force.
EnableExplicit
Declare main()
Declare.i mi(a.i, b.i)
If OpenConsole("MODULAR-INVERSE")
main() : Input() : End
EndIf
Macro ModularInverse(a, b)
PrintN(~"\tMODULAR-INVERSE(" + RSet(Str(a),5) + "," +
RSet(Str(b),5)+") = " +
RSet(Str(mi(a, b)),5))
EndMacro
Procedure main()
ModularInverse(42, 2017) ; = 1969
ModularInverse(40, 1) ; = 0
ModularInverse(52, -217) ; = 96
ModularInverse(-486, 217) ; = 121
ModularInverse(40, 2018) ; = -1
EndProcedure
Procedure.i mi(a.i, b.i)
Define x.i = 1,
y.i = Int(Abs(b)),
r.i = 0
If y = 1 : ProcedureReturn 0 : EndIf
While x < y
r = (a * x) % b
If r = 1 Or (y + r) = 1
Break
EndIf
x + 1
Wend
If x > y - 1 : x = -1 : EndIf
ProcedureReturn x
EndProcedure
- Output:
MODULAR-INVERSE( 42, 2017) = 1969 MODULAR-INVERSE( 40, 1) = 0 MODULAR-INVERSE( 52, -217) = 96 MODULAR-INVERSE( -486, 217) = 121 MODULAR-INVERSE( 40, 2018) = -1
Python
Builtin function
Since python3.8, builtin function "pow" can be used directly to compute modular inverses by specifying an exponent of -1:
>>> pow(42, -1, 2017)
1969
Iteration and error-handling
Implementation of this pseudocode with this.
>>> def extended_gcd(aa, bb):
lastremainder, remainder = abs(aa), abs(bb)
x, lastx, y, lasty = 0, 1, 1, 0
while remainder:
lastremainder, (quotient, remainder) = remainder, divmod(lastremainder, remainder)
x, lastx = lastx - quotient*x, x
y, lasty = lasty - quotient*y, y
return lastremainder, lastx * (-1 if aa < 0 else 1), lasty * (-1 if bb < 0 else 1)
>>> def modinv(a, m):
g, x, y = extended_gcd(a, m)
if g != 1:
raise ValueError
return x % m
>>> modinv(42, 2017)
1969
>>>
Recursion and an option type
Or, using functional composition as an alternative to iterative mutation, and wrapping the resulting value in an option type, to allow for the expression of computations which establish the absence of a modular inverse:
from functools import (reduce)
from itertools import (chain)
# modInv :: Int -> Int -> Maybe Int
def modInv(a):
return lambda m: (
lambda ig=gcdExt(a)(m): (
lambda i=ig[0]: (
Just(i + m if 0 > i else i) if 1 == ig[2] else (
Nothing()
)
)
)()
)()
# gcdExt :: Int -> Int -> (Int, Int, Int)
def gcdExt(x):
def go(a, b):
if 0 == b:
return (1, 0, a)
else:
(q, r) = divmod(a, b)
(s, t, g) = go(b, r)
return (t, s - q * t, g)
return lambda y: go(x, y)
# TEST ---------------------------------------------------
# Numbers between 2010 and 2015 which do yield modular inverses for 42:
# main :: IO ()
def main():
print (
mapMaybe(
lambda y: bindMay(modInv(42)(y))(
lambda mInv: Just((y, mInv))
)
)(
enumFromTo(2010)(2025)
)
)
# -> [(2011, 814), (2015, 48), (2017, 1969), (2021, 1203)]
# GENERIC ABSTRACTIONS ------------------------------------
# enumFromTo :: Int -> Int -> [Int]
def enumFromTo(m):
return lambda n: list(range(m, 1 + n))
# bindMay (>>=) :: Maybe a -> (a -> Maybe b) -> Maybe b
def bindMay(m):
return lambda mf: (
m if m.get('Nothing') else mf(m.get('Just'))
)
# Just :: a -> Maybe a
def Just(x):
return {'type': 'Maybe', 'Nothing': False, 'Just': x}
# mapMaybe :: (a -> Maybe b) -> [a] -> [b]
def mapMaybe(mf):
return lambda xs: reduce(
lambda a, x: maybe(a)(lambda j: a + [j])(mf(x)),
xs,
[]
)
# maybe :: b -> (a -> b) -> Maybe a -> b
def maybe(v):
return lambda f: lambda m: v if m.get('Nothing') else (
f(m.get('Just'))
)
# Nothing :: Maybe a
def Nothing():
return {'type': 'Maybe', 'Nothing': True}
# MAIN ---
main()
- Output:
[(2011, 814), (2015, 48), (2017, 1969), (2021, 1203)]
Quackery
[ dup 1 != if
[ tuck 1 0
[ swap temp put
temp put
over 1 > while
tuck /mod swap
temp take tuck *
temp take swap -
again ]
2drop
temp release
temp take
dup 0 < if
[ over + ] ]
nip ] is modinv ( n n --> n )
42 2017 modinv echo
- Output:
1969
Using Extended Euclidean Algorithm
Handles negative args. Returns -1 for non-coprime args.
[ dup 0 = iff
[ 2drop 1 0 ]
done
dup unrot /mod
dip swap recurse
tuck 2swap *
dip swap - ] is egcd ( n n --> n n )
[ dup 0 < if negate
over 0 < if
[ swap negate
over tuck mod
- swap ]
dup rot 2dup egcd
2swap 2over rot *
unrot * + 1 != iff
[ drop 2drop -1 ]
done
nip swap mod ] is modinv ( n n --> n )
say " 42 2017 modinv --> " 42 2017 modinv echo cr ( --> 1969 )
say " 40 1 modinv --> " 40 1 modinv echo cr ( --> 0 )
say " 52 -217 modinv --> " 52 -217 modinv echo cr ( --> 96 )
say "-486 217 modinv --> " -486 217 modinv echo cr ( --> 121 )
say " 40 2018 modinv --> " 40 2018 modinv echo cr ( --> -1 )
- Output:
42 2017 modinv --> 1969 40 1 modinv --> 0 52 -217 modinv --> 96 -486 217 modinv --> 121 40 2018 modinv --> -1
Racket
(require math)
(modular-inverse 42 2017)
- Output:
1969
Raku
(formerly Perl 6)
sub inverse($n, :$modulo) {
my ($c, $d, $uc, $vc, $ud, $vd) = ($n % $modulo, $modulo, 1, 0, 0, 1);
my $q;
while $c != 0 {
($q, $c, $d) = ($d div $c, $d % $c, $c);
($uc, $vc, $ud, $vd) = ($ud - $q*$uc, $vd - $q*$vc, $uc, $vc);
}
return $ud % $modulo;
}
say inverse 42, :modulo(2017)
REXX
/*REXX program calculates and displays the modular inverse of an integer X modulo Y.*/
parse arg x y . /*obtain two integers from the C.L. */
if x=='' | x=="," then x= 42 /*Not specified? Then use the default.*/
if y=='' | y=="," then y= 2017 /* " " " " " " */
say 'modular inverse of ' x " by " y ' ───► ' modInv(x,y)
exit /*stick a fork in it, we're all done. */
/*──────────────────────────────────────────────────────────────────────────────────────*/
modInv: parse arg a,b 1 ob; z= 0 /*B & OB are obtained from the 2nd arg.*/
$= 1 /*initialize modular inverse to unity. */
if b\=1 then do while a>1
parse value a/b a//b b z with q b a t
z= $ - q * z
$= trunc(t)
end /*while*/
if $<0 then $= $ + ob /*Negative? Then add the original B. */
return $
- output when using the default inputs of: 42 2017
modular inverse of 42 by 2017 ───► 1969
Ring
see "42 %! 2017 = " + multInv(42, 2017) + nl
func multInv a,b
b0 = b
x0 = 0
multInv = 1
if b = 1 return 0 ok
while a > 1
q = floor(a / b)
t = b
b = a % b
a = t
t = x0
x0 = multInv - q * x0
multInv = t
end
if multInv < 0 multInv = multInv + b0 ok
return multInv
Output:
42 %! 2017 = 1969
RPL
Using complex numbers allows to ‘parallelize’ calculations and keeps the stack depth low: never more than 4 levels despite the simultaneous use of 6 variables: r, r’, u, u’, q - and b for the final touch.
RPL code | Comment |
---|---|
≪ DUP ROT 1 R→C ROT 0 R→C WHILE DUP RE REPEAT OVER RE OVER RE / FLOOR OVER * NEG ROT + END DROP C→R ROT MOD SWAP 1 == SWAP 0 IFTE ≫ ‘MODINV’ STO |
MODINV ( a b -- x ) with ax = 1 mod b 3: b 2: (r,u)←(a,1) 1:(r',u')←(b,0) While r' ≠ 0 q ← r // r' (r - q*r', u - q*u') Forget (r',u') and calculate u mod b Test r and return zero if a and b are not co-prime |
- Input:
123 456 MODINV 42 2017 MODINV
- Output:
2: 1969 1: 0
Ruby
#based on pseudo code from http://en.wikipedia.org/wiki/Extended_Euclidean_algorithm#Iterative_method_2 and from translating the python implementation.
def extended_gcd(a, b)
last_remainder, remainder = a.abs, b.abs
x, last_x, y, last_y = 0, 1, 1, 0
while remainder != 0
last_remainder, (quotient, remainder) = remainder, last_remainder.divmod(remainder)
x, last_x = last_x - quotient*x, x
y, last_y = last_y - quotient*y, y
end
return last_remainder, last_x * (a < 0 ? -1 : 1)
end
def invmod(e, et)
g, x = extended_gcd(e, et)
if g != 1
raise 'The maths are broken!'
end
x % et
end
> invmod(42,2017) => 1969
Simplified equivalent implementation
def modinv(a, m) # compute a^-1 mod m if possible
raise "NO INVERSE - #{a} and #{m} not coprime" unless a.gcd(m) == 1
return m if m == 1
m0, inv, x0 = m, 1, 0
while a > 1
inv -= (a / m) * x0
a, m = m, a % m
inv, x0 = x0, inv
end
inv += m0 if inv < 0
inv
end
> modinv(42,2017) => 1969
The OpenSSL module has modular inverse built in:
require 'openssl'
p OpenSSL::BN.new(42).mod_inverse(2017).to_i
Run BASIC
print multInv(42, 2017)
end
function multInv(a,b)
b0 = b
multInv = 1
if b = 1 then goto [endFun]
while a > 1
q = a / b
t = b
b = a mod b
a = t
t = x0
x0 = multInv - q * x0
multInv = int(t)
wend
if multInv < 0 then multInv = multInv + b0
[endFun]
end function
- Output:
1969
Rust
fn mod_inv(a: isize, module: isize) -> isize {
let mut mn = (module, a);
let mut xy = (0, 1);
while mn.1 != 0 {
xy = (xy.1, xy.0 - (mn.0 / mn.1) * xy.1);
mn = (mn.1, mn.0 % mn.1);
}
while xy.0 < 0 {
xy.0 += module;
}
xy.0
}
fn main() {
println!("{}", mod_inv(42, 2017))
}
- Output:
1969
Alternative implementation
fn modinv(a0: isize, m0: isize) -> isize {
if m0 == 1 { return 1 }
let (mut a, mut m, mut x0, mut inv) = (a0, m0, 0, 1);
while a > 1 {
inv -= (a / m) * x0;
a = a % m;
std::mem::swap(&mut a, &mut m);
std::mem::swap(&mut x0, &mut inv);
}
if inv < 0 { inv += m0 }
inv
}
fn main() {
println!("{}", modinv(42, 2017))
}
- Output:
1969
Scala
Based on the Handbook of Applied Cryptography, Chapter 2. See http://cacr.uwaterloo.ca/hac/ .
def gcdExt(u: Int, v: Int): (Int, Int, Int) = {
@tailrec
def aux(a: Int, b: Int, x: Int, y: Int, x1: Int, x2: Int, y1: Int, y2: Int): (Int, Int, Int) = {
if(b == 0) (x, y, a) else {
val (q, r) = (a / b, a % b)
aux(b, r, x2 - q * x1, y2 - q * y1, x, x1, y, y1)
}
}
aux(u, v, 1, 0, 0, 1, 1, 0)
}
def modInv(a: Int, m: Int): Option[Int] = {
val (i, j, g) = gcdExt(a, m)
if (g == 1) Option(if (i < 0) i + m else i) else Option.empty
}
Translated from C++ (on this page)
def modInv(a: Int, m: Int, x:Int = 1, y:Int = 0) : Int = if (m == 0) x else modInv(m, a%m, y, x - y*(a/m))
- Output:
scala> modInv(2,4) res1: Option[Int] = None scala> modInv(42, 2017) res2: Option[Int] = Some(1976)
Seed7
The library bigint.s7i defines the bigInteger function modInverse. It returns the modular multiplicative inverse of a modulo b when a and b are coprime (gcd(a, b) = 1). If a and b are not coprime (gcd(a, b) <> 1) the exception RANGE_ERROR is raised.
const func bigInteger: modInverse (in var bigInteger: a,
in var bigInteger: b) is func
result
var bigInteger: modularInverse is 0_;
local
var bigInteger: b_bak is 0_;
var bigInteger: x is 0_;
var bigInteger: y is 1_;
var bigInteger: lastx is 1_;
var bigInteger: lasty is 0_;
var bigInteger: temp is 0_;
var bigInteger: quotient is 0_;
begin
if b < 0_ then
raise RANGE_ERROR;
end if;
if a < 0_ and b <> 0_ then
a := a mod b;
end if;
b_bak := b;
while b <> 0_ do
temp := b;
quotient := a div b;
b := a rem b;
a := temp;
temp := x;
x := lastx - quotient * x;
lastx := temp;
temp := y;
y := lasty - quotient * y;
lasty := temp;
end while;
if a = 1_ then
modularInverse := lastx;
if modularInverse < 0_ then
modularInverse +:= b_bak;
end if;
else
raise RANGE_ERROR;
end if;
end func;
Original source: [1]
Sidef
Built-in:
say 42.modinv(2017)
Algorithm implementation:
func invmod(a, n) {
var (t, nt, r, nr) = (0, 1, n, a % n)
while (nr != 0) {
var quot = int((r - (r % nr)) / nr);
(nt, t) = (t - quot*nt, nt);
(nr, r) = (r - quot*nr, nr);
}
r > 1 && return()
t < 0 && (t += n)
t
}
say invmod(42, 2017)
- Output:
1969
Swift
extension BinaryInteger {
@inlinable
public func modInv(_ mod: Self) -> Self {
var (m, n) = (mod, self)
var (x, y) = (Self(0), Self(1))
while n != 0 {
(x, y) = (y, x - (m / n) * y)
(m, n) = (n, m % n)
}
while x < 0 {
x += mod
}
return x
}
}
print(42.modInv(2017))
- Output:
1969
Tcl
proc gcdExt {a b} {
if {$b == 0} {
return [list 1 0 $a]
}
set q [expr {$a / $b}]
set r [expr {$a % $b}]
lassign [gcdExt $b $r] s t g
return [list $t [expr {$s - $q*$t}] $g]
}
proc modInv {a m} {
lassign [gcdExt $a $m] i -> g
if {$g != 1} {
return -code error "no inverse exists of $a %! $m"
}
while {$i < 0} {incr i $m}
return $i
}
Demonstrating
puts "42 %! 2017 = [modInv 42 2017]"
catch {
puts "2 %! 4 = [modInv 2 4]"
} msg; puts $msg
- Output:
42 %! 2017 = 1969 no inverse exists of 2 %! 4
Tiny BASIC
2017 causes integer overflow, so I'll do the inverse of 42 modulo 331 instead.
PRINT "Modular inverse."
PRINT "What is the modulus?"
INPUT M
PRINT "What number is to be inverted?"
INPUT X
PRINT "Solution is:"
10 LET A = A + 1
GOTO 20
15 IF B = 1 THEN PRINT A
IF B = 1 THEN END
IF A = M-1 THEN PRINT "nonexistent"
IF A = M-1 THEN END
GOTO 10
20 LET B = A*X
30 IF B < M THEN GOTO 15
LET B = B - M
GOTO 30
- Output:
Modular inverse. What is the modulus? 331 What number is to be inverted? 42 Solution is: 134
Another version:
REM Modular inverse
LET E=42
LET T=2017
GOSUB 100
PRINT M
END
REM Increments E S (step) times until B is greater than T
REM Repeats until B = 1 (MOD = 1) and C (count)
REM B will not be greater than T + E
100 LET D=0
IF E>=T THEN GOTO 130
LET B=E
REM At least one iteration is necessary
LET S=((T-B)/E)+1
LET B=B+S*E
LET C=1+S
LET B=B-T
110 IF B=1 THEN GOTO 120
LET S=((T-B)/E)+1
LET B=B+S*E
LET C=C+S
LET B=B-T
GOTO 110
120 LET D=C
130 LET M=D
RETURN
- Output:
1969
tsql
;WITH Iterate(N,A,B,X0,X1)
AS (
SELECT
1
,CASE WHEN @a < 0 THEN @b-(-@a % @b) ELSE @a END
,CASE WHEN @b < 0 THEN -@b ELSE @b END
,0
,1
UNION ALL
SELECT
N+1
,B
,A%B
,X1-((A/B)*X0)
,X0
FROM Iterate
WHERE A != 1 AND B != 0
),
ModularInverse(Result)
AS (
SELECT
-1
FROM Iterate
WHERE A != 1 AND B = 0
UNION ALL
SELECT
TOP(1)
CASE WHEN X1 < 0 THEN X1+@b ELSE X1 END AS Result
FROM Iterate
WHERE (SELECT COUNT(*) FROM Iterate WHERE A != 1 AND B = 0) = 0
ORDER BY N DESC
)
SELECT *
FROM ModularInverse
TypeScript
// Modular inverse
function modInv(e: number, t: number): number {
var d = 0;
if (e < t) {
var count = 1;
var bal = e;
do {
var step = Math.floor((t - bal) / e) + 1;
bal += step * e;
count += step;
bal -= t;
} while (bal != 1);
d = count;
}
return d;
}
console.log(`${modInv(42, 2017)}`); // 1969
- Output:
1969
uBasic/4tH
Print FUNC(_MulInv(42, 2017))
End
_MulInv Param(2)
Local(5)
c@ = b@
f@ = 0
g@ = 1
If b@ = 1 Then Return
Do While a@ > 1
e@ = a@ / b@
d@ = b@
b@ = a@ % b@
a@ = d@
d@ = f@
f@ = g@ - e@ * f@
g@ = d@
Loop
If g@ < 0 Then g@ = g@ + c@
Return (g@)
Print FUNC(_mul_inv(42, 2017))
Print FUNC(_mul_inv(40, 1))
Print FUNC(_mul_inv(52, -217))
Print FUNC(_mul_inv(-486, 217))
Print FUNC(_mul_inv(40, 2018))
End
_mul_inv Param(2)
Local(6)
If (b@ < 0) b@ = -b@
If (a@ < 0) a@ = b@ - (-a@ % b@)
c@ = 0 : d@ = 1 : e@ = b@ : f@ = a@ % b@
Do Until (f@ = 0)
g@ = e@/f@
h@ = d@ : d@ = c@ - g@*d@ : c@ = h@
h@ = f@ : f@ = e@ - g@*f@ : e@ = h@
Loop
If (e@ > 1) Return (-1) ' No inverse'
If (c@ < 0) c@ = c@ + b@
Return (c@)
- Output:
1969 0 96 121 -1 0 OK, 0:156
UNIX Shell
function invmod {
typeset -i a=$1 n=$2
if (( n < 0 )); then (( n = -n )); fi
if (( a < 0 )); then (( a = n - (-a) % n )); fi
typeset -i t=0 nt=1 r=n nr q tmp
(( nr = a % n ))
while (( nr )); do
(( q = r/nr ))
(( tmp = nt ))
(( nt = t - q*nt ))
(( t = tmp ))
(( tmp = nr ))
(( nr = r - q*nr ))
(( r = tmp ))
done
if (( r > 1 )); then
return 1
fi
while (( t < 0 )); do (( t += n )); done
printf '%s\n' "$t"
}
invmod 42 2017
- Output:
1969
VBA
Private Function mul_inv(a As Long, n As Long) As Variant
If n < 0 Then n = -n
If a < 0 Then a = n - ((-a) Mod n)
Dim t As Long: t = 0
Dim nt As Long: nt = 1
Dim r As Long: r = n
Dim nr As Long: nr = a
Dim q As Long
Do While nr <> 0
q = r \ nr
tmp = t
t = nt
nt = tmp - q * nt
tmp = r
r = nr
nr = tmp - q * nr
Loop
If r > 1 Then
mul_inv = "a is not invertible"
Else
If t < 0 Then t = t + n
mul_inv = t
End If
End Function
Public Sub mi()
Debug.Print mul_inv(42, 2017)
Debug.Print mul_inv(40, 1)
Debug.Print mul_inv(52, -217) '/* Pari semantics for negative modulus */
Debug.Print mul_inv(-486, 217)
Debug.Print mul_inv(40, 2018)
End Sub
- Output:
1969 0 96 121 a is not invertible
V (Vlang)
fn main() {
println("42 %! 2017 = ${mult_inv(42, 2017)}")
}
fn mult_inv(aa int, bb int) int {
mut a, mut b := aa, bb
mut x0, mut t := 0, 0
mut b0 := b
mut x1 := 1
if b == 1 {return 1}
for a > 1 {
q := a / b
t = b
b = a % b
a = t
t = x0
x0 = x1 - q * x0
x1 = t
}
if x1 < 0 {x1 += b0}
return x1
}
- Output:
42 %! 2017 = 1969
Wren
import "./big" for BigInt
var a = BigInt.new(42)
var b = BigInt.new(2017)
System.print(a.modInv(b))
- Output:
1969
XPL0
code IntOut=11, Text=12;
int X;
def A=42, M=2017;
[for X:= 2 to M-1 do
if rem(A*X/M) = 1 then [IntOut(0, X); exit];
Text(0, "Does not exist");
]
- Output:
1969
zkl
fcn gcdExt(a,b){
if(b==0) return(1,0,a);
q,r:=a.divr(b); s,t,g:=gcdExt(b,r); return(t,s-q*t,g);
}
fcn modInv(a,m){i,_,g:=gcdExt(a,m); if(g==1) {if(i<0)i+m} else Void}
divr(a,b) is [integer] (a/b,remainder)
- Output:
modInv(2,4) //-->Void modInv(42,2017) //-->1969
- Programming Tasks
- Solutions by Programming Task
- 11l
- 8th
- Action!
- Ada
- ALGOL 68
- Arturo
- ATS
- AutoHotkey
- AWK
- BASIC
- ASIC
- BASIC256
- Chipmunk Basic
- Minimal BASIC
- QBasic
- True BASIC
- Batch File
- BCPL
- Bracmat
- C
- C sharp
- C++
- Clojure
- CLU
- Comal
- Common Lisp
- Cowgol
- Craft Basic
- Crystal
- D
- Dc
- Delphi
- Draco
- EasyLang
- EchoLisp
- Elixir
- ERRE
- F Sharp
- Factor
- Forth
- Fortran
- FreeBASIC
- Frink
- FunL
- Go
- GW-BASIC
- Haskell
- Icon
- Unicon
- IS-BASIC
- J
- Java
- JavaScript
- Jq
- Julia
- Kotlin
- Lambdatalk
- M4
- MAD
- Maple
- Mathematica
- Wolfram Language
- Maxima
- Mercury
- МК-61/52
- Modula-2
- NewLISP
- Nim
- Oberon-2
- ObjectIcon
- OCaml
- Oforth
- Owl Lisp
- PARI/GP
- Pascal
- Perl
- Phix
- PHP
- PicoLisp
- PL/I
- PowerShell
- Prolog
- PureBasic
- Python
- Quackery
- Racket
- Raku
- REXX
- Ring
- RPL
- Ruby
- Run BASIC
- Rust
- Scala
- Seed7
- Sidef
- Swift
- Tcl
- Tiny BASIC
- Tsql
- TypeScript
- UBasic/4tH
- UNIX Shell
- VBA
- V (Vlang)
- Wren
- Wren-big
- XPL0
- Zkl