function invModGcdEx(x, domain)
{
    if(x === 1)
    {
        return 1;
    }
    else
    {
        //В случае 0 или делителя нуля возвращается 0, означающий некий "некорректный результат"
        if(x === 0 || domain % x === 0)
        {
            return 0;
        }
        else
        {
            //Расширенный алгоритм Евклида, вычисляющий такое число tCurr, что tCurr * x + rCurr * N = 1
            //Другими словами, существует такое число rCurr, при котором tCurr * x mod N = 1
            let tCurr = 0;
            let rCurr = domain;
            let tNext = 1;
            let rNext = x;
            while(rNext !== 0)
            {
                let quotR = Math.floor(rCurr / rNext);
                let tPrev = tCurr;
                let rPrev = rCurr;

                tCurr = tNext;
                rCurr = rNext;

                tNext = Math.floor(tPrev - quotR * tCurr);
                rNext = Math.floor(rPrev - quotR * rCurr);
            }

            tCurr = (tCurr + domain) % domain;
            return tCurr;
        }
    }
}

function wholeMod(x, domain)
{
    return ((x % domain) + domain) % domain;
}

function mulSubRow(rowLeft, rowRight, mulValue, domain)
{
    for(let i = 0; i < rowLeft.length; i++)
    {
        rowLeft[i] = wholeMod(rowLeft[i] - mulValue * rowRight[i], domain);
    }
}

function mulRow(row, mulValue, domain)
{
    for(let i = 0; i < row.length; i++)
    {
        row[i] = (row[i] * mulValue) % domain;
    }
}

function mulSubRowCached(rowLeft, rowRight, mulValue, wholeModCache, cacheIndexOffset)
{
    for(let i = 0; i < rowLeft.length; i++)
    {
        rowLeft[i] = wholeModCache[rowLeft[i] - mulValue * rowRight[i] + cacheIndexOffset];
    }
}

function mulRowCached(row, mulValue, domain, wholeModCache, cacheIndexOffset)
{
    for(let i = 0; i < row.length; i++)
    {
        row[i] = domain - wholeModCache[cacheIndexOffset - row[i] * mulValue];
    }
}

function invertMatrix(matrix, domain)
{
    let matrixSize = matrix.length;

    //Инициализируем обратную матрицу единичной
    let invMatrix = [];
    for(let i = 0; i < matrixSize; i++)
    {
        let matrixRow = new Uint8Array(matrixSize);
        matrixRow.fill(0);

        matrixRow[i] = 1;
        invMatrix.push(matrixRow);
    }

    //Прямой ход: приведение матрицы к ступенчатому виду
    for(let i = 0; i < matrixSize; i++)
    {
        let thisRowFirst = matrix[i][i];
        if(thisRowFirst === 0 || (thisRowFirst !== 1 && domain % thisRowFirst === 0)) //Первый элемент строки 0 или делитель нуля, меняем строку местами со следующей строкой, у которой первый элемент не 0
        {
            for(let j = i + 1; j < matrixSize; j++)
            {
                let otherRowFirst = matrix[j][i];
                if(otherRowFirst !== 0 && (otherRowFirst === 1 || domain % otherRowFirst !== 0)) //Нашли строку с ненулевым первым элементом
                {
                    thisRowFirst = otherRowFirst;

                    let tmpMatrixRow = matrix[i];
                    matrix[i]        = matrix[j];
                    matrix[j]        = tmpMatrixRow;

                    let tmpInvMatrixRow = invMatrix[i];
                    invMatrix[i]        = invMatrix[j];
                    invMatrix[j]        = tmpInvMatrixRow;

                    break;
                }
            }
        }

        //Обнуляем первые элементы всех строк после первой, отнимая от них (otherRowFirst / thisRowFirst) * x mod N
        let invThisRowFirst = invModGcdEx(thisRowFirst, domain);
        for(let j = i + 1; j < matrixSize; j++)
        {
            let otherRowFirst = matrix[j][i];
            let mulValue      = invThisRowFirst * otherRowFirst;

            if(otherRowFirst !== 0 && (otherRowFirst === 1 || domain % otherRowFirst !== 0))
            {
                mulSubRow(matrix[j],    matrix[i],    mulValue, domain);
                mulSubRow(invMatrix[j], invMatrix[i], mulValue, domain);
            }
        }
    }

    //Обратный ход - обнуление всех элементов выше главной диагонали
    let matrixRank = matrixSize;
    for(let i = matrixSize - 1; i >= 0; i--)
    {
        let thisRowLast    = matrix[i][i];
        let invThisRowLast = invModGcdEx(thisRowLast, domain);
        for(let j = i - 1; j >= 0; j--)
        {
            let otherRowLast = matrix[j][i];
            let mulValue     = invThisRowLast * otherRowLast;

            if(otherRowLast !== 0 && (otherRowLast === 1 || domain % otherRowLast !== 0))
            {
                mulSubRow(matrix[j],    matrix[i],    mulValue, domain);
                mulSubRow(invMatrix[j], invMatrix[i], mulValue, domain);
            }
        }

        if(thisRowLast !== 0 && domain % thisRowLast !== 0)
        {
            mulRow(matrix[i],    invThisRowLast, domain);
            mulRow(invMatrix[i], invThisRowLast, domain);
        }

        if(matrix[i].every(val => val === 0))
        {
            matrixRank -= 1;
        }
    }

    return {inverse: invMatrix, rank: matrixRank};
}

function invertMatrixCachedInverses(matrix, domain)
{
    let matrixSize = matrix.length;

    //Инициализируем обратную матрицу единичной
    let invMatrix = [];
    for(let i = 0; i < matrixSize; i++)
    {
        let matrixRow = new Uint8Array(matrixSize);
        matrixRow.fill(0);

        matrixRow[i] = 1;
        invMatrix.push(matrixRow);
    }

    //Вычисляем все обратные элементы заранее
    let domainInvs = [];
    for(let d = 0; d < domain; d++)
    {
        domainInvs.push(invModGcdEx(d, domain));
    }

    //Прямой ход: приведение матрицы к ступенчатому виду
    for(let i = 0; i < matrixSize; i++)
    {
        let thisRowFirst = matrix[i][i];
        if(domainInvs[thisRowFirst] === 0) // <--- Первый элемент строки 0 или делитель нуля, меняем строку местами со следующей строкой, у которой первый элемент не 0
        {
            for(let j = i + 1; j < matrixSize; j++)
            {
                let otherRowFirst = matrix[j][i];
                if(domainInvs[otherRowFirst] !== 0) // <--- Нашли строку с ненулевым первым элементом
                {
                    thisRowFirst = otherRowFirst;

                    let tmpMatrixRow = matrix[i];
                    matrix[i]        = matrix[j];
                    matrix[j]        = tmpMatrixRow;

                    let tmpInvMatrixRow = invMatrix[i];
                    invMatrix[i]        = invMatrix[j];
                    invMatrix[j]        = tmpInvMatrixRow;

                    break;
                }
            }
        }

        //Обнуляем первые элементы всех строк после первой, отнимая от них (otherRowFirst / thisRowFirst) * x mod N
        let invThisRowFirst = domainInvs[thisRowFirst]; // <---
        for(let j = i + 1; j < matrixSize; j++)
        {
            let otherRowFirst = matrix[j][i];
            let mulValue      = invThisRowFirst * otherRowFirst;

            if(domainInvs[otherRowFirst] !== 0) // <---
            {
                mulSubRow(matrix[j],    matrix[i],    mulValue, domain);
                mulSubRow(invMatrix[j], invMatrix[i], mulValue, domain);
            }
        }
    }

    //Обратный ход - обнуление всех элементов выше главной диагонали
    let matrixRank = matrixSize;
    for(let i = matrixSize - 1; i >= 0; i--)
    {
        let thisRowLast    = matrix[i][i];
        let invThisRowLast = domainInvs[thisRowLast]; // <---
        for(let j = i - 1; j >= 0; j--)
        {
            let otherRowLast = matrix[j][i];
            let mulValue     = invThisRowLast * otherRowLast;

            if(domainInvs[otherRowLast] !== 0) // <---
            {
                mulSubRow(matrix[j],    matrix[i],    mulValue, domain);
                mulSubRow(invMatrix[j], invMatrix[i], mulValue, domain);
            }
        }

        if(domainInvs[thisRowLast] !== 0) // <---
        {
            mulRow(matrix[i],    invThisRowLast, domain);
            mulRow(invMatrix[i], invThisRowLast, domain);
        }

        if(matrix[i].every(val => val === 0))
        {
            matrixRank -= 1;
        }
    }

    return {inverse: invMatrix, rank: matrixRank};
}

function invertMatrixCachedMod(matrix, domain)
{
    let matrixSize = matrix.length;

    //Инициализируем обратную матрицу единичной
    let invMatrix = [];
    for(let i = 0; i < matrixSize; i++)
    {
        let matrixRow = new Uint8Array(matrixSize);
        matrixRow.fill(0);

        matrixRow[i] = 1;
        invMatrix.push(matrixRow);
    }

    //Вычисляем все обратные элементы заранее
    let domainInvs = [];
    for(let d = 0; d < domain; d++)
    {
        domainInvs.push(invModGcdEx(d, domain));
    }

    //Вычисляем кэш деления по модулю
    const cacheIndexOffset = (domain - 1) * (domain - 1);

    let wholeModCache = new Uint8Array((domain - 1) * (domain - 1) + domain);
    for(let i = 0; i < wholeModCache.length; i++)
    {
        let divisor      = i - cacheIndexOffset; //[-cacheIndexOffset, domain - 1]
        wholeModCache[i] = wholeMod(divisor, domain);
    }

    //Прямой ход: приведение матрицы к ступенчатому виду
    for(let i = 0; i < matrixSize; i++)
    {
        let thisRowFirst = matrix[i][i];
        if(domainInvs[thisRowFirst] === 0) //Первый элемент строки 0 или делитель нуля, меняем строку местами со следующей строкой, у которой первый элемент не 0
        {
            for(let j = i + 1; j < matrixSize; j++)
            {
                let otherRowFirst = matrix[j][i];
                if(domainInvs[otherRowFirst] !== 0) //Нашли строку с ненулевым первым элементом
                {
                    thisRowFirst = otherRowFirst;

                    let tmpMatrixRow = matrix[i];
                    matrix[i]        = matrix[j];
                    matrix[j]        = tmpMatrixRow;

                    let tmpInvMatrixRow = invMatrix[i];
                    invMatrix[i]        = invMatrix[j];
                    invMatrix[j]        = tmpInvMatrixRow;

                    break;
                }
            }
        }

        //Обнуляем первые элементы всех строк после первой, отнимая от них (otherRowFirst / thisRowFirst) * x mod N
        let invThisRowFirst = domainInvs[thisRowFirst];
        for(let j = i + 1; j < matrixSize; j++)
        {
            let otherRowFirst = matrix[j][i];
            let mulValue      = domain - wholeModCache[cacheIndexOffset - invThisRowFirst * otherRowFirst]; // <---

            if(domainInvs[otherRowFirst] !== 0)
            {
                mulSubRowCached(matrix[j],    matrix[i],    mulValue, wholeModCache, cacheIndexOffset); // <---
                mulSubRowCached(invMatrix[j], invMatrix[i], mulValue, wholeModCache, cacheIndexOffset); // <---
            }
        }
    }

    //Обратный ход - обнуление всех элементов выше главной диагонали
    let matrixRank = matrixSize;
    for(let i = matrixSize - 1; i >= 0; i--)
    {
        let thisRowLast    = matrix[i][i];
        let invThisRowLast = domainInvs[thisRowLast];
        for(let j = i - 1; j >= 0; j--)
        {
            let otherRowLast = matrix[j][i];
            let mulValue     = domain - wholeModCache[cacheIndexOffset - invThisRowLast * otherRowLast]; // <---

            if(domainInvs[otherRowLast] !== 0)
            {
                mulSubRowCached(matrix[j],    matrix[i],    mulValue, wholeModCache, cacheIndexOffset); // <---
                mulSubRowCached(invMatrix[j], invMatrix[i], mulValue, wholeModCache, cacheIndexOffset); // <---
            }
        }

        if(domainInvs[thisRowLast] !== 0)
        {
            mulRowCached(matrix[i],    invThisRowLast, domain, wholeModCache, cacheIndexOffset); // <---
            mulRowCached(invMatrix[i], invThisRowLast, domain, wholeModCache, cacheIndexOffset); // <---
        }

        if(matrix[i].every(val => val === 0))
        {
            matrixRank -= 1;
        }
    }

    return {inverse: invMatrix, rank: matrixRank};
}

function createRandomMatrix(matrixSize, domain)
{
    let matrix = [];
    for(let i = 0; i < matrixSize; i++)
    {
        let matrixRow = new Uint8Array(matrixSize);
        for(let j = 0; j < matrixSize; j++)
        {
            matrixRow[j] = Math.floor(Math.random() * domain);
        }

        matrix.push(matrixRow);
    }

    return matrix;
}

function main()
{
    let matrixSize = 500;
    let domains    = [2, 3, 5, 7, 11, 13, 17, 19, 23, 29, 31, 37, 41, 43, 47, 53, 59, 61, 67, 71, 73, 79, 83, 89, 97]

    const test_count = 5;
    for(let i = 0; i < domains.length; i++)
    {
        let domain = domains[i];
        console.log("");
        console.log("DOMAIN: " + domain);

        for(let i = 0; i < test_count; i++)
        {
            let matrix = createRandomMatrix(matrixSize, domain);

            console.time("Matrix regular inversion " + i);
            invMatrixData = invertMatrix(matrix, domain);
            console.timeEnd("Matrix regular inversion " + i);
        }

        for(let i = 0; i < test_count; i++)
        {
            let matrix = createRandomMatrix(matrixSize, domain);

            console.time("Matrix domain inverses " + i);
            invMatrixData = invertMatrixCachedInverses(matrix, domain);
            console.timeEnd("Matrix domain inverses "  + i);
        }

        for(let i = 0; i < test_count; i++)
        {
            let matrix = createRandomMatrix(matrixSize, domain);

            console.time("Matrix cached all " + i);
            invMatrixData = invertMatrixCachedMod(matrix, domain);
            console.timeEnd("Matrix cached all " + i);
        }
    }
}