Mdtryhht

I have a matrix $A$ for which I want to compute the quantity $Tlambda_j = Pi_lambda_ine lambda_j fracA - lambda_i Ilambda_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]

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]

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

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]

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_]]