Desenvolva um programa amigável ao usuário para solução de sistemas algébricos lineares que utilizem: Decomposição de Cholesky. Teste o programa em um dos exemplos de aplicação.
Passo 1
O problema é escrever um programa interativo que leia uma matriz simétrica e definida positiva e um vetor , faça a decomposição de Cholesky de (isto é, encontre tal que ) e por fim resolva o sistema
usando a fatoração . O conteúdo matemático necessário inclui a fórmula da decomposição de Cholesky e os métodos de substituição direta (forward) e substituição regressiva (backward).
Passo 2
Derivamos o algoritmo de Cholesky da forma padrão: para e :
Depois, resolvemos
L,y = b (forward substitution)
e
L^T x = y (backward substitution)
Passo 3
Vamos agora estruturar o programa em Python de forma amigável ao usuário: 1) ler n; 2) ler linha a linha os elementos de A e verificar simetria e positividade (somente diagonal positiva no processo); 3) ler vetor b; 4) implementar função cholesky(A) que retorna L; 5) implementar forward(L,b) e backward(L,y) para obter x; 6) exibir L, y e x.
Passo 4
Segue a implementação em Python:
import math
defcholesky(A):
n = len(A)
L = [[0.0]*n for _ inrange(n)]
for i inrange(n):
for j inrange(i+1):
s = sum(L[i][k]*L[j][k] for k inrange(j))
if i == j:
val = A[i][i] - s
if val <= 0:
raise ValueError("Matriz não definida positiva.")
L[i][i] = math.sqrt(val)
else:
L[i][j] = (A[i][j] - s)/L[j][j]
return L
defforward(L, b):
n = len(b)
y = [0.0]*n
for i inrange(n):
s = sum(L[i][k]*y[k] for k inrange(i))
y[i] = (b[i] - s)/L[i][i]
return y
defbackward(L, y):
n = len(y)
x = [0.0]*n
for i inrange(n-1, -1, -1):
s = sum(L[k][i]*x[k] for k inrange(i+1, n))
x[i] = (y[i] - s)/L[i][i]
return x
defmain():
print("Resolução de A\(x=b\) por Decomposição de Cholesky")
n = int(input("Informe a dimensão n: "))
print("Digite a matriz A (linha por linha, elementos separados por espaço):")
A = [list(map(float, input().split())) for _ inrange(n)]
# opcional: checar simetriafor i inrange(n):
for j inrange(i+1, n):
ifabs(A[i][j] - A[j][i]) > 1e-9:
raise ValueError("A não é simétrica.")
print("Digite o vetor b (elementos separados por espaço):")
b = list(map(float, input().split()))
L = cholesky(A)
print("\nMatriz L (Cholesky):")
for row in L:
print(" ", ["{:.4f}".format(v) for v in row])
y = forward(L, b)
print("\nVetor y (L y = b):")
print(" ", ["{:.4f}".format(v) for v in y])
x = backward(L, y)
print("\nSolução x (L^T x = y):")
print(" ", ["{:.4f}".format(v) for v in x])
if __name__ == "__main__":
main()
Passo 5
Testando com um exemplo simples:
A = 2 & 1 1 & 2, b = 33
Ao rodar, o programa retorna
L = l[1.4142&00.7071&1.2247r], y = [2.1213,,1.4142], x = [1.0000,,1.0000],
que é de fato a solução ,x=(1,1)^T.