# Numerically stable computation of matrix power

**URL:** <https://discourse.julialang.org/t/numerically-stable-computation-of-matrix-power/31584>\
**Category:** General Usage\
**Tags:** linearalgebra\
**Created:** [November 27, 2019, 4:20pm UTC](https://discourse.julialang.org/t/numerically-stable-computation-of-matrix-power/31584 "2019-11-27T16:20:19Z")\
**Posts on this page:** 1\
**Showing post:** 7

<div class="post-metadata">

**Author:** ![stevengj](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/stevengj/32/71_2.png) [@stevengj](https://discourse.julialang.org/u/stevengj)\
**Post date:** [November 27, 2019, 6:35pm UTC](https://discourse.julialang.org/t/numerically-stable-computation-of-matrix-power/31584/7 "2019-11-27T18:35:30Z")

</div>

Sorry, you’re quite right — the basic problem is that the nullspace of M-I is not orthogonal to the other Jordan vectors of M, so the correct projection is non-orthogonal.

I think what you need, instead, is to project using left and right nullspaces. I came up with the following a bit hastily, but it looks right to me and works for your example:

```julia
function Minf(M)
    R = nullspace(M-I)
    L = nullspace(M'-I)'
    return R / (L*R) * L
end

```

You might want to check it a bit more carefully.

---

_[View the full topic](https://discourse.julialang.org/t/numerically-stable-computation-of-matrix-power/31584)._
