2
$\begingroup$

I have a matrix $A$ for which I want to compute the quantity $T\lambda_j = \Pi_{\lambda_i\ne \lambda_j} \frac{A - \lambda_i I}{\lambda_j-\lambda_i}$, where $\lambda_i$ ($\lambda_j$) denote the eigenvalues of $A$. How can this be implemented in Mathematica? Just gave a try here:

A = {{1, 0, 0, 1},{0, 1, 2, 0},{1, 1, 0, 2},{0, 0, 0, 1}};
Eigenvalues[A]

{2, -1, 1, 1}

Tj = Product[(A - Eigenvalues[A][[i]] IdentityMatrix[4])/(
  Eigenvalues[A][[j]] - Eigenvalues[A][[i]]), {i, 1, 4}]
$\endgroup$
  • $\begingroup$ in which part exactly you want to exclude it in Tj !?? $\endgroup$ – Alrubaie Mar 29 at 18:02
  • $\begingroup$ do you want it to be skipped put not Zero right !? $\endgroup$ – Alrubaie Mar 29 at 18:07
  • $\begingroup$ @Alrubaie, there was a typo in my post. Just edited it. I want the denominator to be non-zero and hence avoid the case for which $i=j$. $\endgroup$ – H. Kenan Mar 29 at 18:08
  • $\begingroup$ @Alrubaie, my $i$ and $j$ are not the indices in my question. They are the eigenvalues. I should have used something like $\lambda_i$ and $\lambda_j$. $\endgroup$ – H. Kenan Mar 29 at 18:12
  • 2
    $\begingroup$ That product is presumably a matrix multiplication? $\endgroup$ – J. M. is away Mar 29 at 18:34
4
$\begingroup$

Here is my pedestrian implementation of your formula:

a = {{1, 0, 0, 1}, {0, 1, 2, 0}, {1, 1, 0, 2}, {0, 0, 0, 1}};

ClearAll[t]
t[amat_, j_] := Module[
  {evals, usable},
  evals = Eigenvalues[amat];
  usable = DeleteDuplicates@Cases[evals, Except@evals[[j]] ];
  Dot @@ 
   Table[
     (amat - i IdentityMatrix[Length[amat]])/(evals[[j]] - i),
     {i, usable}
   ]
]

t[a, 4]

Mathematica graphics

You do not provide an example of desired output, so I will let you check whether this is what you expect.

$\endgroup$
  • $\begingroup$ Thanks, @MarcoB. It leads precisely to the expected result. However, it looks too complicated. Nevertheless, it is fine as it works. $\endgroup$ – H. Kenan Mar 29 at 18:46
3
$\begingroup$

Something like this?

Clear[A, evals, T]
A = {{1, 0, 0, 1}, {0, 1, 2, 0}, {1, 1, 0, 2}, {0, 0, 0, 1}};
T[A_?MatrixQ, j_Integer] := With[
  {evals = Eigenvalues[A], id = IdentityMatrix@Length@A},
  Dot @@ Table[
    If[evals[[j]] - evals[[i]] == 0, id, (A - evals[[i]] id)/(evals[[j]] - evals[[i]])],
    {i, Length@A}
    ]
  ]

MatrixForm /@ Array[T[A, #] &, 4]

enter image description here

$\endgroup$
2
$\begingroup$

This

A = {{1, 0, 0, 1},{0, 1, 2, 0},{1, 1, 0, 2},{0, 0, 0, 1}};
e=Eigenvalues[A];
Map[(A-e[[#[[1]]]]*IdentityMatrix[4])/(e[[#[[2]]]]-e[[#[[1]]]])&,
  DeleteCases[Tuples[Range[4],2],{i_,i_}]]

generates your twelve matricies with i not equal to j.

Put Dot@@ in front of that Map to form the dot product of the 12 matricies.

That works by forming every possible distinct i,j pair and then using those in the Map

If it might be easier to read you can also write it this way

Map[(ei=e[[#[[1]]]];ej=e[[#[[2]]]];
  (A-ei*IdentityMatrix[4])/(ej-ei))&,
  DeleteCases[Tuples[Range[4],2],{i_,i_}]]
$\endgroup$
  • $\begingroup$ Should e[[#[[2]]]]-e[[[[1]]]] be e[[#[[2]]]]-e[[#[[1]]]]? $\endgroup$ – That Gravity Guy Mar 29 at 19:12
  • $\begingroup$ @ThatGravityGuy Yes! Good catch. Thank you! Corrected. $\endgroup$ – Bill Mar 29 at 19:14
0
$\begingroup$

Another way:

ClearAll[t];
t[j_Integer, A_?SquareMatrixQ] := t[j, A, Eigenvalues@A];   (* add the eigenvalues *)
t[j_Integer, A_?SquareMatrixQ, evs_?VectorQ] /; Length@A == Length@evs := (* arg checks *)=
  Fold[
   #1.(A - #2 IdentityMatrix[Length@A])/(evs[[j]] - #2) &,
   IdentityMatrix[Length@A],
   Pick[evs, Unitize[evs - evs[[j]]], 1]                    (* Pick nonzero differences *)
   ];

Performance tuning: One can use DeleteCases[evs, e_ /; e == evs[[j]]] to pick the eigenvalues that give a nonzero difference. It makes no consistent difference to the timing on a 101 x 101 machine real matrix. One can save a little time by computing the identity matrix once and using With[] to inject it in the two places it occurs. One can also save time using dot = (dot = Dot; #2) & instead of Dot to skip the multiplication by the identity matrix (the first step of Fold[]). The differences evs - evs[[j]] appear twice, so they can be replaced by a single computation like the identity matrix. It can make up to a 10% improvement.

$\endgroup$

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.