數(shù)值計算方法編程作業(yè)C語言版講解

1、1:第二章(1)二分法求解非線性方程： #include#include #define f(x) (x*x-1)*x-1)void main() float a,b,x,eps;int k=0;printf(intput epsn);/* 容許誤差 */scanf(%f,&eps);printf(a,b=n);for(;)scanf(%f, %f,&a ,&b);if(f(a)*f(b)=0) /* 判斷是否符合二分法使用的條件 */printf( 二分法不可使用 ,請重新輸入 :n);else break;do x=(a+b)/2;k+;if(f(a)*f(x)0) /* 如果 f(a)*

2、f(x)0) /* 否則根在區(qū)間的右半部分 */a=x;else break;while(fabs(b-a)eps);/* 判斷是否達(dá)到精度要求 ,若沒有達(dá)到 ,繼續(xù)循環(huán) */ x=(a+b)/2; /* 取最后的小區(qū)間中點(diǎn)作為根的近似值 */ printf(n The root is x=%f, k=%dn,x,k);運(yùn)行結(jié)果：intput eps0.00001a,b=2,-5The root is x=1.324721, k=20Press any key to continue 總結(jié)：本題關(guān)鍵在于兩個端點(diǎn)的取值和誤差的判斷，此程序較容易。二分法收斂速度較快， 但缺點(diǎn)是只能求解單根。(2)

3、牛頓法求解非線性方程：#include #include float f(float x) /* 定義函數(shù) f(x) */ return(-3*x+4)*x-5)*x+6; float f1(float x) /* 定義函數(shù) f(x) 的導(dǎo)數(shù) */ return (-9*x+8)*x-5; void main() float eps,x0,x1=1.0;printf(input eps:n);scanf(%f,&eps); /* 輸入容許誤差 */do x0=x1; /* 準(zhǔn)備下一次迭代的初值 */x1=x0-f(x0)/f1(x0); /* 牛頓迭代 */while(fabs(x1-x0)e

4、ps); /* 當(dāng)滿足精度 ,輸出近似根 */ printf(x=%fn,x1);程序運(yùn)行結(jié)果：x=1.265328總結(jié)：關(guān)鍵是牛頓迭代的應(yīng)用， 程序中最大缺點(diǎn)是函數(shù)及其導(dǎo)數(shù)已唯一給出確定不可求的隨 意函數(shù)的根，牛頓法比二分法收斂快，可以求重根。2:第三章(1)列主元素消去法求解線性方程： #include#include#define N 20using namespace std;void load();float aNN;int m;int main()int i,j;int c,k,n,p,r;float xN,lNN,s,d;coutm;coutendl;cout 請按順序輸入增廣矩

5、陣 a:endl;load();for(i=0;im;i+)for(j=i;jfabs(aii)?j:i;/* 找列最大元素 */for(n=0;nm+1;n+)s=ain; ain=acn; acn=s;/* 將列最大數(shù)防在對角線上 */for(p=0;pm+1;p+)coutaipt;coutendl;for(k=i+1;km;k+)lki=aki/aii;for(r=i;r=0;i-)d=0;for(j=i+1;jm;j+)d=d+aij*xj; xi=(aim-d)/aii;cout 該方程組的解為 :endl; for(i=0;im;i+) coutxi=xit; /system(p

6、ause); return 0;void load()int i,j;for(i=0;im;i+)for(j=0;jaij;運(yùn)行結(jié)果：/* 化成三角陣 */*求解 */1 2 3 45 1 0 84 6 9 246920-6.5-11.255.50-1.86265e-008-0.115385該方程組的解為請按順序輸入增廣矩陣 a:3.92308下面請輸入未知數(shù)的個數(shù) m=3 x0=-9.99999 x1=58 x2=-34 Press any key to continue 總結(jié)：列主元素消去法的目的是為了防止減去一個較小的數(shù)時大數(shù)淹沒小數(shù)， 而使結(jié)果產(chǎn)生 較大誤差， 本程序關(guān)鍵在每次消元時找

7、到相應(yīng)列中的最大項， 然后交換兩行位置， 在進(jìn)行計 算。(2)LU 分解法求解線性方程：#includevoid solve(float l100,float u100,float b,float x,int n)int i,j;float t,s1,s2;float y100; for(i=1;i=n;i+) /* 第一次回代過程開始 */s1=0;for(j=1;j=1;i-) /* 第二次回代過程開始 */ s2=0;for(j=n;ji;j-)t=-uij; s2=s2+t*xj; xi=(yi+s2)/uii;void main()float a100100,l100100,u100

8、100,x100,b100;int i,j,n,r,k;float s1,s2;1*/for(i=1;i=99;i+)/* 將所有的數(shù)組置零，同時將 L 矩陣的對角值設(shè)為 for(j=1;j=99;j+)lij=0,uij=0;if(j=i) lij=1;printf (input n:n);/* 輸入方程組的個數(shù) */ scanf(%d,&n);printf (input array A:n);/* 讀取原矩陣 A*/ for(i=1;i=n;i+)for(j=1;j=n;j+) scanf(%f,&aij);printf (input array B:n);/* 讀取列矩陣 B*/for(

9、i=1;i=n;i+)scanf(%f,&bi);for(r=1;r=n;r+)/* 求解矩陣 L for(i=r;i=n;i+)s1=0;for(k=1;k=r-1;k+)s1=s1+lrk*uki;uri=ari-s1; for(i=r+1;i=n;i+)s2=0;for(k=1;k=r-1;k+)s2=s2+lik*ukr;lir=(air-s2)/urr;printf(array L:n);/* 輸出矩陣for(i=1;i=n;i+)for(j=1;j=n;j+)printf(%7.3f ,lij);printf(n);printf(array U:n);/* 輸出矩陣for(i=1;

10、i=n;i+)for(j=1;j=n;j+)printf(%7.3f ,uij);printf(n);solve(l,u,b,x,n);printf( 解為 :n); for(i=1;i=n;i+) printf(x%d=%fn,i,xi);運(yùn)行結(jié)果：input n:3input array A:2 2 34 7 7-2 4 5和 U*/L*/U*/input array B:3 1 -7array L:1.0000.0000.0002.0001.0000.000-1.0002.0001.000array U:2.0002.0003.0000.0003.0001.0000.0000.0006.

11、000解為 :x1=2.000000 x2=-2.000000x3=1.000000Press any key to continue總結(jié)：關(guān)鍵是把矩陣分解為 L、 U 兩個三角矩陣，回代過程比較簡單。3:第四章(1) 拉格朗日差值多項式 ; #include#include#define MAX 100void main() int i,j,k,m,n,N,mi;float tmp,mx;float XMAXMAX,YMAX,xMAX,yMAX,aMAX;printf(n 輸入擬合多項式的次數(shù) :n); scanf(%d,&m);printf(n 輸入給定點(diǎn)的個數(shù) n 及坐標(biāo) (x,y):n

12、); scanf(%d,&N);printf(n); for(i=0;iN;i+)scanf(%f,%f,&xi,&yi); for(i=0;i=m;i+)for(j=i;j=m;j+)tmp=0; for(k=0;kN;k+) tmp=tmp+pow(xk,(i+j);Xij=tmp;Xji=Xij;for(i=0;i=m;i+)tmp=0;for(k=0;kN;k+)tmp=tmp+yk*pow(xk,i);Yi=tmp;for(j=0;jm;j+)for(i=j+1,mi=j,mx=fabs(Xjj);imx)mi=i;mx=fabs(Xij);if(jmi)tmp=Yj;Yj=Ymi;

13、Ymi=tmp;for(k=j;k=m;k+)tmp=Xjk;Xjk=Xmik;Xmik=tmp;for(i=j+1;i=m;i+)tmp=-Xij/Xjj;Yi+=Yj*tmp;for(k=j;k=0;i-)ai=Yi;for(j=i+1;j=m;j+)ai-=Xij*aj;ai/=Xii;printf(n 所求的二次多項式為 :n);printf(P(x)=%f,a0); for(i=1;i=m;i+) printf(+(%f)*x%d,ai,i);運(yùn)行結(jié)果：輸入擬合多項式的次數(shù) :5輸入給定點(diǎn)的個數(shù) n 及坐標(biāo) (x,y):31,25,34,2所求的二次多項式為 :P(x)=1.9804

14、17+(0.282759)*x1+(-0.299937)*x2+(0.022071)*x3+(0.016624)*x4+(-0.0 01934)*x5Press any key to continue 總結(jié)：拉格朗日計算公式中，只需要知道各個點(diǎn)即可4:第五章(1)曲線擬合：#include#include#define MAX 100void main() int i,j,k,m,n,N,mi;float tmp,mx;float XMAXMAX,YMAX,xMAX,yMAX,aMAX;printf(n 輸入擬合多項式的次數(shù) :n);scanf(%d,&m);printf(n 輸入給定點(diǎn)的個數(shù)

15、 n 及坐標(biāo) (x,y):n);scanf(%d,&N);printf(n);for(i=0;iN;i+)scanf(%f,%f,&xi,&yi);for(i=0;i=m;i+) for(j=i;j=m;j+)tmp=0; for(k=0;kN;k+) tmp=tmp+pow(xk,(i+j); Xij=tmp;Xji=Xij;for(i=0;i=m;i+)tmp=0;for(k=0;kN;k+)tmp=tmp+yk*pow(xk,i);Yi=tmp;for(j=0;jm;j+)for(i=j+1,mi=j,mx=fabs(Xjj);imx)mi=i;mx=fabs(Xij);if(jmi)t

16、mp=Yj;Yj=Ymi;Ymi=tmp;for(k=j;k=m;k+)tmp=Xjk;Xjk=Xmik;Xmik=tmp;for(i=j+1;i=m;i+)tmp=-Xij/Xjj;Yi+=Yj*tmp;for(k=j;k=0;i-)ai=Yi;for(j=i+1;j=m;j+)ai-=Xij*aj;ai/=Xii;printf(n 所求的二次多項式為 :n); printf(P(x)=%f,a0);for(i=1;i=m;i+) printf(+(%f)*x%d,ai,i);輸入擬合多項式的次數(shù) :2輸入給定點(diǎn)的個數(shù) n 及坐標(biāo) (x,y):51,25,32,48,3-1,5所求的二次多項

17、式為 :P(x)=3.952280+(-0.506315)*x1+(0.050877)*x2Press any key to continue5:第六章(1)辛普生求積方法：#include #define N 16/* 等分?jǐn)?shù) */float func(float x) float y;y=4.0/(1+x*x);return(y);void gedianzhi(float y,float a,float h) int i;for(i=0;i=N;i+) yi=func(a+i*h);float simpson(float y,float h) float s,s1,s2;int i;s1=

18、y1;s2=0.0;for(i=2;i=N-2;i=i+2) s1+=yi+1; /* 計算奇數(shù)項的函數(shù)值之和 */s2+=yi; /* 計算偶數(shù)項的函數(shù)值之和 */ s=y0+yN+4.0*s1+2.0*s2;return(s*h/3.0);main() float a,b,h,s,fN+1; scanf(%f,%f,&a,&b);h=(b-a)/( float)N;gedianzhi(f,a,h);s=simpson(f,h);printf(s=%fn,s);運(yùn)行結(jié)果：1,3 s=1.854590Press any key to continue總結(jié)：辛普生算法是一種積分方法，采用三點(diǎn)法插

19、值，如果 h 較小的話，誤差很小，因為它 b a h 4的插值余項 R( f ) b a h f (4)( ) ，辛普生算法比較精確，程序關(guān)鍵是對所取的點(diǎn) 180 2的取和，注意6:第七章(1)改進(jìn)歐拉法求解常微分方程的初值問題#include float func(float x,float y) return(y-x);float euler(float x0,float xn,float y0,int N) float x,y,yp,yc,h;int i;x=x0;y=y0;h=(xn-x0)/(float)N;for(i=1;i=N;i+) yp=y+h*func(x,y);x=x0+

20、i*h; yc=y+h*func(x,yp);y=(yp+yc)/2.0;return(y);main() float x0,xn,y0,e;int n;printf(ninput n:n);scanf(%d,&n);printf(input x0,xn:n );scanf(%f,%f,&x0,&xn);printf(input y0:n ); scanf(%f,&y0);e=euler(x0,xn,y0,n); printf(y(%f)=%6.4f,y0,e);input n:20 input x0,xn:1,6 input y0:2 y(2.000000)=7.0000Press any

21、key to continue (2)四階龍格 庫塔法 #include float func(float x,float y) return(x-y);float runge_kutta(float x0,float xn,float y0,int N) float x,y,y1,y2,h,xh;float d1,d2,d3,d4;int i;x=x0;y=y0;h=(xn-x0)/(float)N;for(i=1;i=N;i+) xh=x+h/2;d1=func(x,y);d2=func(xh,y+h*d1/2.0);d3=func(xh,y+h*d2/2.0); d4=func(xh,y+h*d3);y=y+h*(d1+2*d2+2*d3+d4)/6.0;x=x0+i*h; return(y);main() float x0,xn