How to write the Taylor expansion of the sum of several matrices?
dotk[ab_, k_] := Dot @@ ConstantArray[ab, k]
dotk[ab_, 0] := IdentityMatrix[Length@ab];
mInvSeries[a_, b_, n_] := Sum[(-1)^k Inverse[a].dotk[(b.Inverse[a]), k], {k, 0, n}]
Usage
a = 10 {{6, 3}, {2, 8}};
ai = mInvSeries[a, IdentityMatrix[2], 150];
(a.ai ) // N
(* {{0.981341, 0.00691085}, {0.00460723, 0.985948}} *)
How about using NonCommutativeMultiply
as matrix multiplication?
Then, we can let -1 escape the noncommutative multiply
Unprotect[NonCommutativeMultiply];
NonCommutativeMultiply[H___, Times[-1, M_], T___] := - NonCommutativeMultiply[H, M, T]
And do the Taylor series about the big
matrix by the small
matrix by hand:
inv[n_][big_, small_] := Total@FoldList[
#1 ** #2 &,
Inverse[big],
-(small ** Inverse[big]) & /@ Range[n]
]
Then, inv[3][A, B]
gives
Inverse[A] - Inverse[A] ** B ** Inverse[A] +
Inverse[A] ** B ** Inverse[A] ** B ** Inverse[A] -
Inverse[A] ** B ** Inverse[A] ** B ** Inverse[A] ** B ** Inverse[A]
and inv[3][A, B + C]
gives
Inverse[A] - Inverse[A] ** (B + C) ** Inverse[A] +
Inverse[A] ** (B + C) ** Inverse[A] ** (B + C) ** Inverse[A] -
Inverse[A] ** (B + C) ** Inverse[A] ** (B + C) **
Inverse[A] ** (B + C) ** Inverse[A]
In the same spirit to Dr. belisarius' answer:
abSeries[A_, B_, a_, n_] := Series[Inverse[A + a B], {a, 0, n}]
After obtaining the series you can get the matrix without the a
by using
ReplaceAll[Normal[abSeries[A, B, a, n]],a->1]
There is definitely a way to generate the expression you quoted for "general" matrices in mathematica, but since I don't know the math in detail I don't know how to teach mathematica to derive it either. Of course programming a known expression is easy, similarly to what Dr. belisarius did.
Addition
This might be more akin to what was originally intended:
abSeries[A_, B_, n_] := Plus @@ NestList[-Dot[Inverse[A].B, #] &, Inverse[A], n]
It gives an expression for general expressions A
and B
meaning you do not have to input a literal matrix for the expressions for this to work. Of course to get the final literal result you still need to put in the matrices using ReplaceAll
after the calculation.
If you want to use three matrices you can simply put in the addition of two matrices as A
or B
. i.e. B=C+D
. You can then use Distribute
and Factor
to put the resulting expression into different forms.
Example:
In[1]:=ReplaceAll[abSeries[a, ee b, 5].(a + ee b), {a -> PauliMatrix[1], b -> PauliMatrix[2]}] // Simplify
Out[1]:={{1 + ee^6, 0}, {0, 1 + ee^6}}