数值计算方法上机试验报告

华北电力大学试验报告||试验名称 数值计算方法》上机试验课程名称 数值计算方法 专业班级：电力实08 指导教师：郝育黔教师学号：202301001008 学生姓名：李超然成绩： 试验日期：2023年04月数值计算方法上机试验报告一、各算法的算法原理及计算机程序框图1、牛顿法求解非线性方程(1)算法原理：* k kf(x)0X的一个近似值X，将f(x)X处开放成一阶泰勒公* k k式” f“() 2f(X)f(Xk)f(Xk)(XXk)----------------------(xXk)2!无视高次项，有f(X)f(Xk)f”(Xk)(XXk)* 右端是直线方程，用这个直线方程来近似非线性方程 f(X)。将非线性方程f(X)0的根X代入f(X) 0* f(Xk)f”(Xk)(x*Xk)0解出将右端取为Xk1，则Xk 1是比Xk更接近于X的近似值，即X X f(Xk)kk1k这就是牛顿迭代公式。计算机程序框图：(见)

f”(Xk)输入变量、输出变量说明：o输入变量：X迭代初值，迭代精度，N迭代最大次数输出变量：k当前迭代次数，X1当前迭代值o具体算例及求解结果：完毕例：导出计算匸〔c0〕的牛顿迭代公式，并计算115。〔P392-16〕求解结果：10.75000010.72383710.72380510.7238052、列主元素消去法求解线性方程组算法原理：高斯消去法是利用现行方程组初等变换中的一种变换，即用一个不为零的数乘一个方程后加只另一个方程，使方程组变成同解的上三角方程组，然后再自下而上对上三角方程组求解。kk列选主元是当高斯消元到第k步时，从k列的a 以下〔包括a 〕的各元素中选出绝对kkkk值最大的，然后通过行交换将其交换到akk的位置上。交换系数矩阵中的两行〔包括常数项〕，只相当于两个方程的位置交换了，因此，列选主元不影响求解的结果。计算机程序框图：〔见下页〕输入变量、输出变量说明：输入变量：aj系数矩阵元素，b常向量元素输出变量：bi,b2bn解向量元素具体算例及求解结果：例：用列选主元法求解以下线性方程组〔课本P653-3〕0.50x11.10x23.10X36.002.00x14.50x20.36X30.0205.00x10.96X26.50X30.96求解结果：x12.600000x21.000000x32.0000003、LU分解法求解线性方程组算法原理：AxbALULUxb，这时可归用递推计算相继求解两个三角形〔系数矩阵为三角矩阵〕方程组，用顺代，由LybyUxyx。

结为利计算机程序框图：〔见下页〕输入变量、输出变量说明：输入变量：aj系数矩阵元素，b常向量元素从主程序来开头i,j1,2,..,n

读入数据,a j a U1idi,i1,2,...,nU1idi,i1,2,...,nli1i,-U2,3,...,n1Ur1riarilU,irr kir,r1,..,nlir(ak1r1.irlU)/U,iik krrrr 1,r2,...,nk1yi1yib,ybik1n

2,3,...,nX y/U ,X (% UXj/uJ n1,,...,2,1n n nn

ikki1完毕n,输出变量：b1,b2,...,b解向量元素n,(4)具体算例及求解结果:例：用杜里特尔分解法求解方程组(课本P743-8)2 2 3 x1 34 7 7 x2 12 4 5 x3 7求解结果：x12.000000X22.000000x31.0000004、拉格朗日插值法(1)算法原理：|构造基函数|k

(x)

i0xkik

,可以证明基函数满足以下条件:xxlk(x)

0 ik1ik，1对于给定(n1)个节点, n次拉格朗日插值多项式由下式给出:XXkXXL(x) nn ykXXk0i0 k iklk(x)xnL(x)xn次的代数多项式。xXi时，L(xJy，满足插值条件。计算机程序框图：(见下页)输入变量、输出变量说明:输入变量：(x，%)插值节点输出变量：y插值所得到被插函数在插值点的近似值具体算例及求解结果:f(x)sinx的值如下表所示。f(x)sinx的值XX06丄2432sinx0占212sin—的估量值12求解结果：0.2585885、最小二乘法的曲线拟合(1)算法原理：对于给定的一组数据(x，f(xj),i 1,2,...,m，要在给定的函数空间Span{ , }0 1,…，

n*(x)a；o(x)a；i(x)...a；n(x)

na；i(x)n使(x)满足m2 *(X)f(X)]mi2

m(X)f(X)] 22 i i i i的

〔X〕的方法称为曲线拟合的最小二乘法，〔X〕称为最小二乘法最小二乘解。计算机程序框图:输入变量、输出变量说明:输入变量：〔x，yj数据点输出变量：ai拟合多项式的系数(4)具体算例及求解结果：yf(x)的实例数据表，试用最小二乘法求二次拟合多项式。P1863)Xi0123456yi15141414141516求解结果:a014.928572

(课本a,0.892857a,a2 0.178571y1.3181713.431811x0.386363X26、变步长梯形求积分算法原理:设将积分区间［a,b］nnXk步长

akh,k0,1,...,n，其中XkXn对于子区间［Xk,Xk1］，利用体型求其积分近似值对于子区间［a,b］有

2【f(Xk) f(Xk1)]T n1h[f(X)f(X)]n 2k0 k k1对于子区间［Xk,Xk1］再取其中点1(X X )rk k k1作节点，此时区间数增加了一倍为 ，2n对子区间［Xk,Xk1］，其积分近似值4[f(x

) 2f(x )mk k1k2对区间［a,b］有

n1h山

4[5)

2f(x)g)]0 k2h()n1 f(X()

h)] n1 f(x )计算机程序框图:

[fX ki k124ko k k2bah,£[f(a) f(b)] T>0S,ahx2T1T22 2Sf(x)T1T22 2完毕输入变量、输出变量说明：输入变量：［a,b］积分区间，精度输出变量：T2积分结果具体算例及求解结果:例：用变步长梯形公式求积法计算sin例：用变步长梯形公式求积法计算sinx10x求解结果：0.94608277、改进欧拉法〔1〕算法原理：h取值较小时，让梯形法的迭代公式只迭代一次就完毕。这样先用欧拉公式求得初步近似值n〔〕，称之为预报值，预报值的精度不高，用它替代梯形法右端的ni，yn1，并称之为校正值，这时得到预报-校正公式。将预报-校正公式Yn1Yn1〔0〕 nYhf〔X,y〕n nYn1 nhYf〔X,Y〕fX ,Y〕2nn1n1〔0〕称为改进欧拉公式。〔2〕计算机程序框图：〔见下页〕〔3〕输入变量、输出变量说明：输入变量：〔X,Y〕处置点，h区间长度，N计算次数0 0输出变量：〔论,力〕初值问题的数值解法结果〔4〕具体算例及求解结果：例：求解初值问题〔课本P2427-2〕2xyy——,0x1YY〔0〕1求解结果:XnYny(Xn)XnYny(Xn)0.1P1.095909:1.0959090.61.4859561.4859550.21.1840971.1840970.71.5625141.5525140.3[1.2662011.2662010.81.6164751.6164740.4r1.343360:1.3433600.91.6783201.6781660.51.4164021.4164021.01.7378671.737867C始)1 1读入C始)1 1读入X,y,h,N011ony。hf(Xo,y°)yy。hf] p c(X,y)y] p cp (yp ]]/Xy7]]r / -------------1 n1n完毕8四阶龙格-库塔法求解常微分方程的初值问题算法原理：Yn1Ynh(X)4k4kYn1Ynh(X)4k4k1f(X,y)2k23k：3n nk2f(Xn1h,ynk11k1h)3f(Xn2h,yn21 ikhk22k2h)4f(Xn3h,yn31k1h32k2h33k3h)k ki其中，i,i,j(i1,2,3,j1,2,...,i)均为待定系数。k2,k3,k4Xnhyni并进展花间，然后与y(X )在X点上的泰勒开放式比较，使其两式比较，使其两式右端直到h的系数相等，ni n 4经过简洁的数学演算可得到关于i,i,ij的一组特解1211221221313203121122122131320333143211613从而得到以下常用的经典公式yiyn2k3k)412k2kif(X,y)n nk2f(x,y1 nn2k3f(x』n1 nk4f(X ,yn1 nhk)35

截断误计算机程序框图：(见下页)输入变量、输出变量说明：o输入变量：(x,y。)处置点，h区间长度，N计算次数输出变量：(xi,yi)初值问题的数值解法结果o具体算例及求解结果：h0.2x0x1，用经典公式求解初值问题2xy(0)iXnyXnyny(X)n0.20.40.60.81.01.1832291.3416671.4832811.6125141.7321421.1832291.3416671.4832811.6125141.732142二、上机体验与收获本次上机内容为牛顿法求解非线性方程、列主元素消去法求解线性方程组、LU分解法求解线性方程组、拉格朗日插值、最小二乘法的曲线拟合、变步长梯形求积分、改进欧拉方法求 5解常微分方程的初值问题、四阶龙格—库塔法求解常微分方程的初值问题在各个算法程序的编制中，我生疏到从一个数值计算方法到一个具体程序的简洁过程，使我对各种数值计算方法有了更深的了解。同时，补充了C++吾言的应用，把握了C++吾言中文件的操作，这是原先