The restriction on the eigenvalues is usually given in the form of a condition number. This is the $\kappa$ that you see in all the runtimes in your table. $\kappa = |\lambda_{\rm{max}}/\lambda_{\rm{min}}|$ where $\lambda_{\rm{max}}$ and $\lambda_{\rm{min}}$ are the maximum and minimum eigenvalues respectively.
In all runtimes listed in your table, it is assumed that the condition number is known. One does not usually think of "calculating the condition number" as part of the algorithm for solving $Ax=b$, for example. If the condition number is larger, the system is harder to solve, and if it is smaller the system is easier to solve (assuming all other parameters, including the maximum desired error $\epsilon$ are held fixed).
In terms of needing to know that $\lambda_{\rm{max}} < M$ and $\lambda_{\rm{min}}>L$, there are lots of examples where we can know the bounds on the eigenvalues without actually going through the effort of calculating the eigenvalues. In this way, HHL can be a great way to find the state you're looking for, without the cost of calculating the condition number or any eigenvalues.
Let me give just one real-world example. Let's say I want to find the molecular vibrational state $|\psi\rangle$ such that after $t=10$ps of evolving under its Hamiltonian $H$, the molecule ends up in state $|b\rangle $. This can be described by the equation:
$$
e^{-iHt}|\psi\rangle = |b\rangle
$$
where the $|\psi\rangle$ satisfying this equation is what you want to know. You can find your desired $|\psi\rangle$ by using the HHL algorithm with $A = e^{-\frac{i}{\hbar}Ht}$ and $|\psi\rangle = |x\rangle$.
Obtaining the smallest and largest eigenvalues of a molecular Hamiltonian to arbitrary precision is extremely costly on a classical compter, but knowing that they lie within the range $(L,M)$ can be determined with no cost at all. For example if the molecule is the nitrogen dimer we know the lowest and highest vibrational states have energies (eigenvalues) between 0 and 10 eV and since $e^{0}=1$ we have $L=1$ and $M = e^{-\frac{i}{\hbar} 10 \rm{eV} \cdot 10 \rm{ps}}$. You can convert eV to Hz, and ps to seconds to evaluate $M$ numerically, and then you can obtain the lower and upper bounds that you need to use when scaling your matrix the way you described in your previous question. At no point did I need to calculate the eigenvalues of a 14-electron molecular Hamiltonian (which would be extremely hard and would defeat the purpose of using HHL, because if I could calculate the eigenvalues I could just calculate $A$ and invert it to get $|\psi\rangle$). I just used the dissociation energy of the molecule to come up with the bounds on its vibrational energies. I could have come up with even better bounds by using the semi-classical WKB approximation, also with much less cost than actually calculating the eigenvalues, but the first example is already enough.
So now let's address all your individual questions:
First group of questions: I read plenty of papers on HHL and none of
them even mentioned this restriction. Why? Is this restriction known
but considered weak (i.e. it's easy to have this kind of information)?
Or the restriction was not known? Is there any research paper that
mention this restriction?
Out of the 539 papers that have (at present) cited the original HHL paper, many of them will not know the finer details like the dependence of its performance on the condition number or eigenvalues. Some of the papers will certainly know that the performance of the algorithm will depend on the condition number or eigenvalues of the matrix, namely, the papers listed in your table on improvements to the HHL algorithm. Robin Kothari also mentioned it, for example, at the very beginning of his talk in 2016 on the CKS algorithm (which is mentioned in your table).
Second group of questions: Is there a better algorithm in term of
complexity? If not, then why is the HHL algorithm still presented as
an exponential improvement over classical algorithms?
The algorithm you mention, suggested by DaftWulie, to estimate the bounds on the eigenvalues, is not going to be improved over $\mathcal{O}(\sqrt{N})$ because the dominant cost in that algorithm is in searching through all $N$ rows for the maximum and minimum values. The cost of everything else is small because the matrix is assumed to have a sparsity of $s \lll N$. There is no way to do this search faster in faster than $\mathcal{O}(\sqrt{N})$ time (unless you have some other extra knowledge of the system) because Grover's algorithm has been proven to be optimal.
You are right, people should mention the caveats of algorithms more often in their papers. In terms of your specific question "why is the HHL algorithm still presented as an exponential improvement over classical algorithms," I think the original authors HHL did do their due diligence in explaining the algorithm and its caveats, in that they said that there's an exponential scaling but the cost grows quadratically with the condition number and sparsity and inversely with the size of the error you are willing to tolerate. Why do most other people after HHL not mention all the caveats? Well many of them don't know the caveats, and those that do might have felt it wasn't necessary because calculating the condition number is not part of the algorithm. Knowing the condition number will tell you how well the algorithm will work, but it is assumed you already know this like in the molecular vibrations example I gave above!
