In the second-order Tikhonov regularization approach
$$ \min \left\| G m - d \right\|_2^2 + \alpha^2 \left\| \Gamma x \right\|_2^2 \tag{1} $$
given that $\Gamma\ $ contains second order derivatives, $f_i $ are the filter factors, $G=U \Lambda\ X^T $ and $\Gamma\ =V M X^T $ and $Y=X^{-T} $ the solution for model parameter $m$ is
$$m_{\alpha\ \Gamma\ }=\sum_{i=1}^n f_i \dfrac{U_{.,i+k}^T d}{\lambda\ _i} Y_{.,I} \tag{2}$$
where $k=n-m$ when $m \leq n$, and $m$ and $n$ are dimensions of $G$.
In my case I have a $m\times n=16\times84$ matrix of $G$, so $U$ should be a $16\times 16$ matrix.
Now, the question is how to calculate $m_{\alpha\ \Gamma\ }$ from equation (2) while the index $i+k$ is always greater than the dimension of $U^T$?!
Are you sure that the solution for model parameter is given with that form? I'd suppose it is only if $\Gamma = \mathbb{1}$, so it applies to the zeroth-order Tikhonov regularisation.
In fact, you should rather solve the normal equation, so rather $m_{\alpha} = (G^T G + \alpha^2 \Gamma^T \Gamma)^{-1} G^T d$. If your approach is correct, it is new to me.