مصطلحات
العربية: طريقة تشوليسكي
الإنجليزية: Cholesky Method
الفرنسية: Méthode de Cholesky
المبدأ
-
نعتبر A مصفوفة متماثلة وموجبة.
توجد مصفوفة R مثلثية سفلية بحيث A=RRT (تفكيك شوليسكي).
-
طريقة تشوليسكي تتجلى في تحديد (بتحقيق الذاتية) المصفوفة R، ثم حل النظمتين المثلثيتين على التوالي Ry=b و RTx=y.
الخوارزم
المرحلة 1: حساب R بحيث A = RRT، مع R مثلثية سفلية
المرحلة 1: حل النظمة المثلثية R y = b
المرحلة 1: حل النظمة المثلثية RT x = y
مثال
لنحل النظمة Ax=b بحيث:
و
حل يدوي (باستعمال الكسر والجذور)
المرحلة 1: حساب R
وبالتالي
المرحلة 3: حل النظمة المثلثية R y = b
المرحلة 4: حل النظمة المثلثية RT x = y
البرمجة
#include <math.h>
#include <conio.h>
#define NMAX 5
int sl_decomp_cholesky(double a[NMAX][NMAX],double r[NMAX][NMAX],int n);
void sl_resol_cholesky(double r[NMAX][NMAX],double b[NMAX],int n);
main()
{
double a[NMAX][NMAX],r[NMAX][NMAX],b[NMAX];
int i,j,n,err;
clrscr();
n=4;
printf("Méthode de Cholesky\n");
b[1]=1;b[2]=0;b[3]=0;b[4]=-1;
for(i=1;i<=n;i++)
for(j=1;j<=n;j++)
{
a[i][j]=0;
if (j==i+1) a[i][i+1]=1;
if (j==i-1) a[i][i-1]=1;
if (j==i) a[i][j]=2;
}
printf(" A b\n");
for(i=1;i<=n;i++)
{
for(j=1;j<=n;j++)
printf("%12.6e ",a[i][j]);
printf("I %12.6e \n",b[i]);
}
err=sl_decomp_cholesky(a,r,n);
if(err==1)
{
printf(" R \n");
for(i=1;i<=n;i++)
{
for(j=1;j<=i;j++)
printf("%12.6e ",r[i][j]);
for(j=i+1;j<=n;j++)
printf("%12.6e ",0.0);
printf("\n");
}
sl_resol_cholesky(r,b,n);
printf("Solution:\n");
for(i=1;i<=n;i++)
printf("x%d =%23.16e \n",i,b[i]);
}
else printf("Erreur");
}
نتائج البرنامج
يظهر البرنامج ما يلي :
-
المصفوفة A و العنصر الثاني b للنظمة؛
-
المصفوفة R؛
A b
2.00000e+00 1.00000e+00 0.00000e+00 0.00000e+00 I 1.00000e+00
1.00000e+00 2.00000e+00 1.00000e+00 0.00000e+00 I 0.00000e+00
0.00000e+00 1.00000e+00 2.00000e+00 1.00000e+00 I 0.00000e+00
0.00000e+00 0.00000e+00 1.00000e+00 2.00000e+00 I -1.00000e+00
R
1.41421e+00 0.00000e+00 0.00000e+00 0.00000e+00
7.07107e-01 1.22474e+00 0.00000e+00 0.00000e+00
0.00000e+00 8.16497e-01 1.15470e+00 0.00000e+00
0.00000e+00 0.00000e+00 8.66025e-01 1.11803e+00
Solution:
x1 = 1.000000000000000e+00
x2 = -1.000000000000000e+00
x3 = 1.000000000000000e+00
x4 = -1.000000000000000e+00
الدالة الرئيسية sl_decomp_cholesky
- معايير الدالة
a : المصفوفة الرئيسية
r : المصفوفة R
n : درجة المصفوفة
القيمة الرجعية
ترجع القيمة 1 عندما تكون المصفوفة متماثلة معرفة موجبة والقيمة 0 إذا كان عكس ذلك.
إذا كانت قيمة الرجوع هي 1 فسيكون الحل في الجدول r.
خلال كل كرة يتم إظهار المصفوفة Ak والعنصر الثاني bk.
معايير الدالة sl_resol_cholesky
r : المصفوفة R
b : المصفوفة الرئيسية
n : درجة المصفوفة
ترجع الحل ممثلا في الجدول b.
تساوي الثابتة NMAX البعد القصوي للمصفوفة زائد 1، وتساوي الثابتة N2MAX القيمة NMAX+1.
int sl_decomp_cholesky(double a[NMAX][NMAX],double r[NMAX][NMAX],int n)
{
int i,j,k,err;
double s;
err=1;
k=1;
while (err==1 && k<=n)
{
s=a[k][k];
for(i=1;i<k;i++)
s-=r[k][i]*r[k][i];
if(s>=0)
{
r[k][k]=sqrt(s);
for(j=k+1;j<=n;j++)
{
s=a[k][j];
for(i=1;i<k;i++)
s-=r[k][i]*r[j][i];
r[j][k]=s/r[k][k];
}
}
else err=0;
k++;
}
return(err);
}
void sl_resol_cholesky(double r[NMAX][NMAX],double b[NMAX],int n)
{
int i,j;
double s;
double y[NMAX];
y[1]=b[1]/r[1][1];
for(i=2;i<=n;i++)
{
s=b[i];
for(j=1;j<i;j++)
s-=r[i][j]*y[j];
y[i]=s/r[i][i];
}
b[n]=y[n]/r[n][n];
for(i=n-1;i>=1;i--)
{
s=y[i];
for(j=i+1;j<=n;j++)
s-=r[j][i]*b[j];
b[i]=s/r[i][i];
}
}