I had written an algorithm that searches for the optimal weight parameter to be implemented in the successive-over relaxation (SOR) method which worked cleanly by vectorizing the interval and for each ω the spectral radius of the iteration matrix is computed.
However, I was advised not to use this approach for large sparse matrices as it is expensive to compute (the same way computing condition number of a large matrix is unfeasible) and rather use it as a demonstration tool. Therefore, I was wondering what strategy is the best to approximate the optimal weight parameter for large sparse systems () that would allow the best convergence of the SOR which is given as follow starting by decomposing $A$ to lower, diagonal, and upper triangular matrix $A=D-L-U$. Thus, for a weight parameter $\omega$, the SOR is given by : $$ \left(\frac{1}{\omega}D+L\right)x_{k+1}=\left(1-\frac{1}{\omega}\right)Dx_{k}+Ux_{k} $$
Furthermore, as a result of my question I was wondering if classical iterative stationary methods such as Jacobi, Gauss-Seidel, and the SOR are worthy to be used nowadays in dealing with large sparse systems or are they outdated and the default preference now are Krylov methods?