Resolución numérica de sistemas lineales y problemas propios

Cada fila de Ax=b es una ecuación y cada columna corresponde a la misma incógnita en todas ellas. Una matriz cuadrada invertible da solución única. Ser cuadrada por sí solo no basta: si es singular, el sistema puede ser incompatible o tener múltiples soluciones.

Resolver no exige construir A⁻¹. La eliminación reorganiza ecuaciones equivalentes y termina en una forma fácil de despejar. En Gauss basta llegar a triangular y sustituir; Gauss-Jordan continúa eliminando también por encima de los pivotes. En cualquier caso conservamos el término independiente en cada operación.

Usamos un pivote no nulo para hacer ceros bajo él. Reemplazar una fila por ella menos un múltiplo de otra conserva las soluciones porque la operación se puede deshacer. Al terminar, despejamos de abajo arriba en el sistema triangular. Registramos la fila objetivo y calculamos todas sus entradas, incluido b.

Si el pivote es cero, intercambiamos filas cuando haya una entrada adecuada, o estudiamos el rango. En punto flotante también preocupa que sea muy pequeño: genera multiplicadores grandes. El pivoteo parcial elige el mayor módulo de la columna activa entre las filas restantes, e intercambia sus filas completas.

Con la convención PA=LU, P registra permutaciones, L es triangular inferior con diagonal uno y U triangular superior.

Definimos y=Ux, resolvemos Ly=Pb por sustitución hacia delante y después Ux=y hacia atrás. y es un vector auxiliar concreto.

Para una matriz densa general, factorizar cuesta O(n³) y cada nuevo b requiere O(n²) en sustituciones. Si la matriz no cambia podemos reutilizar los factores; si cambia, esa reutilización exacta ya no está justificada. Dispersión y estructura de banda cambian almacenamiento, orden de eliminación y coste.

El pivoteo parcial evita pivotes innecesariamente pequeños en su columna; sus multiplicadores tienen módulo ≤1. Aun así, las entradas intermedias pueden crecer mucho en algunas matrices, por lo que no proporciona crecimiento pequeño universal. Las garantías de redondeo dependen de ese crecimiento y de la estructura del cálculo.

Con A invertible, x−x̃=A⁻¹r; un residuo pequeño necesita escalas y condición para convertirse en error pequeño. Si A es simétrica definida positiva, Cholesky usa A=CCᵀ con C triangular y diagonal positiva: aprovecha la estructura, sin hacer aplicable el método a cualquier matriz simétrica.

La matriz devuelve un múltiplo del vector inicial, de modo que conserva su recta. La escala es libre: cualquier múltiplo no nulo de v es otro autovector del mismo λ. El vector cero no sirve, pues cumpliría la igualdad para cualquier λ sin identificar una dirección.

En una diagonal, los ejes son autovectores y las entradas diagonales sus autovalores. En una matriz real general pueden aparecer autovalores complejos. Si la matriz real es simétrica, el teorema espectral garantiza autovalores reales y una base ortonormal de autovectores. Esa estructura permite interpretar mejor los residuos propios.

Normalizar evita crecimiento de escala; no cambia la recta. Si wₖ=0, no se puede dividir por su norma: la semilla o el paso cae en el núcleo y hay que interpretarlo.

Si |λ₁|>|λⱼ| para todos los demás y c₁≠0, dividir por λ₁ᵏ deja c₁u₁ más términos cⱼ(λⱼ/λ₁)ᵏ que tienden a cero. Esto prueba convergencia de dirección bajo esas hipótesis. La razón de módulos controla la tasa básica; una semilla sin componente dominante no la crea en aritmética exacta.

Si el autovalor dominante es negativo, su potencia alterna signo. El vector normalizado puede alternar entre dos orientaciones mientras la recta que genera converge. Si lo que buscamos es una dirección, debemos reconocer esa equivalencia de signos en vez de declarar fallo solo por comparar vₖ₊₁−vₖ.

Si dos autovalores de módulo máximo empatan, la prueba de dominancia estricta deja de funcionar. Puede haber alternancias, rotaciones o una componente que nunca se hace despreciable. Observar dos pasos parecidos tampoco certifica haber hallado un autovector: se necesita comprobar la ecuación propia.

Para v≠0, el cociente de Rayleigh es ρ=v*Av/(v*v), donde * es traspuesta conjugada. El residuo propio es r=Av−ρv; para comparar entre escalas usamos ‖r‖/‖v‖. En matrices hermíticas, ρ es real y la expansión en una base ortonormal propia muestra dist(ρ, σ(A))≤‖r‖/‖v‖: la norma al cuadrado del residuo es Σ|λⱼ−ρ|²|cⱼ|².

La cota controla cercanía de un valor al espectro hermítico, no automáticamente el ángulo con un autovector concreto sin separación espectral. En matrices no normales la interpretación cambia. Para cualquier matriz, si tomamos un λ candidato y r=Av−λv, la perturbación E=−rv*/‖v‖² satisface (A+E)v=λv y ‖E‖₂=‖r‖₂/‖v‖₂: esto es un control hacia atrás del par.

LU atiende sistemas; potencia aproxima una dirección dominante, no todo el espectro. Iteración inversa con desplazamiento μ resuelve (A−μI)w=v y normaliza: requiere matriz desplazada invertible. En el caso diagonalizable con un autovalor más próximo a μ que todos los demás y componente inicial adecuada, esa dirección domina porque los factores son 1/(λⱼ−μ).

QR y métodos de Krylov atienden otras necesidades y tamaños.

Para sistemas grandes también hay iteraciones. Jacobi separa la diagonal: si A=D+R con diagonal no nula, xₖ₊₁=D⁻¹(b−Rxₖ). Con una norma subordinada que cumpla ‖D⁻¹R‖<1 se obtiene una contracción y convergencia al único punto fijo. No todo sistema converge con esta iteración. En cada método declaramos hipótesis, coste y la comprobación apropiada: Ax≈b o Av≈λv.