2014 dxdy logo

Научный форум dxdy

Математика, Физика, Computer Science, Machine Learning, LaTeX, Механика и Техника, Химия,
Биология и Медицина, Экономика и Финансовая Математика, Гуманитарные науки




Начать новую тему Ответить на тему
 
 Нахождение собственных ВЕКТОРОВ
Сообщение30.10.2006, 03:04 


30/10/06
33
Имеется симметричная действительная неотрицательно определенная матрица с габаритами 512х512. Нужен самый быстрый алгоритм для нахождения первого десятка ее собственных векторов, соответствующих ее самым СТАРШИМ собственным значениям.
Какие будут советы по алгоритму?
Что пробовал делать: приводил к трехдиагональной форме по Хаусхолдеру (включая вычисление матрицы перехода), затем QL-алгоритмом находил ВСЕ собственные значения (соб. вектора при этом не вычислялись), полученные собственные значения сортировались в порядке убывания, и из списка выбрался первый десяток. Затем повторял QL-алгоритм, но уже с вычислением векторов, используя ранее найденные соб. значения в качестве сдвигов (в монотонно убывающем порядке). Ничего хорошего не получилось – на первом же расщеплении на приводимые части (ноль в середине поддиагонали) алгоритм срывался.
Проводить QL-итерации сразу с векторами очень накладно – именно на это тратится основное время вычислений. Вычисление даже всех собственных значений без поворотов векторов на каждой итерации по скорости вполне приемлемо. Но повороты векторов требуют много умножений, в то время как из 512 мне нужны только 10.
Согласно теории, использование соб. значения в качестве сдвига должно приводить к сходимости QL (или QR) алгоритма за одну (!) итерацию. На этом я и хотел выиграть в скорости. Однако данное утверждение относится только к неразложимой матрице.
На практике не вышло ничего даже от "подсказки" сдвигов, сделанных без сортировки. Сходимость наблюдается, но часто не к тому значению, к которому ожидал, особенно это касается малых собственных чисел. К сожалению, моя матрица 512x512 далеко не полного ранга. И те 10 значений представляют максимальный ранг, на который я могу надеяться. Остальные соб.значения – следствие "шумов" в данных.

 Профиль  
                  
 
 
Сообщение30.10.2006, 07:31 


20/12/05
31
Я бы вам посоветовал почитать книжки по вычислительным методам линейной алгебры, та проблема о которой вы говорите называется частичной проблемой собственных значений.
Посмотрите, например, Д.К. Фадеев Численные методы линейной алгебры. Или любую подобную книгу.

 Профиль  
                  
 
 
Сообщение30.10.2006, 10:59 


30/10/06
33
zhe
Книжки читаем :-). Однако моя задача вовсе не "частичная проблема собственных значений", поскольку собственные значения я знаю ВСЕ. Моя проблема именно в том, чтобы построить собственные ВЕКТОРА (а не значения!) для определенных соб.значений, которые я сам выберу из полного набора.
В книжках заморачиваются только поиском соб. значений, забывая о том, что каждая итерация требует синхронной перестройки матрицы вращения, которая в итоге должна превратиться в набор собственных векторов.
Вести поиск (итерации) соб. значений, одновременно отслеживая результирующую матрицу поворота, с вычислительной точки зрения самоубийственно. Это все равно, что искать грибы, таская за собой железнодорожный вагон для их складывая :-). Потому-то и хочется сначала найти грибы налегке, а уж только потом подогнать в нужное место транспорт. Причем, кратчайшим путем. И минуя найденные мухоморы и поганки :-). Тем более, что у меня и лес такой, что на 10 съедобных грибов приходится 502 поганых. Если искать соб. значения налегке, то почти без разницы каким алгоритмом пользоваться – "энергозатраты" при этом малы. Поэтому проще не мучиться с бисекциями, а вычислить сразу полный набор, т.к. такой алгоритм наиболее устойчив и гарантирует, что кратные соб. значения не будут приняты за одно.
Очень хотелось бы получить прямую наводку на алгоритм (без его описания) на способ получения собственного вектора для отдельного (строго указанного!) соб. значения. Ответы "знатоков" типа "книжки читай – в них все написано" мне не нужны.

 Профиль  
                  
 
 
Сообщение31.10.2006, 05:16 
Заслуженный участник


15/05/05
3445
USA
На этой странице: http://www.srcc.msu.su/num_anal/lib_na/cat/cat571.htm есть ссылка на алгоритм вычисления собственных значений, принадлежащих заданному интервалу, их номеров и соответствующих собственных векторов вещественной симметричной матрицы. Там еще упоминаются алгоритмы вычисления а) нескольких экстремальных СЗ и СВ и б)нескольких крайних СЗ и соответствующих СВ методом Ланцоша.

Это не совсем то, что Вы ищете, но может быть натолкнет Вас на какие-нибудь мысли.

 Профиль  
                  
 
 
Сообщение31.10.2006, 14:35 


30/10/06
33
Yuri Gendelman
За участие спасибо. С того я и начинал, что первым делом попытался ознакомиться с процедурами, включенными в известные библиотечные пакеты. На форум я гораздо позже написал, когда уже понял, что во всех библиотеках имеется, по существу, один и тот же стандартный набор возможностей. И на определенные мысли это меня, конечно, натолкнуло. А мысли вот такие.
Обширные пакеты программ создают в "общественном сознании" почти непреодолимое чувство, что на все вопросы уже найдены ответы. Кстати, точно то же самое впечатление возникает в головах обывателей по отношении к науке вообще. И лишь личная углубленность в какую-либо проблему позволяет осознать, насколько поверхностно то знание, которое считается непреложным.
Из учебника в учебник, из книги в книгу, из статьи в статью почти в буквальном пересказе перекочевывает материал, ставший классическим. Несомненно, что образовательным целям такая практика целиком отвечает. Однако, имеет вредный побочный эффект, выражающийся в замазывании подводных камней и рифов в алгоритмах. Доверяясь громкому имени капитана, пассажир считает все форватеры проходимыми. А зря!
Обращение к первоисточникам (в первую очередь к "СПРАВОЧНИК алгоритмов на языке Алгол. Линейная алгебра" под ред. Уилкинсона и Райнша, вышедшую в русском переводе еще в 1960 г.) показывает, что все нынешние алгоритмы взяты именно оттуда.
Математики придумывали алгоритмы на бумаге, а программисты затем перекладывали их на алгоритмические языки. Персональных компьютеров, я полагаю, тогда не было. А если и были ЭВМ, то авторы алгоритмов, по-видимому, к ним не обращались :-). Сейчас, когда за шагами итерации можно удобно наблюдать на экране, многие утверждения авторов алгоритмов уже не кажутся такими непреложными.
Например, по утверждению Парлетта (Parlett B.N., The Symmetric Eigenvalue Problem, 1980) QL- алгоритм со сдвигом, равным собственному значению, сходится за одну итерацию (при точной арифметике). И он подсознательно полагал, что сойдется он именно к этому сдвигу. А вот наблюдение за процессом итерации на персоналке показывает, что процесс может сойтись вовсе не к тому СЗ, которое было выбрано в качестве сдвига, а к соседнему, если оно встретится на пути итерации раньше.
Величина первоначального сдвига (с которого начинается итерация) выбирается по оценке Уилкинсона, которая гарантирует сходимость. Но то, что итерация сойдется к СЗ, ближайшему к этой оценке, в общем случае неверно. Из-за этого мне не удается направить сходимость к желаемому мною СЗ. И все из-за того, что по дороге от краевого диагонального элемента к нужному мне значению СЗ могут встречаться (и почти всегда встречаются!) другие, не нужные мне СЗ. Такое почти всегда встречается, когда краевой элемент достаточно мал. И тогда на пути к максимальному СЗ будет суждено пройти почти через все остальные СЗ. И тут при каждом перешагивании может случиться так, что на кого-то наступишь. Тогда элемент на поддиагонали сразу же обнуляется и итерация вынужденно прекращается на том СЗ, на которое нечаянно наступил. Вот и получай тогда собственный вектор, которого не заказывал.
С бисекциями (алгоритм Штурма) тоже большая закавыка. Тут надо обязательно указывать ИНТЕРВАЛ, на котором ищешь СЗ. Но получается так, что невозможно указать такого интервала, чтобы в нем гарантированно было ровно 10 штук СЗ и не одной ни больше, и ни меньше. Ведь расстояния между соседними СЗ могут быть сколь угодно малы, не говоря уже о неточности их определения.
Хотя я все сильнее склоняюсь к тому, что использовать надо этот алгоритм, а не первый. Т.е. уповать на процедуру типа TINVIT и действовать в духе статьи http://nadav.harel.org.il/papers/eigen.ps.gz
В общем, вопросов больше, чем ответов. А главное не видно никого, кто бы в эти вопросы захотел вникнуть – настолько велик всеобщий гипноз готовых алгоритмических решений.

 Профиль  
                  
 
 
Сообщение02.11.2006, 02:01 
Заслуженный участник


15/05/05
3445
USA
Oam писал(а):
Обширные пакеты программ создают в "общественном сознании" почти непреодолимое чувство, что на все вопросы уже найдены ответы.

По-моему Вы несколько преувеличиваете трагизмситуации. Стандартные пакеты покрывают 99% всех потребностей. Для нестандартных задач создаются специализированные алгоритмы. Но они обычно не оформляются как стандартные и "рассеяны" по приложениям.

Oam писал(а):
Из учебника в учебник ... почти в буквальном пересказе перекочевывает материал, ставший классическим. ... такая практика ... имеет вредный побочный эффект, выражающийся в замазывании подводных камней и рифов в алгоритмах.

Согласен насчет "перекочевывания" и считаю это естественным. Но не согласен насчет замазывания. Другое дело, что не всегда документацию читают, особенно к статистическим пакетам.

Oam писал(а):
А если и были ЭВМ, то авторы алгоритмов, по-видимому, к ним не обращались :-).

Вы видимо младше, чем я думал :-)

Oam писал(а):
В общем, вопросов больше, чем ответов. А главное не видно никого, кто бы в эти вопросы захотел вникнуть – настолько велик всеобщий гипноз готовых алгоритмических решений.

Дело не в гипнозе. Вам нужен алгоритм - кто больше Вас должен вникать? Придумайте алгоритм и опубликуйте его.

Техническое замечание:
"СПРАВОЧНИК алгоритмов на языке Алгол. Линейная алгебра" Уилкинсона и Райнша вышел на русском в 1976 (оригинал – в 1971). Т.е. переводу не 46 лет, а всего 30 :-) .
Более ранняя монография: Уилкинсон Дж.Х. Алгебраическая проблема собственных значений вышла по-русски в 1970 (оригинал – в 1965).
Уилкинсон в этом деле (СЗ-СВ) голова. Потому и алгоритмы его встречаются в всех книгах.

 Профиль  
                  
 
 
Сообщение07.11.2006, 20:51 
Заслуженный участник
Аватара пользователя


31/10/06
371
РФ, РК, г.Симферополь
Может я поздновато зашел.
Из всего, что прочел, я понял только первые три предложения...
А разве нельзя использовать какие-нибудь модификации степенного метода, которые будут одновременно давать и несколько первых собственных значений и соответствующие им собственные векторы? Ведь это самые быстрые методы.

 Профиль  
                  
 
 QL - алгоритм
Сообщение27.01.2009, 17:33 


27/01/09
4
Ивановская область
Если есть у кого нибудь QL-алгоритм с преобразованием Хаусхолдера очень прошу помочь.Желательно исходник с main().

 Профиль  
                  
 
 Re: Нахождение собственных ВЕКТОРОВ
Сообщение05.02.2009, 11:18 
Заблокирован
Аватара пользователя


13/01/09

335
Oam писал(а):
Имеется симметричная ...

Все ответы на ваши вопросы есть в пакете Lapack (в свободном доступе) или в его быстрой реализации Intel MKL.

 Профиль  
                  
 
 
Сообщение05.02.2009, 17:43 


11/07/06
201
allsolovey в сообщении #181702 писал(а):
Если есть у кого нибудь QL-алгоритм с преобразованием Хаусхолдера очень прошу помочь.Желательно исходник с main().


Посмотрите "Numerical recipies in C". Там точно есть исходники
QR-алгоритма с двойным неявным сдвигом, возможно и с одним сдвигом есть.

 Профиль  
                  
Показать сообщения за:  Поле сортировки  
Начать новую тему Ответить на тему  [ Сообщений: 10 ] 

Модераторы: Karan, Toucan, PAV, maxal, Супермодераторы



Кто сейчас на конференции

Сейчас этот форум просматривают: нет зарегистрированных пользователей


Вы не можете начинать темы
Вы не можете отвечать на сообщения
Вы не можете редактировать свои сообщения
Вы не можете удалять свои сообщения
Вы не можете добавлять вложения

Найти:
Powered by phpBB © 2000, 2002, 2005, 2007 phpBB Group