DECOMPOSITION DE CHOLESKY

/*****Resolution du systeme lineaire (triangle inferieur)*****/
void resyst_tri_inf(float mat[dim][dim],float b1[dim],float x1[dim],int n)
{
 	int i,k;
	float S;
 	x1[1]=b1[1]/mat[1][1];
 	for(i=2;i<=n;i++)
 	{
		S=0;
  		for(k=1;k<=i;k++)
			S+=mat[i][k]*x1[k];
		x1[i]=(b1[i]-S)/mat[i][i];
 	}
}
/*****Resolution du systeme lineaire (triangle superieur)*****/
void resyst_tri_sup(float mat[dim][dim],float b1[dim],float x1[dim],int n)
{
 	int i,k;
	float S;
 	x1[n]=b1[n]/mat[n][n];
 	for(i=n-1;i>=1;i--)
 	{
		S=0;
  		for(k=n;k>=i+2;k--)
			S+=mat[i][k]*x1[k];
		x1[i]=(b1[i]-S)/mat[i][i];
 	}
}
/*****cholesky : résolution d'un systeme de la forme Ax=b,*/
/*Avec A = L.LT, A étant une matrice sysmétrique définie positive de dimension n*/
/*L:matrice triangulaire inférieur, et LT:transposée de L*/
void cholesky(float A[dim][dim],float L[dim][dim],int n)
{
	int i,j,k;
	float S;
	for(j=1;j<=n;j++)
	{
		S=0;
		for(k=1;k<j;k++)
		S+=pow(L[j][k],2);
		L[j][j]=sqrt(A[j][j]-S);
		for(i=j;i<=n;i++)
		{
			S=0;
			for(k=1;k<j;k++)
			S+=L[i][k]*L[j][k];
			L[i][j]=(A[i][j]-S)/L[j][j];
		}
	}
}