Giáo trình Turbo C nâng cao và C++ - Chương 10: Giải hệ phương trình đại số tuyến tính

145 Ch−ơng 10 : Giải hệ ph−ơng trình đại số tuyến tính Đ1.Ph−ơng pháp Gauss Có nhiều ph−ơng pháp để giải một hệ ph−ơng trình tuyến tính dạng AX = B. Ph−ơng pháp giải sẽ đơn giản hơn nếu ma trận A có dạng tam giác nghĩa là có dạng : ⎟⎟ ⎟ ⎠ ⎞ ⎜⎜ ⎜ ⎝ ⎛ 333231 2221 11 aaa 0aa 00a hay ⎟⎟ ⎟ ⎠ ⎞ ⎜⎜ ⎜ ⎝ ⎛ 33 2322 131211 a00 aa0 aaa Trong tr−ờng hợp đầu tiên,ma trận đ−ợc gọi là ma trận tam giác d−ới và tr−ờng hợp thứ hai ma trận đ−ợc gọi là ma trận tam giác trên.Ph−ơng trình t−ơng ứng với ma trận tam giác d−ới có dạng t−ờng minh là : 11 1 2 3 1 21 1 2 3 2 31 1 2 3 3 0 0 022 32 33 a x x x b a x a x x b a x a x a x b + + = + + = + + = ⎧ ⎨ ⎪⎪ ⎩ ⎪⎪ Với ph−ơng trình dạng này chúng ta sẽ giải ph−ơng trình từ trên xuống. Ch−ơng trình giải ph−ơng trình ma trận tam giác d−ới là : Ch−ơng trình 10-1 #include #include #include #include #include #define max 10 void main() { float a[max][max]; float b[max],x[max]; int i,j,k,n,t; float s,c; char tl; clrscr(); printf("Cho so phuong trinh n = "); scanf("%d",&n); printf("Cho cac phan tu cua ma tran a\n"); for (i=1;i<=n;i++) for (j=1;j<=n;j++) { printf("a[%d][%d] = ",i,j); scanf("%f",&a[i][j]); } printf("\n"); printf("Ma tran a ma ban da nhap\n"); printf("\n"); 146 for (i=1;i<=n;i++) { for (j=1;j<=n;j++) printf("%15.5f",a[i][j]); printf("\n"); } printf("\n"); t=1; flushall(); while (t) { printf("Co sua ma tran a khong(c/k)?"); scanf("%c",&tl); if (toupper(tl)=='C') { printf("Cho chi so hang can sua : "); scanf("%d",&i); printf("Cho chi so cot can sua : "); scanf("%d",&j); printf("a[%d][%d] = ",i,j); scanf("%f",&a[i][j]); } if (toupper(tl)=='K') t=0; } printf("Ma tran a ban dau\n"); printf("\n"); for (i=1;i<=n;i++) { for (j=1;j<=n;j++) printf("%15.5f",a[i][j]); printf("\n"); } printf("\n"); printf("Cho cac phan tu cua ma tran b\n"); for (i=1;i<=n;i++) { printf("b[%d] = ",i); scanf("%f",&b[i]); } printf("\n"); printf("Ma tran b ma ban da nhap"); printf("\n"); for (i=1;i<=n;i++) printf("b[%d] = %10.5f\n",i,b[i]); printf("\n"); flushall(); t=1; while (t) { 147 printf("Co sua ma tran b khong(c/k)?"); scanf("%c",&tl); if (toupper(tl)=='C') { printf("Cho chi so hang can sua : "); scanf("%d",&i); printf("b[%d] = ",i); scanf("%f",&b[i]); } if (toupper(tl)=='K') t=0; } printf("\n"); printf("Ma tran b ban dau"); printf("\n"); for (i=1;i<=n;i++) printf("%15.5f\n",b[i]); { if (a[1][1]==0) if (b[1]!=0) printf("He da cho vo nghiem\n"); else { printf("He da cho co vo so nghiem"); x[n]=c; } else x[1]=b[1]/a[1][1]; for (i=2;i<=n;i++) { s=0; for (k=1;k<=i-1;k++) s=s+a[i][k]*x[k]; x[i]=(b[i]-s)/a[i][i]; } printf("\n"); printf("Nghiem cua he da cho la"); printf("\n"); for (i=1;i<=n;i++) printf("x[%d] = %10.5f\n",i,x[i]); getch(); } } Ph−ơng trình t−ơng ứng với ma trận tam giác trên có dạng t−ờng minh là : 11 1 2 3 1 1 2 3 2 1 2 3 3 12 13 22 23 33 0 0 0 a x a x a x b x a x a x b x x a x b + + = + + = + + = ⎧ ⎨ ⎪⎪ ⎩ ⎪⎪ Với ph−ơng trình này chúng ta giải từ d−ới lên. 148 Ch−ơng trình giải ph−ơng trình ma trận tam giác trên là : Ch−ơng trình 10-2 #include #include #include #include #include #define max 10 void main() { float a[max][max]; float b[max],x[max]; int i,j,k,n,t; float s,c; char tl; clrscr(); printf("Cho so phuong trinh n = "); scanf("%d",&n); printf("Cho cac phan tu cua ma tran a :\n"); for (i=1;i<=n;i++) for (j=1;j<=n;j++) { printf("a[%d][%d] = ",i,j); scanf("%f",&a[i][j]); } printf("\n"); printf("Ma tran a ma ban da nhap\n"); printf("\n"); for (i=1;i<=n;i++) { for (j=1;j<=n;j++) printf("%15.5f",a[i][j]); printf("\n"); } printf("\n"); t=1; flushall(); while (t) { printf("Co sua ma tran a khong(c/k)?"); scanf("%c",&tl); if (toupper(tl)=='C') { printf("Cho chi so hang can sua : "); scanf("%d",&i); 149 printf("Cho chi so cot can sua : "); scanf("%d",&j); printf("a[%d][%d] = ",i,j); scanf("%f",&a[i][j]); } if (toupper(tl)=='K') t=0; } printf("Ma tran a ban dau"); printf("\n"); for (i=1;i<=n;i++) { for (j=1;j<=n;j++) printf("%15.5f",a[i][j]); printf("\n"); } printf("\n"); printf("Cho cac phan tu cua ma tran b : \n"); for (i=1;i<=n;i++) { printf("b[%d] = ",i); scanf("%f",&b[i]); } printf("\n"); printf("Ma tran b ma ban da nhap"); printf("\n"); for (i=1;i<=n;i++) printf("b[%d] = %10.5f\n",i,b[i]); printf("\n"); flushall(); t=1; while (t) { printf("Co sua ma tran b khong(c/k)?"); scanf("%c",&tl); if (toupper(tl)=='C') { printf("Cho chi so hang can sua : "); scanf("%d",&i); printf("b[%d] = ",i); scanf("%f",&b[i]); } if (toupper(tl)=='K') t=0; } printf("\n"); printf("Ma tran b ban dau\n"); printf("\n"); for (i=1;i<=n;i++) printf("b[%d] = %10.5f\n",i,b[i]); 150 printf("\n"); { if (a[n][n]==0) if (b[n]!=0) printf("He da cho vo nghiem"); else { printf("He da cho co vo so nghiem"); x[n]=c; } else x[n]=b[n]/a[n][n]; for (i=n-1;i>=1;i--) { s=0; for (k=i+1;k<=n;k++) s=s+a[i][k]*x[k]; x[i]=(b[i]-s)/a[i][i]; } printf("\n"); printf("Nghiem cua he da cho la\n"); printf("\n"); for (i=1;i<=n;i++) printf("x[%d] = %10.5f\n",i,x[i]); getch(); } } Tuy nhiên, các hệ ph−ơng trình đơn giản hiếm khi gặp trong thực tế. Các hệ ph−ơng trình tuyến tính có thể biểu diễn d−ới dạng tam giác nếu định thức của nó khác không, nghĩa là ph−ơng trình có nghiệm.Chúng ta biết rằng các nghiệm của hệ không đổi nếu ta thay một hàng bằng tổ hợp tuyến tính của các hàng khác.Nh− vậy bằng một loạt các biến đổi ta có thể đ−ahệ ban đầu về dạng tam giác. Đó chính là nội dung của ph−ơng pháp loại trừ Gauss. Chúng ta hãy xét hệ ph−ơng trình : 11 1 2 3 1 21 1 2 3 2 31 1 2 3 3 12 13 22 23 32 33 a x a x a x b a x a x a x b a x a x a x b + + = + + = + + = ⎧ ⎨ ⎪⎪ ⎩ ⎪⎪ Nhân hàng thứ nhất với a21/a11 ta có : 21 1 21 11 12 2 21 11 13 3 21 11 1a x a a a x a a a x a a b + + = Số hạng đầu của ph−ơng trình bằng số hạng đầu của hàng thứ hai trong hệ ph−ơng trình ban đầu.Khi trừ hàng một đã đ−ợc biến đổi cho hàng 2 ta nhận đ−ợc hàng 2 mới 0 1 22 21 11 12 2 23 21 11 13 3 2 21 11 1x a a a a x a a a a x b a a b + − + − = −( ) ( ) Ta tiếp tục cách này để loại trừ x1 ra khỏi hàng thứ 3.Ph−ơng trình trở thành : 11 12 13 22 23 32 33 1 2 3 1 2 3 0 0 , , , , , , , , ' , a a a a a a a x x x b b b ⎛ ⎝ ⎜⎜⎜⎜⎜⎜ ⎞ ⎠ ⎟⎟⎟⎟⎟⎟ ⎛ ⎝ ⎜⎜⎜⎜⎜ ⎞ ⎠ ⎟⎟⎟⎟⎟ ⎛ ⎝ ⎜⎜⎜⎜⎜⎜ ⎞ ⎠ ⎟⎟⎟⎟⎟⎟ ì = 151 với a,11 = a11 ; a , 12 = a12 ; a , 13 = a13 ; a , 13 = a13 ; b , 1 = b1 22 22 21 11 12 , a a a a a = − 23 23 21 11 13 , a a a a a = − 32 32 31 11 12 , a a a a a = − 33 33 31 11 13 , a a a a a = − 2 2 21 11 1 , b b a a b = − 3 3 31 11 1 , b b a a b = − Ta loại trừ số hạng chứa x3 trong dòng thứ 3 bằng cách t−ơng tự.Ta nhân hàng thứ 2 trong hệ A'X = B' với a,32/a,22 và đem trừ đi hàng thứ 3 trong hệ mới.Nh− vậy số hạng chứa x3 biến mất và ta nhận đ−ợc ma trận tam giác trên. 11 12 13 22 23 33 1 2 3 1 2 3 0 0 0 ,, ,, ,, ,, ,, ,, ,, ,, ,, a a a a a a x x x b b b ⎛ ⎝ ⎜⎜⎜⎜⎜⎜ ⎞ ⎠ ⎟⎟⎟⎟⎟⎟ ⎛ ⎝ ⎜⎜⎜⎜⎜ ⎞ ⎠ ⎟⎟⎟⎟⎟ ⎛ ⎝ ⎜⎜⎜⎜⎜⎜ ⎞ ⎠ ⎟⎟⎟⎟⎟⎟ ì = với 11 11 ,, , a a= 12 12,, ,a a= 13 13,, ,a a= 1 1,, ,b b= 22 22,, ,a a= 23 23,, ,a a= 2 2,, ,b b= 33 33 32 22 23 ,, , , , , a a a a a= − 3 3 33 22 2 ,, , , , , b b a a b= − Các phép tính này chỉ thực hiện đ−ợc khi a11 ≠ 0 và a,11 ≠ 0. Với một hệ có n ph−ơng trình,thuật tính hoàn toàn t−ơng tự.Sau đây là ch−ơng trình giải hệ ph−ơng trình n ẩn số bằng ph−ơng pháp loại trừ Gauss. Ch−ơng trình 10-3 #include #include #include #include #include #define max 10 void main() { float b[max],x[max]; float a[max][max]; int i,j,k,n,t; float c,s,d; char tl; clrscr(); printf("Cho so phuong trinh n = "); scanf("%d",&n); printf("Cho cac phan tu cua ma tran a :\n"); for (i=1;i<=n;i++) for (j=1;j<=n;j++) { printf("a[%d][%d] = ",i,j); scanf("%f",&a[i][j]); } printf("\n"); printf("Ma tran a ma ban da nhap\n"); printf("\n"); 152 for (i=1;i<=n;i++) { for (j=1;j<=n;j++) printf("%15.5f",a[i][j]); printf("\n"); } printf("\n"); t=1; flushall(); while (t) { printf("Co sua ma tran a khong(c/k)?"); scanf("%c",&tl); if (toupper(tl)=='C') { printf("Cho chi so hang can sua : "); scanf("%d",&i); printf("Cho chi so cot can sua : "); scanf("%d",&j); printf("a[%d][%d] = ",i,j); scanf("%f",&a[i][j]); } if (toupper(tl)=='K') t=0; } printf("Ma tran a ban dau\n"); printf("\n"); for (i=1;i<=n;i++) { for (j=1;j<=n;j++) printf("%15.5f",a[i][j]); printf("\n"); } printf("\n"); printf("Cho cac phan tu cua ma tran b : \n"); for (i=1;i<=n;i++) { printf("b[%d] = ",i); scanf("%f",&b[i]); } printf("\n"); printf("Ma tran b ma ban da nhap\n"); printf("\n"); for (i=1;i<=n;i++) printf("b[%d] = %15.5f\n",i,b[i]); printf("\n"); flushall(); t=1; while (t) { 153 printf("Co sua ma tran b khong(c/k)?"); scanf("%c",&tl); if (toupper(tl)=='C') { printf("Cho chi so hang can sua : "); scanf("%f",&i); printf("b[%d] = ",i); scanf("%f",&b[i]); } if (toupper(tl)=='K') t=0; } printf("\n"); printf("Ma tran b\n"); for (i=1;i<=n;i++) printf("b[%d] = %15.5f\n",i,b[i]); printf("\n"); for (k=1;k<=n-1;k++) { for (i=k+1;i<=n;i++) { b[i]=b[i]-b[k]*a[i][k]/a[k][k]; for (j=k+1;j<=n;j++) a[i][j]=a[i][j]-a[k][j]*a[i][k]/a[k][k]; } } { if (a[n][n]==0) if (b[n]==0) printf("He da cho vo nghiem"); else { printf("He da cho co vo so nghiem"); x[n]=c; } else x[n]=b[n]/a[n][n]; for (i=n-1;i>=1;i--) { s=0; for (k=i+1;k<=n;k++) s=s+a[i][k]*x[k]; x[i]=(b[i]-s)/a[i][i]; } printf("\n"); printf("Nghiem cua he da cho la\n"); printf("\n"); for (i=1;i<=n;i++) printf("x[%d] = %15.5f\n",i,x[i]); 154 getch(); } } Đ2.Ph−ơng pháp Gauss-Jordan Xét hệ ph−ơng trình AX=B. Khi giải hệ bằng ph−ơng pháp Gauss ta đ−a nó về dạng ma trận tam giác sau một loạt biến đổi. Ph−ơng pháp khử Gauss-Jordan cải tiến khử Gauss bằng cách đ−a hệ về dạng : EX = B* và khi đó nghiệm của hệ chính là B*.Trong ph−ơng pháp Gauss-Jordan mỗi b−ớc tính phải tính nhiều hơn ph−ơng pháp Gauss nh−ng lại không phải tính nghiệm.Để đ−a ma trận A về dạng ma trận E tại b−ớc thứ i ta phải có aii = 1 và aij=0.Nh− vậy tại lần khử thứ i ta biến đổi : 1.aij = aij/aii (j=i+1,i+2,...,n) 2.k=1,2,...,n akj = akj - aijaki (j=i+1,i+2,...,n) bk = bk - biaki Ví dụ : Cho hệ ⎟⎟ ⎟ ⎠ ⎞ ⎜⎜ ⎜ ⎝ ⎛ = ⎟⎟ ⎟⎟ ⎠ ⎞ ⎜⎜ ⎜⎜ ⎝ ⎛ ì ⎟⎟ ⎟ ⎠ ⎞ ⎜⎜ ⎜ ⎝ ⎛ 21 26 32 24 x x x x 9440 45.652 45104 0248 4 3 2 1 Biến đổi lần 1 : ta chia hàng 1 cho a11 = 8;nhân hàng 1 vừa nhận đ−ợc với 4 và lấy hàng 2 trừ đi; nhân hàng 1 vừa nhận đ−ợc với 2 và lấy hàng 3 trừ đi;giữ nguyên hàng 4 vì phần tử đầu tiên đã bằng 0 ta có ⎟⎟ ⎟ ⎠ ⎞ ⎜⎜ ⎜ ⎝ ⎛ = ⎟⎟ ⎟⎟ ⎠ ⎞ ⎜⎜ ⎜⎜ ⎝ ⎛ ì ⎟⎟ ⎟ ⎠ ⎞ ⎜⎜ ⎜ ⎝ ⎛ 21 20 20 3 x x x x 9440 4640 4480 025.05.01 4 3 2 1 Biến đổi lần 2 : ta chia hàng 2 cho a22 = 8;nhân hàng 2 vừa nhận đ−ợc với 0.5 và lấy hàng 1 trừ đi; nhân hàng 2 vừa nhận đ−ợc với 4 và lấy hàng 3 trừ đi; nhân hàng 2 vừa nhận đ−ợc với 4 và lấy hàng 4 trừ đi ta có : ⎟⎟ ⎟ ⎠ ⎞ ⎜⎜ ⎜ ⎝ ⎛ = ⎟⎟ ⎟⎟ ⎠ ⎞ ⎜⎜ ⎜⎜ ⎝ ⎛ ì ⎟⎟ ⎟ ⎠ ⎞ ⎜⎜ ⎜ ⎝ ⎛ − 11 10 5.2 75.1 x x x x 7200 2400 5.05.010 25.0001 4 3 2 1 Biến đổi lần 3 : ta chia hàng 3 cho a33 = 4;giữ nguyên hàng 1;nhân hàng 3 vừa nhận đ−ợc với 0.5 và lấy hàng 2 trừ đi; nhân hàng 3 vừa nhận đ−ợc với 2 và lấy hàng 4 trừ đi ta có : ⎟⎟ ⎟ ⎠ ⎞ ⎜⎜ ⎜ ⎝ ⎛ = ⎟⎟ ⎟⎟ ⎠ ⎞ ⎜⎜ ⎜⎜ ⎝ ⎛ ì ⎟⎟ ⎟ ⎠ ⎞ ⎜⎜ ⎜ ⎝ ⎛ − 6 5.2 25.1 75.1 x x x x 6000 5.0100 25.0010 25.0001 4 3 2 1 155 Biến đổi lần 4 : ta chia hàng 4 cho a44 = 6;nhân hàng 4 vừa nhận đ−ợc với -0.25 và lấy hàng 1 trừ đi; nhân hàng 4 vừa nhận đ−ợc với 0.25 và lấy hàng 2 trừ đi; nhân hàng 4 vừa nhận đ−ợc với 0.5 và lấy hàng 3 trừ đi ta có : ⎟⎟ ⎟ ⎠ ⎞ ⎜⎜ ⎜ ⎝ ⎛ = ⎟⎟ ⎟⎟ ⎠ ⎞ ⎜⎜ ⎜⎜ ⎝ ⎛ ì ⎟⎟ ⎟ ⎠ ⎞ ⎜⎜ ⎜ ⎝ ⎛ 1 2 1 2 x x x x 1000 0100 0010 0001 4 3 2 1 và ta có ngay vec tơ nghiệm. Ch−ơng trình10-4 #include #include #include #include #include #define spt 10 void main() { float a[spt][2*spt]; float b[spt]; int i,j,k,n,m,t; float max,c; char tl; clrscr(); printf("Cho so phuong trinh n = "); scanf("%d",&n); printf("Cho cac phan tu cua ma tran a :\n"); for (i=1;i<=n;i++) for (j=1;j<=n;j++) { printf("a[%d][%d] = ",i,j); scanf("%f",&a[i][j]); } printf("\n"); printf("Ma tran a ma ban da nhap"); printf("\n"); for (i=1;i<=n;i++) { for (j=1;j<=n;j++) printf("%15.5f",a[i][j]); printf("\n"); } printf("\n"); t=1; flushall(); 156 while (t) { printf("Co sua ma tran a khong(c/k)?"); scanf("%c",&tl); if (toupper(tl)=='C') { printf("Cho chi so hang can sua : "); scanf("%d",&i); printf("Cho chi so cot can sua : "); scanf("%d",&j); printf("a[%d][%d] = ",i,j); scanf("%f",&a[i][j]); } if (toupper(tl)=='K') t=0; } printf("Ma tran a\n"); printf("\n"); for (i=1;i<=n;i++) { for (j=1;j<=n;j++) printf("%15.5f",a[i][j]); printf("\n"); } printf("\n"); printf("Cho cac phan tu cua ma tran b : \n"); for (i=1;i<=n;i++) { printf("b[%d] = ",i); scanf("%f",&b[i]); } printf("\n"); printf("Ma tran b ma ban da nhap\n"); printf("\n"); for (i=1;i<=n;i++) printf("b[%d] = %15.5f\n",i,b[i]); printf("\n"); t=1; flushall(); while (t) { printf("Co sua ma tran b khong(c/k)?"); scanf("%c",&tl); if (toupper(tl)=='C') { printf("Cho chi so hang can sua : "); scanf("%d",&i); printf("b[%d] = ",i); scanf("%f",&b[i]); } 157 if (toupper(tl)=='K') t=0; } printf("\n"); printf("Ma tran b\n"); printf("\n"); for (i=1;i<=n;i++) printf("%15.5f\n",b[i]); printf("\n"); t=1; flushall(); i=1; while (t) { if (a[i][i]==0) { max=0; m=i; for (k=i+1;k<=n;k++) if (max<fabs(a[k][i])) { m=k; max=fabs(a[i][i]); } if (m!=i) { for (j=i;j<=n;j++) { c=a[i][j]; a[i][j]=a[m][j]; a[m][j]=c; } c=b[i]; b[i]=b[m]; b[m]=c; } if (m==i) { t=0; printf("MA TRAN SUY BIEN"); } } if (a[i][i]!=0) { c=1/a[i][i]; for (j=i;j<=n;j++) a[i][j]=a[i][j]*c; b[i]=b[i]*c; for (k=1;k<=n;k++) if (k!=i) 158 { c=a[k][i]; for (j=i;j<=n;j++) a[k][j]=a[k][j]-a[i][j]*c; b[k]=b[k]-b[i]*c; } } i=i+1; if (i==(n+1)) t=0; } if (i==(n+1)) { printf("NGHIEM CUA HE"); printf("\n"); for (i=1;i<=n;i++) printf("x[%d] = %15.5f\n",i,b[i]); } getch(); } Đ3.Ph−ơng pháp Cholesky Trong ph−ơng pháp Cholesky một ma trận đối xứng A đ−ợc phân tích thành dạng A = RTR trong đó R là một ma trận tam giác trên.Hệ ph−ơng trình lúc đó chuyển thành AX = RTRX = B.Nh− vậy tr−ớc hết ta phân tích ma trận A thành tích hai ma trận.Sau đó giải hệ ph−ơng trình RTY = B và cuối cùng là hệ RX = Y.Ch−ơng trình mô tả thuật toán này đ−ợc cho d−ới đây : Ch−ơng trình 10-5 #include #include #include #include #include #define max 6 void main() { float a[max][max],r[max][max]; float b[max],x[max],y[max]; int i,j,k,l,n,t; float s; char tl; clrscr(); printf("Cho so phuong trinh n = "); 159 scanf("%d",&n); printf("Cho cac phan tu cua ma tran a : \n"); for (i=1;i<=n;i++) for (j=1;j<=n;j++) { printf("a[%d][%d] = ",i,j); scanf("%f",&a[i][j]); } printf("\n"); printf("Ma tran a ma ban da nhap\n"); printf("\n"); for (i=1;i<=n;i++) { for (j=1;j<=n;j++) printf("%15.5f",a[i][j]); printf("\n"); } printf("\n"); flushall(); t=1; while (t) { printf("Co sua ma tran a khong(c/k)?"); scanf("%c",&tl); if (toupper(tl)=='C') { printf("Cho chi so hang can sua : "); scanf("%d",&i); printf("Cho chi so cot can sua : "); scanf("%d",&j); printf("a[",i,",",j,"] = "); scanf("%f",&a[i][j]); } if (toupper(tl)=='K') t=0; } printf("Ma tran a\n"); printf("\n"); for (i=1;i<=n;i++) { for (j=1;j<=n;j++) printf("%15.5f",a[i][j]); printf("\n"); } printf("\n"); printf("Cho cac phan tu cua ma tran b : \n"); for (i=1;i<=n;i++) { printf("b[%d] = ",i); scanf("%f",&b[i]); 160 } printf("\n"); printf("Ma tran b ma ban da nhap\n"); printf("\n"); for (i=1;i<=n;i++) printf("b[%d] = %15.5f\n",i,b[i]); printf("\n"); flushall(); t=1; while (t) { printf("Co sua ma tran b khong(c/k)?"); scanf("%c",&tl); if (toupper(tl)=='C') { printf("Cho chi so hang can sua : "); scanf("%d",&i); printf("b[%d] = ",i); scanf("%f",&b[i]); } if (toupper(tl)=='K') t=0; } printf("\n"); printf("Ma tran b\n"); printf("\n"); for (i=1;i<=n;i++) printf("b[%d] = %15.5f\n",i,b[i]); printf("\n"); for(i=1;i<=n;i++) for(j=1;j<=n;j++) r[i][j]=0.0; for (i=1;i<=n;i++) { if (a[i][i]>=0) { r[i][i]=sqrt(a[i][i]); for (j=1+i;j<=n;j++) r[i][j]=a[i][j]/r[i][i]; for (k=i+1;k<=n;k++) for (l=k;l<=n;l++) a[k][l]=a[k][l]-r[i][k]*r[i][l]; } } for (k=1;k<=n;k++) { s=b[k]; if (k!=1) for (i=1;i<=k-1;i++) s=s+r[i][k]*y[i]; 161 y[k]=-s/r[k][k]; } for (i=n;i>=1;i--) { s=-y[i]; if (i!=n) for (k=i+1;k<=n;k++) s=s-r[i][k]*x[k]; x[i]=s/r[i][i]; } printf("Nghiem cua he phuong trinh la\n "); for (i=1;i<=n;i++) printf("x[%d] = %10.5f\n",i,x[i]); getch(); } Đ4.Ph−ơng pháp Crout Ph−ơng pháp Crout là một dạng của ph−ơng pháp Gauss.Với ph−ơng pháp Gauss,chúng ta biến đổi ma trận A thành một ma trận tam giác thì ở ph−ơng pháp Crout chúng ta phân tích ma trận này thành tích của ma trận tam giác trên R và ma trận tam giác d−ới L.Trong ma trận L,các hệ số trên đ−ờng chéo chính bằng 1.Nh− vậy ph−ơng trình AX = B đ−ợc viết thành : A.X = L.R.X = B Chúng ta đặt RX = Y nên : LY = B Nh− vậy tr−ớc hết chúng ta phân tích ma trận thành tích của L.R.Tiếp theo ta giải ph−ơng trình LY = B và sau đó giải ph−ơng trình RX = A để tìm nghiệm X. Ch−ơng trình 10-6 #include #include #include #include #include #define max 6 void main() { float b[max],x[max],y[max]; float a[max][max],r[max][max],l[max][max]; int i,j,k,n,t; float c,tr,tl,s; char tloi; clrscr(); printf("Cho so phuong trinh n = "); 162 scanf("%d",&n); printf("Cho cac phan tu cua ma tran a : \n"); for (i=1;i<=n;i++) for (j=1;j<=n;j++) { printf("a[%d][%d] = ",i,j); scanf("%f",&a[i][j]); } printf("\n"); printf("Ma tran a ma ban da nhap"); printf("\n"); for (i=1;i<=n;i++) { for (j=1;j<=n;j++) printf("%10.5f",a[i][j]); printf("\n"); } printf("\n"); t=1; flushall(); while (t) { printf("Co sua ma tran a khong(c/k)?"); scanf("%c",&tloi); if (toupper(tloi)=='C') { printf("Cho chi so hang can sua : "); scanf("%d",&i); printf("Cho chi so cot can sua : "); scanf("%d",&j); printf("a[%d][%d] = ",i,j); scanf("%f",&a[i][j]); flushall(); } if (toupper(tloi)=='K') t=0; } printf("Ma tran a\n"); printf("\n"); for (i=1;i<=n;i++) { for (j=1;j<=n;j++) printf("%10.5f",a[i][j]); printf("\n"); } printf("\n"); printf("Cho cac phan tu cua ma tran b : \n"); for (i=1;i<=n;i++) { printf("b[%d] = ",i); 163 scanf("%f",&b[i]); } printf("\n"); printf("Ma tran b ma ban da nhap"); printf("\n"); for (i=1;i<=n;i++) printf("b[%d] = %10.5f\n",i,b[i]); printf("\n"); t=1; flushall(); while (t) { printf("Co sua ma tran b khong(c/k)?"); scanf("%c",&tloi); if (toupper(tloi)=='C') { printf("Cho chi so hang can sua : "); scanf("%d",&i); printf("b[%d] = ",i); scanf("%f",&b[i]); flushall(); } if (toupper(tloi)=='K') t=0; } printf("\n"); printf("Ma tran b\n"); printf("\n"); for (i=1;i<=n;i++) printf("b[%d] = %10.5f\n",i,b[i]); printf("\n"); for(i=1;i<=n;i++) for(j=1;j<=n;j++) { r[i][j]=0.0; l[i][j]=0.0; } for (i=1;i<=n;i++) { r[1][i]=a[1][i]; l[i][i]=1.0; l[i][1]=a[i][1]/a[1][1]; } for (k=2;k<=n;k++) { for (j=k;j<=n;j++) { tr=0.0; for (i=1;i<=k;i++) tr=tr+l[k][i]*r[i][j]; 164 r[k][j]=a[k][j]-tr; } if (k!=n) { for (i=1;i<=n;i++) { tl=0.0; for (j=1;j<=k-1;j++) tl=tl+l[i][j]*r[j][k]; l[i][k]=(a[i][k]-tl)/r[k][k]; } } else printf("\n"); } if (l[1][1]==0.0) if (b[1]==0.0) printf("He da cho vo nghiem\n"); else { printf("He da cho co vo so nghiem\n"); y[n]=c; } else y[1]=b[1]/l[1][1]; for (i=2;i<=n;i++) { s=0.0; for (k=1;k<=i-1;k++) s=s+l[i][k]*y[k]; y[i]=(b[i]-s)/l[i][i]; } if (r[n][n]==0.0) if (y[n]==0.0) printf("He da cho vo nghiem\n"); else { printf("He da cho co vo so nghiem\n"); x[n]=c; } else x[n]=y[n]/r[n][n]; for (i=n-1;i>=1;i--) { s=0.0; for (k=i+1;k<=n;k++) s+=r[i][k]*x[k]; x[i]=(y[i]-s)/r[i][i]; } printf("\n"); 165 printf("Nghiem cua he da cho la\n"); printf("\n"); for (i=1;i<=n;i++) printf("x[%d] = %15.5f\n",i,x[i]); getch(); } Đ5.Ph−ơng pháp lặp đơn Xét hệ ph−ơng trình AX = F.Bằng cách nào đó ta đ−a hệ ph−ơng trình về dạng X = BX + G trong đó B = (bij)n,n G = (g1,g2,...,gn) T Chọn vectơ X = ( x1 (o),x2 (o),....,xn (o) )T làm xấp xỉ thứ 0 của nghiệm đúng và xây dựng xấp xỉ X(m+1) = BX(m) + G ( m = 0,1,....) Ng−ời ta chứng minh rằng nếu ph−ơng trình ban đầu có nghiệm duy nhất và một trong ba chuẩn của ma trận B nhỏ hơn 1 thì dãy xấp xỉ hội tụ về nghiệm duy nhất đó.(Cho một ma trận B,chuẩn của ma trận B,kí hiệu là || B || là một trong 3 số : 1 1 B b i ij j n = = ∑max 2 1 B b ji ij i n = = ∑max 3 1 2 1 2 1 B b i n ij j n = = = ∑ ∑ / ( ) (Chuẩn của ma trận quan hệ tới sự hội tụ của ph−ơng pháp lặp) Ví dụ chúng ta có ph−ơng trình 1 2 3 1 2 3 1 2 3 10 2 10 10 2 12 10 8 x x x x x x x x x + + = + + = + + = ⎧ ⎨ ⎪⎪ ⎩ ⎪⎪ Chúng ta đ−a ph−ơng trình về dạng : 1 2 3 2 1 3 3 1 2 1 5 1 10 1 1 10 1 5 6 5 1 10 1 10 4 5 x x x x x x x x x = − − + = − − + = − − + ⎧ ⎨ ⎪⎪⎪⎪⎪ ⎩ ⎪⎪⎪⎪⎪ Nh− vậy : B = − − − − − − ⎡ ⎣ ⎢⎢⎢⎢⎢⎢⎢ ⎤ ⎦ ⎥⎥⎥⎥⎥⎥⎥ 0 1 5 1 10 1 10 0 1 5 1 10 1 10 0 và G = ⎡ ⎣ ⎢⎢⎢⎢⎢⎢ ⎤ ⎦ ⎥⎥⎥⎥⎥⎥ 1 5 6 4 5 Dễ thấy || B ||1 = 3/10 ; || B ||2 = 3/10 và || B ||3 = 12/100 nên phép lặp hội tụ .Ch−ơng trình lặp đơn là : 166 Ch−ơng trình 10-7 #include #include #include #include #include #define max 10 void main() { float a[max][max]; float f[max],x0[max],x1[max]; int i,j,k,n,l,t; float s,c,epsi; char tl; clrscr(); printf("Cho so phuong trinh n = "); scanf("%d",&n); printf("Cho cac phan tu cua ma tran a : \n"); for (i=1;i<=n;i++) for (j=1;j<=n;j++) { printf("a[%d][%d] = ",i,j); scanf("%f",&a[i][j]); } printf("\n"); printf("Ma tran a ma ban da nhap"); printf("\n"); for (i=1;i<=n;i++) { for (j=1;j<=n;j++) printf("%15.5f",a[i][j]); printf("\n"); } printf("\n"); t=1; flushall(); while (t) { printf("Co sua ma tran a khong(c/k)?"); scanf("%c",&tl); if (toupper(tl)=='C') { printf("Cho chi so hang can sua : "); scanf("%d",&i); printf("Cho chi so cot can sua : "); scanf("%d",&j); printf("a[%d][%d] = ",i,j); 167 scanf("%f",&a[i][j]); flushall(); } if (toupper(tl)=='K') t=0; } printf("Ma tran a\n"); printf("\n"); for (i=1;i<=n;i++) { for (j=1;j<=n;j++) printf("%10.5f",a[i][j]); printf("\n"); } printf("\n"); printf("Cho cac phan tu cua ma tran f : \n"); for (i=1;i<=n;i++) { printf("f[%d] = ",i); scanf("%f",&f[i]); } printf("\n"); printf("Ma tran f ma ban da nhap"); printf("\n"); for (i=1;i<=n;i++) printf("f[%d] = %10.5f\n",i,f[i]); printf("\n"); t=1; flushall(); while (t) { printf("Co sua ma tran f khong(c/k)?"); scanf("%c",&tl); if (toupper(tl)=='C') { printf("Cho chi so hang can sua : "); scanf("%d",&i); printf("f[%d] = ",i); scanf("%f",&f[i]); flushall(); } if (toupper(tl)=='K') t=0; } printf("\n"); printf("Ma tran f"); printf("\n"); for (i=1;i<=n;i++) printf("f[%d] = %10.5f\n",i,f[i]); printf("\n"); 168 for (i=1;i<=n;i++) x0[i]=0.0; x0[1]=1.0; printf("Cho so lan lap l = "); scanf("%d",&l); epsi=1e-5; for (i=1;i<=n;i++) { c=1.0/a[i][i]; for (j=1;j<=n;j++) if (j!=i) a[i][j]*=c; f[i]*=c; a[i][i]=0.0; } k=1; t=0; do { for (i=1;i<=n;i++) { x1[i]=f[i]; for (j=1;j<=n;j++) x1[i]=x1[i]-a[i][j]*x0[j]; } s=0.0; for (i=1;i<=n;i++) s=s+fabs(x1[i]-x0[i]); if (s>=epsi) for (i=1;i<=n;i++) x0[i]=x1[i]; if (s<epsi) { t=1; printf("\n"); printf("Phep lap hoi tu sau %d buoc tinh",k); printf("\n"); printf("NGHIEM CUA HE"); printf("\n"); for (i=1;i<=n;i++) printf("x[%d] = %10.5f\n",i,x1[i]); } k=k+1; if (k>l) { t=1; printf("Phep lap khong hoi tu sau %d buoc tinh",k-1); } } while (t==0); 169 getch(); } Đ6. Ph−ơng pháp lặp Gauss-Seidel Ph−ơng pháp lặp Gauss-Seidel đ−ợc cải tiến từ ph−ơng pháp lặp đơn . Nội dung cơ bản của ph−ơng pháp là ở chỗ khi tính nghiệm xấp xỉ thứ (k+1) của ẩn xi ta sở dụng các xấp xỉ thứ (k+1) đã tính của các ẩn x1,...,xi-1.Giả sử đã cho hệ : AX = B và ta có nghiệm : i i ij j n jx x i n= + =∑ =β α1 1( ,...., ) Lấy xấp xỉ ban đầu tuỳ ý x1 (o) , x2 (o) ,...., xn (o) và tất nhiên ta cố gắng lấy chúng t−ơng ứng với x1,x2 ,...,xn (càng gần càng tốt). Tiếp theo ta giả sử rằng đã biết xấp xỉ thứ k xi (k) của nghiệm . Theo Seidel ta sẽ tìm xấp xỉ thứ (k+1) của nghiệm theo các công thức sau : 1 1 1 1 1 ( ) ( )k j j n j k x x + = = + ∑β α 2 1 2 1 1 2 2 21 ( ) ( ) ( )k k j j n j k x x x + + = = + + ∑β α α ...... i k i ij j i j k ij j n j k x x x ( ) ( ) ( )+ = − + = = + +∑ ∑1 1 1 1 1 β α α ...... ni k n ij j n j k nn n k x x x ( ) ( ) ( )+ = − += + +∑1 1 1 1β α α Thông th−ờng ph−ơng pháp Gauss - Seidel hội tụ nhanh hơn ph−ơng pháp lặp đon nh−ng tính toán phức tạp hơn.Dể dễ hiểu ph−ơng pháp này chúng ta xét một ví dụ cụ thể : Cho hệ ph−ơng trình : 1 2 3 1 2 3 1 2 3 10 12 2 10 13 2 2 10 14 x x x x x x x x x + + = + + = + + = ⎧ ⎨ ⎪⎪ ⎩ ⎪⎪ nghiệm đúng của hệ là (1 , 1, 1) Ta đ−a về dạng thuận tiện cho phép lặp : 1 2 3 2 1 3 3 2 2 1 2 01 01 1 3 0 2 01 1 4 0 2 0 2 x x x x x x x x x = − − = − − = − − ⎧ ⎨ ⎪⎪ ⎩ ⎪⎪ , . . , . . , . . Lấy x1 (o) = 1.2 ; x2 (o) = 0 ; x3 (o) = 0; Sử dụng ph−ơng pháp lặp Seidel ta có : ⎪⎩ ⎪⎨ ⎧ =ì−ì−= =ì−ì−= =ì−ì−= 948.006.12.02.12.04.1x 06.101.02.12.03.1x 2.101.001.02.1x 1 3 1 2 1 1 1 2 2 2 3 2 1 2 01 106 01 0 948 0 9992 1 3 0 2 0 9992 01 0 948 100536 1 4 0 2 0 9992 0 2 100536 0 999098 x x x = − − = = − − = = − − = ⎧ ⎨ ⎪⎪⎪ ⎩ ⎪⎪⎪ ì ì ì ì ì ì , . . . . . , . . . . . , . . . . . và cứ thế tiếp tục . Ch−ơng trình mô tả thuật toán Gauss - Seidel nh− sau : 170 Ch−ơng trình 10-8 #include #include #include #include #include #define max 6 void main() { float b[max],x[max]; float a[max][max]; int i,j,k,n,dem,t1; float t,s,d,w,epsi; char tl; clrscr(); printf("Cho so an so n = "); scanf("%d",&n); printf("Cho cac phan tu cua ma tran a : \n"); for (i=1;i<=n;i++) for (j=1;j<=n;j++) { printf("a[%d][%d] = ",i,j); scanf("%f",&a[i][j]); } printf("\n"); printf("Ma tran a ma ban da nhap\n"); printf("\n"); for (i=1;i<=n;i++) { for (j=1;j<=n;j++) printf("%10.5f",a[i][j]); printf("\n"); } printf("\n"); t1=1; flushall(); while (t1) { printf("Co sua ma tran a khong(c/k)?"); scanf("%c",&tl); if (toupper(tl)=='C') { printf("Cho chi so hang can sua : "); scanf("%d",&i); printf("Cho chi so cot can sua : "); scanf("%d",&j); printf("a[%d][%d] = ",i,j); 171 scanf("%f",&a[i][j]); flushall(); } if (toupper(tl)=='K') t1=0; } printf("Ma tran a\n"); printf("\n"); for (i=1;i<=n;i++) { for (j=1;j<=n;j++) printf("%15.5f",a[i][j]); printf("\n"); } printf("\n"); printf("Cho cac phan tu cua ma tran b : \n"); for (i=1;i<=n;i++) { printf("b[%d] = ",i); scanf("%f",&b[i]); } printf("\n"); printf("Ma tran b"); printf("\n"); for (i=1;i<=n;i++) printf("b[%d] = %10.5f\n",i,b[i]); printf("\n"); t1=1; flushall(); while (t1) { printf("Co sua ma tran b khong(c/k)?"); scanf("%c",&tl); if (toupper(tl)=='C') { printf("Cho chi so hang can sua : "); scanf("%d",&i); printf("b[%d] = ",i); scanf("%f",&b[i]); flushall(); } if (toupper(tl)=='K') t1=0; } printf("\n"); printf("Ma tran b"); printf("\n"); for (i=1;i<=n;i++) printf("b[%d] = %10.5f\n",i,b[i]); printf("\n"); 172 printf("Cho so lan lap k : "); scanf("%d",&k); printf("\n"); w=1; epsi=1e-8; for (i=1;i<=n;i++) x[i]=0.0; dem = 0; do { dem=dem+1; for (i=1;i<=n;i++) { s=0.0; for (j=1;j<=n;j++) s=s+a[i][j]*x[j]; d=x[i]; x[i]=(1-w)*d+w*(-s+a[i][i]*d+b[i])/a[i][i]; t=fabs(d-x[i]); } } while ((demepsi*fabs(x[n]))); if (t<epsi*fabs(x[n])) { printf("Nghiem sau %d lan lap la :\n",dem); for (i=1;i<=n;i++) printf("x[%d] = %12.8f\n",i,x[i]); } else { printf("Khong dat do chinh xac sau %d lan lap\n",k); printf("Nghiem cua lan lap cuoi cung la : \n"); for (i=1;i<=n;i++) printf("x[%d] = %12.8f\n",i,x[i]); } getch(); } Đ7. Ph−ơng pháp Cramer Một tr−ờng hợp riêng của hệ ph−ơng trình , trong đó số ph−ơng trình bằng số ẩn , nghĩa là hệ có dạng : 11 1 12 2 13 3 1n 1 21 1 22 2 23 3 2 2 1 1 2 2 3 3 a x a x a x a x b a x a x a x a x b a x a x a x a x b n n n n n n nn n n + + + + = + + + + = + + + + = ⎧ ⎨ ⎪⎪⎪ ⎩ ⎪⎪⎪ ... ... ... ...... trong đó A= (aij) là một ma trận không suy biến . một hệ ph−ơng trình nh− vậy có tên là hệ thống Cramer 173 Định lí Crame : Hệ thống Crame có nghiệm duy nhất đ−ợc cho bởi công thức : i i x A A n= = ( ) (i ,..., )1 trong đó A(i) là ma trận nhận đ−ợc từ A bằng cách thay cột thứ i bởi cột B=[ b1 ,...,bn]T Nh− vậy để giải hệ bằng ph−ơng pháp Cramer chúng ta lần l−ợt tính các định thức của ma trận và ma trận thay thế rồi tìm nghiệm theo công thức Cramer. Ch−ơng trình sau mô tả thuật toán này, Ch−ơng trình 10-.9 // Cramer; #include #include #include #define max 50 void main() { float r[max][max],a[max][max]; float b[max],x[max]; float delta[max]; int i,j,k,l,t,n,ok1,ok2,t1; float c,d; char tl; clrscr(); printf("Cho so an cua phuong trinh n = "); scanf("%d",&n); printf("Cho cac phan tu cua ma tran a : \n"); for (i=1;i<=n;i++) for (j=1;j<=n;j++) { printf("a[%d][%d] = ",i,j); scanf("%f",&a[i][j]); } printf("\n"); printf("Ma tran a ma ban da nhap\n"); printf("\n"); for (i=1;i<=n;i++) { for (j=1;j<=n;j++) printf("%10.5f",a[i][j]); printf("\n"); } printf("\n"); t1=1; flushall(); while (t1) { 174 printf("Co sua ma tran a khong(c/k)?"); scanf("%c",&tl); if (toupper(tl)=='C') { printf("Cho chi so hang can sua : "); scanf("%d",&i); printf("Cho chi so cot can sua : "); scanf("%d",&j); printf("a[",i,",",j,"] = "); scanf("%f",&a[i][j]); flushall(); } if (toupper(tl)=='K') t1=0; } printf("Ma tran a \n"); printf("\n"); for (i=1;i<=n;i++) { for (j=1;j<=n;j++) printf("%10.5f",a[i][j]); printf("\n"); } printf("\n"); printf("Cho cac phan tu cua ma tran b : \n"); for (i=1;i<=n;i++) { printf("b[%d] = ",i); scanf("%f",&b[i]); } printf("\n"); printf("Ma tran b ma ban da nhap\n"); printf("\n"); for (i=1;i<=n;i++) printf("b[%d] = %10.5f\n",i,b[i]); printf("\n"); t1=1; flushall(); while (t1) { printf("Co sua ma tran b khong(c/k)?"); scanf("%c",&tl); if (toupper(tl)=='C') { printf("Cho chi so hang can sua : "); scanf("%d",&i); printf("b[%d] = ",i); scanf("%f",&b[i]); flushall(); } 175 if (toupper(tl)=='K') t1=0; } printf("\n"); printf("Ma tran b\n"); printf("\n"); for (i=1;i<=n;i++) printf("%10.5f",b[i]); printf("\n"); //thay b vao tung cot cua a de tinh cac dinh thuc for (k=0;k<=n;k++) { for (i=1;i<=n;i++) for (j=1;j<=n;j++)//thay cot b vao a if (i==k) r[j][i]=b[j]; else r[j][i]=a[j][i]; //tinh dinh thuc d=1.0; i=1; ok2=1; while (ok2&&(i<=n)) { if (r[i][i]==0.0) { ok1=1; t=t+1; while (ok1&&(t<=n)) if (r[t][i]!=0) { for (j=i;j<=n;j++) { c=r[i][j]; r[i][j]=r[t][j]; r[k][j]=c; } d=-d; ok1=0; } else t=t+1; if (k>n) { printf("\n"); printf("** MA TRAN SUY BIEN **"); ok2=0; d=0.0; } 176 } if (r[i][i]!=0) { c=r[i][i]; for (j=i+1;j<=n;j++) r[i][j]=r[i][j]/c; for (t=i+1;t<=n;t++) { c=r[t][i]; for (j=i+1;j<=n;j++) r[t][j]=r[t][j]-r[i][j]*c; } } i=i+1; } if (ok2) for (i=1;i<=n;i++) d=d*r[i][i]; delta[k]=d; } if (delta[0]==0.0) printf("He da cho vo nghiem\n"); else { printf("\nNGHIEM CUA HE :\n"); printf("\n"); for (i=1;i<=n;i++) { x[i]=delta[i]/delta[0]; printf("x[%d] = %10.5f\n",i,x[i]); } } getch(); } Đ8. Hệ ph−ơng trình số phức Giả sử ta có một hệ ph−ơng trình dạng số phức dạng AX = B trong đó A = C + jD , B = E +jF và X = Y + jZ . Ta viết lại ph−ơng trình d−ới dạng : (C + jD)(Y + jZ) = (E +jF) Nhân biểu thức về trái và cân bằng phần thực với phần thực và phần ảo với phần ảo ta nhận đ−ợc hệ mới : ⎩⎨ ⎧ + F = CZ DY E =DZ - CY Nh− vậy chúng ta nhận đ−ợc một hệ gồm 2n ph−ơng trình số thực . Giải hệ này và kết hợp các phần thực à phần ảo ta nhận đ−ợc nghiệm của hệ ph−ơng trình ban đầu . Ch−ơng trình giải hệ ph−ơng trình nh− vậy cho ở d−ới đây : 177 Ch−ơng trình 10-10 #include #include #include #include #include #define max 20 void main() { int i,j,k,l,n,m; float s,t,a[max][max],b[max][max],x[max]; clrscr(); printf("Cho so an so cua phuong trinh n = "); scanf("%d",&n); printf("Cho phan thuc cua cac he so,ke ca ve phai\n"); for (i=1;i<=n;i++) for (j=1;j<=n+1;j++) { printf("a[%d][%d] = ",i,j); scanf("%f",&a[i][j]); } printf("\n"); printf("Cho phan ao cua cac he so,ke ca ve phai\n"); for (i=1;i<=n;i++) for (j=1;j<=n+1;j++) { printf("b[%d][%d] = ",i,j); scanf("%f",&b[i][j]); } for (i=1;i<=n;i++) a[i][2*n+1]=a[i][n+1]; for (i=n+1;i<=2*n;i++) a[i][2*n+1]=b[i-n][n+1]; for (i=n+1;i<=2*n;i++) for (j=1;j<=n;j++) a[i][j]=b[i-n][j]; for (i=1;i<=n;i++) for (j=n+1;j<=2*n;j++) a[i][j]=-b[i][j-n]; for (i=n+1;i<=2*n;i++) for (j=n+1;j<=2*n;j++) a[i][j]=a[i-n][j-n]; m=2*n; for (k=1;k<=m-1;k++) { s=0.0; 178 for (i=k;i<=m;i++) { t=fabs(a[i][k]); if (s<=t) { s=t; l=i; } } for (j=k;j<=m+1;j++) { s=a[k][j]; a[k][j]=a[l][j]; a[l][j]=s; } if (fabs(a[k][k]/a[1][1])<=1e-08) { printf("Ma tran suy bien\n"); getch(); exit(1); } for (i=k+1;i<=m;i++) { s=a[i][k]/a[k][k]; a[i][k]=0.0; for (j=k+1;j<=m+1;j++) a[i][j]-=s*a[k][j]; } } x[m]=a[m][m+1]/a[m][m]; for (i=m-1;i>=1;i--) { s=a[i][m+1]; for (j=i+1;j<=m;j++) s-=a[i][j]*x[j]; x[i]=s/a[i][i]; } printf("\n"); printf("Nghiem phuc cua he\n"); for (i=1;i<=n;i++) if (x[i+n]<0) printf("%10.5f-%10.5fj\n",x[i],fabs(x[i+n])); else printf("%10.5f+%10.5fj\n",x[i],x[i+n]); getch(); } Dùng ch−ơng trình này giải hệ ph−ơng trình : 179 ( ) ( ) ( ) ( ) ( ) ( ) ( ) ( ) ( ) ( ) ( ) ( ) ( ) 3 7 2 4 1 3 4 2 8 36 5 6 2 5 3 4 10 4 5 1 2 5 6 13 3 2 4 1 2 3 10 6 + + − + + − + + = + − + + + − + = + + + + + − − + = − + + − + − = − + ⎧ ⎨ ⎪⎪ ⎩ ⎪⎪ j x j y j z j r j j x j z j r j j x j y j z r j j x j y j r j Ta nhận đ−ợc các nghiệm x = 2 + 3j ; y = 1 - 2j ; z = -1 + 4j và r = 1- j Ngoài các ph−ơng pháp nêu trên ta thấy rằng từ hệ ph−ơng trình AX = B ta có thể tìm nghiệm X của hệ bằng cách viết lại ph−ơng trình d−ới dạng X = B/A =A-1B với A-1 là ma trận nghịch đảo của A . Do vậy tr−ớc hết ta cần tìm A-1 và sau đó tính tích A-1B.