# Optimising subdiagonal shift matrix generation time

I want a square matrix $$A$$ such that $$A_{i,i-1} = 1 \ \forall i\\ A_{i,j} = 0 \text{ for } \ j \ne i-1$$

I have tried

tminus[n_] := Table[If[i == j + 1, 1, 0], {i, 1, n}, {j, 1, n}]


This is much better than

tminus[n_] := Module[{m = ConstantArray[0, {n, n}], i, j},
For[i = 1, i <= n, i++,
For[j = 1, j <= n, j++,
If[i == j + 1, m[[i, j]] = 1, m[[i, j]] = 0]
];
];
Return[m];
]


However, I wanted to know a better and optimised version so that I can create matrices of order 10000 and higher.

• try not to use For loops in Mathematica. There are faster, efficient functional ways of implementing stuff without using loops. For loops are error-prone Jun 18, 2019 at 19:12

Since your matrix is large and has very few nonzero elements it is much more efficient to create and store it as a SparseArray. In addition many matrix operations are optimized to work faster with SparseArrays.

In order to efficiently generate sparse subdiagonal matrix you can use following function:

tminus7[n_] := DiagonalMatrix[ConstantArray[1, n-1, SparseArray], -1]


Both DiagonalMatrix and ConstantArray can output SparseArrays, but you have to explicitly specify SparseArray option, otherwise a large non-sparse array will be generated at intermediate step and the performance will be much lower (as shown by tminus2 function in kglr's answer).

Let's compare this method with the fastest function from kglr's answer (tminus1):

tminus1[n_] := ArrayPad[IdentityMatrix[n - 1, SparseArray], {{1, 0}, {0, 1}}]

n=10000;
t1=First[RepeatedTiming[r1=tminus1[n];]];
t7=First[RepeatedTiming[r7=tminus7[n];]];
Equal[r1,r7]


True

{t1,t7}


{0.000072, 0.000025}

Thus, for large matrices this function is about 3 times faster than the fastest function (tminus1) from the top answer.

ClearAll[tminus1, tminus2, tminus3, tminus4, tminus5]
tminus1[n_] := ArrayPad[IdentityMatrix[n - 1, SparseArray], {{1, 0}, {0, 1}}]
tminus2[n_] := DiagonalMatrix[ConstantArray[1, n-1], -1]
tminus3[n_] := SparseArray[{i_, j_} /; j == (i - 1) -> 1, {n, n}]
tminus4[n_] := IdentityMatrix[n+1, SparseArray][[;;-2, 2;;]]
tminus5[n_] := Drop[IdentityMatrix[n+1, SparseArray], {-1}, {1}]

tminus1[5] // MatrixForm // TeXForm


$$\left( \begin{array}{ccccc} 0 & 0 & 0 & 0 & 0 \\ 1 & 0 & 0 & 0 & 0 \\ 0 & 1 & 0 & 0 & 0 \\ 0 & 0 & 1 & 0 & 0 \\ 0 & 0 & 0 & 1 & 0 \\ \end{array} \right)$$

Timings: Including tminus0 (from OP) , tminus6 (from Roman's answer) and tminus7 (from Shadowray's answer):

tminus0[n_] := Table[If[i == j + 1, 1, 0], {i, 1, n}, {j, 1, n}]
tminus6[n_] := SparseArray[Band[{2, 1}] -> 1, {n, n}]
tminus7[n_] := DiagonalMatrix[ConstantArray[1, n-1, SparseArray], -1]

funcs ={"tminus1 - arraypad", "tminus3 - sparsearray", "tminus4 - part",
"tminus5 - drop","tminus6 - sparsearray-band", "tminus7 - diagonalmatrix-sparseArray"};


Version/OS:

\$Version (on Wolfram Cloud)


12.0.0 for Linux x86 (64-bit) (April 7, 2019)

n = 100;
t0 = First[AbsoluteTiming[r0 = tminus0[n];]];
t1 = First[AbsoluteTiming[r1 = tminus1[n];]];
t2 = First[AbsoluteTiming[r2 = tminus2[n];]];
t3 = First[AbsoluteTiming[r3 = tminus3[n];]];
t4 = First[AbsoluteTiming[r4 = tminus4[n];]];
t5 = First[AbsoluteTiming[r5 = tminus5[n];]];
t6 = First[AbsoluteTiming[r6 = tminus6[n];]];
t7 = First[AbsoluteTiming[r7 = tminus7[n];]];
Equal[r0, r1, r2, r3, r4, r5, r6, r7]


True

timings = {t0, t1, t2, t3, t4, t5, t6, t7};
Column[{"n = 100",
Grid[Prepend[SortBy[Transpose[{funcs, timings}], Last],
{"function", "timing"}], Dividers -> All, Alignment -> "."]}]


I exclude tminus0 and tminus2 (due limited cloud credit) for n = 10000 and for n = 100000:

• Astonishing timings! Thanks, this comparison is very useful. Jun 12, 2019 at 15:01
• @kglr, could you add a column to your tables to show specifically which functions tminusX go with which times and function cores? I keep scrolling up and down to mentally correlate. Super useful data, thanks! Jun 12, 2019 at 17:29
• Also, which version of Mathematica used. Jun 12, 2019 at 17:31
• @MikeY, please see the update.
– kglr
Jun 12, 2019 at 18:41
• Worth adding Shadowray's fully SparseArray implementation of tminus2 to your table Jun 14, 2019 at 20:43
tminus[n_Integer /; n >= 2] := SparseArray[Band[{2, 1}] -> 1, {n, n}]