Suppose that I have a linear system $AX=b$ with $A\in\mathbb{Z}^{n\times m}$ and $b\in\mathbb{Z}^{n}$. Assume that $AX=b$ has exactly one rational solution. Then how can I obtain this solution efficiently?

I know LU decomposition is a choice for this purpose. But it seems that the denominators may expand in the process of computing a LU decomposition, which I hope to avoid if the rational solution has small denominators.

If I already have a numerical solution of $AX=b$, does it help?

  • $\begingroup$ It is not necessary that $m=n$ though there is exactly one solution. If the rational solution has no large denominators, then I hope no large denominators appear in the computation. $\endgroup$
    – Jie Wang
    Jan 1 '20 at 23:28
  • 2
    $\begingroup$ Maybe this is relevant: cs6505.wordpress.com/2018/10/02/linear-equations-i. One can solve in $O(n^\omega L)$ time an $n\times n$ system with integer entries which take at most $L$ bits to represent. $\endgroup$ Jan 2 '20 at 15:57
  • 7
    $\begingroup$ I've voted to reopen. It is a perfectly legitimate question, concerned with the well known phenomenon of "intermediate expression swell" in computer algebra. Such linear algebra problems may easily arise in applied mathematical research, where practitioners are not expected to be proficient in linear algebra over rings (like the integers). Hence this question could be useful to them. $\endgroup$ Jan 2 '20 at 20:32

The following answer is only theoretical. Suppose you first reduce your matrix to Smith normal form $D = S A T$, where $S$ and $T$ are integral matrices with integral inverses (the determinants of $S$ and $T$ are each either $+1$ or $-1$) and $D$ is diagonal, whose diagonal elements are the so-called elementary divisors of $A$. Then you reduce your problem to solving $DY=c$, where $c = S b$ and $X = T Y$. Since by $c$ is still integral and $D$ is now diagonal, the denominators in $Y$ will be no larger than the largest entry in $D$ (the largest elementary divisor of $A$). Then, since $X$ is obtained from $Y$ by multiplication with another integral matrix, the denominators in $X$ will be no larger than the lcm of the denominators in $Y$.

Unfortunately, I don't think that reduction to Smith normal form is a cheap operation, so whether it is worth going through that step in your problem will depend on other factors.

EDIT: Actually, the Smith normal form is overkill for your problem. The Hermite normal form is basically the integer version of the LU decomposition that you already had in mind: $U = S A$, where $U$ is integral upper triangular, and $S$ is integral with integral inverse. The new system you need to solve is $U X = S b$, which can be done by back substitution, where at each step you only divide by the corresponding diagonal element of $U$. Hence, the final denominators in $X$ will be no larger than the lcm of the diagonal elements of $U$.

The linked Wikipedia article on the Hermite normal form mentions that modern algorithms for its computation have polynomially bounded time and space complexity (respectively, in terms of the matrix size and logarithmic size of the matrix entries). If you (or your computer algebra system) use such an algorithm, it might be your best bet for provably avoiding intermediate expression swell. The practicality of this method for your problem of course would still be subject to experimentation.

EDIT 2: I've also found a very well written survey of several algorithms for solving rational linear systems, with practical timing comparisons in the case of sparse systems.

Cook, William; Steffy, Daniel E., Solving very sparse rational systems of equations, ACM Trans. Math. Softw. 37, No. 4, Paper No. 39, 21 p. (2011). ZBL1365.65121. (author's copy)

  • 1
    $\begingroup$ Part of the hypothesis of OP is that the system has exactly one rational solution. Is there anything to be gained through use of this hypothesis? $\endgroup$ Jan 17 at 5:44
  • $\begingroup$ @GerryMyerson I don't see how, at least in this factorization approach. $\endgroup$ Jan 17 at 23:17

This is standard and can be done using $p$-adic lifting, see for example the following article:

  • J.D. Dixon, Exact solution of linear equations using $P$-adic expansions, Numer. Math. 40 (1982) pp. 137–141, doi:10.1007/BF01459082.

This also exploits that there is a unique solution.


Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.