数值分析上机第四次作业.doc
《数值分析上机第四次作业.doc》由会员分享,可在线阅读,更多相关《数值分析上机第四次作业.doc(8页珍藏版)》请在咨信网上搜索。
数值分析上机第四次作业 实验项目:共轭梯度法求解对称正定得线性方程组 实验内容:用共轭梯度法求解下面方程组 (1) 迭代20次或满足时停止计算。 (2) ,A就是1000阶得Hilbert矩阵或如下得三对角矩阵, A[i,i]=4,A[i,i-1]=A[i-1,i]=-1,i=2,3,、、,n b[1]=3, b[n]=3, b[i]=2,i=2,3,…,n-1 迭代10000次或满足时停止计算。 (3)*考虑模型问题,方程为 用正方形网格离散化,若取,得到得线性方程组,并用共轭梯度法(CG法)求解,并对解作图。 实验要求:迭代初值可以取,计算到停止.本题有精确解,这里表示以为分量得向量,表示在相应点上取值作为分量得向量. 实验一: (1) 编制函数子程序CGmethod。 function [x,k]=CGmethod(A,b) n=length(A);x=zeros(n,1);r=b-A*x;rho=r'*r; k=0; while rho>10^(-12) & k<20 k=k+1; if k==1 p=r; else beta=rho/rho1; p=r+beta*p; end w=A*p; alpha=rho/(p'*w); x=x+alpha*p; r=r-alpha*w; rho1=rho; rho=r'*r; end 编制主程序shiyan1_1: clear,clc A=[2,-1,0,0;-1,3,-1,0;0,-1,4,1;0,0,-1,5]; b=[3,-2,1,5]'; [x,k]=CGmethod(A,b) 运行结果为: x = 1、3882 -0、2855 -0、0222 0、9367 k = 20 (2) 编制函数子程序CGmethod_1 function [x,k]=CGmethod_1(A,b) n=length(A);x(1:n,1)=0;r=b-A*x;r1=r; k=0; while norm(r1,1)>=10^(-7)&k<10^4 k=k+1; if k==1 p=r; else beta=(r1'*r1)/(r'*r);p=r1+beta*p; end r=r1; w=A*p; alpha=(r'*r)/(p'*w); x=x+alpha*p; r1=r-alpha*w; end 编制主程序shiyan1_2: clear,clc n=1000; A=hilb(n); b=sum(A')'; [x,k]=CGmethod_1(A,b) 运行结果为:x得值,均接近1,迭代次数k=32 实验二 实验目得:用复化Simpson方法、自适应复化梯形方法与Romberg方法求数值积分。 实验内容:计算下列定积分 (1) (2) (3) 实验要求: (1)分别用复化Simpson公式、自适应复化梯形公式计算要求绝对误差限为,输出每种方法所需得节点数与积分近似值,对于自适应方法,显示实际计算节点上离散函数值得分布图; (2)分析比较计算结果。 2、实验目得:高斯数值积分方法用于积分方程求解。 实验内容:线性得积分方程得数值求解,可以被转化为线性代数方程组得求解问题。而线性代数方程组所含未知数得个数,与用来离散积分得数值方法得节点个数相同。在节点数相同得前提下,高斯数值积分方法有较高得代数精度,用它通常会得到较好得结果。对第二类Fredholm积分方程 首先将积分区间[a,b]等分成n份,在每个子区间上离散方程中得积分就得到线性代数方程组。 实验要求:分别使用如下方法,离散积分方程中得积分 1、复化梯形方法;2、复化辛甫森方法;3、复化高斯方法。求解如下得积分方程 ,方程得准确解为, 并比较各算法得优劣。 实验二 1、复化Simpson方法) 输入积分区间下限0 输入积分区间上限2 输入等分份数20 输入被积函数(以x为自变量)x^6/10-x^2+x S = 1、1619 输入积分区间下限0 输入积分区间上限1 输入等分份数20 输入被积函数(以x为自变量)x*sqrt(x) S = 0、4000 输入积分区间下限5 输入积分区间上限200 输入等分份数20 输入被积函数(以x为自变量)1/sqrt(x) S = 23、8218 2、自动变步长Simpson方法 函数1: 输入积分区间下限0 输入积分区间上限2 输入为课本得第几个函数(第一个这输入1):1 S =1、619(过程省略) i = 19 函数2: 输入积分区间下限0 输入积分区间上限1 输入为课本得第几个函数(第一个这输入1):2 S =0、4(过程省略) i = 17 函数3: 输入积分区间下限5 输入积分区间上限200 输入为课本得第几个函数(第一个这输入1):3 S=23、8121(过程省略) i = 111 编制程序如下: Clear,clc syms x a=input('输入积分区间下限'); b=input('输入积分区间上限'); n=input('输入等分份数'); ff=input('输入被积函数(以x为自变量)'); h=(b-a)/n; f=inline(ff,'x'); sum1=0; sum2=0; for i=0:n-1 sum1=sum1+f(a+i*h+0、5*h); end for i=1:n-1 sum2=sum2+f(a+i*h); end for i=0:2*n X(i+1,1)=f((b-a)*i/(n*2)+a); end S=h/6*(f(a)+4*sum1+2*sum2+f(b)) function S = zdsps( n ) a=0; b=1; h=(b-a)/4; f=inline('x^(3/2)','x'); sum1=0; sum2=0; for i=0:n-1 sum1=sum1+f(a+i*h+0、5*h); end for i=1:n-1 sum2=sum2+f(a+i*h); end for i=0:2*n x(i+1,1)=f((b-a)*i/(n*2)+a); end S=h/6*(f(a)+4*sum1+2*sum2+f(b)); end function S = zpsgs(a,b,n,ff ) h=(b-a)/n; sum1=0; sum2=0; sum3=0; sum4=0; if ff==1 f=inline('x^6/10-x^2+x','x'); end if ff==2 f=inline('x^(3/2)','x'); end if ff==3 f=inline('1/sqrt(x)','x'); end for i=0:n-1 sum1=sum1+f(a+i*h+0、25*h); sum2=sum2+f(a+i*h+0、75*h); sum4=sum4+f(a+i*h+0、5*h); end for i=1:n-1 sum3=sum3+f(a+i*h); end for i=0:4*n x(i+1,1)=f((b-a)*i/(n*4)+a); end S=h/(6*2)*(f(a)+4*sum1+4*sum2+2*(sum3+sum4)+f(b)); end clear, clc a=input('输入积分区间下限'); b=input('输入积分区间上限'); ff=input('输入为课本得第几个函数(第一个这输入1):'); for i=2:300 S(i)=zpsgs(a,b,(i),ff); S(i+1)=zpsgs(a,b,(i+1),ff); if abs(S(i+1)-S(i))<0、5*10^(-7) break end end S %所求积分值 i %所分份数 实验三 1、对常微分方程初值问题 分别使用Euler显示方法、改进得Euler方法与经典RK法与四阶Adams预测-校正算法,求解常微分方程数值解,并与其精确解进行作图比较。其中多步法需要得初值由经典RK法提供。 2、实验目得:Lorenz问题与混沌 实验内容:考虑著名得Lorenz方程 其中s, r, b为变化区域有一定限制得实参数。该方程形式简单,表面上瞧并无惊人之处,但由该方程揭示出得许多现象,促使“混沌”成为数学研究得崭新领域,在实际应用中也产生了巨大得影响。 实验方法:先取定初值Y0=(x, y, z)=(0, 0, 0),参数s=10, r=28, b=8/3,用MATLAB编程数值求解,并与MATLAB函数ods45得计算结果进行对比。 实验要求: (1)对目前取定得参数值s, r与b,选取不同得初值Y0进行运算,观察计算得结果有什么特点?解得曲线就是否有界?解得曲线就是不就是周期得或趋于某个固定点? (2)在问题允许得范围内适当改变其中得参数值s, r, b,再选取不同得初始值Y0进行试算,观察并记录计算得结果有什么特点?就是否发现什么不同得现象? 3、 定义函数子程序为: function z=f(x,y) z=-y+2*cos(x); return 主程序为: clear,clc b=pi; a=0; n=100; y(1)=1; h=(b-a)/n; x=a:h:b; for i=1:n y(i+1)=y(i)+h*f(x(i),y(i)); end t1=plot(x,y,'r-') hold on for i=1:n K1=f(x(i),y(i)); K2=f(x(i+1),y(i)+h*K1); y(i+1)=y(i)+h*(K1+K2)/2; end t2=plot(x,y,'b+') for i=1:n K1=f(x(i),y(i)); K2=f(x(i)+0、5*h,y(i)+0、5*h*K1); K3=f(x(i)+0、5*h,y(i)+0、5*h*K2); K4=f(x(i),y(i)+h*K3); y(i+1)=y(i)+h*(K1+2*K2+2*K3+K4)/6; end t3=plot(x,y,'ko') for i=1:3 K1=f(x(i),y(i)); K2=f(x(i)+0、5*h,y(i)+0、5*h*K1); K3=f(x(i)+0、5*h,y(i)+0、5*h*K2); K4=f(x(i),y(i)+h*K3); y(i+1)=y(i)+h*(K1+2*K2+2*K3+K4)/6; end for i=4:n z=y(i)+h/24*(55*f(x(i),y(i))-59*f(x(i-1),y(i-1))、、、 +37*f(x(i-2),y(i-2))-9*f(x(i-3),y(i-3))); y(i+1)=y(i)+h/24*(9*f(x(i+1),z)+19*f(x(i),y(i))、、、 -5*f(x(i-1),y(i-1))+f(x(i-2),y(i-2))); end t4=plot(x,y,'g*') t5=ezplot('sin(x)+cos(x)',[0,pi]) xlabel('x轴','FontWeight','bold') ylabel('y轴', 'FontWeight','bold') legend([t1,t2,t3,t4,t5],'向前Euler法','¸改进得Euler法','经典四阶Runge-Kutta法','四阶Adams公式','精确解') 原图为: 局部放大图为: 由图可得:四阶Adams公式及改进得欧拉法将为精确。- 配套讲稿:
如PPT文件的首页显示word图标,表示该PPT已包含配套word讲稿。双击word图标可打开word文档。
- 特殊限制:
部分文档作品中含有的国旗、国徽等图片,仅作为作品整体效果示例展示,禁止商用。设计者仅对作品中独创性部分享有著作权。
- 关 键 词:
- 数值 分析 上机 第四 作业
咨信网温馨提示:
1、咨信平台为文档C2C交易模式,即用户上传的文档直接被用户下载,收益归上传人(含作者)所有;本站仅是提供信息存储空间和展示预览,仅对用户上传内容的表现方式做保护处理,对上载内容不做任何修改或编辑。所展示的作品文档包括内容和图片全部来源于网络用户和作者上传投稿,我们不确定上传用户享有完全著作权,根据《信息网络传播权保护条例》,如果侵犯了您的版权、权益或隐私,请联系我们,核实后会尽快下架及时删除,并可随时和客服了解处理情况,尊重保护知识产权我们共同努力。
2、文档的总页数、文档格式和文档大小以系统显示为准(内容中显示的页数不一定正确),网站客服只以系统显示的页数、文件格式、文档大小作为仲裁依据,个别因单元格分列造成显示页码不一将协商解决,平台无法对文档的真实性、完整性、权威性、准确性、专业性及其观点立场做任何保证或承诺,下载前须认真查看,确认无误后再购买,务必慎重购买;若有违法违纪将进行移交司法处理,若涉侵权平台将进行基本处罚并下架。
3、本站所有内容均由用户上传,付费前请自行鉴别,如您付费,意味着您已接受本站规则且自行承担风险,本站不进行额外附加服务,虚拟产品一经售出概不退款(未进行购买下载可退充值款),文档一经付费(服务费)、不意味着购买了该文档的版权,仅供个人/单位学习、研究之用,不得用于商业用途,未经授权,严禁复制、发行、汇编、翻译或者网络传播等,侵权必究。
4、如你看到网页展示的文档有www.zixin.com.cn水印,是因预览和防盗链等技术需要对页面进行转换压缩成图而已,我们并不对上传的文档进行任何编辑或修改,文档下载后都不会有水印标识(原文档上传前个别存留的除外),下载后原文更清晰;试题试卷类文档,如果标题没有明确说明有答案则都视为没有答案,请知晓;PPT和DOC文档可被视为“模板”,允许上传人保留章节、目录结构的情况下删减部份的内容;PDF文档不管是原文档转换或图片扫描而得,本站不作要求视为允许,下载前自行私信或留言给上传者【天****】。
5、本文档所展示的图片、画像、字体、音乐的版权可能需版权方额外授权,请谨慎使用;网站提供的党政主题相关内容(国旗、国徽、党徽--等)目的在于配合国家政策宣传,仅限个人学习分享使用,禁止用于任何广告和商用目的。
6、文档遇到问题,请及时私信或留言给本站上传会员【天****】,需本站解决可联系【 微信客服】、【 QQ客服】,若有其他问题请点击或扫码反馈【 服务填表】;文档侵犯商业秘密、侵犯著作权、侵犯人身权等,请点击“【 版权申诉】”(推荐),意见反馈和侵权处理邮箱:1219186828@qq.com;也可以拔打客服电话:4008-655-100;投诉/维权电话:4009-655-100。
1、咨信平台为文档C2C交易模式,即用户上传的文档直接被用户下载,收益归上传人(含作者)所有;本站仅是提供信息存储空间和展示预览,仅对用户上传内容的表现方式做保护处理,对上载内容不做任何修改或编辑。所展示的作品文档包括内容和图片全部来源于网络用户和作者上传投稿,我们不确定上传用户享有完全著作权,根据《信息网络传播权保护条例》,如果侵犯了您的版权、权益或隐私,请联系我们,核实后会尽快下架及时删除,并可随时和客服了解处理情况,尊重保护知识产权我们共同努力。
2、文档的总页数、文档格式和文档大小以系统显示为准(内容中显示的页数不一定正确),网站客服只以系统显示的页数、文件格式、文档大小作为仲裁依据,个别因单元格分列造成显示页码不一将协商解决,平台无法对文档的真实性、完整性、权威性、准确性、专业性及其观点立场做任何保证或承诺,下载前须认真查看,确认无误后再购买,务必慎重购买;若有违法违纪将进行移交司法处理,若涉侵权平台将进行基本处罚并下架。
3、本站所有内容均由用户上传,付费前请自行鉴别,如您付费,意味着您已接受本站规则且自行承担风险,本站不进行额外附加服务,虚拟产品一经售出概不退款(未进行购买下载可退充值款),文档一经付费(服务费)、不意味着购买了该文档的版权,仅供个人/单位学习、研究之用,不得用于商业用途,未经授权,严禁复制、发行、汇编、翻译或者网络传播等,侵权必究。
4、如你看到网页展示的文档有www.zixin.com.cn水印,是因预览和防盗链等技术需要对页面进行转换压缩成图而已,我们并不对上传的文档进行任何编辑或修改,文档下载后都不会有水印标识(原文档上传前个别存留的除外),下载后原文更清晰;试题试卷类文档,如果标题没有明确说明有答案则都视为没有答案,请知晓;PPT和DOC文档可被视为“模板”,允许上传人保留章节、目录结构的情况下删减部份的内容;PDF文档不管是原文档转换或图片扫描而得,本站不作要求视为允许,下载前自行私信或留言给上传者【天****】。
5、本文档所展示的图片、画像、字体、音乐的版权可能需版权方额外授权,请谨慎使用;网站提供的党政主题相关内容(国旗、国徽、党徽--等)目的在于配合国家政策宣传,仅限个人学习分享使用,禁止用于任何广告和商用目的。
6、文档遇到问题,请及时私信或留言给本站上传会员【天****】,需本站解决可联系【 微信客服】、【 QQ客服】,若有其他问题请点击或扫码反馈【 服务填表】;文档侵犯商业秘密、侵犯著作权、侵犯人身权等,请点击“【 版权申诉】”(推荐),意见反馈和侵权处理邮箱:1219186828@qq.com;也可以拔打客服电话:4008-655-100;投诉/维权电话:4009-655-100。
关于本文