7
$\begingroup$

Let $A$ be a $n\times n$ matrix over a commutative ring. I'm looking for a good method to compute its adjugate matrix.

My current approach is to use the Cayley-Hamilton theorem: $$\text{adj}(A) = -(A^{n-1} + c_1A^{n-2}+\ldots +c_{n-2}A+c_{n-1}\text{I}) $$ where $$\lambda^n + c_1\lambda^{n-1} + \ldots + c_{n-1}\lambda +c_n = \det(\lambda\text{I}-A).$$ i.e. the $c_i$ are the coefficients of the characteristic polynomial, which I can determine using the Berkowitz algorithm. In particular, I never have to make divisions, so all works over an arbitrary commutative ring.

However, this method needs $\mathcal{O}(n^4)$ arithmetical ring operations, so there might be better alternatives. I hope for something like $\mathcal{O}(n^3)$, as for matrix inversion over a field.

$\endgroup$
  • $\begingroup$ If you take the generic matrix (the one with indeterminate entries $a_{11}\dots$), you can compute its inverse (inside the field of fractions), multiply by the determinant to eventually get the generic adjugate matrix. Anyway I'd expect the practical complexity to depend on the entries of the given matrix and not only on the size of the matrix...? $\endgroup$ – YCor Dec 5 '14 at 11:25
  • 1
    $\begingroup$ I think that you can run Gaussian elimination keeping the denominators written formally instead of computing them. To compute sums you just reduce everything to common denominator. In the end you'll end up with a "formal denominator" equal to the product of the pivots, that is, the determinant. You just strike it out and it's done. $\endgroup$ – Federico Poloni Dec 5 '14 at 12:51
  • $\begingroup$ @YCor: This means that you use the formula for the entries of the adjugate matrix. But this is a polynomial of degree $n-1$ in $n^2$ variables for each entry, which is very ineffective. $\endgroup$ – Nicolas Malebranche Dec 5 '14 at 12:54
  • $\begingroup$ @Federico Poloni: this will work only over integral domains. But if I do so, I will need to simplify ratios, which requires a good knowledge of the ring I'm using and is maybe hard to carry out. But I'm hoping to find a method that is independent of the chosen base ring. $\endgroup$ – Nicolas Malebranche Dec 5 '14 at 13:03
  • $\begingroup$ @NicolasMalebranche: but the case of the generic matrix is a particular case of what you ask, if you consider general commutative rings (or fields), since one particular case is the function ring on $n^2$ indeterminates. If you have something else in mind, you should specify which kind of (computable) commutative rings you want to consider. $\endgroup$ – YCor Dec 5 '14 at 15:39
2
$\begingroup$

A reference dated 2001-2003 is http://www4.ncsu.edu/~kaltofen/bibliography/01/KaVi01.pdf

Over a commutative ring, the complexities of the calculation of $\det(A)$ and $adj(A)$ are the same. With the standard multiplication, the authors obtain a bit-complexity in $n^{10/3+\epsilon}N^{1+\epsilon}$ where $N$ is the length of the entries of $A$. Yet, the programming seems complicated ; there is a simpler method with exponent $3.5$.

Conclusion: If you have a robust and simple method with exponent $4$, then it is perfect.

$\endgroup$
1
$\begingroup$

(EDIT: there is a flaw when using $\det A$ below, as pointed out by @NicolasMalebranche.)

This answer expands on my comment in the hope to convince you that Gaussian elimination can be run formally without the need to ever simplify ratios. Your observation that this only applies to an integral domain is correct, though.

Consider the $3\times 3$ case for simplicity. Let us call $p_1,p_2,p_3$ the pivots encountered during Gaussian elimination, and suppose for simplicity that no row exchanges are needed. Let us suppose at first that all the pivots are nonzero.

At the first step, we multiply our matrix $A$ on the left by $$ L_1 = \begin{bmatrix} 1 & 0 & 0\\ -a_{21} & p_1 & 0\\ -a_{31} & 0 & p_1 \end{bmatrix}, $$ which zeroes out the entries in position (2,1) and (3,1). Then we premultiply by a matrix with form $$ L_2 = \begin{bmatrix} 1 & 0 & 0\\ 0 & 1 & 0\\ 0 & * & p_2 \end{bmatrix}, $$ to zero out the (3,2) entry. Let us call the resulting upper triangular matrix $$ U= \begin{bmatrix} p_1 & * & *\\ 0 & p_2 & *\\ 0 & 0 & p_3 \end{bmatrix}. $$ We have written $L_2L_1A=U$, so $\operatorname{Adj} A = (p_1p_2p_3)A^{-1} = (p_1p_2p_3) U^{-1}L_2L_1$.

Now I'll show that the result of this computation can be determined by a sort of back-substitution. Let $b$ be a generic column of $L_2L_1$. We write $Ux=b$ as $$ \begin{bmatrix} p_1 & * & *\\ 0 & p_2 & *\\ 0 & 0 & p_3 \end{bmatrix} \begin{bmatrix} x_1\\x_2\\x_3 \end{bmatrix} = \begin{bmatrix} b_1\\b_2\\b_3 \end{bmatrix} $$ We aim to solve this linear system not for $x$ but for $p_1p_2p_3x$.

From the third row, we can determine $p_3x_3$ without divisions. Multiplying the second equation by $p_3$, we can get $p_2p_3x_2$ without divisions. Multiplying the first equation by $p_2p_3$, we can get $p_1p_2p_3x_1$ without divisions.

So we just need to multiply by the remaining factors to obtain the entries of $p_1p_2p_3x$.

For a generic $n\times n$ matrix, everything should be perfectly analogous. One needs to precompute the products $p_1p_2\dots p_i$ and $p_i p_{i+1}\dots p_n$ for each $i$; once you do that, I think that all the computations can be performed in $O(n^3)$.

What if some pivots are zero? If we use total pivoting, the only breakdown case is when we arrive at a point when the last $n-k$ rows of $L_{k}\dots L_2L_1U$ are zero. We can restrict to the case $k=n-1$, otherwise $A$ has rank $n-2$ or lower and $\operatorname{Adj} A$ is trivially zero. So $p_n=0$ and all the previous pivots are nonzero. Set $L_n=I$, and factor as above to get $LA=U$, with $L$ lower triangular and invertible. We have $\operatorname{Adj} A = \operatorname{Adj} U \operatorname{Adj} L^{-1}$. The term $\operatorname{Adj} L^{-1}$ is a multiple of $L$, so it should pose little problems. The term $\operatorname{Adj} U$ has nonzero only in its last row (since the last row of $U$ is zeros), and it should be possible to compute these determinants in $O(n^2)$ each, since each is the determinant of an invertible triangular matrix + a rank-1 matrix (failing everything, use the matrix determinant lemma and the trick explained previously to compute $\det B B^{-1}u$, where $B$ is the top $(n-1)\times(n-1)$ block of $U$).

All of this is absolutely inelegant, I concur.

$\endgroup$
  • 1
    $\begingroup$ Thank you very much for your answer. But from $L_2L_1A = U$ I read off $p_1^2p_2\det A = p_1p_2p_3$, so I can't see how $\det A$ should equal $p_1p_2p_3$. $\endgroup$ – Nicolas Malebranche Dec 7 '14 at 9:02
  • $\begingroup$ The example I had in mind when I was remarking that Gaussian elimination will need ring divisions was some unimodular matrix like $\begin{bmatrix} 2 & 3 \\ 3 & 5 \end{bmatrix}$ where $\det A=1$, so $\text{Adj}(A) = A^{-1}$. That means that any formal divisions are forbidden, since the common denominator, the determinant is $1$, right? $\endgroup$ – Nicolas Malebranche Dec 7 '14 at 9:05
  • $\begingroup$ @NicolasMalebranche You are correct about the determinant. The algorithm I suggested won't work then. I am still not convinced that there cannot be a simple trick to make those simplifications 'automatically', though. $\endgroup$ – Federico Poloni Dec 7 '14 at 9:36

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.