2021-03-24 13:04:30 +00:00
|
|
|
/**
|
|
|
|
* Absolute value. abs(a)==a if a>=0. abs(a)==-a if a<0
|
|
|
|
*
|
|
|
|
* @param a
|
|
|
|
*
|
|
|
|
* @returns The absolute value of a
|
|
|
|
*/
|
|
|
|
function abs(a) {
|
|
|
|
return (a >= 0) ? a : -a;
|
|
|
|
}
|
|
|
|
|
|
|
|
/**
|
|
|
|
* Returns the bitlength of a number
|
|
|
|
*
|
|
|
|
* @param a
|
|
|
|
* @returns The bit length
|
|
|
|
*/
|
|
|
|
function bitLength(a) {
|
2021-03-24 14:11:07 +00:00
|
|
|
if (typeof a === 'number')
|
|
|
|
a = BigInt(a);
|
2021-03-24 13:04:30 +00:00
|
|
|
if (a === 1n) {
|
|
|
|
return 1;
|
|
|
|
}
|
|
|
|
let bits = 1;
|
|
|
|
do {
|
|
|
|
bits++;
|
|
|
|
} while ((a >>= 1n) > 1n);
|
|
|
|
return bits;
|
|
|
|
}
|
|
|
|
|
|
|
|
/**
|
|
|
|
* An iterative implementation of the extended euclidean algorithm or extended greatest common divisor algorithm.
|
|
|
|
* Take positive integers a, b as input, and return a triple (g, x, y), such that ax + by = g = gcd(a, b).
|
|
|
|
*
|
|
|
|
* @param a
|
|
|
|
* @param b
|
|
|
|
*
|
2022-09-12 08:59:51 +00:00
|
|
|
* @throws {@link RangeError} if a or b are <= 0
|
2021-03-24 17:30:45 +00:00
|
|
|
*
|
2021-03-24 13:04:30 +00:00
|
|
|
* @returns A triple (g, x, y), such that ax + by = g = gcd(a, b).
|
|
|
|
*/
|
|
|
|
function eGcd(a, b) {
|
2021-03-24 14:11:07 +00:00
|
|
|
if (typeof a === 'number')
|
|
|
|
a = BigInt(a);
|
|
|
|
if (typeof b === 'number')
|
|
|
|
b = BigInt(b);
|
|
|
|
if (a <= 0n || b <= 0n)
|
2021-03-24 13:04:30 +00:00
|
|
|
throw new RangeError('a and b MUST be > 0'); // a and b MUST be positive
|
|
|
|
let x = 0n;
|
|
|
|
let y = 1n;
|
|
|
|
let u = 1n;
|
|
|
|
let v = 0n;
|
2021-03-24 14:11:07 +00:00
|
|
|
while (a !== 0n) {
|
|
|
|
const q = b / a;
|
|
|
|
const r = b % a;
|
2021-03-24 13:04:30 +00:00
|
|
|
const m = x - (u * q);
|
|
|
|
const n = y - (v * q);
|
2021-03-24 14:11:07 +00:00
|
|
|
b = a;
|
|
|
|
a = r;
|
2021-03-24 13:04:30 +00:00
|
|
|
x = u;
|
|
|
|
y = v;
|
|
|
|
u = m;
|
|
|
|
v = n;
|
|
|
|
}
|
|
|
|
return {
|
2021-03-24 14:11:07 +00:00
|
|
|
g: b,
|
2021-03-24 13:04:30 +00:00
|
|
|
x: x,
|
|
|
|
y: y
|
|
|
|
};
|
|
|
|
}
|
|
|
|
|
|
|
|
/**
|
2022-09-12 08:59:51 +00:00
|
|
|
* Greatest common divisor of two integers based on the iterative binary algorithm.
|
2021-03-24 13:04:30 +00:00
|
|
|
*
|
|
|
|
* @param a
|
|
|
|
* @param b
|
|
|
|
*
|
|
|
|
* @returns The greatest common divisor of a and b
|
|
|
|
*/
|
|
|
|
function gcd(a, b) {
|
2021-03-24 14:11:07 +00:00
|
|
|
let aAbs = (typeof a === 'number') ? BigInt(abs(a)) : abs(a);
|
|
|
|
let bAbs = (typeof b === 'number') ? BigInt(abs(b)) : abs(b);
|
2021-03-24 13:04:30 +00:00
|
|
|
if (aAbs === 0n) {
|
|
|
|
return bAbs;
|
|
|
|
}
|
|
|
|
else if (bAbs === 0n) {
|
|
|
|
return aAbs;
|
|
|
|
}
|
|
|
|
let shift = 0n;
|
|
|
|
while (((aAbs | bAbs) & 1n) === 0n) {
|
|
|
|
aAbs >>= 1n;
|
|
|
|
bAbs >>= 1n;
|
|
|
|
shift++;
|
|
|
|
}
|
|
|
|
while ((aAbs & 1n) === 0n)
|
|
|
|
aAbs >>= 1n;
|
|
|
|
do {
|
|
|
|
while ((bAbs & 1n) === 0n)
|
|
|
|
bAbs >>= 1n;
|
|
|
|
if (aAbs > bAbs) {
|
|
|
|
const x = aAbs;
|
|
|
|
aAbs = bAbs;
|
|
|
|
bAbs = x;
|
|
|
|
}
|
|
|
|
bAbs -= aAbs;
|
|
|
|
} while (bAbs !== 0n);
|
|
|
|
// rescale
|
|
|
|
return aAbs << shift;
|
|
|
|
}
|
|
|
|
|
|
|
|
/**
|
|
|
|
* The least common multiple computed as abs(a*b)/gcd(a,b)
|
|
|
|
* @param a
|
|
|
|
* @param b
|
|
|
|
*
|
|
|
|
* @returns The least common multiple of a and b
|
|
|
|
*/
|
|
|
|
function lcm(a, b) {
|
2021-03-24 14:11:07 +00:00
|
|
|
if (typeof a === 'number')
|
|
|
|
a = BigInt(a);
|
|
|
|
if (typeof b === 'number')
|
|
|
|
b = BigInt(b);
|
|
|
|
if (a === 0n && b === 0n)
|
2021-03-24 13:04:30 +00:00
|
|
|
return BigInt(0);
|
2022-01-17 10:04:55 +00:00
|
|
|
// return abs(a * b) as bigint / gcd(a, b)
|
|
|
|
return abs((a / gcd(a, b)) * b);
|
2021-03-24 13:04:30 +00:00
|
|
|
}
|
|
|
|
|
|
|
|
/**
|
2022-09-12 08:59:51 +00:00
|
|
|
* Maximum. max(a,b)==a if a>=b. max(a,b)==b if a<b
|
2021-03-24 13:04:30 +00:00
|
|
|
*
|
|
|
|
* @param a
|
|
|
|
* @param b
|
|
|
|
*
|
|
|
|
* @returns Maximum of numbers a and b
|
|
|
|
*/
|
|
|
|
function max(a, b) {
|
|
|
|
return (a >= b) ? a : b;
|
|
|
|
}
|
|
|
|
|
|
|
|
/**
|
2022-09-12 08:59:51 +00:00
|
|
|
* Minimum. min(a,b)==b if a>=b. min(a,b)==a if a<b
|
2021-03-24 13:04:30 +00:00
|
|
|
*
|
|
|
|
* @param a
|
|
|
|
* @param b
|
|
|
|
*
|
|
|
|
* @returns Minimum of numbers a and b
|
|
|
|
*/
|
|
|
|
function min(a, b) {
|
|
|
|
return (a >= b) ? b : a;
|
|
|
|
}
|
|
|
|
|
|
|
|
/**
|
|
|
|
* Finds the smallest positive element that is congruent to a in modulo n
|
2021-03-24 14:11:07 +00:00
|
|
|
*
|
|
|
|
* @remarks
|
|
|
|
* a and b must be the same type, either number or bigint
|
|
|
|
*
|
2021-03-24 17:30:45 +00:00
|
|
|
* @param a - An integer
|
|
|
|
* @param n - The modulo
|
|
|
|
*
|
2022-09-12 08:59:51 +00:00
|
|
|
* @throws {@link RangeError} if n <= 0
|
2021-03-24 13:04:30 +00:00
|
|
|
*
|
2021-03-24 17:30:45 +00:00
|
|
|
* @returns A bigint with the smallest positive representation of a modulo n
|
2021-03-24 13:04:30 +00:00
|
|
|
*/
|
|
|
|
function toZn(a, n) {
|
2021-03-24 14:11:07 +00:00
|
|
|
if (typeof a === 'number')
|
|
|
|
a = BigInt(a);
|
|
|
|
if (typeof n === 'number')
|
|
|
|
n = BigInt(n);
|
|
|
|
if (n <= 0n) {
|
2021-03-24 17:30:45 +00:00
|
|
|
throw new RangeError('n must be > 0');
|
2021-03-24 13:04:30 +00:00
|
|
|
}
|
2021-03-24 14:11:07 +00:00
|
|
|
const aZn = a % n;
|
|
|
|
return (aZn < 0n) ? aZn + n : aZn;
|
2021-03-24 13:04:30 +00:00
|
|
|
}
|
|
|
|
|
|
|
|
/**
|
|
|
|
* Modular inverse.
|
|
|
|
*
|
|
|
|
* @param a The number to find an inverse for
|
|
|
|
* @param n The modulo
|
|
|
|
*
|
2022-09-12 08:59:51 +00:00
|
|
|
* @throws {@link RangeError} if a does not have inverse modulo n
|
2021-03-24 17:30:45 +00:00
|
|
|
*
|
|
|
|
* @returns The inverse modulo n
|
2021-03-24 13:04:30 +00:00
|
|
|
*/
|
|
|
|
function modInv(a, n) {
|
2021-03-24 17:30:45 +00:00
|
|
|
const egcd = eGcd(toZn(a, n), n);
|
|
|
|
if (egcd.g !== 1n) {
|
|
|
|
throw new RangeError(`${a.toString()} does not have inverse modulo ${n.toString()}`); // modular inverse does not exist
|
2021-03-24 13:04:30 +00:00
|
|
|
}
|
2021-03-24 17:30:45 +00:00
|
|
|
else {
|
|
|
|
return toZn(egcd.x, n);
|
2021-03-24 13:04:30 +00:00
|
|
|
}
|
|
|
|
}
|
|
|
|
|
|
|
|
/**
|
|
|
|
* Modular exponentiation b**e mod n. Currently using the right-to-left binary method
|
|
|
|
*
|
|
|
|
* @param b base
|
|
|
|
* @param e exponent
|
|
|
|
* @param n modulo
|
|
|
|
*
|
2022-09-12 08:59:51 +00:00
|
|
|
* @throws {@link RangeError} if n <= 0
|
2021-03-24 17:30:45 +00:00
|
|
|
*
|
|
|
|
* @returns b**e mod n
|
2021-03-24 13:04:30 +00:00
|
|
|
*/
|
|
|
|
function modPow(b, e, n) {
|
2021-03-24 14:11:07 +00:00
|
|
|
if (typeof b === 'number')
|
|
|
|
b = BigInt(b);
|
|
|
|
if (typeof e === 'number')
|
|
|
|
e = BigInt(e);
|
|
|
|
if (typeof n === 'number')
|
|
|
|
n = BigInt(n);
|
|
|
|
if (n <= 0n) {
|
2021-03-24 17:30:45 +00:00
|
|
|
throw new RangeError('n must be > 0');
|
2021-03-24 13:04:30 +00:00
|
|
|
}
|
2021-03-24 14:11:07 +00:00
|
|
|
else if (n === 1n) {
|
2021-03-24 17:30:45 +00:00
|
|
|
return 0n;
|
2021-03-24 13:04:30 +00:00
|
|
|
}
|
2021-03-24 14:11:07 +00:00
|
|
|
b = toZn(b, n);
|
2021-03-24 13:04:30 +00:00
|
|
|
if (e < 0n) {
|
2021-03-24 14:11:07 +00:00
|
|
|
return modInv(modPow(b, abs(e), n), n);
|
2021-03-24 13:04:30 +00:00
|
|
|
}
|
|
|
|
let r = 1n;
|
|
|
|
while (e > 0) {
|
|
|
|
if ((e % 2n) === 1n) {
|
2021-03-24 14:11:07 +00:00
|
|
|
r = r * b % n;
|
2021-03-24 13:04:30 +00:00
|
|
|
}
|
|
|
|
e = e / 2n;
|
2021-03-24 14:11:07 +00:00
|
|
|
b = b ** 2n % n;
|
2021-03-24 13:04:30 +00:00
|
|
|
}
|
|
|
|
return r;
|
|
|
|
}
|
|
|
|
|
|
|
|
export { abs, bitLength, eGcd, gcd, lcm, max, min, modInv, modPow, toZn };
|
2022-09-12 08:59:51 +00:00
|
|
|
//# sourceMappingURL=data:application/json;charset=utf-8;base64,eyJ2ZXJzaW9uIjozLCJmaWxlIjoiaW5kZXgubm9kZS5qcyIsInNvdXJjZXMiOlsiLi4vLi4vc3JjL3RzL2Ficy50cyIsIi4uLy4uL3NyYy90cy9iaXRMZW5ndGgudHMiLCIuLi8uLi9zcmMvdHMvZWdjZC50cyIsIi4uLy4uL3NyYy90cy9nY2QudHMiLCIuLi8uLi9zcmMvdHMvbGNtLnRzIiwiLi4vLi4vc3JjL3RzL21heC50cyIsIi4uLy4uL3NyYy90cy9taW4udHMiLCIuLi8uLi9zcmMvdHMvdG9abi50cyIsIi4uLy4uL3NyYy90cy9tb2RJbnYudHMiLCIuLi8uLi9zcmMvdHMvbW9kUG93LnRzIl0sInNvdXJjZXNDb250ZW50IjpudWxsLCJuYW1lcyI6W10sIm1hcHBpbmdzIjoiQUFBQTs7Ozs7O0FBTUc7QUFDRyxTQUFVLEdBQUcsQ0FBRSxDQUFnQixFQUFBO0FBQ25DLElBQUEsT0FBTyxDQUFDLENBQUMsSUFBSSxDQUFDLElBQUksQ0FBQyxHQUFHLENBQUMsQ0FBQyxDQUFBO0FBQzFCOztBQ1RBOzs7OztBQUtHO0FBQ0csU0FBVSxTQUFTLENBQUUsQ0FBZ0IsRUFBQTtJQUN6QyxJQUFJLE9BQU8sQ0FBQyxLQUFLLFFBQVE7QUFBRSxRQUFBLENBQUMsR0FBRyxNQUFNLENBQUMsQ0FBQyxDQUFDLENBQUE7SUFFeEMsSUFBSSxDQUFDLEtBQUssRUFBRSxFQUFFO0FBQUUsUUFBQSxPQUFPLENBQUMsQ0FBQTtBQUFFLEtBQUE7SUFDMUIsSUFBSSxJQUFJLEdBQUcsQ0FBQyxDQUFBO0lBQ1osR0FBRztBQUNELFFBQUEsSUFBSSxFQUFFLENBQUE7QUFDUCxLQUFBLFFBQVEsQ0FBQyxDQUFDLEtBQUssRUFBRSxJQUFJLEVBQUUsRUFBQztBQUN6QixJQUFBLE9BQU8sSUFBSSxDQUFBO0FBQ2I7O0FDVkE7Ozs7Ozs7Ozs7QUFVRztBQUNhLFNBQUEsSUFBSSxDQUFFLENBQWdCLEVBQUUsQ0FBZ0IsRUFBQTtJQUN0RCxJQUFJLE9BQU8sQ0FBQyxLQUFLLFFBQVE7QUFBRSxRQUFBLENBQUMsR0FBRyxNQUFNLENBQUMsQ0FBQyxDQUFDLENBQUE7SUFDeEMsSUFBSSxPQUFPLENBQUMsS0FBSyxRQUFRO0FBQUUsUUFBQSxDQUFDLEdBQUcsTUFBTSxDQUFDLENBQUMsQ0FBQyxDQUFBO0FBRXhDLElBQUEsSUFBSSxDQUFDLElBQUksRUFBRSxJQUFJLENBQUMsSUFBSSxFQUFFO0FBQUUsUUFBQSxNQUFNLElBQUksVUFBVSxDQUFDLHFCQUFxQixDQUFDLENBQUE7SUFFbkUsSUFBSSxDQUFDLEdBQUcsRUFBRSxDQUFBO0lBQ1YsSUFBSSxDQUFDLEdBQUcsRUFBRSxDQUFBO0lBQ1YsSUFBSSxDQUFDLEdBQUcsRUFBRSxDQUFBO0lBQ1YsSUFBSSxDQUFDLEdBQUcsRUFBRSxDQUFBO0lBRVYsT0FBTyxDQUFDLEtBQUssRUFBRSxFQUFFO0FBQ2YsUUFBQSxNQUFNLENBQUMsR0FBRyxDQUFDLEdBQUcsQ0FBQyxDQUFBO0FBQ2YsUUFBQSxNQUFNLENBQUMsR0FBVyxDQUFDLEdBQUcsQ0FBQyxDQUFBO1FBQ3ZCLE1BQU0sQ0FBQyxHQUFHLENBQUMsSUFBSSxDQUFDLEdBQUcsQ0FBQyxDQUFDLENBQUE7UUFDckIsTUFBTSxDQUFDLEdBQUcsQ0FBQyxJQUFJLENBQUMsR0FBRyxDQUFDLENBQUMsQ0FBQTtRQUNyQixDQUFDLEdBQUcsQ0FBQyxDQUFBO1FBQ0wsQ0FBQyxHQUFHLENBQUMsQ0FBQTtRQUNMLENBQUMsR0FBRyxDQUFDLENBQUE7UUFDTCxDQUFDLEdBQUcsQ0FBQyxDQUFBO1FBQ0wsQ0FBQyxHQUFHLENBQUMsQ0FBQTtRQUNMLENBQUMsR0FBRyxDQUFDLENBQUE7QUFDTixLQUFBO0lBQ0QsT0FBTztBQUNMLFFBQUEsQ0FBQyxFQUFFLENBQUM7QUFDSixRQUFBLENBQUMsRUFBRSxDQUFDO0FBQ0osUUFBQSxDQUFDLEVBQUUsQ0FBQztLQUNMLENBQUE7QUFDSDs7QUMzQ0E7Ozs7Ozs7QUFPRztBQUNhLFNBQUEsR0FBRyxDQUFFLENBQWdCLEVBQUUsQ0FBZ0IsRUFBQTtJQUNyRCxJQUFJLElBQUksR0FBRyxDQUFDLE9BQU8sQ0FBQyxLQUFLLFFBQVEsSUFBSSxNQUFNLENBQUMsR0FBRyxDQUFDLENBQUMsQ0FBQyxDQUFDLEdBQUcsR0FBRyxDQUFDLENBQUMsQ0FBVyxDQUFBO0lBQ3RFLElBQUksSUFBSSxHQUFHLENBQUMsT0FBTyxDQUFDLEtBQUssUUFBUSxJQUFJLE1BQU0sQ0FBQyxHQUFHLENBQUMsQ0FBQyxDQUFDLENBQUMsR0FBRyxHQUFHLENBQUMsQ0FBQyxDQUFXLENBQUE7SUFFdEUsSUFBSSxJQUFJLEtBQUssRUFBRSxFQUFFO0FBQ2YsUUFBQSxPQUFPLElBQUksQ0FBQTtBQUNaLEtBQUE7U0FBTSxJQUFJLElBQUksS0FBSyxFQUFFLEVBQUU7QUFDdEIsUUFBQSxPQUFPLElBQUksQ0FBQTtBQUNaLEtBQUE7SUFFRCxJQUFJLEtBQUssR0FBRyxFQUFFLENBQUE7SUFDZCxPQUFPLENBQUMsQ0FBQyxJQUFJLEdBQUcsSUFBSSxJQUFJLEVBQUUsTUFBTSxFQUFFLEVBQUU7UUFDbEMsSUFBSSxLQUFLLEVBQUUsQ0FBQTtRQUNYLElBQUksS0FBSyxFQUFFLENBQUE7QUFDWCxRQUFBLEtBQUssRUFBRSxDQUFBO0FBQ1IsS0FBQTtBQUNELElBQUEsT0FBTyxDQUFDLElBQUksR0FBRyxFQUFFLE1BQU0sRUFBRTtRQUFFLElBQUksS0FBSyxFQUFFLENBQUE7SUFDdEMsR0FBRztBQUNELFFBQUEsT0FBTyxDQUFDLElBQUksR0FBRyxFQUFFLE1BQU0sRUFBRTtZQUFFLElBQUksS0FBSyxFQUFFLENBQUE7UUFDdEMsSUFBSSxJQUFJLEdBQUcsSUFBSSxFQUFFO1lBQ2YsTUFBTSxDQUFDLEdBQUcsSUFBSSxDQUFBO1lBQ2QsSUFBSSxHQUFHLElBQUksQ0FBQTtZQUNYLElBQUksR0FBRyxDQUFDLENBQUE7QUFDVCxTQUFBO1FBQ0QsSUFBSSxJQUFJLElBQUksQ0FBQTtLQUNiLFFBQVEsSUFBSSxLQUFLLEVBQUUsRUFBQzs7SUFHckIsT0FBTyxJQUFJLElBQUksS0FBSyxDQUFBO0FBQ3RCOztBQ3BDQTs7Ozs7O0FBTUc7QUFDYSxTQUFBLEdBQUcsQ0FBRSxDQUFnQixFQUFFLENBQWdCLEVBQUE7SUFDckQsSUFBSSxPQUFPLENBQUMsS0FBSyxRQUFRO0FBQUUsUUFBQSxDQUFDLEdBQUcsTUFBTSxDQUFDLENBQUMsQ0FBQyxDQUFBO0lBQ3hDLElBQUksT0FBTyxDQUFDLEtBQUssUUFBUTtBQUFFLFFBQUEsQ0FBQyxHQUFHLE1BQU0sQ0FBQyxDQUFDLENBQUMsQ0FBQTtBQUV4QyxJQUFBLElBQUksQ0FBQyxLQUFLLEVBQUUsSUFBSSxDQUFDLEtBQUssRUFBRTtBQUFFLFFBQUEsT0FBTyxNQUFNLENBQUMsQ0FBQyxDQUFDLENBQUE7O0FBRTFDLElBQUEsT0FBTyxHQUFHLENBQUMsQ0FBQyxDQUFDLEdBQUcsR0FBRyxDQUF
|