Modular multiplicative inverse

23

1

Your task is to given two integer numbers, a and b calculate the modular multiplicative inverse of a modulo b, if it exists.

The modular inverse of a modulo b is a number c such that ac ≡ 1 (mod b). This number is unique modulo b for any pair of a and b. It exists only if the greatest common divisor of a and b is 1.

The Wikipedia page for modular multiplicative inverse can be consulted if you require more information about the topic.

Input and Output

Input is given as either two integers or a list of two integers. Your program should output either a single number, the modular multiplicative inverse that is in the interval 0 < c < b, or a value indicating there is no inverse. The value can be anything, except a number in the range (0,b), and may also be an exception. The value should however be the same for cases in which there is no inverse.

0 < a < b can be assumed

Rules

  • The program should finish at some point, and should solve each test case in less than 60 seconds
  • Standard loopholes apply

Test cases

Test cases below are given in the format, a, b -> output

1, 2 -> 1
3, 6 -> Does not exist
7, 87 -> 25
25, 87 -> 7
2, 91 -> 46
13, 91 -> Does not exist
19, 1212393831 -> 701912218
31, 73714876143 -> 45180085378
3, 73714876143 -> Does not exist

Scoring

This is code golf, so the shortest code for each language wins.

This and this are similar questions, but both ask for specific situations.

Halvard Hummel

Posted 2017-08-24T15:41:12.173

Reputation: 3 131

Sandbox – Halvard Hummel – 2017-08-24T15:41:36.377

Can our output be empty in case there is no such value? – Mr. Xcoder – 2017-08-24T15:51:09.367

@Mr.Xcoder Yes, but in that case it would need to be empty for all cases in which there is no such value – Halvard Hummel – 2017-08-24T15:53:50.570

6It follows from Fermat's Little Theorem that the multiplicative inverse of a, if it exists, can be computed efficiently as a^(phi(b)-1) mod b, where phi is Euler's totient function: phi(p0^k0 * p1^k1 * ...) = (p0-1) * p0^(k0-1) * (p1-1) * p1^(k1-1) * ... Not saying it leads to shorter code :) – ngn – 2017-08-24T16:08:29.737

Can we output a optional type? I'm assuming we can but I just wanted to check.

– Post Rock Garf Hunter – 2017-08-24T16:22:08.873

can I take as input [a,-1,b] ? – J42161217 – 2017-08-24T16:31:38.363

1@Jenny_mathy Taking additional input is generally disallowed. – Mr. Xcoder – 2017-08-24T16:32:49.707

@Mr.Xcoder yes but I would beat mathematica's built in which is generally great! – J42161217 – 2017-08-24T16:34:35.047

5I count six answers that seem to be brute forcing, and unlikely to run all test cases in 60 seconds (some of them give a stack or memory error first). – Ørjan Johansen – 2017-08-24T23:28:40.677

1@ngn : You've conflated Fermat's Little Theorem (FLT) with Euler's improvement to it. Fermat did not know about the Euler phi function. Further, FLT and Euler's improvement only apply if gcd(a,b) = 1. Finally, in the form you have written it, "a^(\phi(b)-1) mod b" is congruent to 1, not a^(-1). To get a^(-1), use a^(\phi(b)-2) mod b. – Eric Towers – 2017-08-25T04:49:59.090

1@EricTowers Euler's is a consequence. Regarding "gcd(a,b)=1" - I did say "if it [the inverse] exists". Are you sure about phi(b)-2? – ngn – 2017-08-25T06:45:29.833

@ngn : Hmm... Seems you got me turned around. It's a^(b -2) = a^(-1) mod p in Fermat's theorem and a^(phi(b) -1) = a^(-1) mod p in Euler's. – Eric Towers – 2017-08-25T12:56:57.097

@EricTowers Looks correct now, assuming by "p" you mean "b=p, where p is prime". I admit it would have been more proper, though just as true, to quote Euler's theorem (which is the generalisation to multiple prime factors) instead of Fermat's Little. – ngn – 2017-08-25T13:19:37.337

60 seconds on which machine? – totallyhuman – 2017-08-25T21:10:21.740

@totallyhuman Reasonable normal machine – Halvard Hummel – 2017-08-25T21:12:28.633

So at least @Mr.Xcoder claims that his brute force solution (in Python) does run in less than 60 seconds on a normal machine - barely. Which is exactly what you'd want to avoid in choosing a non-scoring time limit. – Ørjan Johansen – 2017-08-26T17:14:19.013

Why limit the algorithm for just positive results and arguments, when during the computation one see values can be negative, the same the result,and that algorithm it is for all integers (negative positive zero)? – RosLuP – 2017-08-29T09:01:01.193

Answers

11

Mathematica, 14 bytes

Obligatory Mathematica builtin:

ModularInverse

It's a function that takes two arguments (a and b), and returns the inverse of a mod b if it exists. If not, it returns the error ModularInverse: a is not invertible modulo b..

Tutleman

Posted 2017-08-24T15:41:12.173

Reputation: 571

7

JavaScript (ES6), 79 73 62 61 bytes

Returns false if the inverse does not exist.

It uses the extended Euclidean algorithm and solves all test cases almost instantly.

f=(a,b,c=!(n=b),d=1)=>a?f(b%a,a,d,c-(b-b%a)/a*d):b<2&&(c+n)%n

Test cases

f=(a,b,c=!(n=b),d=1)=>a?f(b%a,a,d,c-(b-b%a)/a*d):b<2&&(c+n)%n

console.log(f(1, 2)) // -> 1
console.log(f(3, 6)) // -> Does not exist
console.log(f(7, 87)) // -> 25
console.log(f(25, 87)) // -> 7
console.log(f(2, 91)) // -> 46
console.log(f(13, 91)) // -> Does not exist
console.log(f(19, 1212393831)) // -> 701912218
console.log(f(31, 73714876143)) // -> 45180085378
console.log(f(3, 73714876143)) // -> Does not exist

Arnauld

Posted 2017-08-24T15:41:12.173

Reputation: 111 334

Why is it not possible to write the name of function f, as in f(c,a,b=0,d=1,n=a)=>c?f(a%c,c,d,b-(a-a%c)/c*d,n):a<2&&(b+n)%n ? – RosLuP – 2017-08-24T17:47:06.917

@RosLup f(x,y) is always parsed as a function call, except if it is explicitly preceded by the function keyword. An anonymous arrow function, on the other hand, is declared as (x,y)=>something and f=(x,y)=>something assigns the function to the f variable. – Arnauld – 2017-08-24T18:01:11.523

5

Python 2, 34 bytes

f=lambda a,b:a==1or-~b*f(-b%a,a)/a

Try it online!

Recursive function that gives True for print f(1,2), which I believe to be acceptable, and errors for invalid inputs.

We are trying to find \$x\$ in \$a\cdot x\equiv 1\pmod{b}\$.

This can be written as \$a\cdot x-1=k\cdot b\$ where \$k\$ is an integer.

Taking \$\mod{a}\$ of this gives \$-1\equiv k\cdot b\pmod{a}\$. Moving the minus gives \$-k\cdot b\equiv1\pmod{a}\$, where we have to solve for \$k\$.

Seeing how it resembles the initial scenario, allow us to recurse to solve for \$k\$ by calling the function with \$f(-b\%a,a)\$ (works because Python gives positive values for modulo with a negative argument).

The program recurses for until \$a\$ becomes 1, which only happens if the original \$a\$ and \$b\$ are coprime to each other (ie there exists a multiplicative inverse), or ends in an error caused by division by 0.

This value of \$k\$ can be substituted in the equation \$a\cdot x-1=k\cdot b\$ to give \$x\$ as \$\frac{k\cdot b+1}{a}\$.

user41805

Posted 2017-08-24T15:41:12.173

Reputation: 16 320

4

Jelly, 2 bytes

æi

Try it online!

This uses a builtin for modular inverse, and returns 0 for no modular inverse.

Jelly, 7 bytes

R×%⁸’¬T

Try it online!

Outputs empty set (represented as empty string) on no modular inverse. Runs out of memory on TIO for the largest test-cases, but should work given enough memory.

How it Works

R×%⁸’¬T  
R        Generate range of b
 ×       Multiply each by a
  %⁸     Mod each by b
    ’    Decrement (Map 1 to 0 and all else to truthy)
     ¬   Logical NOT
      T  Get the index of the truthy element.

If you want to work for larger test-cases, try this (relatively ungolfed) version, which requires much time rather than memory:

Jelly, 9 bytes

×⁴%³’¬ø1#

Try it online!

How it Works

×⁴%³’¬ø1#
        #   Get the first
      ø1      one integer
            which meets:
×⁴            When multiplied by a
  %³          And modulo-d by b
    ’         Decrement
     ¬        Is falsy

fireflame241

Posted 2017-08-24T15:41:12.173

Reputation: 7 021

4

Mathematica, 18 bytes

PowerMod[#,-1,#2]&

input

[31, 73714876143]

J42161217

Posted 2017-08-24T15:41:12.173

Reputation: 15 931

3

R + numbers, 15 bytes

numbers::modinv

returns NA for those a without inverses mod b.

R-Fiddle to try it!

R, 33 bytes (non-competing)

This will fail on very large b since it actually creates a vector of size 32*b bits.

function(a,b)which((1:b*a)%%b==1)

Try it online!

Returns integer(0) (an empty list) for those a without inverses mod b.

Giuseppe

Posted 2017-08-24T15:41:12.173

Reputation: 21 077

3

Japt, 9 8 bytes

Takes the inputs in reverse order. Outputs -1 for no match. Craps out as the bigger integer gets larger.

Ç*V%UÃb1

Test it

  • Saved 1 byte thanks to ETH pointing out an errant, and very obvious, space.

Shaggy

Posted 2017-08-24T15:41:12.173

Reputation: 24 623

The test input 73714876143,31 seems to produce an out-of-memory error on Firefox (and to crash Chromium). I don't think this is a valid answer. – Ilmari Karonen – 2017-08-25T00:54:47.077

@IlmariKaronen: I clearly pointed out that fact in my solution. We can assume infinite memory for the purposes of code golf so the memory issues and crashes do not invalidate this solution. – Shaggy – 2017-08-25T08:51:55.840

3Unfortunately the memory issues also make it impossible to tell whether your code would actually solve the test cases in 60 seconds as stipulated by the challenge. I suspect it would not, even if there was sufficient memory available to make it not crash, but without a computer that can actually run the program for that long there's no way to tell for sure. – Ilmari Karonen – 2017-08-25T10:09:09.580

2

Python 2, 51 49 54 53 51 49 bytes

-1 byte thanks to officialaimm
-1 byte thanks to Shaggy

a,b=input()
i=a<2
while(a*i%b-1)*b%a:i+=1
print+i

Try it online!

Prints 0 when there is no solution.

Rod

Posted 2017-08-24T15:41:12.173

Reputation: 17 588

1Outputs 0 for a=1 and b=2; from the test cases, it should output 1. – Shaggy – 2017-08-24T16:03:36.913

As a recursive algo – Mr. Xcoder – 2017-08-24T16:03:43.827

1As Shaggy pointed out, fails for 2, 1 – Mr. Xcoder – 2017-08-24T16:06:10.683

@Shaggy should work now – Rod – 2017-08-24T16:14:59.770

1This fails to return an answer in 60 seconds (on TIO) for the input 31,73714876143. – Ilmari Karonen – 2017-08-25T00:49:37.770

1This goes into infinite loop for the input 4, 6, since there is no inverse, but b%a != 0. (For inverse to exist, it's not enough that b is not divisible by a, you need them to be coprime.) – Arthur – 2017-08-25T07:41:38.047

2

Pyth, 10 bytes

3 bytes saved thanks to @Jakube.

xm%*szdQQ1

Try it here!

Returns -1 for no multiplicative inverse.

Code Breakdown

xm%*szdQQ1      Let Q be the first input.
 m      Q       This maps over [0 ... Q) with a variable d.
   *szd         Now d is multiplied by the evaluated second input.
  %    Q        Now the remained modulo Q is retrieved.
x        1      Then, the first index of 1 is retrieved from that mapping.

Pyth, 15 13 bytes

KEhfq1%*QTKSK

Throws an exception in case no multiplicative inverse exists.

Try it here!

Pyth, 15 bytes

Iq1iQKEfq1%*QTK

This adds lots of bytes for handling the case where no such number exists. The program can be shortened significantly if that case would not need to be handled:

fq1%*QTK

Try it here!

Mr. Xcoder

Posted 2017-08-24T15:41:12.173

Reputation: 39 774

2 bytes saved with KExm%*QdKK1 – Jakube – 2017-08-25T06:51:07.423

Or 3 bytes if you swap the order of inputs: xm%*szdQQ1 – Jakube – 2017-08-25T06:53:03.790

@Jakube Thanks a lot, editing! – Mr. Xcoder – 2017-08-25T07:26:12.650

How does this work? – user41805 – 2018-11-20T20:04:41.893

@Cowsquack I've added a completely primitive code breakdown but rn I don't have time to include a complete explanation. Hopefully it is clear enough for now but I'll try to add a more complete explanation soon. – Mr. Xcoder – 2018-11-21T14:42:21.250

@Mr.Xcoder no rush, I wondered how the algorithm worked, thanks – user41805 – 2018-11-21T15:00:04.407

2

Python 3 + gmpy, 23 bytes

I don't think it can get any shorter in Python.

gmpy.invert
import gmpy

Try it online! (won't work if you do not have gmpy installed)

Mr. Xcoder

Posted 2017-08-24T15:41:12.173

Reputation: 39 774

2

Python 3, 49 bytes

lambda a,b:[c for c in range(b)if-~c*a%b==1][0]+1

Try it online!

Python 3, 50 bytes

lambda a,b:[c for c in range(1,b+1)if c*a%b==1][0]

Try it online!

This throws IndexError: list index out of range in case there is no modular multiplicative inverse, as it is allowed by the rules.

Mr. Xcoder

Posted 2017-08-24T15:41:12.173

Reputation: 39 774

This fails to return a result for the input 31,73714876143 in 60 seconds (on TIO). – Ilmari Karonen – 2017-08-25T00:58:07.103

@IlmariKaronen Seems to finish in 56 seconds on my machine (Macbook Pro '15) – Mr. Xcoder – 2017-08-25T04:58:15.667

2

8th, 6 bytes

Code

invmod

Explanation

invmod is a 8th word that calculates the value of the inverse of a, modulo b. It returns null on overflow or other errors.

Usage and test cases

ok> 1 2 invmod .
1
ok> 3 6 invmod .
null
ok> 7 87 invmod .
25
ok> 25 87 invmod .
7
ok> 2 91 invmod .
46
ok> 13 91 invmod .
null
ok> 19 1212393831 invmod .
701912218
ok> 31 73714876143 invmod .
45180085378
ok> 3 73714876143 invmod .
null

Chaos Manor

Posted 2017-08-24T15:41:12.173

Reputation: 521

2

Pari/GP, 22 bytes

a->b->lift(1/Mod(a,b))

Throws an error when there is no inverse.

Try it online!

alephalpha

Posted 2017-08-24T15:41:12.173

Reputation: 23 988

2

J, 28 bytes

4 :'(1=x+.y)*x y&|@^<:5 p:y'

Try it online!

Uses Euler's theorem. Returns 0 if the inverse does not exist.

Explanation

4 :'(1=x+.y)*x y&|@^<:5 p:y'  Input: a (LHS), b (RHS)
4 :'                       '  Define an explicit dyad - this is to use the special
                              form `m&|@^` to perform modular exponentiation
                          y   Get b
                      5 p:    Euler totient
                    <:        Decrement
             x                Get a
                   ^          Exponentiate
               y&|@             Modulo b
       x+.y                   GCD of a and b
     1=                       Equals 1
            *                 Multiply

miles

Posted 2017-08-24T15:41:12.173

Reputation: 15 654

1

C (gcc), 48 110 104 bytes

#define f(a,b)g(a,b,!b,1,b)
long g(a,b,c,d,n)long a,b,c,d,n;{a=a?g(b%a,a,d,c-(b-b%a)/a*d):!--b*(c+n)%n;}

Try it online!

This should work with all inputs (that fit within a long) within 60 seconds.

Edit. I'm already abusing the n variable so I might as well assume that gcc puts the first assignment in %rax.

ceilingcat

Posted 2017-08-24T15:41:12.173

Reputation: 5 503

2Alas, this gives wrong results even for fairly small inputs due to integer overflow inside the loop. For example, f(3,1000001) returns 717, which is obviously nonsense (the correct answer is 333334). Also, even if this bug was fixed by using a wider integer type, this brute-force approach would certainly time out for some of the larger test cases given in the challenge. – Ilmari Karonen – 2017-08-25T01:20:15.813

1

C (gcc), 115 bytes

#define L long long
L g(L a,L b,L c,L d){return a?g(b%a,a,d-b/a*c,c):b-1?0:d;}L f(L a,L b){return(g(a,b,1,0)+b)%b;}

Try it online!

Extended Euclidean algorithm, recursive version

C (gcc), 119 bytes

long long f(a,b,c,d,t,n)long long a,b,c,d,t,n;{for(c=1,d=0,n=b;a;a=t)t=d-b/a*c,d=c,c=t,t=b%a,b=a;return b-1?0:(d+n)%n;}

Try it online!

Extended Euclidean algorithm, iterative version

nwellnhof

Posted 2017-08-24T15:41:12.173

Reputation: 10 037

0

Python 2 + sympy, 74 bytes

from sympy import*
def f(a,m):i,_,g=numbers.igcdex(a,m);return g==1and i%m

Try it online!

Taken from Jelly source code.

totallyhuman

Posted 2017-08-24T15:41:12.173

Reputation: 15 378

0

JavaScript (ES6), 42 41 39 38 bytes

Outputs false for no match. Will throw a overflow error as the second number gets to be too large.

x=>y=>(g=z=>x*z%y==1?z:z<y&&g(++z))(1)

Shaggy

Posted 2017-08-24T15:41:12.173

Reputation: 24 623

0

Python 2, 52 bytes

-3 bytes thanks to Mr. Xcoder.

f=lambda a,b,i=1:i*a%b==1and i or i<b and f(a,b,i+1)

Try it online!

Outputs False on no solution and errors out as b gets larger.

Embedded TIO

I'm just testing out iframes in Stack Snippets and they work absolutely fantastic.

html,body{height:100%;}iframe{height:100%;width:100%;border:none;}
<iframe src="https://tio.run/##PY3BCsIwEETv/Yq5CEndy6ZotbhfIh4SanBB09Lm4tdHU8HTvIHhzfzOjym5UqI8/SuMHp4CqfCgrd8FEfZphGJaoJeAWqLZJnu2Jd/XvEJwNUxwlmA6wrFmTzj1FdzhT4QzV@DuR7emidULTdhMA@ZFU/4@tGrLBw"></iframe>

totallyhuman

Posted 2017-08-24T15:41:12.173

Reputation: 15 378

I'm not certain this works, can't i*a%b be 0? – Post Rock Garf Hunter – 2017-08-24T16:29:21.340

1Fails with "maximum recursion depth exceeded" error for input (31,73714876143). – Ilmari Karonen – 2017-08-25T01:09:34.713

0

Axiom, 45 bytes

f(x:PI,y:PI):NNI==(gcd(x,y)=1=>invmod(x,y);0)

0 for error else return z with x*z Mod y =1

RosLuP

Posted 2017-08-24T15:41:12.173

Reputation: 3 036

0

Jelly, 27 bytes

²%³
⁴Ç⁹Сx⁸
ÆṪ’BṚçL$P%³×gỊ¥

Try it online!

Uses Euler's theorem with modular exponentiation. Since Jelly doesn't have a builtin to perform modular exponentiation, it had to be implemented, and it took most of the bytes.

miles

Posted 2017-08-24T15:41:12.173

Reputation: 15 654

0

Axiom, 99 bytes

w(a,b,x,u)==(a=0=>(b*b=1=>b*x;0);w(b rem a,a,u,x-u*(b quo a)));h(a,b)==(b=0=>0;(b+w(a,b,0,1))rem b)

it use the function h(); h(a,b) return 0 if error [not exist inverse] otherwise it return the z such that a*z mod b = 1 This would be ok even if arguments are negative...

this would be the general egcd() function that retunr a list of int (so they can be negative too)

egcd(aa:INT,bb:INT):List INT==
   x:=u:=-1   -- because the type is INT
   (a,b,x,u):=(aa,bb,0,1)
   repeat
      a=0=>break
      (q,r):=(b quo a, b rem a)
      (b,a,x,u):=(a,r,u,x-u*q)
   [b,x, (b-x*aa)quo bb]

this is how use it

(7) -> h(31,73714876143)
   (7)  45180085378
                                                    Type: PositiveInteger

i find the base algo in internet from https://pastebin.com/A13ybryc

RosLuP

Posted 2017-08-24T15:41:12.173

Reputation: 3 036

0

APL (Dyalog Unicode), 36 38 bytes

{⌈/i×1=⍵|(i←⍳⍵)×⍵|⍺}

Try it online!

Explanation:

                   ⍵|⍺} ⍝ Get ⍺ mod ⍵
             (i←⍳⍵)×     ⍝ Multiply the result by all numbers up to ⍵
          ⍵|            ⍝ Take result mod ⍵
     i×1=                ⍝ Find all numbers (1,⍵) where the mod is 1
{⌈/                      ⍝ And take the largest

Much thanks to Adam in the APL Orchard chatroom for the help with this one!

Formula obtained from this site

First iteration:

{((⍵|⍺),⍵){+/⌈/⍵×1=(¯1↑⍺)|⍵×⊃⍺}⍳⍵}

JPeroutek

Posted 2017-08-24T15:41:12.173

Reputation: 734

1what about larger tests? the challenge says: "The program should finish at some point, and should solve each test case in less than 60 seconds" – ngn – 2020-02-17T16:15:18.160

1we usually count 1 char = 1 byte in apl – ngn – 2020-02-17T16:15:50.997

Asked: 2017-08-24T15:41:12.173

Viewed: 2 564 times

Active: 2020-02-18T16:10:02.307