Solution 1:

Yes, the method is completely rigorous, because here $t$ is, algebraically, an indeterminate.

If the coefficients of the matrix are supposed to be in a field $F$, then the computation you make takes place in the field $F(t)$ of rational functions in the indeterminate $t$. No problem here in considering $t^{-1}$, because $t$ is a nonzero element of the field.

To be picky, induction should be needed, but laying out the argument is very easy.

Solution 2:

The argument can be made rigorous by the following identity: let $U$ be the final triangular matrix and let $D$ be the matrix containing only a subdiagonal of $1, \cdots, 1$. Then one has

$$ \left(I - \frac{1}{t} D\right) U = (A + t I) $$ hence $\det(A+tI)= \det(I - \frac{1}{t} D) \det (U) = \det(U)$

It means that this process is actually an $LU$ decomposition of $A+tI$.

For a complete calculation, let

$$\renewcommand{\arraystretch}{2} \begin{array}{rcl}{P}_{0}&=&0\\ {P}_{i}&=&\displaystyle \sum _{p = 0}^{i-1} {a}_{p} {t}^{p} \quad i \geqslant 1\\ L&=&{\left({{\delta}}_{i}^{j}-{t}^{{-1}} {{\delta}}_{i}^{j+1}\right)}_{i , j}\\ U&=&{\left(t {{\delta}}_{i}^{j}+{t}^{1-i} {P}_{i} {{\delta}}_{j}^{n}\right)}_{i , j} \end{array}$$

Then, using ${P}_{i}-{P}_{i-1} = {a}_{i-1} {t}^{i-1}$ for $i \geqslant 1$,

$$\renewcommand{\arraystretch}{2} \begin{array}{rcl}{\left(L U\right)}_{i , j}&=&\displaystyle \sum _{k = 1}^{n} \left({{\delta}}_{i}^{k}-{t}^{{-1}} {{\delta}}_{i}^{k+1}\right) \left(t {{\delta}}_{k}^{i}+{t}^{1-k} {P}_{k} {{\delta}}_{j}^{n}\right)\\ &=&t {{\delta}}_{i}^{j}+{t}^{1-i} {P}_{i} {{\delta}}_{j}^{n}-{{\delta}}_{i}^{j+1}-{t}^{1-i} {P}_{i-1} {{\delta}}_{j}^{n}\\ &=&t {{\delta}}_{i}^{j}-{{\delta}}_{i}^{j+1}+{a}_{j-1} {{\delta}}_{j}^{n}\\ &=&{\left(A+t I\right)}_{i , j} \end{array}$$

where $\delta_i^j$ is Kronecker's delta.