求解约束线性方程组

Question

我有一个y = Ax + b形式的方程组,其中y,x和b是n×1向量,A是×n(对称)矩阵.

所以这里是皱纹.并非所有的x都是未知的.指定了某些x行,并且y的相应行未知.以下是一个例子

| 10  |   |  5  -2  1 | | * |   | -1 |
|  *  | = | -2   2  0 | | 1 | + |  1 |
|  1  |   |  1   0  1 | | * |   |  2 |

其中,*指定未知数.

我已经为Fortran中的上述问题构建了一个求解器,但是我想知道是否有一个不错的鲁棒求解器作为Lapack或MLK的一部分用于这些类型的问题？

---

我的求解器基于一个排序矩阵,称为根据已知和未知pivot = [1,3,2]重新排列x和y向量

| 10  |   |  5  1 -2 | | * |   | -1 |
|  1  |   |  1  1  0 | | * | + |  2 |
|  *  |   | -2  0  2 | | 1 |   |  1 |

并使用块矩阵解和LU分解求解

! solves a n×n system of equations where k values are known from the 'x' vector
function solve_linear_system(A,b,x_known,y_known,pivot,n,k) result(x)
use lu
integer(c_int),intent(in) :: n, k, pivot(n)
real(c_double),intent(in) :: A(n,n), b(n), x_known(k), y_known(n-k)
real(c_double) :: x(n), y(n), r(n-k), A1(n-k,n-k), A3(n-k,k), b1(n-k)
integer(c_int) :: i, j, u, code, d, indx(n-k)

    u = n-k
    !store known `x` and `y` values
    x(pivot(u+1:n)) = x_known
    y(pivot(1:u)) = y_known

    !define block matrices
    ! |y_known| = | A1  A3 | |    *    | + |b1|
    | |   *   | = | A3` A2 | | x_known |   |b2|

    A1 = A(pivot(1:u), pivot(1:u))
    A3 = A(pivot(1:u), pivot(u+1:n))
    b1 = b(pivot(1:u))

    !define new rhs vector
    r = y_known -matmul(A3, x_known)-b1

    % solve `A1*x=r` with LU decomposition from NR book for 'x'
    call ludcmp(A1,u,indx,d,code)
    call lubksb(A1,u,indx,r)

    % store unknown 'x' values (stored into 'r' by 'lubksb')
    x(pivot(1:u)) = r

end function

对于上面的例子,解决方案是

    | 10.0 |        |  3.5 | 
y = | -4.0 |    x = |  1.0 |
    |  1.0 |        | -4.5 |

PS.线性系统通常具有n<=20方程式.

Answer 1

只有未知数的问题是线性最小二乘问题。

您的先验知识可以通过等式约束（固定一些变量）引入，将其转换为线性等式约束最小二乘问题。

lapack 中确实有一个算法可以解决后者，称为xGGLSE。

以下是一些概述。

（看来，在您的情况下，您需要将 b 与 -1 相乘才能与定义兼容）

编辑：经过进一步检查，我错过了y. 哎哟。这不好。