Вызов ретрансляционной программы ARPACK с обратным сообщением

Question

Я пытаюсь написать драйвер на C++ для вычисления собственных значений для асимметричной разрешающей матрицы с реальными значениями, используя функции fortran, предлагаемые ARPACK, но у меня есть немного проблема с обратным коммуникационным подходом.

Вообще, я пытаюсь решить нормальное собственное уравнение:

A*v = lambda*v 

и любое взаимодействие с матрицей А делается в ARPACK с помощью функции «ау»:

av(n, workd[ipntr[0]], workd[ipntr[1]]) 

, который умножает вектор, удерживаемый в массиве «workd», начинающийся с местоположения «ipntr [0]», и вставляет результат в массив «workd», начинающийся с местоположения «ipntr [1]». Примеры такого подхода приведены в руководстве по адресу http://www.caam.rice.edu/software/ARPACK/, а также в коде ARPACK/ПРИМЕР/ПРОСТО/dnsimp.f.

Что я хотел бы знать, так это как я могу связать матрицу A? Если он не передается в подпрограмму, то как можно найти его действие на предоставленном векторе?

В примере кода dnsimp.f их матрица A вычисляется внутри функции 'av' и является «полученной из стандартной центральной дискретизации развязки двумерного оператора конвекции-диффузии». Однако, я считаю, что это проблема конкретной? Также не кажется слишком полезным кодировать вывод матрицы A в функцию. Я не могу найти много информации об этом из руководства.

Это, кажется, не слишком большая проблема, поскольку поскольку это определенная пользователем функция, я могу просто изменить определение «av», чтобы включить матрицу A в качестве параметра. Однако я хотел бы знать, как это делается правильно в случае каких-либо потенциальных проблем совместимости.

Спасибо!

Answer 1

Вам не нужно поставлять матрицу ARPACK.

Все, что вам нужно сделать, это умножить матрицу на возвращаемые векторы (таким образом, обратную связь) до достижения желаемой конвергенции.

Для получения информации об алгоритмах, вы должны взглянуть на users guide и особенно на chapter about the underlying algorithms.

Ответ на комментарий: Основной алгоритм является формой Арнольди Итерация. Основной алгоритм показан в wikipedia и показывает, что матрица A не будет доступна. Ни прямо, ни косвенно.

В частности, алгоритмы начинаются с произвольного нормированного вектора q_1. Этот вектор возвращается пользователю. Пользователь умножает этот вектор на матрицу A, используя свою любимую процедуру (обычно некоторое эффективное разреженное преобразование матрицы-вектора), и возвращает результат в Arnoldi Iteration для вычисления части матрицы Гессенберга H (собственные значения которой обычно сходятся к экстремальные собственные значения А) и следующий вектор q_2. Это нужно повторить, пока ваши результаты не будут сведены.