Есть ли эффективный способ вычислить 2-норму обратной матрицы (Matlab)?

Контекст: у меня есть некоторые теоретические апостериорные границы погрешности для двухсеточной конечно-элементной схемы, которые я использую для решения проблемы потери устойчивости. Однако вычисление одного из членов непомерно дорого, и мне интересно, не делаю ли я это наивно.

Проблема: Обозначив A мою (разреженную) матрицу жесткости, я вычисляю факторизацию Холецкого:

L = chol(A,'lower');

Затем мне нужно вычислить 2-норму inv(L). В настоящее время я использую «inv» и svds:

Linv = inv(L);
LinvNorm = svds(Linv,1);

Обратите внимание, что я использую inv, поскольку на самом деле это рекомендуемый синтаксис на веб-сайте Mathworks: "Для разреженных входных данных inv(X) создает разреженную идентификационную матрицу и использует обратную косую черту, X\speye(size(X)). "

Вопрос: Это, конечно, очень медленно (особенно вычисление обратного, хотя svds тоже недешев). Я пропустил трюк здесь?

Редактировать: я пытался использовать svds(L,1,0) (чтобы получить обратное значение), но это не сходится. У меня 2015a, поэтому я не могу увеличить размер Крылова, что является рекомендуемым исправлением в более новых версиях.


person Kino    schedule 22.08.2017    source источник
comment
Ничего не пропало. Вы выполняете дорогостоящую математическую операцию, она требует времени   -  person Ander Biguri    schedule 22.08.2017
comment
@AnderBiguri Я так и думал, я просто был оптимистичен - лучше спросить, может быть, Mathematics SE. В таком случае мне придется попытаться связать эти термины.   -  person Kino    schedule 22.08.2017
comment
@Kino Вы также можете задать вопрос в Computational Science SE: scicomp.stackexchange.com   -  person rayryeng    schedule 22.08.2017


Ответы (1)


Что о

LinvNorm = 1 / svds(L, 1, 0)

то есть обратное значение наименьшего единственного числа L? (Сингулярные значения обратной матрицы являются величинами, обратными сингулярным значениям исходной матрицы).

Альтернативное выражение

LinvNorm = 1 / sqrt(eigs(A, 1, 0))

or

LinvNorm = 1 / sqrt(svds(A, 1, 0))

которые происходят от A == L*L'.

person Stefano M    schedule 22.08.2017
comment
Извините, да, я пробовал это изначально, но в этом случае svds не сходятся. Я на 2015a, поэтому не могу увеличить размер пространства Крылова, что рекомендуется для этого в более новых выпусках. - person Kino; 23.08.2017
comment
@kino Вы пытались вычислить eigs(A,1,0)? Может быть, это более стабильно (и наверняка более редко). - person Stefano M; 23.08.2017
comment
@kino Если приведенного выше предложения, основанного на элементарной линейной алгебре, недостаточно, вы можете попробовать задать новый вопрос на scicomp.stackexchange.com о лучшем способе оценки ваших границ ошибки. - person Stefano M; 23.08.2017