diff --git a/src/__tests__/index.spec.ts b/src/__tests__/index.spec.ts index f34bd7e..6670681 100644 --- a/src/__tests__/index.spec.ts +++ b/src/__tests__/index.spec.ts @@ -22,6 +22,8 @@ import { monzoToCents, tenneyHeight, wilsonHeight, + modInv, + mmod, } from '../index'; const FUZZ = 'FUZZ' in process.env; @@ -621,3 +623,59 @@ describe('Wilson complexity measure', () => { } }); }); + +describe('Modular inverse calculator', () => { + it('finds modular inverses when they exist', () => { + for (let a = -30; a < 30; ++a) { + for (let b = 2; b < 20; ++b) { + const c = Math.abs(gcd(a, b)); + if (c !== 1) { + expect(() => modInv(a, b)).toThrow(); + } else { + expect(mmod(a * modInv(a, b), b)).toBe(1); + } + } + } + }); + + it("finds the next best alternative even if a true modular inverse doesn't exist", () => { + for (let a = -30; a < 30; ++a) { + for (let b = 1; b < 20; ++b) { + let c = Math.abs(gcd(a, b)); + if (c === b) { + c = 0; + } + expect(mmod(a * modInv(a, b, false), b)).toBe(c); + } + } + }); + + it('can be used to break 15edo into 3 interleaved circles of fifths', () => { + const edo = 15; + const gen = 9; + const genInv = modInv(gen, edo, false); + const numCycles = mmod(gen * genInv, edo); + function f(step: number) { + const c = mmod(step, numCycles); + return mmod((step - c) * genInv, edo) + c; + } + expect([...Array(16).keys()].map(f)).toEqual([ + 0, // Root + 1, // Root + 1 + 2, // Root + 2 + 6, // Gen * 2 + 7, // +1 + 8, // +2 + 12, // Gen * 4 + 13, // +1 + 14, // +2 + 3, // Gen + 4, // Gen + 1 + 5, // Gen + 1 + 9, // Gen * 3 + 10, // +1 + 11, // +2 + 0, // Circle closes + ]); + }); +}); diff --git a/src/index.ts b/src/index.ts index a633706..df58177 100644 --- a/src/index.ts +++ b/src/index.ts @@ -128,6 +128,23 @@ export function iteratedEuclid(params: Iterable) { return coefs; } +/** + * Find modular inverse of a (mod b). + * @param a Number to find modular inverse of. + * @param b Modulus. + * @param strict Ensure that a * modInv(a, b) = 1 (mod b). If `strict = false` we have a * modInv(a, b) = gdc(a, b) (mod b) instead. + * @returns The modular inverse in the range {0, 1, ..., b - 1}. + */ +export function modInv(a: number, b: number, strict = true) { + const {gcd, coefA} = extendedEuclid(a, b); + if (strict && gcd !== 1) { + throw new Error( + `${a} does not have a modular inverse modulo ${b} since they're not coprime` + ); + } + return mmod(coefA, b); +} + /** * Collection of unique fractions. */