Fatoração LU, valor de y não está certo e eu não sei porque

[Cairo Alberto]

Tenho que faze um codigo que resolva Matrizes A*x=b, pelo metodo de fatoração LU, aonde o programa recebe A e b e tenho que resolver x (sem a necessidade de trocar linhas)
eu estou usando esses valores:
Ann => n=3 => A{3, 2, 4; 1, 1, 2;4, 3, -2}   bn =>  b{1, 2, 3}   é esperado o retorno de yn => n=3 => y{1, 1.6666..., 0} e xn=> x{-3, 5, 0}

porém, não sei porque, y, está com o valor errado, e como se usa y para conseguir x, x acaba ficando errado também.

o algoritmo é é suposto eu implementar é esse:
https://imgur.com/a/QB1pU7g

OBs: o algoritmo da foto, as posições do vetor começa a contar a partir de 1 (1, 2, 3,..), e em c++ o vetor começa a contar a partir de 0 (0, 1, 2, 3,...)
então se você ver por exemplo, na foto: "y1=b1", e no meu codigo está: "y[0] = b[0];", esse é o porque.

 

 

e meu cogigo está assim: 

 

 

#include <iomanip>
#include <iostream>

using std::cout; using std::cin;
using std::setprecision;

int main(){
	cout<<setprecision(2);
	setlocale(LC_ALL,"");
	//definir tamanho da matriz A
	int n,e;
	cout<<"Defina tamanho da matriz A nxn:\n";
	cout<<"n = ";	cin>>n;
	
	//criar matriz A e b, e definir os valores
	float A[n][n];
	float b[n];
	
	for(int i=0;i<n;i++){
	  for(int j=0;j<n;j++){
	  	cout<<"Defina o valor da matriz na posição "<<i+1<<", "<<j+1<<"\n";
	  	cin>>e;
	    A[i][j] = e;
	  }
	}
	for(int j=0;j<n;j++){
	  	cout<<"Defina o valor do vetor b na posição , "<<j+1<<"\n";
	  	cin>>e;
	    b[j] = e;
	  }
	  
	//imprimir matriz A original aumentada(b)
	cout<<"\n\n   ";
	for(int i=0;i<n+1;i++){
		if(i!=n){
			cout<<"c"<<i+1<<"	";
		}else{
			cout<<"b";	
		}
	}
		cout<<"\n";
	for(int i=0;i<n;i++){
		
		cout<<"l"<<i+1<<"|";
		for(int j=0;j<n+1;j++){
			if(j!=n){
			
			cout<<A[i][j]<<"	";
			}else{
			cout<<b[i]<<"	";
			}
		}
		cout<<"\n";
	}
  //_______________________________________________________________________________________________________________
	//Criar e organizar matriz A=LU	
	float U[n][n];
	float L[n][n];

	//copiar matriz A em U
	for(int i=0;i<n;i++){
	  for(int j=0;j<n;j++){
	  	U[i][j]=A[i][j];
	  }
	}
	//Matriz L = Matriz Identidade
	for(int i=0;i<n;i++){
	  for(int j=0;j<n;j++){
	  	if(i==j){
		   L[i][j] = 1;
		}else{L[i][j] = 0; }
	  }
	}

	//_______________________________________________________________________________________________________________
  
	// Resolver Matrizes L U | U = triangular superior | L = triangular inferior
	for(int j=0;j<n;j++){
		for(int i=j+1;i<n;i++){ 
			if(U[i][j]!=0){
				L[i][j] = U[i][j]/U[j][j]; 
				for(int c=j;c<n;c++){
					U[i][c]=U[i][c]+U[j][c]*(-1*(L[i][j]));
				}
			}
		}
	}
	float soma=0;
	float y[n];
	float x[n];
	
	//resolver L*y=b
	//valores de texte A = 					|3, 2, 4; 1, 2, 2; 4, 3, -2|     b = |1; 2; 3|
  	//valores esperados:
  	//Matriz 	U=					|3, 2, 4; 0, 0.33... ,0.66... ; 0, 0, -8|	(correto)
  	//matriz 	L=					|1, 0, 0; 0.33..., 1, 0; 1.33.. , 1, 1|		(correto)
	//vetor  	y=					|1, 1.66..., 0|					(errado)
  	//vetor		x=					|-3, 5, 0|					(errado)					

	y[0]=b[0]; //A nxn = A 3x3,  l=3
	for(int i=1;i<=n-1;i++){				//i[1]						|i[2]			|i[2]
          soma=0;						//soma=0					|			|	
		for(int j=0;j<=i-1;j++){			//j[0]						|j[0]			|j[1]
			soma+= L[i][j]*y[j];			//soma=(0,333..)				|soma+=(1,33..)		|soma=(1,33..+1*1,66..)=2,99..
		}						//						|soma=1,333		|
		y[i] = b[i]-soma; 				//y[1]=(2-0,33)=1,66...				|y[2]=			|y[2]=3-3=0
	}
	
	// resolver U*y=x
	
	x[n-1]= y[n-1]/A[n-1][n-1]; //Xn = Yn/Ann
	
	for(int i=n-2;i>-1;i--){
		soma=0;
		for(int j=i+1;j<n;j++){
			soma = soma+A[i][j]*x[j];	
		}
		x[i]=(y[i]-soma)/A[i][i];
	}
	

	
	
	
	
	//__________________________________________________________________________________________________________________
	cout<<"\n\n		Matriz U   \n   ";
	for(int i=0;i<n+1;i++){
		if(i!=n){
			cout<<"c"<<i+1<<"		";
			
		}
	}
		cout<<"\n";
	for(int i=0;i<n;i++){
		
		cout<<"l"<<i+1<<"|";
		for(int j=0;j<n+1;j++){
			if(j!=n){
			
			cout<<U[i][j]<<"		";
		
			}
		}
		cout<<"\n";
	}
	
	
	cout<<"\n\n		Matriz L   \n   ";
	for(int i=0;i<n+1;i++){
		if(i!=n){
			cout<<"c"<<i+1<<"		";
			
		}
	}
		cout<<"\n";
	for(int i=0;i<n;i++){
		
		cout<<"l"<<i+1<<"|";
		for(int j=0;j<n+1;j++){
			if(j!=n){
			
			cout<<L[i][j]<<"		";
		
			}
		}
		cout<<"\n";
	}
	cout<<"y = \n";
	for(int i=0;i<n;i++){
		
	cout<<"y"<<i+1<<" = "<<y[i]<<"		";
		
	}
	//cout<<setprecision(5);
	cout<<"\n\n";
	for(int i=0;i<n;i++){
		
	cout<<"x"<<i+1<<" = "<<x[i]<<"		\n";
		
	}
	
}

 

Editado Abril 7, 2022 por Cairo Alberto