›
›

Giáo trình Phương pháp tính - Chương 2, Phần 2: Giải gần đúng phương trình đại số và siêu việt

Download
với  
Xi = { si, pi}T và Xi+1 = { si+1, pi+1}T  
f(s ,p )  
g(si ,pi )  
i
i
F(Xi )   
f f  
s p  
g g  
J(Xi )   
s p  
Quan hệ : J(Xi)X = -F(Xi) với X = {si+1 - si,pi+1 - pi}T tương ứng với một hệ phương trình tuyến tính hai ẩn số  
s = si+1 - si và p = pi+1 - pi :  
f  
f  
s   
s   
p  f(si ,pi )  
p  g(si ,pi )  
s  
g  
s  
p  
g  
p  
Theo công thức Cramer ta có :  
g  
p p  
f  
f g  
s  
g  
f  
g f  
s s  
p  
g  
g  
f  
f  
  
s p p s  
g g  
f f  
Để dùng được công thức này ta cần tính được các đạo hàm  
s , p , s , p . Các đạo hàm này được tính  
theo công thức truy hồi.  
Do bo = ao nên  
b0  
s  
b0  
p  
0  
0  
0  
b1 = a1 + sbo nên  
b1  
s  
b1  
p  
b0  
b2 = a2 + sb1- pbo nên  
(pb0 )  
s  
b2 a2 (sb1 )  
s  
s  
s  
Mặt khác :  
(pb0 )  
0  
a2  
(sb1 )  
s  
(b1 )  
s  
0  
s  
b1  
s  
s  
b2  
nên:  
b1 sb0  
s  
b3 = a3 + sb2- pb1 nên:  
b3  
b2  
s  
b1  
s  
b2 s  
p  
s  
Nếu chúng ta đặt :  
35  
bk  
s  
ck1  
thì :  
co = bo  
(2)  
c1 = b1 + sbo = b1 + sco  
c2 = b2 + sc1 - pco  
....................  
ck = bk + sck-1 - pck-2  
cn-1 = bn-1 + scn-2 - pcn-3  
Như vậy các hệ số cũng được tính theo cách như các hệ số bk. Cuối cùng với f = bn-1 và g = bn ta được:  
f  
s  
f  
s  
f  
s  
f  
s  
cn2  
cn3  
cn1  
cn2  
bn1cn2 bncn3  
cn1cn3 c2n2  
s   
p   
(3)  
(4)  
bn1cn1 bncn2  
cn1cn3 c2n2  
Sau khi phân tích xong Pn(x) ta tiếp tục phân tích Pn-2(x) theo phương pháp trên. Các bước tính  
toán gồm:  
- Chọn các giá trị ban đầu bất kì s0 và p0  
- Tính các giá trị bo,..,bn theo (1)  
- Tính các giá trị co,...,cn theo (2)  
- Tính so và po theo (3) và (4)  
- Tính s1 = s0 + so và p1 = po+ po  
- Lặp lại bước 1 cho đến khi pi+1 = pi = p và si+1 = si = s  
- Giải phương trình x2 - sx + p để tìm 2 nghiệm của đa thức  
- Bắt đầu quá trình trên cho đa thức Pn-2(x)  
Ví dụ: Tìm nghiệm của đa thức P4(x) = x4 - 1.1x3 + 2.3x2 + 0.5x2 + 3.3.  
Với lần lặp ban đầu ta chọn s = -1 và p =1, nghĩa là tam thức có dạng: x2+x+1  
a0  
1
a1  
-1.1  
-1  
a2  
a3  
a4  
2.3  
2.1  
-1  
0.5  
-3.4  
2.1  
3.3  
0.8  
-3.4  
sbi  
-pbi-1  
bi 1  
sbi  
-2.1  
3.4  
-0.8 = bn-1  
0.7=bn  
-1.0  
-3.1  
3.1  
-5.5  
-pbi-1  
ci  
-1.0  
5.5  
3.1  
-3.2  
1
0.8 3.1  
0.7 5.5  
5.5 3.1  
3.2 5.5  
s   
0.11  
36  
5.5 0.8  
3.2 0.7  
5.5 3.1  
3.2 5.5  
p   
0.06  
s* = -1 + 0.11 = -0.89  
p* = 1 + 0.06 = 1.06  
Tiếp tục lặp lần 2 với s1 = s* và p1 = p* ta có :  
a0  
1
a1  
-1.1  
-0.89  
a2  
2.3  
1.77  
-1.06  
a3  
0.5  
-2.68  
2.11  
a4  
3.3  
0.06  
-3.17  
sbi  
-pbi-1  
bi 1  
sbi  
-pbi-1  
ci  
-1.99  
3.01  
-1.0  
-0.07 = bn-1  
0.17=bn  
-0.89  
-2.88  
2.56  
4.51  
-4.01  
-1.03  
3.1  
1
0.07 2.88  
0.7 5.5  
4.51 2.88  
s   
 0.01  
0.04  
1.03 4.51  
4.51  
0.07  
1.03 0.17  
4.51 2.88  
1.03 4.51  
p   
s* = -0.89 - 0.01 = -0.9  
p* = 1.06 + 0.04 = 1.1  
Như vậy:  
P4(x) = (x2 + 0.9x + 1.1)(x2 + 2x + 3)  
Chương trình sau áp dụng lí thuyết vừa nêu để tìm nghiệm của đa thức.  
Chương trình 2-10  
//phuong phap Bairstow  
#include <conio.h>  
#include <stdio.h>  
#include <math.h>  
#include <stdlib.h>  
#define m 10  
void main()  
{
float a[m],b[m],c[m];  
int i,n,v;  
float s,e1,t,p,q,r,p1,q1;  
clrscr();  
printf("Cho bac cua da thuc n = ");  
scanf("%d",&n);  
printf("Cho cac he so cua da thuc can tim nghiem\n");  
37  
for (i=n;i>=0;i--)  
{
printf("a[%d] = ",n-i);  
scanf("%f",&a[i]);  
}
printf("\n");  
e1=0.0001;  
if (n<=2)  
if (n==1)  
{
printf("Nghiem cua he\n");  
printf("%.8f",(a[0]/(-a[1])));  
getch();  
exit(1);  
}
do  
{
v=0;  
p=1;  
q=-1;  
b[n]=a[n];  
c[n]=a[n];  
do  
{
b[n-1]=b[n]*p+a[n-1];  
c[n-1]=b[n-1]+b[n]*p;  
for (i=n-2;i>=0;i--)  
{
b[i]=b[i+2]*q+b[i+1]*p+a[i];  
c[i]=c[i+2]*q+c[i+1]*p+b[i];  
}
r=c[2]*c[2]-c[1]*c[3];  
p1=p-(b[1]*c[2]-b[0]*c[3])/r;  
q1=q-(b[0]*c[2]-b[1]*c[1])/r;  
if ((fabs(b[0])<e1)&&(fabs(b[1])<e1))  
goto tt;  
v=v+1;  
p=p1;  
q=q1;  
}
while (v<=40);  
if(v>40)  
{
printf("Khong hoi tu sau 40 lan lap");  
getch();  
exit(1);  
}
tt:s=p1/2;  
t=p1*p1+4*q1;  
if(t<0)  
{
printf("Nghiem phuc\n");  
printf("%.8f+%.8fj\n",s,(sqrt(-t)/2));  
printf("%.8f-%.8fj\n",s,(sqrt(-t)/2));  
printf("\n");  
38  
}
else  
{
printf("Nghiem thuc\n");  
printf("%.8f\n",(s+sqrt(t)/2));  
printf("%.8f\n",(s-sqrt(t)/2));  
printf("\n");  
}
for (i=2;i<=n;i++)  
a[i-2]=b[i];  
n=n-2;  
}
while ((n>2)&(r!=0.0));  
s=-a[1]/(2*a[2]);  
t=a[1]*a[1]-4*a[2]*a[0];  
if (t<0)  
{
printf("Nghiem phuc\n");  
printf("%.8f+%.8fj\n",s,(sqrt(-t)/(2*a[2])));  
printf("%.8f-%.8fj\n",s,(sqrt(-t)/(2*a[2])));  
printf("\n");  
}
else  
{
printf("Nghiem thuc\n");  
printf("%.8f\n",(s-sqrt(t)/(2*a[2])));  
printf("%.8f\n",(s-sqrt(t)/(2*a[2])));  
printf("\n");  
}
getch();  
}
Dùng chương trình trên để xác định nghiệm của đa thức :  
x6 - 2x5 - 4x4 + 13x3 - 24x2 + 18x - 4 = 0  
ta nhận được các nghiệm :  
x1 = 2.61903399  
x2 = -2.73205081  
x3 = 0.732050755  
x4 = 0.381966055  
x5 = 0.500011056 + i*1.3228881  
x6 = 0.500011056 - i*1.3228881  
§11. HỆ PHƯƠNG TRÌNH PHI TUYẾN  
Phương pháp Newton có thể được tổng quát hoá để giải hệ phương trình phi tuyến dạng :  
f (x ,x ,x ,...,x ) 0  
1  
1
2
3
n
f2 (x1 ,x2 ,x3 ,...,xn ) 0  
f (x ,x ,x ,...,x ) 0  
  
3
1
2
3
n
fn (x1 ,x2 ,x3 ,...,xn ) 0  
hay viết gọn hơn dưới dạng :  
F(X) = 0  
Trong đó :  
39  
X = (x1,x2,x3,.....,xn)  
Với một phương trình một biến, công thức Newton là :  
f(xi )  
xi1 xi   
f (xi )  
hay : f'(xi).x = -f(xi)  
với x = xi+1 - xi  
Đối với hệ,công thức lặp là :  
J(Xi)x = -F(Xi)  
Trong đó J(Xi) là toán tử Jacobi. Nó là một ma trận bậc n ( n - tương ứng với số thành phần trong vectơ X)  
có dạng :  
f1 f1 f1  
x1 x2 x3  
f1  
xn  
f2  
  
  
f  
f2 f2  
2
x1 x2 x3  
  
x  
n   
J(X )   
i
  
  
  
fn fn fn  
x1 x2 x3  
fn  
xn  
  
và  
X = Xi+1 - Xi  
Phương pháp Newton tuyến tính hoá hệ và như vậy với mỗi bước lặp cần giải một hệ phương trình  
tuyến tính (mà biến là xi) xác định bởi công thức lặp cho tới khi vectơ X(x1,x2,x3,.....,xn) gần với nghiệm.  
Dưới đây là chương trình giải hệ phương trình phi tuyến  
x13 x22 3x1x2x4 8 0  
x1 x2 x3 x4 5 0  
25 x12 8x3 4 0  
2x1x2x3 x4 8 0  
Ma trận đạo hàm riêng J(Xi) là :  
3x12 3x2x4 3x22 3x1x4  
0
1
8
3x1x   
2   
1
1
0
1
0
x1  
25 x12  
2x2x3  
2x1x3  
2x1x2  
1  
Ma trận này được chương trình đọc vào nhờ thủ tục doc.Trong thủ tục này,các hệ số a[i,5] là các  
hàm fi(x).Vectơ nghiệm ban đầu được chọn là { 0,-1,-1,1}T.Kết quả tính cho ta : x = {0.01328676,-1.94647929,-  
1.12499779,8.05819031 }T với độ chính xác 0.000001.Vectơ số dư r = { 0.00000536,-0.00000011,-0.00000001,-  
0.00000006}T.  
Chương trình 2-11  
//giai he pt phi tuyen  
#include <conio.h>  
#include <stdio.h>  
#include <math.h>  
#include <stdlib.h>  
#define n 4  
float a[n+1][n+2];  
40  
float x[n+1],y[n+1];  
int i,j,k,l,z,r;  
float e,s,t;  
void main()  
{
void doc();  
clrscr();  
printf("Cho cac gia tri nghiem ban dau\n");  
for (i=1;i<=n;i++)  
{
printf("x[%d] = ",i);  
scanf("%f",&x[i]);  
}
e=1e-6;  
z=30;  
for (r=1;r<=z;r++)  
{
doc();  
for (k=1;k<=n-1;k++)  
{
s=0 ;  
for (i=k;i<=n;i++)  
{
t=fabs(a[i][k]);  
if (s<=t)  
{
s=t;  
l=i;  
}
}
for (j=k;j<=n+1;j++)  
{
s=a[k][j];  
a[k][j]=a[l][j];  
a[l][j]=s;  
}
if (a[1][1]==0)  
{
printf("Cac phan tu duong cheo cua ma tran bang khong");  
getch();  
exit(1);  
}
else  
{
if (fabs(a[k][k]/a[1][1])<(1e-08))  
{
printf("Ma tran suy bien");  
goto mot;  
}
}
for (i=k+1;i<=n;i++)  
{
if (a[k][k]==0)  
{
41  
printf("Cac phan tu duong cheo cua ma tran bang khong\n");  
goto mot;  
}
s=a[i][k]/a[k][k];  
a[i][k]=0;  
for (j=k+1;j<=n+1;j++)  
a[i][j]=a[i][j]-s*a[k][j];  
}
y[n]=a[n][n+1]/a[n][n];  
for (i=n-1;i>=1;i--)  
{
s=a[i][n+1];  
for (j=i+1;j<=n;j++)  
s=s-a[i][j]*y[j];  
if (a[i][i]==0)  
{
printf("Cac phan tu duong cheo cua ma tran bang khong\n");  
goto mot;  
}
y[i]=s/a[i][i];  
}
}
if (r!=1)  
for (i=1;i<=n;i++)  
{
if (fabs(y[i])<e*fabs(x[i]))  
goto ba;  
}
for (i=1;i<=n;i++)  
x[i]=x[i]-y[i];  
printf("\n");  
}
printf("Khong hoi tu sau %d lan lap\n",z);  
goto mot;  
clrscr();  
ba:printf("Vec to nghiem\n");  
for (i=1;i<=n;i++)  
printf("%.5f\n",(x[i]-y[i]));  
printf("\n");  
printf("Do chinh xac cua nghiem la %.5f: \n", e);  
printf("\n");  
printf("Vec to tri so du :\n");  
for (i=1;i<=n;i++)  
printf("%.5f\n",(a[i][n+1]));  
mot:printf("\n");  
getch();  
}
void doc()  
{
a[1][1]=3*x[1]*x[1]-3*x[2]*x[4];  
a[1][2]=-3*x[2]*x[2]-3*x[1]*x[4];  
a[1][3]=0;  
a[1][4]=-3*x[1]*x[2];  
a[1][5]=x[1]*x[1]*x[1]-x[2]*x[2]*x[2]-3*x[1]*x[2]*x[4]-8;  
42  
a[2][1]=1;  
a[2][2]=1;  
a[2][3]=1;  
a[2][4]=1;  
a[2][5]=x[1]+x[2]+x[3]+x[4]-5;  
a[3][1]=-x[1]/sqrt(25-x[1]*x[1]);  
a[3][2]=0;  
a[3][3]=8;  
a[3][4]=0;  
a[3][5]=sqrt(25-x[1]*x[1])+8*x[3]+4;  
a[4][1]=2*x[2]*x[3];  
a[4][2]=2*x[1]*x[3];  
a[4][3]=2*x[1]*x[2];  
a[4][4]=-1;  
a[4][5]=2*x[1]*x[2]*x[3]-x[4]+8;  
}
43  
doc 9 trang Thùy Anh 10940 Free
Download
Bạn đang xem tài liệu "Giáo trình Phương pháp tính - Chương 2, Phần 2: Giải gần đúng phương trình đại số và siêu việt", để tải tài liệu gốc về máy hãy click vào nút Download ở trên

File đính kèm:

  • docgiao_trinh_phuong_phap_tinh_chuong_2_phan_2_giai_gan_dung_ph.doc