Эффективное решение линейной системы Ax= b при изменении только одного постоянного члена

Как эффективно решить большую систему линейных уравнений, если меняются лишь некоторые постоянные члены? Например:

В настоящее время у меня есть система Ax= b. Я вычисляю инверсию A один раз, сохраняю ее в матрице, и каждый раз, когда любая запись обновляется в b, я выполняю умножение матрицы на вектор A^-1(b) для повторного вычисления x.

Это неэффективно, так как только несколько записей будут обновлены в b. Существуют ли более эффективные способы решения этой системы, когда A-1 остается постоянным, но конкретные известные значения b изменяются?

Я использую uBlas и Eigen, но не знаю решений, которые решат эту проблему выборочного пересчета. Спасибо за любое руководство.


c++ eigen linear-algebra blas ublas
person Atlas    schedule 01.11.2012    source источник

Ответы (4)


arrow_upward
3
arrow_downward

Вычислите A^-1. Если b_i является i-м компонентом b, то исследуйте d/db_i A^-1 b (производную A^-1 по i-му компоненту b) -- он равен столбцу A^-1 (в частности, i-му столбцу). А производные линейных функций постоянны в своей области определения. Итак, если у вас есть b и b', и они отличаются только i-й составляющей, то A^-1 b - A^-1 b' = [d/db_i A^-1] * (b-b')_i. Для нескольких компонентов просто сложите их (поскольку A^-1 является линейным).

Или, короче говоря, вы можете вычислить A^-1 (b'-b) с некоторыми оптимизациями для входных компонентов, равных нулю (которые, если только некоторые компоненты изменятся, будут большинством компонентов). A^-1 b' = A^-1 (b'-b) + A^-1 (b). И если вы знаете, что изменятся только некоторые компоненты, вы можете взять копию соответствующего столбца A^-1, а затем умножить ее на изменение этого компонента b.

person Yakk - Adam Nevraumont    schedule 01.11.2012

arrow_upward
1
arrow_downward

Вы можете воспользоваться линейностью задачи:

x0 = A_(-1)*b0 
x = A_(-1)*b = x0 + A_(-1)*db

где db — это разностная матрица между b и b0, и она должна быть заполнена нулем: вы можете сжать ее до разреженной матрицы.

в Eigen lib есть много интересных функций для разреженных матриц (умножение, инверсия, ...).

person lucasg    schedule 01.11.2012

arrow_upward
1
arrow_downward

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

Библиотека выполнит некоторую декомпозицию, например LU, и использует ее для вычисления x. Если вы выберете итеративный решатель, то он уже делает в значительной степени то, что вы описываете, чтобы сосредоточиться на решении; он начнет с худшего предположения и сгенерирует лучшее, и любая хорошая процедура использует начальное предположение, чтобы ускорить процесс. Во многих случаях у вас в любом случае есть хорошее представление о результате, поэтому имеет смысл использовать это.

Если новый b находится рядом со старым b, то новый x должен быть рядом со старым x, и это послужит хорошей начальной догадкой.

person Phil H    schedule 01.11.2012

arrow_upward
0
arrow_downward

Во-первых, не вычисляйте обратную матрицу, а используйте LU-разложение или QR-разложение (медленнее, чем LU, но стабильнее). Такие разложения масштабируются лучше, чем инверсия, с размером матрицы и обычно более стабильны (особенно QR).

Есть способы обновить QR-разложение, если A немного изменится (например, с помощью матрицы первого ранга), но если B изменится, вам придется снова решать с новым b - вы не можете избежать этого, и это O (n ^ 2).

Однако, если правая часть B изменится только на фиксированный элемент, т.е. B' = B + дБ с заранее известными дБ, вы можете раз и навсегда решить A dx = дБ, и теперь решение x' для Ax' = B' равно x + dX.

Если dB не известен заранее, но всегда является линейной комбинацией нескольких векторов dB_i, вы можете решить для A dx_i = dB_i, но если у вас много таких dB_i, вы получите процесс ^2 (фактически это равносильно вычисление обратного)...

person Alexandre C.    schedule 01.11.2012