当前位置:文档之家› 面向方程的数值积分法仿真

面向方程的数值积分法仿真

面向方程的数值积分法仿真
面向方程的数值积分法仿真

实验一面向方程的数值积分法仿真

1、实验目的

学习用Matlab编写数值积分法仿真程序,深入理解一阶常微分方程组数值积分法的原理和程序结构,加深理解4阶龙格—库塔法的原理以及步长与算法稳定性的关系。

2、实验内容

线性定常系统仿真,二阶系统传递函数为

假设状态初值为零,系统输入为阶跃输入,

=1。

(1)根据所学的原理,将传递函数转换为线性状态空间模型(能控标准型),写出状态空间模型表达式。

(2)对状态空间模型采用定步长 RK4 法进行仿真,取步长 h=0.1,仿真时间长度为20秒,做出阻尼比分别为ξ=0,0.1,0.5,1.0,1.5 时的响应曲线,并叠加显示。

答:由于本实验没有控制作用,所以只需要编写三个函数,即“.m”文件。他们分别是

线性系统数值积分法仿真框架函数

function [t,y]=w_LinearSimu(tstart,tstop,h,x0,u0,A,B,C,D)

t=[tstart:h:tstop];%t数一个行序列

cntt=size(t,2);%返回列数

cnty=size(C,1);%返回y的维数

y=zeros(cnty,cntt);%构造一个空矩阵,用来存储结果

y0=C*x0+D*u0;%

%eval([OutputFile,'(tstart,x0,u0)']);%计算初始输出

y(:,1)=y0 ;%将y0作为输出的第1列

curx=x0; %当前一步的x

curu=u0; %当前一步的u

cury=y0; %当前一步的y

for i=1:1:cntt-1

curx=w_LinearStepIntegral(h,curx,curu,'RK4',A,B);%单步积分运算,针对线性模型

cury=C*curx+D*curu;%计算输出

y(:,i+1)=cury;%将输出加入到输出序列里

end

线性系统单步积分函数

function NewX=w_LinearStepIntegral(h,curx,curu,InteMethod,A,B); % 函数功能:单步积分运算,与模型方程有关

% 输入参数:h是数值积分步长,curx,curu分别是当前的状态和控制向量% InteMethod是积分算法,字符串类型,

% 可以取'EUL','RK2','RK4',由于要用于判断,字符串必须等长

% A,B是线性状态模型的微分方程的系数矩阵

% 输出参数:NewX是这一步计算的新的状态向量

Bu=B*curu;

if InteMethod == 'RK4'

k1=A*curx+Bu; k2=A*(curx+0.5*h*k1)+Bu;

k3=A*(curx+0.5*h*k2)+Bu; k4=A*(curx+h*k3)+Bu;

NewX=curx+h*(k1+2*k2+2*k3+k4)/6;

else if InteMethod == 'RK2'

k1=A*curx+Bu;

k2=A*(curx+h*k1)+Bu;

NewX=curx+0.5*h*(k1+k2);

else %欧拉法EUL

NewX=curx+h*(A*curx+Bu);

end

仿真组织函数

function simu_StateSpace% 函数功能:对一个状态空间模型进行仿真

tstart=0;

tstop=20; h=0.1; B=[0;1]; C=[1,0]; D=[0]; x0=[0;0];

u0=[1];%单位阶跃输入

dampratio=0; A=[0,1;-1,-2*dampratio];

[t,y0]=w_LinearSimu(tstart,tstop,h,x0,u0,A,B,C,D,'RK4');

dampratio=0.1; A(2,2)=-2*dampratio;

[t,y1]=w_LinearSimu(tstart,tstop,h,x0,u0,A,B,C,D,'RK4');

dampratio=0.5; A(2,2)=-2*dampratio;

[t,y5]=w_LinearSimu(tstart,tstop,h,x0,u0,A,B,C,D,'RK4');

dampratio=1; A(2,2)=-2*dampratio;

[t,y10]=w_LinearSimu(tstart,tstop,h,x0,u0,A,B,C,D,'RK4');

dampratio=1.5; A(2,2)=-2*dampratio;

[t,y15]=w_LinearSimu(tstart,tstop,h,x0,u0,A,B,C,D,'RK4');

plot(t,y0,'+k',t,y1,'-r',t,y5,':g',t,y10,'*y',t,y15,'.b');%绘制叠加曲线

xlabel('x'); ylabel('y'); title('二阶系统阶跃输入响应');

在编写完这三个函数后,调用“function simu_StateSpace”,得出仿真图形如下图

(3)取

=1,ξ=0.5,仿真时间长度 20 秒。步长从 h=0.1不断增大,做出不同步长的仿真曲线。验证步长与算法稳定性的关系,找出大概的最大允许步长。

答:本题只需对上一步的“function simu_StateSpace”函数修改为

function simu_StateSpace2% 函数功能:对一个状态空间模型进行仿真

tstart=0;

tstop=20; h=2.25; B=[0;1]; C=[1,0]; D=[0]; x0=[0;0];

u0=[1];%单位阶跃输入

dampratio=0.5; A=[0,1;-1,-2*dampratio];

[t,y5]=w_LinearSimu(tstart,tstop,h,x0,u0,A,B,C,D,'RK4');

plot(t,y5);

然后修改相应的h值,进行仿真即可,得出下列仿真波形

h=0.1 h=0.5

h=1.0 h=1.5

h=2.0 h=2.5

h=3.0

观察这些波形可以看出,随着步长的增大,波形的线性度下降。

当h>2.5时,波形会发散,系统趋于不稳定,所以大概的最大允许步长为2.5。

3、实验总结

本实验主要是针对一个传递函数描述的系统,将其转换为线性状态空间模型,再用matlab根据算法理论进行编程。

刚开始不清楚函数的架构,感觉无从下手,通过看实验讲义,大致弄清了函数的四个架构,即数值积分法仿真函数、单步数值积分法函数、系统控制函数、仿真组织函数,但还是不太确定对不对,最后问了下老师,得到了进一步的理解,由于本次试验没有控制,所以不需要系统控制函数。所以,非常感谢老师的辅

导。在编程的时候大部分借鉴了老师的框架,通过自己的编程也进一步理解了整个函数的组织和函数的具体实现。试验完成后,仔细想想试验并不是很难,关键是对算法的理解还是不透彻,所以,对于仿真这样实践性很强的知识,还是更应该多加练习,勤实践,通过实践更好的理解算法。

matlab数值微积分与方程数值求解

电子一班王申江 实验九数值微积分与方程数值求解 一、实验目的 1、掌握求数值导数和数值积分的方法 2、掌握代数方程数值求解的方法 3、掌握常微分方程数值求解的方法 二、实验内容 1、求函数在指定点的数值导数。 () 23 2 123,1,2,3 026 x x x f x x x x x == >>syms x >>f=[x x^2 x^3;1 2*x 3*x^2;0 2 6*x]; >>F=det(f) F=2*x^3 >>h=0.1 >>x=[0:h:4]; >>f=2*x^3; >>[dy,dx]=diff_ctr(f,h,1); >>y1=dy(dx==1) y1=6.0000 >>y2=dy(dx==2)

y2=24.0000 >>y3=dy(dx==3) y3=54.0000 2、用数值方法求定积分。 (1) 210I π =?的近似值 a=inline('sqrt(cos(t.^2)+4*sin((2*t).^2)+1)'); I=quadl(a,0,2*pi) I = 6.7992 + 3.1526i (2)()1 202ln 11x I dx x +=+? b=inline('log(1+x)./(1+x.^2)'); I=quadl(b,0,1) I = 0.2722 3、分别用3种不同的数值方法解线性方程组。 6525494133422139211 x y z u x y z u x y z u x y u +-+=-??-+-=??++-=??-+=? A=[6,5,-2,5;9,-1,4,-1;3,4,2,-2;3,-9,0,2]; b=[-4,13,1,11]'; x=A\b

数值积分算法与MATLAB实现陈悦5133201讲解

东北大学秦皇岛分校 数值计算课程设计报告 数值积分算法及MATLAB实现 学院数学与统计学院 专业信息与计算科学 学号5133201 姓名陈悦 指导教师姜玉山张建波 成绩 教师评语: 指导教师签字: 2015年07月14日

1 绪论 数值分析是计算数学的一个主要部分,计算数学是数学科学的一个分支,它研究用计算机求解各种数学问题的数值检索方其理论与软件的实现.而数值分析主要研究数值计算. 现科学技术的发展与进步提出了越来越多的复杂的数值计算问题,这些问题的圆满解决已远人工手算所能胜任,必须依靠电子计算机快速准确的数据处理能力.这种用计算机处理数值问题的方法,成为科学计算.今天,科学计算的应用范围非常广泛,天气预报、工程设计、流体计算、经济规划和预测以及国防尖端的一些科研项目,如核武器的研制、导弹和火箭的发射等,始终是科学计算最为活跃的领域. 1.1 数值积分介绍 数值积分是数值分析的重要环节,实际问题当中常常需要计算积分,有些数值方法,如微分方程和积分方程的求解,也都和积分计算相联系. 求某函数的定积分时,在多数情况下,被积函数的原函数很难用初等函数表达出来,因此能够借助微积分学的牛顿-莱布尼兹公式计算定积分的机会是不多的.另外,许多实际问题中的被积函数往往是列表函数或其他形式的非连续函数,对这类函数的定积分,也不能用不定积分方法求解.由于以上原因,数值积分的理论与方法一直是计算数学研究的基本课题.对微积分学做出杰出贡献的数学大师,如I.牛顿、L.欧拉、C.F.高斯、拉格朗日等人都在数值积分这个领域作出了各自的贡献,并奠定了这个分支的理论基础. 构造数值积分公式最通常的方法是用积分区间上的n 次插值多项式代替被积函数,由此导出的求积公式称为插值型求积公式.特别在节点分布等距的情形称为牛顿-科特斯公式,例如梯形公式(Trapezoidal Approximations)与抛物线公式(Approximations Using Parabolas)就是最基本的近似公式.但它们的精度较差.龙贝格算法是在区间逐次分半过程中,对梯形公式的近似值进行加权平均获得准确程度较高的积分近似值的一种方法,它具有公式简练、计算结果准确、使用方便、稳定性好等优点,因此在等距情形宜采用龙贝格求积公式(Rhomberg Integration).当用不等距节点进行计算时,常用高斯型求积公式计算,它在节点数目相同情况下,准确程度较高,稳定性好,而且还可以计算无穷积分.数值积分还是微分方程数值解法的重要依据.许多重要公式都可以用数值积分方程导出.现探讨数值积分算法以及运用MATLAB软件的具体实现

微分方程数值解法

《微分方程数值解法》 【摘要】自然界与工程技术中的很多现象,可以归结为微分方程定解问题。其中,常微分方程求解是微分方程的重要基础内容。但是,对于许多的微分方程,往往很难得到甚至不存在精确的解析表达式,这时候,数值解提供了一个很好的解决思路。,针对于此,本文对常微分方程数值解法进行了简单研究,主要讨论了一些常用的数值解法,如欧拉法、改进的欧拉法、Runge —Kutta 方法、Adams 预估校正法以及勒让德谱方法等,通过具体的算例,结合MA TLAB 求解画图,初步给出了一般常微分方程数值解法的求解过程。同时,通过对各种方法的误差分析,让大家对各种方法的特点和适用范围有一个直观的感受。 【关键词】 常微分方程 数值解法 MA TLAB 误差分析 引言 在我国高校,《微分方程数值解法》作为对数学基础知识要求较高且应用非常广泛的一门课程,不仅 在数学专业,其他的理工科专业的本科及研究生教育中开设这门课程.近四十年来,《微分方程数值解法》不论在理论上还是在方法上都获得了很大的发展.同时,由于微分方程是描述物理、化学和生物现象的数学模型基础,且它的一些最新应用已经扩展到经济、金融预测、图像处理及其他领域 在实际应用中,通过相应的微分方程模型解决具体问题,采用数值方法求得方程的近似解,使具体问题迎刃而解。 2 欧拉法和改进的欧拉法 2.1 欧拉法 2.1.1 欧拉法介绍 首先,我们考虑如下的一阶常微分方程初值问题 ???==0 0)() ,('y x y y x f y (2--1) 事实上,对于更复杂的常微分方程组或者高阶常微分方程,只需要将x 看做向量,(2--1)就成了一个一阶常微分方程组,而高阶常微分方程也可以通过降阶化成一个一阶常微分方程组。 欧拉方法是解常微分方程初值问题最简单最古老的一种数值方法,其基本思路就是把(2--1)中的导数项'y 用差商逼近,从而将一个微分方程转化为一个代数方程,以便求解。 设在[]b a ,中取等距节点h ,因为在节点n x 点上,由(2--1)可得:

偏微分方程数值解法答案

1. 课本2p 有证明 2. 课本812,p p 有说明 3. 课本1520,p p 有说明 4. Rit2法,设n u 是u 的n 维子空间,12,...n ???是n u 的一组基底,n u 中的任一元素n u 可 表为1n n i i i u c ?==∑ ,则,11 11()(,)(,)(,)(,)22j n n n n n n i j i j j i j j J u a u u f u a c c c f ???=== -=-∑∑是12,...n c c c 的二次函数,(,)(,)i j j i a a ????=,令 () 0n j J u c ?=?,从而得到12,...n c c c 满足1 (,)(,),1,2...n i j i j i a c f j n ???===∑,通过解线性方程组,求的i c ,代入1 n n i i i u c ?==∑, 从而得到近似解n u 的过程称为Rit2法 简而言之,Rit2法:为得到偏微分方程的有穷维解,构造了一个近似解,1 n n i i i u c ?== ∑, 利用,11 11()(,)(,)(,)(,)22j n n n n n n i j i j j i j j J u a u u f u a c c c f ???===-=-∑∑确定i c ,求得近似解n u 的过程 Galerkin 法:为求得1 n n i i i u c ? == ∑形式的近似解,在系数i c 使n u 关于n V u ∈,满足(,)(,) n a u V f V =,对任 意 n V u ∈或(取 ,1j V j n ?=≤≤) 1 (,)(,),1,2...n i j i j i a c f j n ???===∑的情况下确定i c ,从而得到近似解1 n n i i i u c ?==∑的过程称 Galerkin 法为 Rit2-Galerkin 法方程: 1 (,)(,)n i j i j i a c f ???==∑ 5. 有限元法:将偏微分方程转化为变分形式,选定单元的形状,对求解域作剖分,进而构 造基函数或单元形状函数,形成有限元空间,将偏微分方程转化成了有限元方程,利用 有效的有限元方程的解法,给出偏微分方程近似解的过程称为有限元法。 6. 解:对求解区间进行网格剖分,节点01......i n a x x x x b =<<<<=得到相邻节点1,i i x x -

几种常用数值积分方法的比较汇总

学科分类号110.3420 州 GUIZHOU NORMAL COLLEGE 本科毕业论文 题目—几种常用数值积分方法的比较_____________ 姓名潘晓祥学号1006020540200 院(系)数学与计算机科学学院 __________________ 专业数学与应用数学年级_____________2010级 指导教师雍进军职称______________________讲师 二O—四年五月

贵州师范学院本科毕业论文(设计)诚信声明本人郑重声明:所呈交的本科毕业论文(设计),是本人在指导老师的指导下,独立进行研究工作所取得的成果,成果不存在知识产权争议,除文中已经注明引用的内容外,本论文不含任何其他个人或集体已经发表或撰写过的作品成果。对本文的研究做出重要贡献的个人和集体均已在文中以明确方式标明。本人完全意识到本声明的法律结果由本人承担。 本科毕业论文作者签名: 年月曰

贵州师范学院本科毕业论文(设计)任务书

研究方法: 本论文主要通过对相关文献和书籍的参考,合自己的见解,复化求积公式,Newton —Cotes求积公式,Romberg求积公式,高斯型求积公式进行讨论并进行上机实验,从代数精度,求积公式误差等角度对这些方法进行分析比较完成期限和采取的主要措施: 本论文计划用6个月的时间完成,阶段的任务如下: (1) 7月份查阅相关书籍和文献; (2) 8月份完成开题报告并交老师批阅; (3) 9月份完成论文初稿并交老师批阅; (4) 10月份完成论文二搞并交老师批阅; (5) 11月份完成论文三搞; (6) 12月份定稿. 主要措施:考相关书籍和文献,合自己的见解,老师的指导下和同学的帮助下完成 主要参考文献及资料名称: [1] 关治?陆金甫?数学分析基础(第二版) [M].北京:等教育出版社.2010.7 [2] 胡祖炽.林源渠.数值分析[M]北京:等教育出版社.1986.3 [3] 薛毅.数学分析与实验[M] 北京:业大学出版社2005.3 [4] 徐士良.数值分析与算法[M].北京:械工业出版社2007.1 [5] 王开荣.杨大地.应用数值分析[M]北京:等教育出版社2010.7 [6] 杨一都.数值计算方法[M].北京:等教育出版社.2008.4 [7] 韩明.王家宝.李林.数学实验(MATLAB版[M].上海:济大学出版社2012.1 [8] 圣宝建.关于数值积分若干问题的研究[J].南京信息工程大学.2009.05.01. : 42 [9] 刘绪军.几种求积公式计算精确度的比较[J].南京职业技术学院.2009. [10] 史万明.吴裕树.孙新.数值分析[M].北京理工大学出版社.2010.4. 指导教师意见: 签名: 年月日

数值积分法

第7章 数值积分法 7.1 实验目的 了解求积公式及代数精度概念,理解并掌握求定积分的求积公式的算法构造和计算,学习用计算机求定积分的一些科学计算方法和简单的编程技术和能用程序实现这些算法。 7.2 概念与结论 1. 求积公式 计算定积分的如下形式的近似公式: 称为求积公式。 2.代数精度 若求积公式 对一切不高于m 次的 多项都准确成立,而对于m+1次多项式等号不成立,则称此求积公式 的代数精度为m 。 代数精度越高,求积公式越好。 3.求积余项 4.Newton-Cotes 求积公式的代数精度 n 点Newton-Cotes 求积公式的代数精度至少可以达到n-1,且当n 为奇数时,可以达到n 。 ?∑=≈b a n k k k x f A dx x f 1 ) ()(?∑=≈b a n k k k x f A dx x f 1 ) ()(?∑=-=b a n k k k x f A dx x f f R 1 ) ()()(

? =b a dx x f I )(5.Richardson 外推定理 设函数F 1(h)逼近量F*的余项为: F*-F 1(h)=a 1h p1 +a 2h p2 +····+a k p k +··· 式中p k >p k-1>···>p 2>p 1>0, F*和a i (i=1,2, ···)都是与h 无关的常数,且k ≥1时,a k ≠0,则由: 定义的函数F 2(h)也逼近F*,且有 F*-F 2(h)= b 2h p2 +····+b k p k +··· 6. 关于复合梯形公式的展开定理 设f(x)在[a,b]区间上无穷次可微,则有如下展开式: T(h)=I+a 1h 2 +a 2h 4 +a 3h 6 +…+a m h 2m +… 式中T(h)是函数f(x)在[a,b]区间上的复化梯形值Tn, 7.3 程序中Mathematica 语句解释 1. 随机函数 Random[] 随机给出闭区间[0,1]内的一个实数 Random[Real, xmax] 随机给出闭区间[0,xmax]内的一个实数 Random[Real, {xmin, xmax}] 随机给出闭区间[xmin,xmax]内的一个实数 Random[Integer] 随机给出整数0或1 Random[Integer, {xmin, xmax}] 随机给出xmin 到xmax 之间的一个整数 Random[Complex] 随机给出单位正方形内的一个复数 2.{a1,a2,…,an} 表示由元素a1,a2,…,an 组成的一个表,元素可以是任何内容。 ) 10(1) ()()(1 1112<<--= q q h F q qh F h F p p

求第一类Fredholm积分方程的离散正则化方法

求第一类Fredholm积分方程的离散正则化方法 【摘要】基于矩阵奇异值分解的离散正则化算法,本文给出了第一类Fredholm积分方程的求解方法。并通过算例验证了此算法的可行性。 【关键词】第一类Fredholm积分方程;矩阵奇异值;正则化方法 0 引言 在实际问题中,有很多数学物理方程反问题的求解最后总要归结为一个第一类算子方程: Kx=y(1) 的求解问题,其中K是从Hilbert空间X到Hilbert空间Y一个有界线性算子,x∈X,y∈Y。通常右端项y是观测数据,因而不可避免的带有一定的误差δ。文中假设方程(1)的右端的扰动数据yδ∈Y满足条件:yδ-y≤δ(C1)。我们需要求解扰动方程Kx=yδ∈Y。(2) 通常境况下,当K为紧算子时,方程(1)的求解时不适定的[1]。即右端数据的小扰动可导致解的巨大变化。消除不稳定性的一个自然的方式是用一族接近适定问题的模型去逼近原问题,比如最著名的Tikhonov正则化方法,用如下适定的算子方程: 去逼近原问题Kx=yδ,其中α>0为一正的“正则参数”,K*表示K的伴随算子。正则化[2-3]是近似求解方程(1)的一种有效方法。Krish应用奇异系统理论提出的正则化子的概念,这给正则化方法的建立提供了新的理论依据。本文利用基于矩阵奇异值分解的离散正则化算法,通过适当选取正则化参数进行不适定问题的求解。 1 基于矩阵奇异值分解的离散正则化算法 矩阵的奇异值分解(SVD)是现代数值线性代数中最重要的基本计算分析工具之一,它具有优良的数值稳定性。其重要应用领域包括矩阵理论以及自动控制理论,力学和物理学等,还有更多的应用方面尚在继续探索中。 对于一般算子方程Kx=y,利用高斯-勒让德求积公式、复化梯形公式或者复化辛普森求积公式等的数值方法将它离散得到一个矩阵方程Ax=y,这样,算子方程Kx=y的求解就转化为矩阵方程: 的求解。 定义设A是m×n实矩阵(m≥n),称n阶方阵ATA的非零特征值的算术平

1面向方程的数值积分方法仿真

实验1面向方程的数值积分方法仿真 1. 实验目的 运用CSS01.C 仿真程序解题,培养阅读及修改仿真程序的能力,学习并了解仿真程序的结构及特点.通过实验,加深理解4阶龙格-库塔法的原理及其稳定域. 2. 实验设备:装有BC 语言的PC 机一台 3. 实验内容 修改CSS01.C 仿真程序,对如下系统进行仿真. (1) 线性定常系统 ??????????321x x x =????????? ?---00600120000120??????????321x x x +??????????60000u y=[]00 1[]T x x x 32 1 ??????????)0()0()0(321x x x =???? ? ?????000, u=1(t) a)、对龙格库塔法进行分析:它是一种数值积分法,也就是微分方程初 值问题数值计算法,是对初值微分方程的离散化求解。对于数值积分法我们常用的是欧拉法以及二阶和四阶龙格库塔法其原理分别如下: 对于形如 dt dy =),(y t f 的微分方程 用欧拉法仿真其迭代公式为 ),(*1k k k k y t f h y y +=+,其中h 为仿真步长(下同);

二阶龙格库塔法其迭代公式为 y k+1= y k +h/2*(k1+k2),其中),(1k k y t f k =, );*,(21h k y h t f k k k ++= 四阶龙格库塔法其迭代公式为:)432221(*6/1k k k k h y y k k ++++=+, 其中),(1k k y t f k =,k2= f(t k +h/2,y k +k1*h/2),k3= f(t k +h/2,y k +k2*h/2), k4= f(t k +h,y k +k3*h); 本程序中大致分以下四次计算: 第一次计算: g[1]*h/2——k[1][2],y[1]——k[1][1],把y[1]的值赋给k[1][1] Y[1]= k[1][1]+ k[1][2]; dt dy =),(y t f ——g[1] , (计算k1), g[1]*h/2——k[1][2],(把k1*h/2存储到k[1][2]) 第二次计算: Y[1]= k[1][1]+ k[1][3]; dt dy =),(y t f ——g[1] , (计算k2), g[1]*h/2——k[1][3],(把k2*h/2存储到k[1][3]) 第三次计算: dt dy =),(y t f ——g[1],(计算k3), Y[1]= k[1][1]+ k[1][4]; k[1][4]= g[1]*h ;(把k3*h 存储到k[1][4]) 第四次计算: dt dy =),(y t f ——g[1],(计算k4), 此时 把k1,k2,k3,k4 的值代入下式: ) 432221(*6/1k k k k h y y k k ++++=+ 即为 6/)4*3**22**21*(1k h k h k h k h y y k k ++++=+ , 由于h*k1=2*k[1][2],2*h*k2=4*k[1][3],2*h*k3=2*k[1][4],h*k4=h*g[1] 所以可得如程序中的计算公式: Y[1]= k[1][1]+(2*k[1][2]+4* k[1][3]+2* k[1][4]+h*g[1]); 可见该仿真程序为四阶龙格库塔法的数值积分仿真程序。

微分方程的分类及其数值解法

微分方程的分类及其数值解法 微分方程的分类: 含有未知函数的导数,如dy/dx=2x 、ds/dt=0.4都是微分方程。 一般的凡是表示未知函数、未知函数的导数与自变量之间的关系的方程,叫做微分方程。未知函数是一元函数的,叫常微分方程;未知函数是多元函数的叫做偏微分方程。微分方程有时也简称方程。 一、常微分方程的数值解法: 1、Euler 法: 00d (,), (1.1)d (), (1.2) y f x y x y x y ?=???=? 001 (),(,),0,1,,1n n n n y y x y y hf x y n N +=??=+=-? (1.4) 其中0,n b a x x nh h N -=+=. 用(1.4)求解(1.1)的方法称为Euler 方法。 后退Euler 公式???+==+++),,(),(111 00n n n n y x hf y y x y y 梯形方法公式 )].,(),([2 111+++++=n n n n n n y x f y x f h y y 改进的Euler 方法11(,),(,),1().2p n n n c n n p n p c y y hf x y y y hf x y y y y ++?=+??=+???=+??? 2、Runge-Kutta 方法: p 阶方法 : 1()O h -=?总体截断误差局部截断误差 二阶Runge-Kutta 方法 ??? ????++==++=+),,(),,(,2212 1211hk y h x f k y x f k k h k h y y n n n n n n

数值积分的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) 一般的,步长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);

1. 积分方程一般概念与弗雷德霍姆方程

第十五章 积分方程 积分方程论是泛函分析的一个重要分支,它是研究数学其他学科(例如偏微分方程边值问题)和各种物理问题的一个重要数学工具。本章叙述线性积分方程,重点介绍弗雷德霍姆积分方程的性质和解法;并简略地介绍了沃尔泰拉积分方程以及一些奇异积分方程;此外,还扼要地叙述积分方程的逐次逼近法和预解核,并举例说明近似解法;最后考察了一个非线性积分方程。 §1 积分方程一般概念与弗雷德霍姆方程 一. 积分方程一般概念 1. 积分方程的定义与分类 [线形积分方程] 在积分号下包含未知函数y (x )的方程 ()()()()(),d b a x y x F x K x y αλξξξ=+? (1) 称为积分方程。式中α(x ),F (x )和K (x,ξ)是已知函数,λ,a,b 是常数,变量x 和ξ可取区间(a,b ) 内的一切值;K (x,ξ)称为积分方程的核,F (x )称为自由项,λ称为方程的参数。如果K (x,ξ)关于x,ξ是对称函数,就称方程(1)是具有对称核的积分方程;如果方程中的未知函数是一次的,就称为线性积分方程,方程(1)就是线性积分方程的一般形式;如果F (x )≡0 ,就称方程(1)为齐次积分方程,否则称为非齐次积分方程。 [一维弗雷德霍姆积分方程(Fr 方程)] 第一类Fr 方程 ()()(),d b a K x y F x ξξξ=? 第二类Fr 方程 ()()()(),d b a y x F x K x y λξξξ=+? 第三类Fr 方程 ()()()()(),d b a x y x F x K x y αλξξξ=+? [n 维弗雷德霍姆积分方程] 111()()()()(),d D P y P F P K P P y P P α=+? 称为n 维弗雷德霍姆积分方程,式中D 是n 维空间中的区域,P ,P 1∈D ,它们的坐标分别是 (x 1,x 2, ,x n )和),,,(21 n x x x ''' ,α(P )=α(x 1,x 2, ,x n ),F (P )=F (x 1,x 2, x n )和K (P ,P 1)=K (x 1,x 2, ,x n , ),,,21 n x x x ''' 是已知函数,f (P )是未知函数。 关于Fr 方程的解法,一维和n (>1)维的情况完全类似,因此在以后的讨论中仅着重考虑一维Fr 方程。 [沃尔泰拉积分方程] 如果积分上限b 改成变动上限,上面三类Fr 方程分别称为第一、第二、第三类沃尔泰拉积分方程。 由于第三类Fr 方程当α(x )在(a ,b )内是正函数时,可以化成

常微分方程数值解法

i.常微分方程初值问题数值解法 常微分方程初值问题的真解可以看成是从给定初始点出发的一条连续曲线。差分法是常微分方程初值问题的主要数值解法,其目的是得到若干个离散点来逼近这条解曲线。有两个基本途径。一个是用离散点上的差商近似替代微商。另一个是先对微分方程积分得到积分方程,再利用离散点作数值积分。 i.1 常微分方程差分法 考虑常微分方程初值问题:求函数()u t 满足 (,), 0du f t u t T dt =<≤ (i.1a ) 0(0)u u = (i.1b) 其中(,)f t u 是定义在区域G : 0t T ≤≤, u <∞上的连续函数,0u 和T 是给定的常数。我们假设(,)f t u 对u 满足Lipschitz 条件,即存在常数L 使得 121212(,)(,), [0,]; ,(,)f t u f t u L u u t T u u -≤-?∈∈-∞∞ (i.2) 这一条件保证了(i.1)的解是适定的,即存在,唯一,而且连续依赖于初值0u 。 通常情况下,(i.1)的精确解不可能用简单的解析表达式给出,只能求近似解。本章讨论常微分方程最常用的近似数值解法-差分方法。先来讨论最简单的Euler 法。为此,首先将求解区域[0,]T 离散化为若干个离散点: 0110N N t t t t T -=<< <<= (i.3) 其中n t hn =,0h >称为步长。 在微积分课程中我们熟知,微商(即导数)是差商的极限。反过来,差商就是微商的近似。在0t t =处,在(i.1a )中用向前差商 10()()u t u t h -代替微商du dt ,便得 10000()()(,())u t u t hf t u t ε=++ 如果忽略误差项0ε,再换个记号,用i u 代替()i u t 便得到 1000(,)u u hf t u -= 一般地,我们有 1Euler (,), 0,1, ,1n n n n u u hf t u n N +=+=-方法: (i.4) 从(i.1b) 给出的初始值0u 出发,由上式可以依次算出1,,N t t 上的差分解1,,N u u 。

实验一 面向微分方程的数值积分法仿真

实验一面向微分方程的数值积分法仿真 一、实验目的 1.掌握数值积分法的基本概念、原理及应用; 2.用龙格-库塔法解算微分方程,增加编写仿真程序的能力; 3.分析数值积分算法的计算步长与计算精度、速度、稳定性的关系; 4. 对数值算法中的“病态问题”进行研究。 二、实验内容 1、已知系统微分方程及初值条件 ,(0)1y t y y =+= 取步长0.1h =,试分别用欧拉方程法和RK4法求2t h =时的y 值,并将求得的值与解析解 ()21t y t e t =--比较(将三个解绘于同一坐标中,且用数值进行比较),说明造成差异的原 因。(①编程完成;②选用MATLAB ode 函数完成。) 程序代码如下: t0=0; tf=2; h=0.1; y1=1; y2=1; y3=1; t1=0; t2=0; t3=0 n=round(tf-t0)/h; for i=1:n y1(i+1)=y1(i)+h*(2*h+y1(i)); t1=[t1,t1(i)+h]; end for i=1:n k1=y2(i)+t2(i); k2=y2(i)+h*k1/2+t2(i)+h/2; k3=y2(i)+h*k2/2+t2(i)+h/2; k4=y2(i)+h*k3+t2(i)+h; y2(i+1)=y2(i)+h*(k1+2*k2+2*k3+k4)/6; t2=[t2,t2(i)+h]; end for i=1:n y3(i+1)=2*exp(t3(i))-t3(i)-1; t3=[t3,t3(i)+h]; end plot(t1,y1,'r',t2,y2,'g',t3,y3,'k') 实验结果如下;

数值积分 (论文)

目录 第一章数值积分计算的重述 (1) 1.1引言 (1) 1.2问题重述 (2) 第二章复化梯形公式 (3) 2.1 复化梯形公式的算法描述 (3) 2.2 复化梯形公式在C语言中的实现 (3) 2.3 测试结果 (4) 第三章复化simpson公式 (6) 3.1 复化simpson公式的算法描述 (6) 3.2 复化simpson公式在C语言中的实现 (6) 3.3 测试结果 (7) 第四章复化cotes公式 (8) 4.1 复化cotes公式的算法描述 (8) 4.2 复化cotes公式在C语言中的实现 (9) 4.3 测试结果 (10) 第五章Romberg积分法 (11) 5.1 Romberg积分法的算法描述 (11) 5.2 Romberg积分法在C中的实现 (12) 5.3 测试结果 (13) 第六章结果对比分析和体会 (144) 参考文献 (16) 附录 (16)

数值积分?-10 2 dx e x (一) 第一章 数值积分计算的重述 1.1引言 数值积分是积分计算的重要方法,是数值逼近的重要内容,是函数插值的最直接应用,也是工程技术计算中常常遇到的一个问题。在应用上,人们常要求算出具体数值,因此数值积分就成了数值分析的一个重要内容。在更为复杂的计算问题中,数值积分也常常是一个基本组成部分。 在微积分理论中,我们知道了牛顿-莱布尼茨(Newton-Leibniz)公式 ()() () b a f x d x F b F a =-? 其中()F x 是被积函数()f x 的某个原函数。但是随着学习的深入,我们发现一个问题: 对很多实际问题,上述公式却无能为力。这主要是因为:它们或是被积函数没有解析形式的原函数,或是只知道被积函数在一些点上的值,而不知道函数的形式,对此,牛顿—莱布尼茨(Newton-Leibniz)公式就无能为力了。此外,即使被积函数存在原函数,但因找原函数很复杂,人们也不愿花费太多的时间在求原函数上,这些都促使人们寻找定积分近似计算方法的研究,特别是有了计算机后,人们希望这种定积分近似计算方法能在计算机上实现,并保证计算结果的精度,具有这种特性的定积分近似计算方法称为数值积分。由定积分知识,定积分只与被积函数和积分区间有关,而在对被积函数做插值逼近时,多项式的次数越高,对被积函数的光滑程度要求也越高,且会出现Runge 现象。如7n >时,Newton-Cotes 公式就是不稳定的。因而,人们把目标转向积分区间,类似分段插值,把积分区间分割成若干小区间,在每个小区间上使用次数较低的Newton-Cotes 公式,然后把每个小区间上的结果加起来作为函数在整个区间上积分的近似,这就是复化的基本思想。本文主要

非线性Volterra积分方程(学习资料)

一类第二种非线性Volterra 积分方程 积分数值解方法 1前言 微分方程和积分方程都是描述物理问题的重要数学工具,各有优点.相对于某种情况来说,对于某种物理数学问题,积分方程对于问题的解决比微分方程更加有优势,使对问题的研究更加趋于简单化,在数学上,利用积分形式讨论存在性、唯一性往往比较方便,结果也比较完美,所以研究积分方程便得越来越有用,日益受到重视. 积分方程的发展,始终是与数学物理问题的研究息息相关.一般认为,从积分发展的源头可以追溯到国外的数学家克莱茵的著作《古今数学思想》,该书是被认为第一个清醒的认为应用积分方程求解的是Abel.Abel 分别于1833年和1826年发表了两篇有关积分方程的文章,但其正式的名称却是由数学家du Bois-Raymond 首次提出的,把该问题的研究正式命名为积分方程。所以最早研究积分方程的是Abel,他在1823年从力学问题时首先引出了积分方程,并用两种方法求出了它的解,第一的积分方程便是以Abel 命名的方程.该方程的形式为:?=-b a a x f dt t x t )()() (?,该方程称为广义Abel 方程,式中a 的值在(0,1)之间.当a=21时,该式子便成为)()(x f dt t x x x a =-??.在此之前,Laplace 于1782年所提出的求Laplace 反变换问题,当时这个问题就要求解一个积分方程.但是Fourier 其实已经求出了一类积分方程的反变换,这就说明在早些时候积分方程就已经在专业性很针对的情况下得到了研究,实际上也说明了Fourier 在研究反变换问题是就相当于解出了一类积分方程.积分方程的形成基础是有两位数学家Fredholm 和V olterra 奠定的,积分方程主要是研究两类相关的方程,由于这两位数学家的突出贡献,所以这两个方程被命名为Fredholm 方程和V olterra 方程。后来又有德国数学家D.Hilbert 进行了重要的研究,并作出了突出的贡献,由于D.Hilbert 领头科学家的研究,所以掀起了一阵研究积分方程的热潮,并出现了很多重要的成果,后来该理论又推广到非线性部分。我国在60年代前,积分方程这部分的理论介绍和相关书本主要靠翻译苏联的相关书籍,那时研究的积分方程基本是一种模式,即用古典的方法来研究相关的积分方程问题,这样使得问题的研究变得繁琐、复杂,在内容方面比较单一、狭隘,甚至有些理论故意把积分方程的研究趋向于复杂化。随着数学研究的高速发展,特别是积分方程近年来的丰富发展,如此单一、刻板的解法已经不能跟上数学研究时代的步伐。在九十年代我国的数学专家路见可、钟寿国出版了《积分方程论》,该书选择2L 空间来讨论古典积分方程,并结合泛函分析的算子理论来分析积分方程的相关问题。最近出版的比较适

微分方程数值解法答案

包括基本概念,差分格式的构造、截断误差和稳定性,这些内容是贯穿整个教材的主线。解答问题关键在过程,能够显示出你已经掌握了书上的内容,知道了解题方法。这次考试题目的类型:20分的选择题,主要是基本概念的理解,后面有五个大题,包括差分格式的构造、截断误差和稳定性。 习题一 1. 略 2. y y x f -=),(,梯形公式:n n n n n n y h h y y y h y y )121(),(2111+-+=+- =+++,所以0122)1(01])121[()121()121(y h h y h h y h h y h h n h h n n n +--+--+-+=+-+==+-+= ,当0→h 时, x n e y -→。 同理可以证明预报-校正法收敛到微分方程的解. 3. 局部截断误差的推导同欧拉公式; 整体截断误差: ? ++++++-++≤1 ),())(,(11111n n x x n n n n n n n dx y x f x y x f R εε 11)(++-++≤n n n y x y Lh R ε,这里R R n ≤ 而111)(+++-=n n n y x y ε,所以 R Lh n n += -+εε1)1(,不妨设1

积分方程

积分方程理论的发展,始终与数学物理问题的研究紧密相联,它在工程、力学等方面有着极其广泛的应用。通常认为,最早自觉应用积分方程并求出解的是阿贝尔(Abel),他在1823年研究质点力学问题时引出阿贝尔方程。此前,拉普拉斯(Laplace)於1782年在数学物理中研究拉普拉斯变换的逆变换以及傅里叶(Fourier)於1811年研究傅里叶变换的反演问题实际上都是解第一类积分方程。随着计算技术的发展,作为工程计算的重要基础之一,积分方程进一步得到了广泛而有效地应用。如今,“物理问题变得越来越复杂,积分方程变得越来越有用”。 积分方程与数学的其他分支,例如,微分方程、泛函分析、复分析、计算数学、位势理论和随机分析等都有着紧密而重要地联系。甚至它的形成和发展是很多重要数学思想和概念的最初来源和模型。例如,对泛函分析中平方可积函数、平均收敛、算子等的形成,对一般线性算子理论的创立,以至於对整个泛函分析的形成都起着重要的推动作用。积分方程论中许多思想和方法,例如,关於第二种弗雷德霍姆(Fredholm)积分方程的弗雷德霍姆理论和奇异积分方程的诺特(Noether)理论以及逐次逼近方法,本身就是数学中经典而优美的理论和方法之一。 编辑本段起源 积分号下含有未知函数的方程。其中未知函数以线性形式出现的,称为线性积分方程;否则称为非线性积分方程。积分方程起源于物理问题。牛顿第二运动定律的出现,促进了微分方程理论的迅速发展,然而对积分方程理论发展的影响却非如此。1823年,N.H.阿贝尔在研究地球引力场中的一个质点下落轨迹问题时提出的一个方程,后人称之为阿贝尔方程,是历史上出现最早的积分方程,但是在较长的时期未引起人们的注意。“积分方程”一词是 P.du B.雷蒙德于1888年首先提出的。19世纪的最后两年,瑞典数学家(E.)I.弗雷德霍姆和意大利数学家V.沃尔泰拉开创了研究线性积分方程理论的先河。从此,积分方程理论逐渐发展成为数学的一个分支。 1899年,弗雷德霍姆在给他的老师(M.)G.米塔-列夫勒的信中,提出如下的方程 公式 , (1) 式中φ(x)是未知函数;λ是参数,K(x,y)是在区域0 ≤x,y≤1上连续的已知函数;ψ(x)是在区间0≤x≤1上连续的已知函数。并认为方程(1)的解可表为关于λ的两个整函数之商。1900年,弗雷德霍姆在

微分方程数值解――

微分方程数值解―― 第二章 习题 1. 设)('x f 为)(x f 的一阶广义导数,试用类似办法定义)(x f 的k 阶广义导数) () (x f k ( ,2,1=k )。 解:对一维情形,函数的广义导数是通过分部积分来定义的。 我们知,)(x f 的一阶广义导数位)(x g ,如果满足 dx x x f dx x x g b a b a )()()()('?? -=?? 类似的,)(x f 的k 阶广义导数为)()() (x f x g k =,如果有 dx x x f dx x x g b a k k b a )()()1()()()(?? -=?? 2. 试建立与边值问题 ?????====<<=+=) 2.1(0)()(,0)()() 1.1(,''44b u b u a u a u b x a f u dx u d Lu 等价的变分问题。 证明: 设}0)()(,0)()(),(|{' '2====∈=b v a v b v a v I H v v V 对方程)1.1(两边同乘以v ,再关于x 在),(b a 上积分)(V v ∈,得 ??=+b a b a fvdx vdx u dx u d )(44 其中 dx dx dv dx u d dx dx dv dx u d dx u d v dx u d d v vdx dx u d b a b a b a b a b a ???? -=-==33 33333344|)( dx dx v d dx u d dx dv dx u d dx u d d dx dv b a b a b a ??+-=-=22222222|)( dx dx v d dx u d b a ? = 2 222 (*) 记dx uv dx v d dx u d v u a b a ?+=)(),(2 222,?=b a fvdx v f ),(。于是我们得到以下等价变分问题的提法:

相关主题
文本预览
相关文档 最新文档