Почему такая подстановка убирает сингулярность?

[10110111]

Есть уравнение:

$$\frac{d^2}{d\rho^2}R(\rho)+\frac5\rho\frac d{d\rho}R(\rho)+\frac{4L}{\rho^2} R(\rho)+2ER(\rho)-\frac2\rho U R(\rho)=0,$$

где $U$ и $L$ — заданные постоянные матрицы, а $R$ — вектор. Краевые условия: $R(0)$ конечно, $R(+\infty)=0$. $E$ — собственное число.
Это задача из вот этой статьи. Около $\rho=0$ функция $R(\rho)$ имеет сингулярности в форме ряда из $\rho^a(\ln\rho)^b$ (разложение Бартлетта-Фока). Поэтому авторы статьи предлагают для нахождения $E$ вычислять вместо $R$ матрицу $F(\rho)$, определённую таким образом:

$$\rho\frac d{d\rho}R(\rho)=F(\rho)R(\rho).$$

После подстановки получаем новое уравнение:

$$\rho\frac{dF(\rho)}{d\rho}+(F(\rho))^2+4F(\rho)+4L+2E\rho^2-2\rho U(\rho)=0,$$

автоматически дающее и начальное условие:

$$(F(0))^2+4F(0)+4L=0.$$

Из этих уравнений для $F(\rho)$ можно вычислить $E$, ища значения, при которых полюс функции $F$ уходит на бесконечность (методом стрельбы). Полюс этой функции фактически соответствуют нулю функции $R$, поскольку определение $F$ очень похоже на определение логарифмической производной (с точностью до множителя $\rho$), поэтому так мы находим фактически $E$, которое позволяет $R$ удовлетворить краевому условию на бесконечности.

Мне вроде даже удалось воспроизвести результаты из статьи (правда, другой, аналогичной — для S-состояний, для простоты) по вычислению $E$, но возникает один большой вопрос: почему это вообще работает? $F$, судя по статье, аналитична около $\rho=0$. Почему использование (аналога) логарифмической производной убирает сингулярность в нуле? Ведь если найти логарифмическую производную от $\rho^a(\ln\rho)^b$, то она определённо не будет аналитичной в нуле, даже после умножения на $\rho$.

И ещё вопрос: а как можно это уравнение решить (численно, конечно) сразу для $R$, минуя вообще всякие $F$? Тут проблема в том, что с какой бы точки ни начинать решение, везде требуется знать мало того, что значение функции, так ещё и её производную (ибо уравнение второго порядка). Я пытался начать решение где-то в далеке от нуля, устанавливая там значение функции $0$, но результат оказывается сильно зависимым от соотношения производных от компонент вектора $R$, так что похоже, что этот вариант плох...