任意拉格朗日-欧拉描述法研究进展

任意拉格朗日-欧拉描述法研究进展
任意拉格朗日-欧拉描述法研究进展

微分方程数值解法

《微分方程数值解法》 【摘要】自然界与工程技术中的很多现象,可以归结为微分方程定解问题。其中,常微分方程求解是微分方程的重要基础内容。但是,对于许多的微分方程,往往很难得到甚至不存在精确的解析表达式,这时候,数值解提供了一个很好的解决思路。,针对于此,本文对常微分方程数值解法进行了简单研究,主要讨论了一些常用的数值解法,如欧拉法、改进的欧拉法、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)可得:

常微分方程作业欧拉法与改进欧拉法

P77 31.利用改进欧拉方法计算下列初值问题,并画出近似解的草图:dy + =t = t y y ≤ ≤ ,2 ;5.0 0,3 )0( )1(= ,1 ? dt 代码: %改进欧拉法 function Euler(t0,y0,inv,h) n=round(inv(2)-inv(1))/h; t(1)=t0; y(1)=y0; for i=1:n y1(i+1)=y(i)+h*fun(t(i),y(i)); t(i+1)=t(i)+h; y(i+1)=y(i)+1/2*h*(fun(t(i),y(i))+ fun(t(i+1),y1(i+1))) end plot(t,y,'*r') function y=fun(t,y); y=y+1; 调用:Euler(0,3,[0,2],0.5) 得到解析解:hold on; y=dsolve('Dy=y+1','(y(0)=3)','t'); ezplot(y,[0,2]) 图像:

dy y =t - t y ;2.0 t = ≤ )0( 0,5.0 ,4 )2(2= ≤ ? ,2 dt 代码: function Euler1(t0,y0,inv,h) n=round(inv(2)-inv(1))/h; t(1)=t0; y(1)=y0; for i=1:n y1(i+1)=y(i)+h*fun(t(i),y(i)); t(i+1)=t(i)+h; y(i+1)=y(i)+1/2*h*(fun(t(i),y(i))+ fun(t(i+1),y1(i+1))) end plot(t,y,'*r') function y=fun(t,y); y=y^2-4*t; 调用: Euler1(0,0.5,[0,2],0.2) 图像:

微分方程常用的两种数值解法:欧拉方法与龙格—库塔法

四川师范大学本科毕业论文 微分方程常用的两种数值解法:欧拉方法与龙 格—库塔法 学生姓名XXX 院系名称数学与软件科学学院 专业名称信息与计算科学 班级2006级 4 班 学号20060640XX 指导教师Xxx 四川师范大学教务处 二○一○年五月

微分方程常用的两种数值解法:欧拉方法与龙格—库塔法 学生姓名:xxx 指导教师:xx 【内容摘要】微分方程是最有生命力的数学分支,在自然科学的许多领域中,都 会遇到常微分方程的求解问题。当前计算机的发展为常微分方程的应用及理论研究提供了非常有力的工具,利用计算机解微分方程主要使用数值方法,欧拉方法和龙格——库塔方法是求解微分方程最典型常用的数值方法。本文详细研究了这两类数值计算方法的构造过程,分析了它们的优缺点,以及它们的收敛性,相容性,及稳定性。讨论了步长的变化对数值方法的影响和系数不同的同阶龙格—库塔方法的差别。通过编制C程序在计算机上实现这两类方法及对一些典型算例的结果分析比较,能更深切体会它们的功能,优缺点及适用场合,从而在实际应用中能对不同类型和不同要求的常微分方程会选取适当的求解方法。 关键词:显式单步法欧拉(Euler)方法龙格—库塔(Runge—Kutta)方法截断误差收敛性 Two commonly used numerical solution of differential equations:Euler method and Runge - Kutta method Student Name: Xiong Shiying Tutor:Zhang Li 【Abstract】The differential equation is the most vitality branch in mathematics. In many domains of natural science, we can meet the ordinary differential equation solution question. Currently, the development of computer has provided the extremely powerful tool for the ordinary differential equation application and the fundamental research, the computer solving differential equation mainly uses value method. The Euler method and the Runge—Kutta method are the most typical commonly value method to solve the differential equation. This article dissects the structure process of these two kinds of values commonly value method to solve the analyses their good and bad points, to their astringency, the compatibility, and the stability has made the proof. At the same time, the article discuss the length of stride to the numerical method changing influence and the difference of the coefficient different same step Runge—kutta method. Through establishing C program on the computer can realize these two kind of methods, Anglicizing some models of calculate example result can sincerely realize their function, the advantage and disadvantage points and the suitable situation, thus the suitable solution method can be selected to solve the different type and the

《图解刚体力学——欧拉运动学方程》

本科生毕业论文 论文题目:图解刚体力学——欧拉运动学方程 学生姓名:罗加宽 学号: 2008021152 专业名称:物理学 论文提交日期: 2012年05月17日 申请学位级别:理学学士 论文评审等级: 指导教师姓名:陈洛恩 职称:教授 工作单位:玉溪师范学院 学位授予单位:玉溪师范学院 玉溪师范学院理学院物理系 2012年05月

图解刚体力学—欧拉运动学方程 罗加宽 (玉溪师范学院理学院物理系 08级物理1班云南玉溪 653100) 指导教师:陈洛恩、杨春艳 摘要:本文阐述了描述刚体定点转动的欧拉角及欧拉运动学方程的图解,以期让复杂的问题转 化得简单清晰而易于学习者的理解,抽象的概念变得直观具体而易于学习者的掌握;并能在一 定程度上对提高学习者的空间思维能力、引导和培养学习者的创新思维能力有一定的帮助。 关键字:图解;刚体;欧拉角;欧拉运动学方程 1.引言 理论力学是研究物体机械运动一般规律的科学;依照牛顿的说法,理论力学“是关于力产生的运动和产生任何运动的力的理论,是精确的论述和证明” [1]。理论力学作为使用数学方法的自然知识的一部分,不仅研究实际物体,而且研究其模型—质点、质点系、刚体和连续介质。从研究次序来看,通常先研究描述机械运动现象的运动学,然后再进一步研究机械运动应当遵循哪些规律的动力学。至于研究平衡问题的静力学,对理科来讲可以作为动力学的一部分来处理,但在工程技术上,静力学却是十分的重要,因此,常把它和动力学分开,自成一个系统[2]。本文图解的内容为刚体力学运动学问题之一的刚体的绕定点的转动。 “图解”的方法,较早见于上海科学技术出版社1988年翻译出版的《图解量子力学》,原书名为The Picture Book of Quantum Mechanics,由Springer-Verlag 出版;类似的书还有Springer-Verlag出版的Visual Quantum Mechanics。其特点是通过将理论物理与数值计算相结合实现可视化来讲解物理知识。国外对物理的可视化教学十分重视,早在1995-1996年间Wiley出版社出版了9本有关物理多媒体教学的丛书,是由大学高等物理软件联盟(The Consortium for Upper-Level Physics Software,CUPS)编写该丛书及其所用的教学软件[3]。如今,图解法已经广泛应用于力学、电磁学、模拟电子技术等方面,理论力学方面同样也有不少人已经采用了图解法。如赵宗杰使用3dsmax建立质点外弹道运动规律的虚拟模型和场景[4];乐山师范学院王峰等利用Matlab分别对质点受力仅为位置、速度或时间的函数进行了图解,并说明了Matlab在理论力学中的应用[5];阜阳师范学院孙美娟、韩修林利用Mathematica进行编程作出了落体的位移—时间图像[6]。通过图解,使很多抽象繁难的物理问题在解析时达到空间立体直观,概念形成清晰,逻辑链路晓畅明朗,数式转换准确易见。 理论力学因理论性较强,与高等数学联系密切,一些概念的形成、公式的推导、逻辑推理等较抽象、繁难、复杂,往往使教授者感到教学很难达到预期的效果,学

第二章 用拉格朗日方程建立系统数学模型

第二章 用拉格朗日方程建立系统的数学模型 §2.1概述 拉格朗日方程——属于能量法,推导中使用标量,直接对整个系统建模 特点:列式简洁、考虑全面、建模容易、过程规范 适合于线性系统也适合于非线性系统,适合于保守系统,也适合于非保守系统。 §2.2拉格朗日方程 1. 哈密尔顿原理 系统总动能 ),,,,,,,(321321N n q q q q q q q q T T = (2-1) 系统总势能 ),,,,(321t q q q q U U N = (2-2) 非保守力的虚功 N N nc q Q q Q q Q W δδδδ ++=2211 (2-3) 哈密尔顿原理的数学描述: 0)(2 1 21 =+-??t t nc t t dt W dt U T δδ (2-4) 2. 拉格朗日方程: 拉格朗日方程的表达式: ),3,2,1()(N i Q q U q T q T dt d i i i i ==??+??-?? (2-5) (推导:) 将系统总动能、总势能和非保守力的虚功的表达式代入哈密尔顿原理式中(变分驻值原理),有 0)( 22112211221122112 1 =+++??-??-??-??++??+??+??+??+??? dt q Q q Q q Q q q T q q U q q U q q T q q T q q T q q T q q T q q T N N N N N N N N t t δδδδδδδδδδδδ (2-6) 利用分步积分

dt q q T dt d q q T dt q q T i t t i t t i i i t t i δδδ?? ??-??=??21212 1 )(][ (2-7) 并注意到端点不变分(端点变分为零) 0)()(21==t q t q i i δδ (2-8) 故 dt q q T dt d dt q q T i i t t i t t i δδ)(212 1 ??-=???? (2-9) 从而有 0)])([2 1 1 =+??-??+??- ?∑=dt q Q q U q T q T dt d i i i t t i i N i δ ( (2-10) 由变分学原理的基本引理: (设 n 维向量函数M(t),在区间],[0f t t 内处处连续,在],[0f t t 内具有二阶连续导 数,在f t t ,0处为零,并对任意选取的n 维向量函数)(t η,有 ? =f t t T dt t M t 0 0)()(η 则在整个区间],[0f t t 内,有 0)(≡t M ) 我们可以得到: 0)(=+??-??+??- i i i i Q q U q T q T dt d (2-11) 即 i i i i Q q U q T q T dt d =??+??-??)( (2-12) 对非保守系统,阻尼力是一种典型的非保守力,如果采用线性粘性阻尼模型, 则阻尼力与广义速度}{q 成正比,在这种情况下,可引入瑞利耗散(耗能)函数D , }]{[}{2 1 q C q D T ≡ (2-13) 阻尼力产生的广义非保守力为:

任意拉格朗日欧拉(ALE)理论基础

暨两岸船舶与海洋工程水动力学研讨会文集 拉格朗日、欧拉和任意拉格朗日-欧拉描述的 有限元分析 孙江龙1杨文玉2 杨侠3 (1 华中科技大学船舶与海洋工程学院,武汉 430074;2 华中科技大学机械科学与工程学院,武汉 430074;3 武 汉工程大学机电工程学院,武汉 430073) 摘要:对拉格朗日、欧拉和任意拉格朗日–欧拉三种描述方法进行了分析,为了便于理解给出了三种描述的参考构形和参考坐标系,在参考坐标系下根据物质导数的定义分别得到相应的速度和加速度,并进行比较,将三种描述方法的区别列于表中,清晰地阐述了三种描述之间的相互关系,并进行了有限元分析。 关键词:拉格朗日;欧拉;任意拉格朗日–欧拉;有限元法 1 引言 自由液面大晃动引起的强非线性往往给问题的求解造成很大困难,对大晃动问题进行数值模拟,要先解决描述方法的选择问题。过去通常采用欧拉法[1-3]和拉格朗日法[4-5]来描述非定常自由面流体流动,它们有着各自的优势和局限性。 采用固定网格的欧拉描述,整个计算过程中计算网格始终保持初始状态,从而可以描述流体质点运动的急剧变化,如碎波等现象。欧拉描述虽然可以有效地分析整个流场内部的运动,但很难精确跟踪流体的自由液面,即很难给出准确的自由面形状和位置。 在拉格朗日描述中,网格结点与流体质点在整个运动过程中始终保持重合,流体质点与网格结点之间不存在相对运动,因此很容易跟踪自由液面,适用于线性小晃动问题。这不仅大大地简化了控制方程地求解,而且还能有效地跟踪流体质点的运动轨迹,准确地描述波动的自由液面。但是,在涉及求解带自由面流体大幅运动时,此时的晃动已经具有很强的非线性特征,如果还采用拉格朗日描述,由于流体质点运动的急剧变化,将导致计算网格的扭曲,会面临网格奇异问题,从而使计算无法继续进行。 拉格朗日描述和欧拉描述虽有各自的优点,但也存在较大的缺陷,如果将它们有机地结合在一起,充分利用各自的优点并克服其缺点,则可以解决各自都难于解决的问题,任意拉格朗日–欧拉描述[6-7](ALE)方法就是基于该思路提出的。在任意拉格朗日–欧拉描述中,网格结点的运动方式比较灵活,网格结点可以跟随流体质点一起运动,也可以固定不变,甚至可以采用网格结点在一个方向上固定而在其他方向上随流体质点一起运动等方式。为了更加清晰地理解这三种描述方法,本研究从以下几个方面进行阐述和比较。 - 164 -

欧拉及改进的欧拉法求解常微分方程

生物信息技术0801 徐聪U200812594 #include #include void f1(double *y,double *x,double *yy) { y[0]=2.0; x[0]=0.0; yy[0]=2.0; for(int i=1;i<=9;i++) { x[i]=x[i-1]+0.2; y[i]=y[i-1]+0.2*(y[i-1]-x[i-1]); yy[i]=x[i]+1+exp(x[i]); printf("若x=%f,计算值是%f,真实值是%f,截断误差是%f\n ",x[i],y[i],yy[i],y[i]-yy[i]); } }; void f2(double *y,double *x,double *yy) { y[0]=1.0; x[0]=0.0; yy[0]=1.0; for(int i=1;i<=9;i++) { x[i]=x[i-1]+0.2; y[i]=y[i-1]+0.2*(2*y[i-1]+x[i-1]*x[i-1]); yy[i]=-0.5*(x[i]*x[i]+x[i]+0.5)+1.25*exp(2*x[i]); printf("若x=%f,计算值是%f,真实值是%f,截断误差是%f\n ",x[i],y[i],yy[i],y[i]-yy[i]); } }; void f3(double *y,double *x,double *yy,double *y0) { y[0]=2.0; x[0]=0.0; yy[0]=2.0; for(int i=1;i<=9;i++) { x[i]=x[i-1]+0.2; y0[i]=y[i-1]+0.2*(y[i-1]-x[i-1]); y[i]=y[i-1]+0.1*(y[i-1]-x[i-1]+y0[i-1]-x[i-1]);

对流体力学欧拉运动方程式的修正

第9节 对流体力学欧拉运动方程式的修正(探讨) 内容提要:本文是探讨性的论文。观念正确如否有待学界审视及实践的检验。流体力学的欧拉运动方程式有修正的必要吗?首先,欧拉运动方程式是在《场论》只具有散度和旋度的数学基础为背景的产物;其次,人们注意到,航天器在飞行运动中存在一未知的莫铭的力。这个莫铭的力应该是欧拉方程尚未虑及的因素造成的。作者在研究《超变函数论》过程中揭示了在三维向量场中除了散度、旋度外尚存在一个为目前所未知的副冲量度【见文献3】。 我们所提出的修正意见就是从这里切入的,即在考虑存在副冲量度这一因素后,欧拉运动方程式应该发生怎样的变化。 关键词:理想流体,时变加速度,位变加速度,欧拉运动方程式,副冲量度,冲量力,压扁的四维空间. 分类号: 一,现在的欧拉运动方程式[见文献4,第77页] 在理想流体场中取出一微小六面体流体微团。微团中心的压力为P ,速度为,,x y z ωωω。微团所受的力有表面力(压力)和体积力(质量力)。六面体各面所受的表面力如下图所示。体积力为,,x y z F F F 。设单位质量的的体积力为X,Y,Z ,则在x 轴方向微团所受的力为 ()()22 ()??+-+???=-?dx dx X dxdydz dydz dydz x x X dxdydz x P P ρP - P P ρ 在x 轴方向微团产生加速度的运动力为 x d dxdydz dt ωρ 【注:其中,总加速度 ???????=+++???????????=+++????y x x x z y x x z x y z d x y z dt t x t y t z t t x y z ωωωωωωωωωωωω 该式右侧第一项称为时变加速度;第二、三、四项总称为位变加速度。】

第8章 常微分方程数值解法 本章主要内容: 1.欧拉法

第8章 常微分方程数值解法 本章主要内容: 1.欧拉法、改进欧拉法. 2.龙格-库塔法。 3.单步法的收敛性与稳定性。 重点、难点 一、微分方程的数值解法 在工程技术或自然科学中,我们会遇到的许多微分方程的问题,而我们只能对其中具有较简单形式的微分方程才能够求出它们的精确解。对于大量的微分方程问题我们需要考虑求它们的满足一定精度要求的近似解的方法,称为微分方程的数值解法。本章我们主要 讨论常微分方程初值问题?????==00 )() ,(y x y y x f dx dy 的数值解法。 数值解法的基本思想是:在常微分方程初值问题解的存在区间[a,b]内,取n+1个节点a=x 0<x 1<…<x N =b (其中差h n = x n –x n-1称为步长,一般取h 为常数,即等步长),在这些节点上把常微分方程的初值问题离散化为差分方程的相应问题,再求出这些点的上的差分方程值作为相应的微分方程的近似值(满足精度要求)。 二、欧拉法与改进欧拉法 欧拉法与改进欧拉法是用数值积分方法对微分方程进行离散化的一种方法。 将常微分方程),(y x f y ='变为() *+=?++1 1))(,()()(n x n x n n dt t y t f x y x y 1.欧拉法(欧拉折线法) 欧拉法是求解常微分方程初值问题的一种最简单的数值解法。 欧拉法的基本思想:用左矩阵公式计算(*)式右端积分,则得欧拉法的计算公式为:N a b h N n y x hf y y n n n n -= -=+=+)1,...,1,0(),(1 欧拉法局部截断误差 11121 )(2 ++++≤≤''=n n n n n x x y h R ξξ或简记为O (h 2)。

微分方程数值解

微 分方程数值解及其应用 绪论 自然界中的许多事物的运动和变化规律都可以用微分方程来描述,因此对工程和科学技术中的实际问题的研究中, 常常需要求解微分方程.但往往只有少数较简单和典型的微分方程可求出其解析解,在大多数情况下,只能用近似法求解,数值解法是一类重要的近似方法.本文主要讨论一阶常微分方程的初值问题的数值解法,探讨这些算法在处理来自生活实际问题中的应用,并结合MATLAB 软件,动手编程予以解决. 1 微分方程的初值问题[1] 1.1 预备知识 在对生活实际问题的研究中,通常需要考虑一阶微分方程的初值问题 00(,)()dy f x y dx y x y ?=???=? (1) 这里(),f x y 是矩形区域R :00,x x a y y b -≤-≤上的连续函数. 对初值问题(1)需要考虑以下问题:方程是否一定有解呢?若有解,有多少个解呢?下面给出相关的概念与定理. 定义1 Lipschitz 条件[1][2]:矩形区域R :00,x x a y y b -≤-≤上的连续函数(),f x y 若满足:存在常数0L >,使得不等式()()1212,,f x y f x y L y y -≤-对所有()()12,,,x y x y R ∈都成立,则称(),f x y 在R 上关于y 满足Lipschitz 条件. 定理 1 解的存在唯一性定理[1][3]:设f 在区域()}{,,D x y a x b y R =≤≤∈上连续,关于y 满足Lipschitz 条件,则对任意的[]00,,∈∈x a b y R ,常微分方程初值问题(1)当[],x a b ∈时存在唯一的连续解()y x . 该定理保证若一个函数(),f x y 关于y 满足Lipschitz 条件,它所对应的微分方程的初值问题就有唯一解.在解的存在唯一性得到保证的前提下,自然要考虑方程的求

牛顿—欧拉方程

牛顿-欧拉方程 欧拉方程(Euler equations),是欧拉运动定律的定量描述,欧拉运动定律是牛顿运动定律的延伸,在牛顿发表牛顿运动定律超过半个世纪后,于1750年,欧拉才成功的用欧拉方程表述了该定律: 该方程是建立在角动量定理的基础上的描述刚体的旋转运动时刚体所受外力矩与角加速度的关系式,大多时候可简写成: 其中,分别为刚体坐标系下三个轴的所受的外力矩,分别为刚体三个坐标轴的转动惯量(刚体坐标系下)。 欧拉方程通常与牛顿的平移运动方程被一起写出,称为牛顿-欧拉方程(Newton-Euler equations): 这里对牛顿的平移运动方程不赘述,只对欧拉方程进行讨论。 1.单质点角动量定理 质点旋转时,有动量定理: 对两边叉乘质点位置矢量:

观察: 因为: 故有: 定义角动量,可以看出为外力矩 故有单质点的角动量定理: 2.刚体的角动量定理 定义刚体的角动量为: 其中:下标G表示该向量为大地坐标系下的,的下标i表示该向量为大地坐标下各个质量元的向量。刚体旋转运动参考的惯性系是大地坐标系,不能把采用刚体的本身坐标系作为参考系,本身坐标系的提出只是方便我们某些量的分析与表述,如角速度、惯性张量。

(这里需要特别说明的是因为刚体质量分布不均匀的原因,角动量的方向往往不与刚体角速度方向一致,这也是无力矩进动的原因,即很多时候刚体角速度不守恒但刚体的角动量守恒了,宏观来看就是因为要保证角动量和动量守恒所以才要产生内力作用使角速度变化达到守恒的效果。) 由牛顿第三定律易知内力矩产生的角动量变化相抵,故有刚体的角动量定理: 其中:为外力矩 把上式展开有: 其中:称为惯性矩阵

MATLAB求解常微分方程数值解

利用MATLAB求解常微分方程数值解

目录 1. 内容简介 (1) 2. Euler Method(欧拉法)求解 (1) 2.1. 显式Euler法和隐式Euler法 (2) 2.2. 梯形公式和改进Euler法 (3) 2.3. Euler法实用性 (4) 3. Runge-Kutta Method(龙格库塔法)求解 (5) 3.1. Runge-Kutta基本原理 (5) 3.2. MATLAB中使用Runge-Kutta法的函数 (7) 4. 使用MATLAB求解常微分方程 (7) 4.1. 使用ode45函数求解非刚性常微分方程 (8) 4.2. 刚性常微分方程 (9) 5. 总结 (9) 参考文献 (11) 附录 (12) 1. 显式Euler法数值求解 (12) 2. 改进Euler法数值求解 (12) 3. 四阶四级Runge-Kutta法数值求解 (13) 4.使用ode45求解 (14)

1.内容简介 把《高等工程数学》看了一遍,增加对数学内容的了解,对其中数值解法比较感兴趣,这大概是因为在其它各方面的学习和研究中经常会遇到数值解法的问题。理解模型然后列出微分方程,却对着方程无从下手,无法得出精确结果实在是让人难受的一件事情。 实际问题中更多遇到的是利用数值法求解偏微分方程问题,但考虑到先从常微分方程下手更为简单有效率,所以本文只研究常微分方程的数值解法。把一个工程实际问题弄出精确结果远比弄清楚各种细枝末节更有意思,因此文章中不追求非常严格地证明,而是偏向如何利用工具实际求解出常微分方程的数值解,力求将课程上所学的知识真正地运用到实际方程的求解中去,在以后遇到微分方程的时候能够熟练运用MATLAB得到能够在工程上运用的结果。 文中求解过程中用到MATLAB进行数值求解,主要目的是弄清楚各个函数本质上是如何对常微分方程进行求解的,对各种方法进行MATLAB编程求解,并将求得的数值解与精确解对比,其中源程序在附录中。最后考察MATLAB中各个函数的适用范围,当遇到实际工程问题时能够正确地得到问题的数值解。 2.Euler Method(欧拉法)求解 Euler法求解常微分方程主要包括3种形式,即显式Euler法、隐式Euler法、梯形公式法,本节内容分别介绍这3种方法的具体内容,并在最后对3种方法精度进行对比,讨论Euler法的实用性。 本节考虑实际初值问题 使用解析法,对方程两边同乘以得到下式

微分方程数值解欧拉法

1.1、求解初值问题()?????=-=-1 0y y xe dx dy x ,已知精确解为 ()()x x x x y -+=22 12 当h=0.1时,解为: n x n y ()n x y ()n n y x y - 0 1 1 0 0.1 0.900000 0.909362 9.3616E-03 0.2 0.819048 0.835105 1.6057E-02 0.3 0.753518 0.774155 2.0637E-02 0.4 0.700391 0.723946 2.3555E-02 0.5 0.657165 0.682347 2.5182E-02 0.6 0.621775 0.647598 2.5823E-02 0.7 0.592526 0.618249 2.5723E-02 0.8 0.568034 0.593114 2.5080E-02 0.9 0.547177 0.571230 2.4053E-02 1.0 0.529051 0.551819 2.2768E-02 0.1 0.2 0.30.40.50.60.70.80.91 当h=0.05时,解为:

n x n y ()n x y ()n n y x y - 0 1 1 0 0.05 0.950000 0.952418 2.4185E-03 0.10 0.904878 0.909362 4.4835E-03 0.15 0.864158 0.870391 6.2326E-03 0.20 0.827406 0.835105 7.6996E-03 0.25 0.794223 0.803138 8.9155E-03 0.30 0.764247 0.774155 9.9084E-03 0.35 0.737147 0.747850 1.0704E-02 0.40 0.712621 0.723946 1.1324E-02 0.45 0.690397 0.702188 1.1791E-02 0.50 0.670223 0.682347 1.2124E-02 0.55 0.651876 0.664213 1.2338E-02 0.60 0.635148 0.647598 1.2450E-02 0.65 0.619855 0.632328 1.2473E-02 0.70 0.605829 0.618249 1.2420E-02 0.75 0.592918 0.605220 1.2302E-02 0.80 0.580985 0.593114 1.2129E-02 0.85 0.569909 0.581819 1.1909E-02 0.90 0.559579 0.571230 1.1651E-02 0.95 0.549896 0.561258 1.1362E-02 1.00 0.540771 0.551819 1.1048E-02 0.1 0.2 0.30.40.50.60.70.80.91

fortran下欧拉法求解常微分方程(实例)

1. Euler 公式 100(,)() i i i i y y hf x y y y x +=+??=? 实例: ,00(,),0,1,01f x y x y x y x =-==≤≤ 精确解为:1x y x e -=+- 程序代码: DIMENSION x(0:20),y(0:20),z(0:20),k(0:21) DOUBLE PRECISION x,y,z,k,h,x0,y0,z0,k0,n f(x,y)=x-y n=20 h=1/n x(0)=0 y(0)=0 DO i=0,n-1 y(i+1)=y(i)+f(x(i),y(i))*h x(i+1)=x(i)+h ENDDO k(0)=0 DO i=0,n z(i)=k(i)+exp(-k(i))-1 k(i+1)=k(i)+h END DO open(10,file='1.txt') WRITE(10,10) (x(i),y(i),z(i),i=0,20) WRITE(*,10) (x(i),y(i),z(i),i=0,20) 10 FORMAT(1x,f10.8,2x,f10.8,2x,f10.8/) END 输出结果: 0.00000000 0.00000000 0.00000000 0.05000000 0.00000000 0.00122942 0.10000000 0.00250000 0.00483742 0.15000000 0.00737500 0.01070798 0.20000000 0.01450625 0.01873075 0.25000000 0.02378094 0.02880078 ???=='00)(),(y x y y x f y ???=='0 0)(),(y x y y x f y

全欧拉法流固耦合问题

摘要 在这项工作中,我们提出一个关于流固耦合的一个完全欧拉结构(FSI)问题即耦合不可压缩Navier-Stokes超弹性固体方程。完全欧拉结构是一个单一变量结构耦合问题。与此相反的任意拉格朗日欧拉(ALE)坐标也是行之有效的,但是完全欧拉结构,包含两个子问题,液体和固体问题。在欧拉坐标,这个概念绕过与ALE各种困难坐标的联系,因为没有其他人的领域成果可以使用,所以该结构是一种首创全新的方法。本方法主要研究该固体的变形,它作为一种固体变形研究的的扩展,建立在初始点的设置,以此检测连接点的位置。由于涉及到大变形,所以本文尽可能的利用了像固体接触边界的变化或其他固体领域的研究。1介绍 我们为流固耦合问题提供了一种完全统一的可变有限单元法。重点强调于大变形结构领域的应用,处于流域和该流域边界接触的结构的自由运动,和其他自接触结构的应用。这项工作中出现的方程是Eulerian-Eulerian型的,并且这种奇特的方法第一次出现是由Dunne引入的[14,15]。 现在有数不清的不同方法来建模和模拟流固耦合问题。在这些方法中我们专注于单片模型,这种模型整个问题被描述为一种包含固体和流体表面的双系统。单片模型允许含糊的描述方案,大的时间步长并且提供使用基于错误假设的可能和最优化方法。他们已经很好的适用于大流体密度问题中如血液动力学。当流体问题被本质的描述为一种混合欧拉或朗格朗日的结构,材料的描述通常是固体问题的基础。所有的描述流固耦合问题的单片模型不管怎么样都要符合这两种结构。 在拉格朗日法或专一的拉格朗日法中,流体问题被定位于一种涉及域的匹配的位置。经典的方法是ALE法,见[29.4.35]或变空间域/稳态空间法(DSD/SST),见实例[53,51]。这些方程有相同之处,运动学和动力学结合很容易嵌入路径空间并有多种技术支持。拉格朗日法的缺点是流体问题的转化可在大变形或大的固体运动中破坏。 欧拉-拉格朗日法为解决流体问题使用一种欧拉混合计算单元,为解决固体 问题使用朗格朗日单元。通过使用力密度法结合这两种结构的应用于大边界[43]或大表面[39]。其它介绍的方法在于额外的表面不同。这些例子都是虚构的域法,其他最近的方法都是基于扩展有限单元法的。表面混合法的调查由Felippa[17] 给出。欧拉-拉格朗日法是一种表面捕捉的方法。这种表面不是修复的欧拉单元

常微分方程欧拉算法

常微分方程欧拉算法 Company Document number:WUUT-WUUY-WBBGB-BWYTT-1982GT

常微分方程欧拉算法 摘要:本文主要论述了常微分方程的欧拉算法的算法原理,误差分析,实例,程序,以及算法比较等内容。 关键词:常微分方程 显式欧拉法 隐式欧拉法 引言:微分方程初值问题模型是常见的一类数学模型。对于一些简单而典型的微分方程模型,譬如线性方程、某些特殊的一阶非线性方程等是可以设法求出其解析解的,并有理论上的结果可资利用。但在数学建模中碰到的常微分方程初值问题模型,通常很难,甚至根本无法求出其解析解,而只能求其近似解。因此,研究其数值方法,以便快速求得数值鳃有其重大意义。 一、欧拉算法原理 对于微分方程初值问题 的解在xy 平面上是一条曲线,称为该微分方程的积分曲线。积分曲线上一点(),x y 的切线斜率等于函数f 在点(),x y 的值,从初始点()000,P x y 出发,向该点的切线方向推进到下一个点()111,P x y ,然后依次做下去,得到后面的未知点。一般地,若知道(),n n n P x y 依上述方法推进到点()111,n n n P x y +++,则两点的坐标关系为: 即 这种方法就是欧拉(Euler )方法(也叫显式欧拉法或向前欧拉法)。当初值0y 已知,则n y 可以逐步算出 对微分方程()=x y dy f dx ,从n x 到1n x +积分,那么有 现在用左矩形公式()(),n n hf x y x 代替()()1 ,n n x x f t y t dt +?,n y 代替()n y x ,1n y +代替() 1n y x +就得到了欧拉方法。如果用右矩形公式()()11,n n hf x y x ++去代替右端积分,则得到另外一 个公式,该方法就称为隐式欧拉法(或后退欧拉法),其公式为 欧拉公式与隐式欧拉公式的区别在于欧拉公式是关于1n y +的一个直接计算公式,然而隐式欧拉公式右端含有1n y +,所以它实际上是关于1n y +的一个函数方程。 二、实例 例 取h=,用Euler 方法解

微分方程数值解欧拉法

dy??x??xey?1.1、求解初值问题,已知精确解为 ????x?2xx?y?2x2当h=0.1时,解为:?dx????01y?1 ????yxy?xyyx nnnnn 1 0 1 9.3616E-03 0.1 0.900000 0.909362 1.6057E-02 0.819048 0.2 0.835105 2.0637E-02 0.774155 0.753518 0.3 2.3555E-02 0.723946 0.4 0.700391 2.5182E-02 0.5 0.682347 0.657165 2.5823E-02 0.621775 0.6 0.647598 2.5723E-02 0.592526 0.618249 0.7 2.5080E-02 0.568034 0.8 0.593114 2.4053E-02 0.547177 0.9 0.571230 2.2768E-02 1.0 0.551819 0.529051

1 0.950.90.850.80.750.70.650.60.550.510.100.20.80.70.90.60.40.30.5时,解为:h=0.05 当. ????x xyy y?yx nnn nn 1 1 0 2.4185E-03 0.952418 0.05 0.950000 4.4835E-03 0.10 0.909362 0.904878 6.2326E-03 0.15 0.864158 0.870391 7.6996E-03 0.827406 0.20 0.835105 8.9155E-03 0.794223 0.25 0.803138 9.9084E-03 0.774155 0.764247 0.30 1.0704E-02 0.737147 0.747850 0.35 1.1324E-02 0.723946 0.40 0.712621 1.1791E-02 0.702188 0.45 0.690397 1.2124E-02 0.50 0.670223 0.682347

MATLABEuler法解常微分方程

Euler法解常微分方程 Euler法解常微分方程算法: Step 1 分别取积分上限、积分下限、步长 Step 2计算判断是否成立,成立转到Step 3,否则继续进行Step 4 Step 3 计算 Step 4 Euler法解常微分方程算程序: function euler2(fun,y0,A,h) %fun--y' %y0---初值 %A----x取值范围 %a----x左区间端点值 %b----x右区间端点值 %h----给定步长 x=min(A); b=max(A); y=y0; while x

Step 3 (1)做显性Euler预测 (2)将带入 Step 4计算判断是否成立,成立返回Step 3,否则继续进行Step 5 Step 5 改进Euler法解常微分方程算程序: function gaijineuler2(fun,y0,A,h) %fun--y' %y0---初值 %A----x取值范围 %a----x左区间端点值 %b----x右区间端点值 %h----给定步长 a=min(A); b=max(A); x=a:h:b; y(1)=y0; for i=1:length(x)-1 w1=feval(fun,x(i),y(i)); y(i+1)=y(i)+h*w1; w2=feval(fun,x(i+1),y(i+1)); y(i+1)=y(i)+h*(w1+w2)/2; end x=x' y=y' 例:用改进Euler法计算下列初值问题(取步长h=0.25) 输入:fun=inline('-x*y^2') gaijineuler2(fun,2,[0 5],0.25) 得到: x = 0.2500 0.5000 0.7500 1.0000 1.2500 1.5000 1.7500 2.0000 2.2500 2.5000 2.7500

欧拉运动微分方程各项的单位

第四章 1 欧 拉 运 动 微 分 方 程 d d u f t p =-?1 ζ 各 项 的 单 位 是: (1) 单 位 质 量 力 (2) 单 位 重 能 量 (3) 单 位 重 的 力 (4) 上 述 回 答 都 不 对 2. 欧 拉 运 动 微 分 方 程 在 每 点 的 数 学 描 述 是: (1)流入的质量流量等于流出的质量流量(2) 单 位 质 量 力 等 于 加 速 度 (3) 能 量 不 随 时 间 而 改 变 (4) 服 从 牛 顿 第 二 定 律 3. 欧 拉 运 动 微 分 方 程: (1) 适 用 于 不 可 压 缩 流 体, 不 适 用 于 可 压 缩 流 体 (2) 适 用 于 恒 定 流, 不 适 用 非 恒 定 流 (3) 适 用 于 无 涡 流, 不 适 用 于 有 涡 流 (4) 适 用 于 上 述 所 提 及 的 各 种 情 况 下 流 体 流 动 4. 水 流 一 定 方 向 应 该 是( ) (1) 从 高 处 向 低 处 流; (2) 从 压 强 大 处 向 压 强 小 处 流; (3) 从 流 速 大 的 地 方 向 流 速 小 的 地 方 流; (4) 从 单 位 重 量 流 体 机 械 能 高 的 地 方 向 低 的 地 方 流。 5. 理 想 流 体 流 经 管 道 突 然 放 大 断 面 时, 其 测 压 管 水 头 线( ) (1) 只 可 能 上 升; (2) 只 可 能 下 降; (3) 只 可 能 水 平; (4) 以 上 三 种 情 况 均 有 可 能。 6 在应用恒定总流的能量方程时,可选用图中的( ) 断 面, 作为计算断面。 (a )1,2,3,4,5 (b )1,3,5 (c )2,4 (d )2,3,4 1 122 3 3 4 4 5 5 7. 设有一恒定汇流,如图所示,Q Q Q 312=+, 根据总流伯努力方程式,则有( ) ()12221111 2 2222 2 3333 2 13 23 z p g V g z p g V g z p g V g h h w w + + ++ + =+ + ++--ραραρα ()()() 22211111 2 22222 2 ρραρραgQ z p g V g gQ z p g V g + + ++ + =++ + ++--ρραρρg Q Q z p g V g gQ h gQ h w w ()()123333 2 12213 23 (3) 上 述 两 式 均 不 成 立, 都 有 错 误;

相关文档
最新文档