matrix - sistemas - Resolviendo un sistema restringido de ecuaciones lineales

resolver sistemas de ecuaciones en r

Question

Tengo un sistema de ecuaciones de la forma y = Ax + b donde y , x y b son n × 1 vectores y A es una matriz × n (simétrica).

Entonces aquí está la arruga. No todo de x es desconocido . Se especifican ciertas filas de x y las filas correspondientes de y son desconocidas. A continuación hay un ejemplo

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

donde * designa cantidades desconocidas.

He desarrollado un solucionador para problemas como el anterior en Fortran, pero quería saber si existe un solucionador robusto decente como parte de Lapack o MLK para este tipo de problemas.

Mi solucionador se basa en una matriz de clasificación llamada pivot = [1,3,2] que reorganiza los vectores xey de acuerdo con lo conocido y lo desconocido.

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

y la solución usando una solución de matriz de bloques y descomposición de 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

Para el ejemplo de arriba, la solución es

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

PD. Los sistemas lineales tienen típicamente n<=20 ecuaciones.

Answer 1

El problema con solo incógnitas es un problema lineal de mínimos cuadrados .

Su conocimiento a priori puede ser introducido con restricciones de igualdad (fijando algunas variables), transformándolo en un problema lineal de mínimos cuadrados restringidos por la igualdad .

De hecho, hay un algoritmo dentro de lapack que soluciona este último, llamado xGGLSE .

Aquí hay una descripción general .

(También parece que necesita multiplicar b por -1 en su caso para que sea compatible con la definición)

Editar: En otra inspección, me perdí las incógnitas dentro de y . Ay. Esto es malo.

Answer 2

Primero, volvería a escribir su sistema en una forma AX = b donde A y b son conocidos. En su ejemplo, y siempre que no haya cometido ningún error, se daría:

    5 0 1     x1         13
A = 2 1 0 X = x2 and b =  3
    1 0 1     x3         -1

Luego puede usar muchos métodos provenientes de varias bibliotecas, como LAPACK o BLAS, dependiendo de las propiedades de su matriz A (positiva definida, ...). Como punto de partida, sugeriría un método simple con una inversión directa de la matriz A, especialmente si su matriz es pequeña. También hay muchos enfoques iterativos (Jacobi, Gradients, Gauss seidel ...) que puede usar para casos más grandes.

Editar : una idea para resolverlo en 2 pasos

Primer paso: puede reescribir su sistema en 2 subsistemas que tienen X e Y como desconocidos, pero las dimensiones son iguales a los números de incógnitas en cada vector. El primer subsistema en X será AX = b que se puede resolver por métodos directos o iterativos.

Segundo paso: el segundo sistema en Y se puede resolver directamente una vez que se conoce X porque Y se expresará en la forma Y = A''X + b ''

Creo que este enfoque es más general.