Todos estamos familiarizados con los muchos métodos computacionales para resolver el sistema lineal estándar.
Sin embargo, tengo curiosidad si hay algún método computacional "estándar" para resolver un sistema lineal más general (dimensión finita) de la forma
donde, por ejemplo, es una matriz , es una matriz , y es un operador lineal que toma matrices matrices , lo que no implica vectorización las matrices , es decir, convertir todo a la forma estándar . m 1 × n 1 B m 2 × n 2 L m 1 × n 1 m 2 × n 2
La razón por la que pregunto es que necesito resolver la siguiente ecuación para :
donde es la transformación de radón 2d, es contigua, y tanto como son matrices (imágenes) 2d. Es posible vectorizar esta ecuación, pero es un dolor, especialmente si vamos a 3D.R ∗ u f
Más en general, ¿qué pasa con las matrices ? Por ejemplo, resolver donde y son matrices 3D (necesitaré hacer esto con Radon transform en algún punto también).L A = B A B
Gracias de antemano, y siéntase libre de enviarme a otro StackExchange si lo necesita.
linear-algebra
icurays1
fuente
fuente
Respuestas:
Una cosa que debe tener cuidado al implementar CG (o enfoques iterativos similares) con operadores lineales generales es implementar el adjunto de su operador lineal correctamente. Es decir, las personas a menudo obtienen correctamente, pero cometen un error al implementar .z = F ∗ ( y )y= F( x ) z= F∗( y)
Recomiendo implementar una prueba simple que aproveche la siguiente identidad: para cualquier e conforme , Entonces, lo que debe hacer es generar valores aleatorios de e , ejecutarlos a través de sus operaciones directas y adjuntas, respectivamente, y calcular los dos productos internos anteriores. Asegúrese de que coincidan con una precisión razonable y repita varias veces.y ⟨ y , F ( x ) ⟩ = ⟨ F * ( y ) , x ⟩ . x yX y
EDITAR: ¿qué haces si se supone que tu operador lineal es simétrico? Bueno, también necesitas verificar esa simetría. Entonces, use la misma prueba, solo observando que --- aplica la misma operación a e . Por supuesto, el OP tiene tanto un operador asimétrico como uno simétrico para tratar ... x yF= F∗ X y
fuente
linop_test.m
rutina por este motivo. Ese paquete también admite matrices, matrices y productos cartesianos en espacios vectoriales.Como resultado, debido a que mi sistema es simétrico y definitivo positivo (dado que mi operador lineal está escrito como ), el gradiente conjugado se puede adaptar para resolver este tipo de ecuación de forma iterativa. La única modificación se produce al calcular productos internos, es decir, un cálculo típico de un producto interno en CG se ve como r T k r k o p T k A p k . En la versión modificada, utilizamos el producto interno Frobenius, que se puede calcular sumando las entradas del producto Hadamard (puntiagudo). Es decirR∗R + λ I rTkrk pagTkA pk
Sospecho que esto pasará bien cuando actualice a las matrices 3D, aunque todavía no he visto el producto interno de Frobenius definido en las matrices 3D (trabajaré bajo el supuesto de que nuevamente puedo sumar el producto puntiagudo).
¡Todavía estaría interesado en métodos más generales si alguien sabe de alguno!
fuente