Matlab中积分函数程序的改进
Ｍａｔｌａｂ中积分函数程序的改进?
丁克会?周建来?倪立学
(淮海工学院?江苏连云港２２２００５)
摘一要:就曲线的长度积分分析了Ｍａｔｌａｂ中三个一维积分函数程序的共性问题?当被积函数可积?函数值在积分区间的端点处无界时?可能会使程序出错?积分失败?分析了程序运行失败的原因?并提出解决方法?提高了Ｍａｔｌａｂ中积分函数的可靠性?增加了用户程序的稳定性?
关键词:Ｍａｔｌａｂ一ｑｕａｄｌ函数一积分一程序
中图分类号:ＴＰ３１１一一一一一文献标识码:Ａ一一一一一文章编号:１００２－６８８６(２０１９)０１－００６３－０３
ＩｍｐｒｏｖｅｍｅｎｔｏｆｉｎｔｅｇｒａｌｆｕｎｃｔｉｏｎｐｒｏｇｒａｍｉｎＭａｔｌａｂ
ＤＩＮＧＫｅｈｕｉ?ＺＨＯＵＪｉａｎｌａｉ?ＮＩｌｉｘｕｅ
Ａｂｓｔｒａｃｔ:Ｔｈｅｃｏｍｍｏｎｐｒｏｂｌｅｍｓｏｆｔｈｒｅｅｏｎｅ－ｄｉｍｅｎｓｉｏｎａｌｉｎｔｅｇｒａｌｆｕｎｃｔｉｏｎｐｒｏｇｒａｍｓｉｎＭａｔｌａｂａｒｅａｎａｌｙｚｅｄｆｏｒｔｈｅｌｅｎｇｔｈｉｎｔｅｇｒａｌｏｆｃｕｒｖｅｓ.Ｗｈｅｎｔｈｅｉｎｔｅｇｒａｎｄｉｓｉｎｔｅｇｒａｂｌｅａｎｄｔｈｅｖａｌｕｅｏｆｔｈｅｆｕｎｃｔｉｏｎｉｓｕｎｂｏｕｎｄｅｄａｔｔｈｅｅｎｄｏｆｔｈｅｉｎｔｅｇｒａｌｉｎ￣ｔｅｒｖａｌ?ｔｈｅｐｒｏｇｒａｍｍａｙｍａｋｅｍｉｓｔａｋｅｓａｎｄｔｈｅｉｎｔｅｇｒａｌｍａｙｆａｉｌ.Ｔｈｅｒｅａｓｏｎｓｆｏｒｔｈｅｆａｉｌｕｒｅｏｆｔｈｅｐｒｏｇｒａｍａｒｅａｎａｌｙｚｅｄ?ａｎｄｔｈｅｓｏｌｕｔｉｏｎｓａｒｅｐｕｔｆｏｒｗａｒｄ.ＩｔｉｍｐｒｏｖｅｓｔｈｅｒｅｌｉａｂｉｌｉｔｙｏｆｉｎｔｅｇｒａｌｆｕｎｃｔｉｏｎｉｎＭａｔｌａｂａｎｄｉｎｃｒｅａｓｅｓｔｈｅｓｔａｂｉｌｉｔｙｏｆｕｓｅｒｐｒｏ￣ｇｒａｍ.
Ｋｅｙｗｏｒｄｓ:Ｍａｔｌａｂ?ｑｕａｄｌｆｕｎｃｔｒｏｎ?ｉｎｔｅｇｒａｌ?ｐｒｏｇｒａｍ
一一时至今日?Ｍａｔｌａｂ在学校已经成为线性代数?自动控制理论?数字信号处理?动态系统仿真等高级课程的基本教学工具?在设计研究单位和工业部门?Ｍａｔｌａｂ被广泛用于科学研究和解决各种具体问题?它除了具备卓越的数值计算能力外?还有符号计算?文字处理?建模仿真和实时控制等功能?
Ｍａｔｌａｂ的主要特点:１)语言简洁紧凑?库函数极其丰富?使用方便灵活?Ｍａｔｌａｂ编程简单?程序书写形式自由?其库函数都由本领域的专家编写?函数的可靠性高?２)源程序的开放性?除内部函数以外?所有Ｍａｔｌａｂ的核心文件和工具箱文件都是可读可改的源文件?用户可通过对源文件的修改以及加入自己的文件构成新的工具箱?本文讨论Ｍａｔｌａｂ中几个积分函数的可靠性问题?
１一Ｍａｔｌａｂ中积分函数简介
Ｍａｔｌａｂ的工具箱中有下列几个数值积分函数:ｑｕａｄ?ｑｕａｄｌ?ｑｕａｄｖ?ｄｂｌｑｕａｄ?ｔｒｉｐｌｅｑｕａｄ?其中ｄｂｌｑｕａｄ和ｔｒｉｐｌｅｑｕａｄ两个函数通过调用ｑｕａｄ?ｑｕａｄｖ?ｑｕａｄｌ三个函数来积分?这三个函数要求被积函数连续?当积分区间端点有界时函数很稳定?不会出现问题?当积分区间端点无界而函数可积时?可能求不到积分而出错?因为?当积分区间端点无界时?被积函数的区间端点的函数值为无穷大?这三个函数用同样的方法:通过将计算端点向区间内偏移一微小的距离?重新计算端点函数值?以避开无穷大值(ｉｎｆ)?本文以Ｍａｔｌａｂ７.１中ｑｕａｄｌ函数为例来讨论这三个函数的这一共性问题?并约定被积函数为无界函数?在给定区间内函数连续可积?经检验?Ｍａｔｌａｂ７.１版后续的版本中积分函数存在同样的问题?
２一积分数学模型
本文以椭圆曲线为例?直角坐标表示的标准椭圆曲线方程如下:
ｘ２
α２＋
ｙ２
β２＝１(１)
３６
其上半部分曲线显式形式为?ｘ?[－α?α]:ｙ＝ｆ(ｘ)＝βαα２－ｘ２(２)
对直角坐标函数曲线ｙ＝ｆ(ｘ)?其在[ａｂ]区间的曲线弧长积分表示为:
Ｌ＝?ｂａ１＋ｙ?ｄｘ(３)
对椭圆上半部分曲线?在积分区间内?(３)式显然收敛?
在运行文献[１]中的程序研究圆弧逼近椭圆曲线时?当逼近到曲线的末端时可能出错?出错的可能性与逼近误差的大小和椭圆曲线的参数有关?报错的原因都是积分函数问题?对这一问题?经过单步深入执行程序?发现问题出现在积分计算末端曲线的长上?在逼近计算到曲线末端时?积分区间很短?ｑｕａｄｌ函数的第７３行试图避开区间的末端点函数值为ｉｎｆ(无穷大)时失败?使得整个逼近过程失败?
３一失败原因分析
３.１一机器数系
电子计算机中数的表示大都采用浮点形式表示?机器数系是一个离散的由有限个有理数所组成的集合?如有一台计算机?ｎ位字长(ｎ＝３)?采用β进制(β＝２)?界码为ｐ?Ｌ?ｐ?Ｕ(硬件常数Ｌ＝－１?Ｕ＝２)?则此计算机表示的数的个数是:
Ｆ＝１＋２(β－１)βｎ－１(Ｕ－Ｌ＋１)＝３３
(４)该机器数系转换为１０进制数为:０??０.２５??０.３１２５??０.３７５??０.４３７５??０.５??０.６２５??０.７５??０.８７５??１??１.２５??１.５??１.７５??２??２.５??３??３.５?可见机器数系是有限不连续?间隔不均匀的?见文献[３]?
３.２一函数ｅｐｓ(ｎ)与ｅｐｓ?ｎ
在Ｍａｔｌａｂ中?ｅｐｓ为浮点相对精度函数?ｅｐｓ为返回从１到下一个最近的较大双精度浮点数的绝对距离?ｅｐｓ＝２(－５２)?ｅｐｓ(ｎ)是返回从绝对值ｎ到下一个最近的较大浮点数的绝对距离?精度同ｎ?ｅｐｓ(ｎ)与ｅｐｓ?ｎ有着不同的意义?后者表示倍数关系?将它们用图表示
?如图１?它们都是非连续的?ｅｐｓ(ｎ)为
一图１一ｅｐｓ(ｎ)与ｅｐｓ?ｎ函数图像阶跃平台函数?
当ｎ＝２ｋ?ｋ取正
整数时?函数值
位于某一平台的
左端点?如ｎ＝
１６或３２时?ｅｐｓ?
ｎ＝ｅｐｓ(ｎ)?当
１６?ｎ<３２时?ｅｐｓ
(ｎ)＝ｅｐｓ(１６)?
所以?此时有ｅｐｓ(ｎ)<ｅｐｓ?ｎ?
３.３一区间端点的取值使端点偏移失败
ｑｕａｄｌ函数的６９和７３行?分别表示对函数区间[ａｂ]的左端点和右端点进行回避无穷大处理?７３行是对被积函数的右端点进行处理?其语句如下:ｙ(１３)＝ｆｅｖａｌ(ｆ?ｂ－ｅｐｓ(ｓｕｐｅｒｉｏｒｆｌｏａｔ(ａ?ｂ))?(ｂ－ａ)?ｖａｒａｒｇｉｎ{:})?
该函数在初始化时?将积分区间分成１２段?共１３个点?ｙ(１３)为区间右端点函数值?当ｆｅｖａｌ(ｆ?ｂ)为ｉｎｆ时?将ｂ用ｂ－ｅｐｓ(ｓｕｐｅｒｉｏｒｆｌｏａｔ(ａ?ｂ))?(ｂ－ａ)替换?其中ｓｕｐｅｒｉｏｒｆｌｏａｔ函数返回单精度或双精度字符串?只要ａ?ｂ有一者为单精度则返回单精度?正常情况下一般都为双精度?而ｅｐｓ( ｄｏｕｂｌｅ )＝ｅｐｓ?所以ｂ－ｅｐｓ(ｓｕｐｅｒｉｏｒｆｌｏａｔ(ａ?ｂ))?(ｂ－ａ)＝ｂ－ｅｐｓ?(ｂ－ａ)?可见函数采用向内偏移ｅｐｓ?(ｂ－ａ)的距离?来试图避开ｉｎｆ?只有当ｂ－ｅｐｓ?(ｂ－ａ)?ｂ时?偏移才成功?才有可能避开ｉｎｆ值?考虑机器数的问题?当ｂ较大?而(ｂ－ａ)较小时?使得ｅｐｓ(ｂ)与ｅｐｓ?(ｂ－ａ)不在一个平台上?ｅｐｓ?(ｂ－ａ)始终在ｅｐｓ(ｂ)平台的左下方?所以ｂ－ｅｐｓ?(ｂ－ａ)＝ｂ?显然?这时不能避开ｉｎｆ值?使得函数运行失败?如２３＝８是函数ｅｐｓ(ｎ)的一个间断点?９和７在间断点两侧有９－ｅｐｓ(７)＝９?增大避开ｉｎｆ的可能性方法可以用ｂ－ｅｐｓ?(ｂ)或ｂ－ｅｐｓ(ｂ)?这两个表达式都有可能偏移成功?而后者较前者偏移的距离一般更小?精度更高?见图１?３.４一端点偏移成功但函数值仍为ｉｎｆ
上述情况并非对所有情况都适用?就是说:即使区间端点偏移成功?但端点函数值计算后仍然为ｉｎｆ?如椭圆的长短半轴半径分别为４０和２０?将椭圆的曲线ｘ２/４０２＋ｙ２/２０２＝１向左平移３８?函数区间[－７８一２]?即ａ＝－７８?ｂ＝２?函数如下:
４６
(ｘ＋３８)２/４０２＋ｙ２/２０２＝１其一阶导数为:
ｙ?＝－(ｘ＋３８)/２/(４０２－(ｘ＋３８)２)１/２
ｙ?函数分母含(４０２－(ｘ＋３８)２)?当ｘ＝ｂ＝２?分母
为０时?其函数值ｙ?为ｉｎｆ?对自变量ｘ＝ｂ－ｅｐｓ?(ｂ－ａ)?当ｂ＝２时有:ｘ＝２－ｅｐｓ?(８０)?ｘ＝２－１.７７６４ｅ－
０１４?可使端点偏移成功?但ｘ值偏移成功后?它与
３８在进行加运算时丢失精度?变为整数?故偏移失败
?
一图２一ｑｕａｄｌ函数程序框图
４一函数区间端点初始函数值的改变对积分结果的影响
一一ｑｕａｄｌ函数初始化计算到的１３点函数值分１１个内点和２个端点?它们的适用范围有所不同?内点的使用仅到判断调整精度结束?见图２第４框?而端点的使用一直持续到积分过程结束?可见?如对初始１３点值进行修改?其影响范围是不同的?
对内点值的修改仅有可能影响到积分精度?而对端点值的修
改肯定会影响到积分结果?两者的影响是两个方面?前者是间接影响?后者是直接影响?若前者有变化?则后者受前者的控制?一般情况下?内点是无需变动的?而端点值的变化是在端点函数值出现ｉｎｆ时才必须改变的?经过实例验证?端点函数值增大时会使得积分值增大?端点函数值减小会使得积分值减小?在积分值是整数时很明显?其积分值的精度是不受影响的?所以调整端点函数值不会改变积分精度?但递归次数会发生变化?数值偏离愈大?函数运行时间会增大?
５一程序改进
考虑以上因素?在函数可积的情况下?当初始化
１３点函数值的端点值为ｉｎｆ?在偏移失败后?可直接指定端点函数值?以增大函数的可靠性?以ｄｕａｑｌ函数为例?在其第７３行下插入下两
行?用于避开右端点是无穷大的情况?ｉｆ一
~ｉｓｆｉｎｉｔｅ(ｙ(１３))一ｙ(１３)＝ｆｅｖａｌ(ｆ?ｂ－ｅｐｓ(ｂ)?ｖａｒａｒｇｉｎ{:})?ｅｎｄ
ｉｆ一~ｉｓｆｉｎｉｔｅ(ｙ(１３))一ｙ(１３)＝１ｅ＋８?ｅｎｄ
在其第６９行下插入下两行?用于避开左端点是无穷大的情况?
ｉｆ一~ｉｓｆｉｎｉｔｅ(ｙ(１))一ｙ(１)＝ｆｅｖａｌ(ｆ?ａ＋ｅｐｓ
(ａ)?ｖａｒａｒｇｉｎ{:})?ｅｎｄ
ｉｆ一
~ｉｓｆｉｎｉｔｅ(ｙ(１))一ｙ(１)＝１ｅ＋８?ｅｎｄ
程序经过改进后不会再出现上述类似的情况?
６一结语
曲线在端点处的切线垂直?导致其求曲线长的被积函数在端点处无界?当曲线长度较短时?在Ｍａｔｌａｂ中用ｑｕａｄｌ等函数积分求曲线的长度时可能失败?本文分析了失败原因?提出了改进积分程序的措施?对可积函数在区间端点无界时都有效?提高了积分函数的可靠性?增加了用户程序的稳定性?
参考文献
[１]一丁克会?席平原?李强化.基于ＭＡＴＬＡＢ单双圆弧混合
逼近曲线的算法及实现[Ｊ].制造技术与机床?２００８(６):１２８－１３１.
[２]一同济大学数学教研室.高等数学[Ｍ].北京:高等教育
出版社?１９９６.
[３]一孙志忠?袁慰平?闻震初.数值分析[Ｍ].南京:东南大
学出版社?２００２.
[４]一李丽?王振领.Ｍａｔｌａｂ工程计算及应用[Ｍ].北京:人民
邮电出版社?２００３.
基金项目:淮海工学院海洋智能装备研究院开放基金:
ＺＢ２０１８００１?
作者简介:丁克会(１９６７－)?男?江苏泰州人?高级实验师?从
事ＣＡＤ/ＣＡＭ的实验教学与研究?
收稿日期:２０１８－０８－２６
５６
几种常见窗函数及其MATLAB程序实现
几种常见窗函数及其MATLAB程序实现 2013-12-16 13:58 2296人阅读评论(0) 收藏举报 分类： Matlab（15） 数字信号处理中通常是取其有限的时间片段进行分析，而不是对无限长的信号进行测量和运算。具体做法是从信号中截取一个时间片段，然后对信号进行傅里叶变换、相关分析等数学处理。信号的截断产生了能量泄漏，而用FFT算法计算频谱又产生了栅栏效应，从原理上讲这两种误差都是不能消除的。在FFT分析中为了减少或消除频谱能量泄漏及栅栏效应，可采用不同的截取函数对信号进行截短，截短函数称为窗函数，简称为窗。 泄漏与窗函数频谱的两侧旁瓣有关，对于窗函数的选用总的原则是，要从保持最大信息和消除旁瓣的综合效果出发来考虑问题，尽可能使窗函数频谱中的主瓣宽度应尽量窄，以获得较陡的过渡带；旁瓣衰减应尽量大，以提高阻带的衰减，但通常都不能同时满足这两个要求。 频谱中的如果两侧瓣的高度趋于零，而使能量相对集中在主瓣，就可以较为接近于真实的频谱。不同的窗函数对信号频谱的影响是不一样的，这主要是因为不同的窗函数，产生泄漏的大小不一样，频率分辨能力也不一样。信号的加窗处理，重要的问题是在于根据信号的性质和研究目的来选用窗函数。图1是几种常用的窗函数的时域和频域波形，其中矩形窗主瓣窄，旁瓣大，频率识别精度最高，幅值识别精度最低，如果仅要求精确读出主瓣频率，而不考虑幅值精度，则可选用矩形窗，例如测量物体的自振频率等；布莱克曼窗主瓣宽，旁瓣小，频率识别精度最低，但幅值识别精度最高；如果分析窄带信号，且有较强的干扰噪声，则应选用旁瓣幅度小的窗函数，如汉宁窗、三角窗等；对于随时间按指数衰减的函数，可采用指数窗来提高信噪比。表1 是几种常用的窗函数的比较。 如果被测信号是随机或者未知的，或者是一般使用者对窗函数不大了解，要求也不是特别高时，可以选择汉宁窗，因为它的泄漏、波动都较小，并且选择性也较高。但在用于校准时选用平顶窗较好，因为它的通带波动非常小，幅度误差也较小。
数值积分的matlab实现
实验10 数值积分 实验目的： 1．了解数值积分的基本原理； 2．熟练掌握数值积分的MATLAB 实现； 3．会用数值积分方法解决一些实际问题。 实验内容： 积分是数学中的一个基本概念，在实际问题中也有很广泛的应用。同微分一样，在《微积分》中，它也是通过极限定义的，由于实际问题中遇到的函数一般都以列表形式给出，所以常常不能用来直接进行积分。此外有些函数虽然有解析式，但其原函数不是初等函数，所以仍然得不到积分的精确值，如不定积分?1 0 d sin x x x 。这时我们一般考虑用数值方法计算其 近似值，称为数值积分。 10.1 数值微分简介 设函数()y f x =在* x 可导，则其导数为 h x f h x f x f h ) ()(lim )(**0* -+='→ （10.1） 如果函数()y f x =以列表形式给出（见表10-1），则其精确值无法求得，但可由下式求得其近似值 h x f h x f x f ) ()()(*** -+≈' （10.2） 表 10-1 一般的，步长h 越小，所得结果越精确。（10.2）式右端项的分子称为函数()y f x =在 *x 的差分，分母称为自变量在*x 的差分，所以右端项又称为差商。数值微分即用差商近似 代替微商。常用的差商公式为： 000()() ()2f x h f x h f x h +--'≈ （10.3） h y y y x f 243)(2 100-+-≈ ' （10.4）
h y y y x f n n n n 234)(12+-≈ '-- （10.5） 其误差均为2 ()O h ，称为统称三点公式。 10.2 数值微分的MATLAB 实现 MATLAB 提供了一个指令求解一阶向前差分，其使用格式为： dx=diff(x) 其中x 是n 维数组，dx 为1n -维数组[]21321,, ,n x x x x x x ---,这样基于两点的数值导 数可通过指令diff(x)/h 实现。对于三点公式，读者可参考例1的M 函数文件diff3.m 。 例1 用三点公式计算()y f x =在=x 1.0,1.2,1.4处的导数值，()f x 的值由下表给 解：建立三点公式的M 函数文件diff3.m 如下： function f=diff3(x,y) n=length(x);h=x(2)-x(1); f(1)=(-3*y(1)+4*y(2)-y(3))/(2*h); for j=2:n-1 f(j)=(y(j+1)-y(j-1))/(2*h); end f(n)=(y(n-2)-4*y(n-1)+3*y(n))/(2*h); 在MATLAB 指令窗中输入指令： x=[1.0,1.1,1.2,1.3,1.4];y=[0.2500,0.2268,0.2066,0.1890,0.1736];diff3(x,y) 运行得各点的导数值为：-0.2470，-0.2170，-0.1890，-0.1650，-0.0014。所以()y f x =在=x 1.0,1.2,1.4处的导数值分别为-0.2470，-0.1890和-0.0014。 对于高阶导数，MATLAB 提供了几个指令借助于样条函数进行求导，详细使用步骤如下： step1：对给定数据点（x,y ），利用指令pp=spline(x,y)，获得三次样条函数数据pp ，供后面ppval 等指令使用。其中，pp 是一个分段多项式所对应的行向量，它包含此多项式的阶数、段数、节点的横坐标值和各段多项式的系数。 step2：对于上面所求的数据向量pp ，利用指令[breaks,coefs,m,n]=unmkpp(pp)进行处理，生成几个有序的分段多项式pp 。 step3：对各个分段多项式pp 的系数，利用函数ppval 生成其相应导数分段多项式的系数，再利用指令mkpp 生成相应的导数分段多项式 step4：将待求点xx 代入此导数多项式，即得样条导数值。 上述过程可建立M 函数文件ppd.m 实现如下： function dy=ppd(pp) [breaks,coefs,m]=unmkpp(pp);
关于MATLAB中分段函数的画法
关于MATLAB中分段函数的画法 最近拿到一题关于MATLAB的分段函数画法的题目,我在网上找了挺久,但没发现很多有用的资料.所以感觉很棘手.但是问题还是要解决,所以我就自己整理了些东西,不怕大家见笑. 我把这些分段函数分为两类: 一.对于y=f(x)这个模型来讲,一类是关于其中一个段是y为常量的一个模型,举例说明. 例 1.y={0,(x<0);1,(x>=0)};在x>-10&x<10区间内的图形 代码如下 : x=-10:0.01:10; y=ones(size(x)); y(x<＝0)=0; plot(x,y); axis([-10 10 -0.5 1.5]); 这样的处理方法就是对于x是变量而Y为常量的而直接定义常数矩阵,再通过判断进行修改,只适合于Y为常量的基础上. ________________________________________________华丽分割线_______________________________________________ 二.第二种是y=f(x),y是关于x的一个变量.需要将x进行赋值的分段函数.这种处理方法比较多. 这里引用一段经典matlab分段画图的例子给大家(代码为蓝色区域): 例 2: x=-3:0.01:3; y1=zeros(size(x)); y2=zeros(size(x)); y3=zeros(size(x)); N=length(x); for k=1:N if x(k)<-1&x(k)>=-3; y1(k)=(-x(k).^2-4*x(k)-3)/2; elseif x(k)>=-1&x(k)<1 ; y2(k)=-x(k).^2+1; else x(k)<=3&x(k)>=1 ; y3(k)=(-x(k).^2+4*x(k)-3)/2; end end y=y1+y2+y3; plot(x,y) 这里运用的是将Y的值设置成三个与x的数量相等的空变量.然后分别依次讲X 的值通过f(x)转换为Y然后画出图形并将三个图形进行组合.
一个简单的Matlab_GUI编程实例
Matlab GUI编程教程（适用于初学者） 1．首先我们新建一个GUI文件：如下图所示； 选择Blank GUI(Default) 2．进入GUI开发环境以后添加两个编辑文本框，6个静态文本框，和一个按钮，布置如下
图所示； 布置好各控件以后，我们就可以来为这些控件编写程序来实现两数相加的功能了。3．我们先为数据1文本框添加代码； 点击上图所示红色方框，选择edit1_Callback，光标便立刻移到下面这段代码的位置。 1. 2. 3.function edit1_Callback(hObject, eventdata, handles) 4.% hObject handle to edit1 (see GCBO) 5.% eventdata reserved - to be defined in a future version of MATLAB
6.% handles structure with handles and user data (see GUIDATA) 7.% Hints: get(hObject,'String') returns contents of edit1 as text 8.% str2double(get(hObject,'String')) returns contents of edit1 as a double 复制代码 然后在上面这段代码的下面插入如下代码： 1. 2.%以字符串的形式来存储数据文本框1的内容. 如果字符串不是数字，则现实空白内容input = str2num(get(hObject,'String')); %检查输入是否为空. 如果为空,则默认显示为0if (isempty(input)) set(hObject,'String','0')endguidata(hObject, handles); 复制代码 这段代码使得输入被严格限制，我们不能试图输入一个非数字。 4．为edit2_Callback添加同样一段代码 5 现在我们为计算按钮添加代码来实现把数据1和数据2相加的目的。 用3中同样的方法在m文件中找到pushbutton1_Callback代码段 如下； 1.function pushbutton1_Callback(hObject, eventdata, handles) 2.% hObject handle to pushbutton1 (see GCBO) 3.% eventdata reserved - to be defined in a future version of MATLAB 4.% handles structure with handles and user data (see GUIDATA) 复制代码
matlab编程实现求解最优解
《现代设计方法》课程 关于黄金分割法和二次插值法的Matlab语言实现在《现代设计方法》的第二章优化设计方法中有关一维搜索的最优化方法的 一节里，我们学习了黄金非分割法和二次插值法。它们都是建立在搜索区间的优先确定基础上实现的。 为了便于方便执行和比较，我将两种方法都写进了一个程序之内，以选择的方式实现执行其中一个。下面以《现代设计方法》课后习题为例。见课本70页，第2—7题。原题如下： 求函数f(x)=3*x^2+6*x+4的最优点，已知单谷区间为[-3,4],一维搜索精度为0.4。 1、先建立函数f(x),f(x)=3*x^2+6*x+4。函数文件保存为：lee.m 源代码为：function y=lee(x) y=3*x^2+6*x+4; 2、程序主代码如下，该函数文件保存为：ll.m clear; a=input('请输入初始点'); b=input('请输入初始步长'); Y1=lee(a);Y2=lee(a+b); if Y1>Y2 %Y1>Y2的情况 k=2; Y3=lee(a+2*b); while Y2>=Y3 %直到满足“大，小，大”为止 k=k+1; Y3=lee(a+k*b); end A=a+b;B=a+k*b; elseif Y1
Matlab考试题库+答案(中北大学)
1、标点符号; _______可以使命令行不显示运算结果， %——用来表示该行为注释行。 2、x为0～4pi，步长为0.1pi的向量，使用命令_______创建。 x=0:0.1*pi:4*pi 3、输入矩阵A=，使用全下标方式用A(2,2) 取出元素“-5”，使用单下标方式用_______取出元素“-5”。 A(5) 4、符号表达式sin(2*a+t)+m中独立的符号变量为_______。 t 5、M脚本文件和M函数文件的主要区别是M脚本文件没有函数定义和M函数文件有函数定义_______。 6. 设x是一维数组，x的倒数第3个元素表示为_______； 设y为二维数组，要删除y的第34行和48列，可使用命令_______; _______； x(_end-2_) y(34,:)=[] y(:,48)=[] 7. 将变量x以Ascii文本格式存储到文件fname.txt，应使用命令_________ _； save _x 8. 在while 表达式, 语句体, End 循环语句中，表达式的值__ __时表示循环条件为真，语句体将被执行，否则跳出该循环语句； 非零 9.要从键盘读入一个字符串并赋值给变量x，且给出提示“Who is she?”，应使用命令_________； x=input(‘Who is she?’,’s’)_ 10．设A=和B=和C=均为m*n矩阵，且存在于WorkSpace中，要产生矩阵D=，可用命令________ _， 计算可用命令________; D=(A-C)/B.^C det(inv(A’*B) 11. 在MATLAB命令窗口中的“>>”标志为MATLAB的_______提示符，“│”标志为_______提示符。 命令行 输入 12.已知A=[1 2 3；4 5 0；7 8 9]；B=[1 0 3；1 5 0；0 1 2]；写出下列各指令运行的结果。 A+B； A.*B； A==B ； ans= [2,2,6;5,10,0;7,9,11] ans= [1,0,9;4,25,0;0,8,18] ans= [1,0,1;0,1,1;0,0,0] 13.已知A是矩阵，求A的对角矩阵函数是_______， 求A的下三角矩阵函数是_______。 diag tril 14.MATLAB的程序文件和Simulink模型文件的扩展名分别是_______、。 .m .mdl 15.MATLAB最基本的绘图函数为_______。 plot() 16. A=[1,2,3;4,5,6]; A(:,[1,3])=[];A=__________________ [2；5] 17. fix(-1.5)=___ ________, round(-1.5)=__ _______________. -1 -2
数字信号处理Matlab实现实例(推荐给学生)
数字信号处理Matlab 实现实例 第１章离散时间信号与系统 例1-1 用MATLAB计算序列{-2 0 1 –1 3}和序列{1 2 0 -1}的离散卷积。 解 MATLAB程序如下： a=[-2 0 1 -1 3]; b=[1 2 0 -1]; c=conv(a,b); M=length(c)-1; n=0:1:M; stem(n,c); xlabel('n'); ylabel('幅度'); 图1.1给出了卷积结果的图形，求得的结果存放在数组c中为：{-2 -4 1 3 1 5 1 -3}。 例1-2 用MATLAB计算差分方程 当输入序列为时的输出结果。 解 MATLAB程序如下： N=41; a=[0.8 -0.44 0.36 0.22]; b=[1 0.7 -0.45 -0.6]; x=[1 zeros(1,N-1)];
k=0:1:N-1; y=filter(a,b,x); stem(k,y) xlabel('n');ylabel('幅度') 图 1.2 给出了该差分方程的前41个样点的输出，即该系统的单位脉冲响应。 例1-3 用MATLAB 计算例1-2差分方程 所对应的系统函数的DTFT 。 解 例1-2差分方程所对应的系统函数为： 123 123 0.80.440.360.02()10.70.450.6z z z H z z z z -------++= +-- 其DTFT 为 23230.80.440.360.02()10.70.450.6j j j j j j j e e e H e e e e ωωωω ωωω--------++= +-- 用MATLAB 计算的程序如下： k=256; num=[0.8 -0.44 0.36 0.02]; den=[1 0.7 -0.45 -0.6]; w=0:pi/k:pi; h=freqz(num,den,w); subplot(2,2,1); plot(w/pi,real(h));grid title('实部') xlabel('\omega/\pi');ylabel('幅度')
CA码生成原理及matlab程序实现
作业：用Matlab写C/A码生成器程序，并画生成码的方波图。 C/A码生成原理 C/A 码是用m 序列优选对组合形成的Gold 码。Gold码是由两个长度相同而互相关极大值为最小的m 序列逐位模2 相加所得到的码序列。它是由两个10 级反馈移位寄存器组合产生的，其产生原理如图1 所示。 图1 C/A码生成原理 发生器的抽头号为3和10，发生器的抽头号为2、3、6、8、9、10；发生器的第10位输出的数字即为码，而码是由的两个抽头的输出结果进行模2相加得到。 卫星的PRN码与延时的量是相关联的,对C/A码来说,每颗卫星都有特别的延时,如第1颗GPS卫星的G2 抽为2、6,第2颗为3、7,第3 颗为4、8,第4 颗为5、9 等，如图2所示。通过G2 相位选择可以产生结构不同的伪随机码,从而可以实现不同卫星之间的码分多址技术与卫星识别。
图2 prn序号与G2抽头、时延对应关系 基于MATLAB的GPS信号实现 编写成“codegen”程序，输入[ca_used]=codegen(svnum)，其中svnum为卫星号，ca_used 为得到的C/A码序列。程序具体实现流程如下： 在程序中定义一个数组，使得卫星号与G2的码片延时一一对应。 gs2=[5;6;7;8;17;18;139;140;141;251;252;254;255;256;257;258;469;470;471;472;473;474;509;512 ;513;514;515;516;859;860;861;862]; 定义两个1×1 023 的数组g1、g2 用来存放生成的Gold 码。定义一个全1 的10 位数组，作为移位寄存器，相当于G1、G2 生成模块的初值均置为全“1”。按原理式
基于matlab程序实现人脸识别
基于m a t l a b程序实现 人脸识别 TYYGROUP system office room 【TYYUA16H-TYY-TYYYUA8Q8-
基于m a t l a b程序实现人脸识别 1.人脸识别流程 基于YCbCr颜色空间的肤色模型进行肤色分割。在YCbCr色彩空间内对肤色进行了建模发现，肤色聚类区域在Cb—Cr子平面上的投影将缩减，与中心区域显着不同。采用这种方法的图像分割已经能够较为精确的将人脸和非人脸分割开来。 人脸识别流程图 2.人脸识别程序 （1）人脸和非人脸区域分割程序 function result = skin(Y,Cb,Cr) %SKIN Summary of this function goes here % Detailed explanation goes here a=; b=; ecx=; ecy=; sita=; cx=; cy=; xishu=[cos(sita) sin(sita);-sin(sita) cos(sita)]; %如果亮度大于230，则将长短轴同时扩大为原来的倍 if(Y>230) a=*a; b=*b; end %根据公式进行计算 Cb=double(Cb); Cr=double(Cr);
t=[(Cb-cx);(Cr-cy)]; temp=xishu*t; value=(temp(1)-ecx)^2/a^2+(temp(2)-ecy)^2/b^2; %大于1则不是肤色，返回0；否则为肤色，返回1 if value>1 result=0; else result=1; end end （2）人脸的确认程序 function eye = findeye(bImage,x,y,w,h) %FINDEYE Summary of this function goes here % Detailed explanation goes here part=zeros(h,w); %二值化 for i=y:(y+h) for j=x:(x+w) if bImage(i,j)==0 part(i-y+1,j-x+1)=255; else part(i-y+1,j-x+1)=0; end end end [L,num]=bwlabel(part,8); %如果区域中有两个以上的矩形则认为有眼睛 if num<2 eye=0;
数值积分用matlab实现
数值积分用m a t l a b实 现
东北大学秦皇岛分校 数值计算课程设计报告 数值积分及Matlab实现 学院数学与统计学院 专业信息与计算科学 学号5133117 姓名楚文玉 指导教师张建波姜玉山 成绩 教师评语： 指导教师签字： 2015年07月14日
1 绪论 在科研计算中，经常会碰到一些很难用公式定理直接求出精确解的积分问题，对于这类问题，我们一般转化为数值积分问题，用计算机来实现求解问题． 1.1 课题的背景 对于定积分()b a f x dx ?在求某函数的定积分时，在一定条件下，虽然有牛顿-莱布里 茨公式()()()b a I f x dx F b F a ==-?可以计算定积分的值，但在很多情况下的原函数() f x 不易求出或非常复杂．被积函数的原函数很难用初等函数表达出来，例如 2 sin (),x x f x e x -= 等；有的函数()f x 的原函数()F x 存在，但其表达式太复杂，计算量太大，有的甚至无法有解析表达式．因此能够借助牛顿-莱布尼兹公式计算定积分的情形是不多的．另外，许多实际问题中的被积函数()f x 往往是列表函数或其他形式的非连续函数，对这类函数的定积分，也不能用不定积分方法求解，只能设法求其近似值．因此，探讨近似计算的数值积分方法是有明显的实际意义的，即有必要研究定积分的数值计算方法，以解决定积分的近似计算．而数值积分就是解决此类问题的一种有效的方法，它的特点是利用被积函数在一些节点上的信息求出定积分的近似值．微积分的发明是人类科学史上一项伟大的成就，在科学技术中，积分是经常遇到的一个重要计算环节数值积分是数学上重要的课题之一，是数值分析中重要的内容之一．随着计算机的出现，近几十年来，对于数值积分问题的研究已经成为一个很活跃的研究领域．现在，数值积分在计算机图形学，积分方程，工程计算，金融数学等应用科学领域都有着相当重要的应用，所以研究数值积分问题有着很重要的意义．国内外众多学者在数值积分应用领域也提出了许多新方法．在很多实际应用中，只能知道积分函数在某些特定点的取值，比如天气测量中的气温、湿度、气压等，医学测量中的血压、浓度等等．通过这个课题的研究，我们将会更好地掌握运用数值积分算法求出特殊积分函数的定积分的一些基本方法、理论基础；并且通过Matlab 软件编程的实现，应用于实际生活中． 1.2 课题的主要内容框架
matlab源代码实例
1.硬币模拟试验 源代码： clear; clc; head_count=0; p1_hist= [0]; p2_hist= [0]; n = 1000; p1 = 0.3; p2=0.03; head = figure(1); rand('seed',sum(100*clock)); fori = 1:n tmp = rand(1); if(tmp<= p1) head_count = head_count + 1; end p1_hist (i) = head_count /i; end figure(head); subplot(2,1,1); plot(p1_hist); grid on; hold on; xlabel('重复试验次数'); ylabel('正面向上的比率'); title('p=0.3试验次数N与正面向上比率的函数图'); head_count=0; fori = 1:n tmp = rand(1); if(tmp<= p2) head_count = head_count + 1; end p2_hist (i) = head_count /i; end figure(head); subplot(2,1,2); plot(p2_hist); grid on; hold on; xlabel('重复试验次数'); ylabel('正面向上的比率'); title('p=0.03试验次数N与正面向上比率的函数图'); 实验结果：
2.不同次数的随机试验均值方差比较 源代码： clear ; clc; close; rand('seed',sum(100*clock)); Titles = ['n=5时' 'n=20时' 'n=25时' 'n=50时' 'n=100时']; Titlestr = cellstr(Titles); X_n_bar=[0]; %the samples of the X_n_bar X_n=[0]; %the samples of X_n N=[5,10,25,50,100]; j=1; num_X_n = 100; num_X_n_bar = 100; h_X_n_bar = figure(1);
遗传算法的原理及MATLAB程序实现
遗传算法的原理及MATLAB程序实现 1 遗传算法的原理 1.1 遗传算法的基本思想 遗传算法(genetic algorithms，GA)是一种基于自然选择和基因遗传学原理，借鉴了生物进化优胜劣汰的自然选择机理和生物界繁衍进化的基因重组、突变的遗传机制的全局自适应概率搜索算法。 遗传算法是从一组随机产生的初始解(种群)开始，这个种群由经过基因编码的一定数量的个体组成，每个个体实际上是染色体带有特征的实体。染色体作为遗传物质的主要载体，其内部表现(即基因型)是某种基因组合，它决定了个体的外部表现。因此，从一开始就需要实现从表现型到基因型的映射，即编码工作。初始种群产生后，按照优胜劣汰的原理，逐代演化产生出越来越好的近似解。在每一代，根据问题域中个体的适应度大小选择个体，并借助于自然遗传学的遗传算子进行组合交叉和变异，产生出代表新的解集的种群。这个过程将导致种群像自然进化一样，后代种群比前代更加适应环境，末代种群中的最优个体经过解码，可以作为问题近似最优解。 计算开始时，将实际问题的变量进行编码形成染色体，随机产生一定数目的个体，即种群，并计算每个个体的适应度值，然后通过终止条件判断该初始解是否是最优解，若是则停止计算输出结果，若不是则通过遗传算子操作产生新的一代种群，回到计算群体中每个个体的适应度值的部分，然后转到终止条件判断。这一过程循环执行，直到满足优化准则，最终产生问题的最优解。图1-1给出了遗传算法的基本过程。 1.2 遗传算法的特点 1.2.1 遗传算法的优点
遗传算法具有十分强的鲁棒性，比起传统优化方法，遗传算法有如下优点: 1. 遗传算法以控制变量的编码作为运算对象。传统的优化算法往往直接利用控制变量的实际值的本身来进行优化运算，但遗传算法不是直接以控制变量的值，而是以控制变量的特定形式的编码为运算对象。这种对控制变量的编码处理方式，可以模仿自然界中生物的遗传和进化等机理，也使得我们可以方便地处理各种变量和应用遗传操作算子。 2. 遗传算法具有内在的本质并行性。它的并行性表现在两个方面，一是遗传 开始 初始化，输入原始参 数及给定参数，gen=1 染色体编码，产生初始群体 计算种群中每个个体的适应值 终止条件的判断, N gen=gen+1 选择 交叉 Y 变异 新种群 输出结果 结束 图1-1 简单遗传算法的基本过程
数值积分算法与MATLAB实现
数值积分算法与MATLAB实现 本文从网络收集而来，上传到平台为了帮到更多的人，如果您需要使用本文档，请点击下载按钮下载本文档（有偿下载），另外祝您生活愉快，工作顺利，万事如意！ 摘要：在求一些函数的定积分时，由于原函数十分复杂难以求出或用初等函数表达，导致积分很难精确求出，只能设法求其近似值，因此能够直接借助牛顿-莱布尼兹公式计算定积分的情形是不多的。数值积分就是解决此类问题的一种行之有效的方法。积分的数值计算是数值分析的一个重要分支；因此，探讨近似计算的数值积分方法是有着明显的实际意义的。本文从数值积分问题的产生出发，详细介绍了一些数值积分的重要方法。 本文较详细地介绍了牛顿-科特斯求积公式，以及为了提高积分计算精度的高精度数值积分公式，即龙贝格求积公式和高斯-勒让德求积公式。除了研究这些数值积分算法的理论外，本文还将这些数值积分算法在计算机上通过MATLAB软件编程实现，并通过实例用各种求积公式进行运算，分析比较了各种求积公式的计算误差。 【关键词】数值积分牛顿-科特斯求积公式高精度求积公式MATLAB软件
前言 对于定积分，在求某函数的定积分时，在一定条件下，虽然有牛顿-莱布里茨公式可以计算定积分的值，但在很多情况下的原函数不易求出或非常复杂。被积函数的原函数很难用初等函数表达出来，例如等；有的函数的原函数存在，但其表达式太复杂，计算量太大，有的甚至无法有解析表达式。因此能够借助牛顿-莱布尼兹公式计算定积分的情形是不多的。另外，许多实际问题中的被积函数往往是列表函数或其他形式的非连续函数，对这类函数的定积分，也不能用不定积分方法求解，只能设法求其近似值。因此，探讨近似计算的数值积分方法是有明显的实际意义的，即有必要研究定积分的数值计算方法，以解决定积分的近似计算。而数值积分就是解决此类问题的一种有效的方法，它的特点是利用被积函数在一些节点上的信息求出定积分的近似值。 微积分的发明是人类科学史上一项伟大的成就，在科学技术中，积分是经常遇到的一个重要计算环节。数值积分是数学上重要的课题之一，是数值分析中重要的内容之一，也是应用数学研究的重点。随着计算机的出现，近几十年来，对于数值积分问题的研究已经成为一个很活跃的研究领域。现在，数值积分在计算
matlab使用指导(2012)
一、基础知识 1.1 常见数学函数 如：输入x=[-4.85 -2.3 -0.2 1.3 4.56 6.75]，则： ceil(x)= -4 -2 0 2 5 7 fix(x) = -4 -2 0 1 4 6 floor(x) = -5 -3 -1 1 4 6 round(x) = -5 -2 0 1 5 7 1.2 系统的在线帮助 1 help 命令： 1.当不知系统有何帮助内容时，可直接输入help以寻求帮助: >>help（回车） 2.当想了解某一主题的内容时，如输入： >> help syntax（了解Matlab的语法规定） 3.当想了解某一具体的函数或命令的帮助信息时，如输入： >> help sqrt （了解函数sqrt的相关信息）
2 lookfor命令 现需要完成某一具体操作，不知有何命令或函数可以完成，如输入： >> lookfor line （查找与直线、线性问题有关的函数） 1.3 常量与变量 系统的变量命名规则：变量名区分字母大小写；变量名必须以字母打头，其后可以是任意字母，数字，或下划线的组合。此外，系统内部预先定义了几个有特殊意 1 数值型向量（矩阵）的输入 1．任何矩阵（向量），可以直接按行方式 ．．．输入每个元素：同一行中的元素用逗号（，）或者用空格符来分隔；行与行之间用分号（；）分隔。所有元素处于一方括号（[ ]）内； 例1： >> Time = [11 12 1 2 3 4 5 6 7 8 9 10] >> X_Data = [2.32 3.43；4.37 5.98] 2
上面函数的具体用法，可以用帮助命令help得到。如help zeros,可查到zeros的具体用法。 例：meshgrid(x,y) 输入x=[1 2 3 4]; y=[1 0 5]; [X,Y]=meshgrid(x, y)，则 X = Y = 1 2 3 4 1 1 1 1 1 2 3 4 0 0 0 0 1 2 3 4 5 5 5 5 目的是将原始数据x，y转化为矩阵数据X，Y。 2 符号向量（矩阵）的输入 1．用函数sym定义符号矩阵： 函数sym实际是在定义一个符号表达式，这时的符号矩阵中的元素可以是任何的符号或者是表达式，而且长度没有限制。只需将方括号置于单引号中。 例2： >> sym_matrix = sym（'[a b c；Jack Help_Me NO_WAY]'） sym_matrix = [ a, b, c] [Jack, Help_Me, NO_WAY] 2．用函数syms定义符号矩阵 先定义矩阵中的每一个元素为一个符号变量，而后像普通矩阵一样输入符号矩阵。 例3： >> syms a b c ； >> M1 = sym（'Classical'）； >> M2 = sym（' Jazz'）； >> M3 = sym（'Blues'）； >> A = [a b c；M1，M2，M3；sym（[2 3 5]）] A = [ a, b, c] [Classical, Jazz, Blues] [ 2, 3, 5]
基于Matlab的动态规划程序实现
动态规划方法的Matlab 实现与应用 动态规划(Dynamic Programming)是求解决策过程最优化的有效数学方法，它是根据“最优决策的任何截断仍是最优的”这最优性原理，通过将多阶段决策过程转化为一系列单段决策问题，然后从最后一段状态开始逆向递推到初始状态为止的一套最优化求解方法。 1．动态规划基本组成 (1) 阶段 整个问题的解决可分为若干个阶段依次进行，描述阶段的变量称为阶段变量，记为k (2) 状态 状态表示每个阶段开始所处的自然状况或客观条件，它描述了研究问题过程的状况。各阶段状态通常用状态变量描述，用k x 表示第k 阶段状态变量，n 个阶段决策过程有n+ 1个状态。 (3) 决策 从一确定的状态作出各种选择从而演变到下一阶段某一状态，这种选择手段称为决策。描述决策的变量称为决策变量，决策变量限制的取值范围称为允许决策集合。用()k k u x 表示第k 阶段处于状态k x 时的决策变量，它是k x 的函数。用()k k D x Dk(xk)表示k x 的允许决策的集合。 (4) 策略 每个阶段的决策按顺序组成的集合称为策略。由第k 阶段的状态k x 开始到终止状态的后部子过程的策略记为{}11(),(),,()k k k k n n u x u x u x ++ 。可供选择的策略的范围称为允许策略集合，允许策略集合中达到最优效果的策略称为最优策略。从初始状态* 11()x x =出发，过程按照最优策略和状态转移方程演变所经历的状态序列{ } **** 121,,,,n n x x x x + 称为最优轨线。 (5) 状态转移方程 如果第k 个阶段状态变量为k x ，作出的决策为k u ，那么第k+ 1阶段的状态变量1k x +也被完全确定。用状态转移方程表示这种演变规律，记为1(,)k k k x T x u +=。 (6) 指标函数 指标函数是系统执行某一策略所产生结果的数量表示，是衡量策略优劣的数量指标，它定义在全过程和所有后部子过程上，用()k k f x 表示。过程在某阶段j 的阶段指标函数是衡量该阶段决策优劣数量指标，取决于状态j x 和决策j u ，用(,)j j j v x u 表示。 2．动态规划基本方程 (){} 11()min ,,(),()k k k k k k k k k k f x g v x u f x u D x ++=∈???? Matlab 实现 (dynprog.m 文件) function [p_opt,fval]=dynprog (x,DecisFun,SubObjFun,TransFun,ObjFun) % x 是状态变量，一列代表一个阶段的所有状态； % M-函数DecisFun(k,x) 由阶段k 的状态变量x 求出相应的允许决策变量; % M-函数SubObjFun(k,x,u) 是阶段指标函数, % M-函数ObjFun(v,f) 是第k 阶段至最后阶段的总指标函数 % M-函数TransFun(k,x,u) 是状态转移函数, 其中x 是阶段k 的某状态变量, u 是相应的决策变量; %输出 p_opt 由4列构成，p_opt=[序号组;最优策略组;最优轨线组;指标函数值组]； %输出 fval 是一个列向量，各元素分别表示p_opt 各最优策略组对应始端状态x 的最优函数值。
利用Matlab实现Romberg数值积分算法----系统建模与仿真结课作业
利用Matlab 实现Romberg 数值积分算法 一、内容摘要 针对于某些多项式积分，利用Newton —Leibniz 积分公式求解时有困难，可以采用数值积分的方法，求解指定精度的近似解，本文利用Matlab 中的.m 文件编写了复化梯形公式与Romberg 的数值积分算法的程序，求解多项式的数值积分，比较两者的收敛速度。 二、数值积分公式 1.复化梯形公式求解数值积分的基础是将区间一等分时的Newton —Cotes 求积公式: I =(x)[f(a)f(b)]2 b a b a f dx -≈ +? 其几何意义是，利用区间端点的函数值、与端点构成的梯形面积来近似(x)f 在区间[a,b]上的积分值，截断误差为： 3" (b a)()12 f η-- (a,b)η∈ 具有一次的代数精度，很明显，这样的近似求解精度很难满足计算的要求，因而，可以采用将积分区间不停地对分，当区间足够小的时候，利用梯形公式求解每一个小区间的积分近似值，然后将所有的区间加起来，作为被求函数的积分，可以根据计算精度的要求，划分对分的区间个数，得到复化梯形公式： I =1 1 (b a)(b a) (x)dx [f(a)f(b)2(a )]2n b a k k f f n n -=--≈+++∑? 其截断误差为：
2" (b a)h ()12 R f η--= (a,b)η∈ 2.Romberg 数值积分算法 使用复化的梯形公式计算的数值积分，其收敛速度比减慢，为此，采用Romberg 数值积分。其思想主要是，根据I 的近似值2n T 加上I 与2n T 的近似误差，作为新的I 的近视，反复迭代，求出满足计算精度的近似解。 用2n T 近似I 所产生的误差可用下式进行估算： 12221 ()3 n n n I T T T -?=-=- 新的I 的近似值： 122 n n j T T -=?+ j =（0 1 2 ….） Romberg 数值积分算法计算顺序 i=0 (1) 002T i=1 (2) 102T (3) 012T i=2 (4) 202T (5) 112T (6) 022T i=3 (7) 302T (8) 212T (9) 122T (10) 032T i=4 (11) 402T (12) 312T (13) 222T (14) 132T … … … … 其中，第一列是二阶收敛的，第二列是四阶收敛的，第三列是六阶收敛的，第四列是八阶收敛的，即Romberg 序列。
用matlab绘制各种数字信号中的函数还有分段函数及翻褶平移
《数字信号处理》 （一）实验目的 使用stem绘图函数分别画出离散时间信号在指定范围内的图形。画图时使用xlabel,ylabel,title,legend等函数进行注释。复习MATLAB的基本应用，如：函数的定义、画图……并巩固理论知识中的多种离散时间信号及其图形，以及延迟与翻褶的函数变换等。 （二）程序的运行与截图 1)用stem绘制单位阶跃序列u（n） clear all;close all;clc;%清除所有变量 n=0:50;%取值范围 y=(n>=0);%n>=0,y=1;n<0,y=0 stem(n,y);%显示出当0<=n<=50 时，函数u(n)的取值范围 xlabel('n');%对横轴进行注释 ylabel('y=u(n)');%对纵轴进行注释 title('y=u(n)的图形');%对图像的标题进行注释 legend('y=u(n)',2);%对图中曲线进行注释，标注在第二象限 2)用stem绘制单位抽样（冲激）序列δ（n） clear all;close all;clc; %清除所有变量
n=0:50; %取值范围 y=(n==0);%n=0,y=1;n!=0,y=1 stem(n,y);%显示出当0<=n<=50 时，函数δ(n)的取值范围xlabel('n');%对横轴进行注释 ylabel('y=δ(n)');%对纵轴进行注释 title('y=δ(n)的图形');%对图像的标题进行注释 legend('y=δ(n)',2);%对图中曲线进行注释，标注在第二象限
3）用stem绘制矩形序列Rn（n）clear all;close all;clc; %清除所有变量 n=0:50; %取值范围 R10=((n>=0)&(n-9)<=0);%0<=n<=10,y=1;n>10,y=0 stem(n,R10);%显示出当0<=n<=50 时，函数Rn（n）的取值范围xlabel('n');%对横轴进行注释 ylabel(' y=R10（n）');%对纵轴进行注释 title('y=R10（n）的图形');%对图像的标题进行注释 legend('y=R10（n）',2);%对图中曲线进行注释，标注在第二象限