Kirill_SalВы сделали много правильных и полезных вещей. Перенумеровали частицы так, что неподвижные крайние имеют номера 

 и 

 (я сделал бы точно так же). Заменили координату на смещение относительно положения равновесия. Ввели обозначение 

. А вот идея использовать здесь цилиндрические функции мне кажется малопродуктивной. Функции Бесселя и Неймана с хитрым нецелым и неполуцелым порядком. Что можно с ними сделать? Разложить в ряд, записать интегральное представление? Самое лучшее — это определить их рекуррентной формулой, что возвращает нас к исходной точке.
Обозначим ещё 

. Амплитуды буду обозначать 

, количество подвижных частиц 

. Тогда уравнения примут вид:

Запишем систему в матричной форме. Чем писать для произвольного 

 c многоточиями, будет нагляднее, если я возьму небольшое конкретное 

, например 

.

Обратите внимание, что в такой записи из системы исчезли какие-либо параметры (вроде массы или жесткости), кроме разве что количества частиц 

.
У нас получилась обобщённая задача о собственных значениях/векторах 
.Это однородная система. Она имеет нетривиальное решение, только если определитель матрицы равен нулю. Это и дает уравнение для определения спектра:

Обозначим такой определитель 

-го порядка через 

. Так как 

 имеет специальный вид, можно продвинуться в его вычислении. Разложив 

 по последней строке, получим (Вы это легко сами проделаете) рекуррентную формулу

При этом 

 и по определению 

.
Итак, определители для различных 

 связаны рекуррентной последовательностью.

 представляет собой полином 

-й степени от 

, поэтому уравнение 

 имеет 

 корней 

 с учетом кратности. Каждый 

 подставляйте в систему и находите соответствующий ему собственный вектор 

.
Хочу подчеркнуть, что собственные векторы, соответствующие различным собственным значениям, будут совсем непохожи друг на друга. Мне показалось, что Вы забыли о том, что здесь будет 

 различных по характеру колебательных мод.