实验一 常微分方程
数值分析常微分方程求解实验
实验报告
一、实验目的
解初值问题各种方法比较。
二、实验题目
给定初值问题
⎪⎩⎪⎨⎧=≤<+=,
0)1(,21y x xe x y dx dy x , 其精确解为)(e e x y x -=,按
(1)改进欧拉法,步长01.0,1.0==h h ;
(2)四阶标准龙格-库塔法,步长1.0=h ;
求在节点)10,...,2,1(1.01=+=k k x k 处的数值解及其误差,比较各个方法的优缺点。
三、实验原理
改进欧拉法程序,四阶标准龙格-库塔法程序。
四、实验内容及结果
五、实验结果分析
实验2中改进欧拉法和四阶标准龙格-库塔法的比较:
结果的第一个ans是x的值与对应的y的值,第二个ans是精确解x的对应值y,第三个ans 是与精确值的误差百分数。
通过误差百分数的比较,可以明显的发现改进欧拉法比四阶标准龙格-库塔法更精确。
2016春作业实验(1)常微分方程
1. 分别用Euler 法和ode45解下列常微分方程并与解析解比较:(1) ,(0)1, 03y x y y x '=+=<<function [ t,y ] = euler(f,ts,y0,h)t=ts(1):h:ts(2);y(1)=y0;for i=1:length(t)-1y(i+1)=y(i)+h*f(t(i),y(i));endt=t';y=y';endf=(t,y)t+y;[t1,y1]=euler(f,[0,3],1,0.05);[t2,y2]=ode45(f,[0,3],1);plot(t1,y1,'.-',t2,y2,'ro')hold ony3=dsolve('Dy=x+y','y(0)=1','x')ezplot(y3,[0,3])hold offlegend('euler','ode45','解析解');(2)22()5()3()45,(0)2,(0)1, 02tx t x t x t e x x t ''''--===<<f=(t,x)[2*x(2);5*x(2)+3*x(1)+45*exp(2*t)];[t1,y1]=ode45(f,[0,2],[2,1]);plot(t1,y1)2. 求一通过原点的曲线,它在(,)x y 处的切线斜率等于22,0 1.57.x y x +<<若x 上限增为1.58,1.60会发生什么?function dy = odefun_2(x,y)dy=2*x+y^2;dy=dy(:);end[t1,y]=ode45('odefun_2',[0,1.58],0) plot(t1,y);[t2,y]=ode45('odefun_2',[0,1.60],0) plot(t2,y);3. 求解刚性方程组:112121221000.25999.750.5,(0)1,050.999.751000.250.5,(0)1,y y y y x y y y y '=-++=⎧<<⎨'=-+=-⎩function Dy=fun(t,y)Dy=zeros(2,1);Dy(1)=-1000.25*y(1)+999.75*y(2)+0.5;Dy(2)=999.75*y(1)-1000.25*y(2)+0.5;[t,y]=ode15s('fun',[0,5],[1,-1]); plot(t,y(:,1),'o',t,y(:,2),'k-','LineWidth',2);4. (广告效应) 某公司生产一种耐用消费品,市场占有率为5%时开始做广告,一段时间的市场跟踪调查后,该公司发现:单位时间购买人口百分比的相对增长率与当时还没有买的百分比成正比,且估得此比例系数为0.5。
常微分方程实验报告
常微分方程实验报告一、实验目的常微分方程是数学分析和实际应用中非常重要的一部分,本次实验的主要目的是通过实际操作和计算,深入理解常微分方程的概念、性质和求解方法,并能够将其应用到实际问题中,提高我们解决数学问题和实际应用问题的能力。
二、实验原理常微分方程是指含有一个自变量和一个未知函数及其导数的等式。
求解常微分方程的方法有很多,常见的有变量分离法、一阶线性方程的求解方法(如常数变易法)、恰当方程的求解方法(通过积分因子)等。
对于一阶常微分方程,形如\(y' + p(x)y = q(x)\)的方程,可以使用积分因子\(e^{\int p(x)dx}\)来求解。
对于可分离变量的方程,形如\(g(y)dy = f(x)dx\),可以通过分别积分求解。
三、实验内容(一)一阶常微分方程的求解1、求解方程\(y' + 2xy = 2x\)首先,计算积分因子\(e^{\int 2xdx} = e^{x^2}\),然后将方程两边乘以积分因子得到:\((ye^{x^2})'= 2xe^{x^2}\)两边积分可得\(ye^{x^2} = e^{x^2} + C\),解得\(y =1 + Ce^{x^2}\)2、求解方程\(xy' y = x^2\)将方程化为\(y' \frac{y}{x} = x\),这里\(p(x) =\frac{1}{x}\),积分因子为\(e^{\int \frac{1}{x}dx} =\frac{1}{x}\)。
方程两边乘以积分因子得到\((\frac{y}{x})'= 1\),积分可得\(\frac{y}{x} = x + C\),即\(y = x^2 + Cx\)(二)二阶常微分方程的求解1、求解方程\(y'' 2y' + y = 0\)特征方程为\(r^2 2r + 1 = 0\),解得\(r = 1\)(二重根),所以通解为\(y =(C_1 + C_2x)e^x\)2、求解方程\(y''+ 4y = 0\)特征方程为\(r^2 + 4 = 0\),解得\(r =\pm 2i\),所以通解为\(y = C_1\cos(2x) + C_2\sin(2x)\)(三)应用常微分方程解决实际问题1、考虑一个物体在受到与速度成正比的阻力作用下的运动,其运动方程为\(m\frac{dv}{dt} = kv\)(其中\(m\)为物体质量,\(k\)为阻力系数),求解速度\(v\)随时间\(t\)的变化。
常微分方程的概念与初值条件
常微分方程的概念与初值条件常微分方程(Ordinary Differential Equation, ODE)是数学建模的一个基本问题。
它描述了一个未知函数在一个自变量上的变化规律。
大部分自然现象都可以通过常微分方程来描述,例如力学、电磁学、光学等。
因此,研究常微分方程不仅在数学上具有重要意义,而且在各个领域的应用中也具有不可替代的地位。
常微分方程的概念十分简单,它是一种描述未知函数变化规律的数学模型。
设未知函数为$y=y(x)$,自变量为$x$,则常微分方程的一般形式为:$$F(x, y, y', y'',\ldots, y^{(n)})=0,$$其中$y', y'',\ldots, y^{(n)}$分别表示未知函数的一阶导数、二阶导数、$\ldots$、$n$阶导数,$F$是一个给定的函数。
通常,我们将$n$称为该方程的阶数。
例如,一阶线性常微分方程的一般形式为:$$y' + p(x)y = q(x),$$其中$p,q$是已知函数。
在解决常微分方程的问题中,初值条件是必不可少的。
初值条件指定了未知函数$y$在某一点$x_0$上的函数值$y(x_0)$,以及该点上的某些导数值。
通常我们将这些导数值称为初始斜率,它们是求解常微分方程的关键。
将初始斜率确定下来,我们就可以唯一地确定未知函数$y(x)$在整个定义域上的函数值,这就是常微分方程解的存在唯一性定理。
因此,初值条件可以直接影响常微分方程的解。
在实际应用中,我们常常需要通过实验或测量来确定初始条件,例如弹簧的振动实验、电路的响应测试等等。
当我们知道了一个系统的初值条件时,就可以用常微分方程来描述该系统的动力学行为,并预测未来的变化趋势。
因此,研究常微分方程的初值问题对于控制、优化和决策等方面有着广泛的应用价值。
解决常微分方程的方法可以分为解析求解和数值求解两种。
解析方法是指通过数学分析得到常微分方程的解析解,即用解析式来表示未知函数$y(x)$。
数值计算实验报告-欧拉法常微分方程
数学与计算科学学院实验报告实验项目名称欧拉法解常微分方程所属课程名称数值计算实验类型验证型实验日期2012-6- 4班级隧道1002班学号201008020233姓名李彬彬成绩一、实验概述:【实验目的】 通过运用相关的数值计算软件,解决最基本的常微分方程的数值计算,并且能够熟练的运用这种方法。
【实验原理】 欧拉法1.对常微分方程初始问题(9.2))((9.1)),(00⎪⎩⎪⎨⎧==y x y y x f dxdy用数值方法求解时,我们总是认为(9.1)、(9.2)的解存在且唯一。
欧拉法是解初值问题的最简单的数值方法。
从(9.2)式由于y (x 0) = y 0已给定,因而可以算出),()('000y x f x y =设x 1 = h 充分小,则近似地有:),()(')()(00001y x f x y hx y x y =≈-(9.3)记 ,n ,,i x y y i i 10 )(== 从而我们可以取),(0001y x hf y y ==作为y (x 1)的近似值。
利用y 1及f (x 1, y 1)又可以算出y (x 2)的近似值:),(1112y x hf y y +=一般地,在任意点x n +1 = (n + 1)h 处y (x )的近似值由下式给出),(1n n n n y x hf y y +=+(9.4)这就是欧拉法的计算公式,h 称为步长。
不难看出,近似解的误差首先是由差商近似代替微商(见(9.3))引起的,这种近似代替所产生的误差称为截断误差。
还有一种误差称为舍入误差,这种误差是由于利用(9.4)进行计算时数值舍入引起的。
【实验环境】Windows XP 环境下运行 NumericalAnalyse 软件二、实验内容:【实验方案】在区间[0,1]上以h=0.1为步长,分别用欧拉法与预估-校正法求初值问题y’=y-2x/y且 y|x=0 =1的数值解。
将上述方程输入到软件NumericalAnalyse中步骤如图选择常微分方程的数值解法。
数学实验——常微分方程数值解
实验4 常微分方程数值解分1 黄浩 43一、实验目的1.掌握用MATLAB软件求微分方程初值问题数值解的方法;2.通过实例学习用微分方程模型解决简化的实际问题;3.了解欧拉方法和龙格-库塔方法的基本思想和计算公式,及稳定性等概念。
二、实验内容1.《数学实验》第一版(问题2)问题叙述:小型火箭初始重量为1400kg,其中包括1080kg燃料。
火箭竖直向上发射时燃料燃烧率为18kg/s,由此产生32000N的推力,火箭引擎在燃烧用尽时关闭。
设火箭上升时空气阻力正比于速度的平方,比例系数为0.4kg/m,求引擎关闭瞬间火箭的高度、速度、加速度,及火箭到达最高点时的高度和加速度,并画出高度、速度、加速度随时间变化的图形。
模型转换及实验过程:(一)从发射到引擎关闭设火箭总质量为m,上升高度为h,瞬时速度为v,瞬时加速度为a,由燃料燃烧时间t=60s,可列如下的方程组:t∈[0,60]−9.8因此,上述方程为二元常微分方程组,选择t为自变量,h和v为因变量进行分析。
初值条件:h(0)=0 ,v(0)=0对上述模型,使用ode45()函数求数值解(程序见四.1、四.2),结果如下:由上表可知,引擎关闭瞬间,火箭的高度为12189.78m,速度为267.26m/s,加速度为0.9170m/s2,火箭至此已飞行60s而高度、速度、加速度随时间的变化曲线如下:(二)从引擎关闭到最高点设引擎关闭时,t1=0,由上一问的结果可知,h(0)=12189.78m,v(0)= 267.26m/s,m=320kg,则可列二元常微分方程组如下:9.8因此,可选择t1为自变量,h、v为因变量进行分析(程序见四.3、四.4),实验结果如下:由上表可知,当t1∈[11,12]时,v(t1)有零点,即该区间内某时刻火箭达到最高点。
再进行更细致的实验(程序略),设步长为0.01,观察该区间内v(t1)的零点,如下表所示:可以看出,当t1=11.30s,即总时间t=71.30s时,火箭达到最高点,高度为13115.36m,加速度为-9.8m/s2。
常微分方程的求解与定性分析实验报告课件.doc
常微分方程的求解与定性分析实验报告一、实验综述1、实验目的及要求归纳和学习求解常微分方程(组)的基本原理和方法;掌握解析、数值解法,并学会用图形观察解的形态和进行解的定性分析;熟悉MATLAB 软件关于微分方程求解的各种命令;通过范例学习建立微分方程方面的数学模型以及求解全过程;通过该实验的学习,使学生掌握微分方程(组)求解方法(解析法、欧拉法、梯度法、改进欧拉法等),对常微分方程的数值解法有一个初步了解,同时学会使用MATLAB 软件求解微分方程的基本命令,学会建立微分方程方面的数学模型。
这对于学生深入理解微分、积分的数学概念,掌握数学的分析思维方法,熟悉处理大量的工程计算问题的方法是十分必要的。
2、实验仪器、设备或软件电脑、matlab7.0二、实验过程(实验步骤、记录、数据、分析)实验内容:根据实验内容和步骤,完成以下实验,要求写出实验报告(实验目的→问题→数学模型→算法与编程→计算结果→分析、检验和结论)1.求微分方程的解析解,并画出它们的图形。
y '= y + 2 x, y (0) = 1, 0< x <1;m=dsolve( 'Dy=y+2*x' , 'y(0)=1' , 'x' )ezplot(m,[0 1])m =3*exp(x) - 2*x –23 exp(x) - 2 x - 243.532.521.510 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1x3u u 0. 1u 0u( 0) 0; u( 0) 01.求微分方程的数值解,要求编写求解程序。
t [010]function dy=vdp1000(t,y)dy=zeros(2,1);dy(1)=y(2);dy(2)=-y(1)+0.1*y(1)^3;[T,Y]=ode15s( 'vdp1000' ,[0 10],[0 0]);plot(T,Y(:,1), '-' )10.80.60.40.2-0.2-0.4-0.6-0.8-10 1 2 3 4 5 6 7 8 9 103.Rossler微分方程组:'x y z'y x ay'z b z( x c)讨论化(如a∈(0,0.65))而方随参数a由小到大变,试当固定参数b=2,c=4时图形,观是否形成混沌状?察空间曲线程解的变化情况,并且画出空间曲线function r=rossler(t,x)global a;global b;global c;r=[-x(2)-x(3);x(1)+a*x(2);b+x(3)*(x(1)-c)];global a;global b;global c;b=2;c=4;t0=[0,200];for a=0:0.1:0.6[t,x]=ode45( 'rossler' ,t0,[0,0,0]);subplot(1,2,1);plot(t,x(:,1), 'r' ,t,x(:,2), 'g' ,t,x(:,3), 'b' );title( 'x(红色),y(绿色),z(蓝色) 随t 的变化情况' );xlabel( 't' );subplot(1,2,2);plot3(x(:,1),x(:,2),x(:,3))title( ' 相图' );xlabel( 'x' );ylabel( 'y' );zlabel( 'z' ); pauseend结果显示:a=0:x(红色),y(绿色),z(篮色)随t变化情况0.6相图0.40.80.20.60.4z-0.20.2-0.4-0.6 00.5-0.8-0.5-10 50 100 150 200t y -1 x-0.5a=0.1:x(红色),y( 绿色),z( 篮色)随t变化情况相图0.60.40.80.20.60.4z-0.20.2-0.4-0.61 -0.8-0.5 0.50 -1y0 50 100 150 200 -1 -0.5xta=0.2:x(红色),y(绿色),z( 篮色)随t变化情况相图43221.511z0.5-1-225-3-2-40 50 100 150 200t y -4 x-5a=0.3:x(红色),y(绿色),z( 篮色)随t变化情况相图54332.52211.5z1-10.5-25-35 -4 0-50 50 100 150 200t y -5 x-5a=0.4:x(红色),y(绿色),z( 篮色)随t变化情况6相图4 65243z21-25-4 105-60 50 100 150 200t y-5-10x-5a=0.5:x(红色),y(绿色),z( 篮色)随t变化情况相图1210128106846z242-25-410-6 5-80 50 100 150 200t y-5-10x-5结果分析:从图像可以看出,当a=0 时,微分方程的解(x,y,z)收敛与(0,0.5,0.5);当a=0.1时,(x,y,z)仍收敛与(0,0.5,0.5),只是收敛速度减慢;当a=0.2时,(x,y,z)已发散,周期性变化;随着a的增大,(x,y,z)接近其极限环的速度加快,空间曲线成混沌状。
《常微分方程指导与实验》第2章:一阶微分方程的解的存在定理
第二章 一阶微分方程的解的存在定理§2.1 一阶微分方程解的基本理论主要内容一 导数已解出方程初值问题解的存在唯一性定理 考虑导数已解出的一阶DE 的初值问题()()00,y f x y y x y '=⎧⎪⎨=⎪⎩(2.1)(2.2)这里()y x f ,是在闭矩形域R : a x x ≤-0,b y y ≤-0上的连续函数。
定义2.1 如果存在常数0>L ,使得对于所有的点()1,y x ,()2,y x R ∈,都有不等式()()2121,,y y L y x f y x f -≤-成立,则称函数()y x f ,在R 上关于y 满足李普希兹(Lipschitz )条件。
1定理2.1 (毕卡存在唯一性定理) 如果()y x f ,在R 上满足条件: 1)连续;2)关于y 满足李普希兹条件,则初值(2.1)和(2.2)在区间h x x ≤-0上存在唯一解()x y y =,其中()M b a h ,m in=,()y x f M R y x ,max ),(∈=。
注1 取数h 的意义。
注意到()y x f M R y x ,max ),(∈=,从而积分曲线()x y y =在任一点()()R x y x ∈,处的切线斜率()M x y ≤'。
于是从点()o y x p ,0引两条斜率分别为M 和M -的直线1l 和2l ,便知过点P 的积分曲线必限制在图2.1和图2.2的阴影区域内。
而直线1l 和2l 相交情形有如下两种可能。
(i )若相交成如图 2.1所示的情况,则a Mb>,积分曲线()x y y =在a x x ≤-0上不越出R ,从而应取a h =。
(ii )若相交成如图 2.2所示的情况,则a Mb >,积分曲线()x y y =在Mb x x ≤-0上不越出R ,从而应取Mb h =。
总之,取()M ba h ,min =,就是为了使初值问题(2.1)和(2.2)的解在h x x ≤-0上总存在。
常微分方程
●凡表示未知函数、未知函数的导数与自变量之间的关系的方程,叫做微分方程。
●未知函数是一元函数的,叫做常微分方程。
未知函数是多元函数的,叫做偏微分方程。
●微分方程中所出现的未知函数的最高阶导数的阶数,叫做微分方程的阶。
●在研究某些实际问题时,首先要建立微分方程,然后找出满足微分方程的函数(解微分方程),就是说,找出这样的函数,把这函数代入微分方程能使该方程成为恒等式。
这个函数就叫该微分方程的解。
●如果微分方程的解中含有任意常数,且任意常数的个数与微分方程的阶数相同,这样的解叫做微分方程的通解。
●设微分方程中的未知函数为y=y(x),如果微分方程是一阶的,通常用来确定任意常数的条件是x=x0时,y=y0,或写成y|x=x0=y0,其中x0,y0都是给定的值;如果微分方程是二阶的,通常用来确定任意常数的条件是:x=x0时,y=y0,y′=y0′,或写成y|x=x0=y0,y′|x=x0=y0′,其中x0,y0和y0′都是给定的值,上述这种条件叫做初始条件。
●确定了通解中的任意常数以后,就得到微分方程的特解。
●求微分方程y′=f(x,y)满足初始条件y|x=x0=y0的特解这样一个问题,叫做一阶微分方程的初值问题,记作{y′=f(x,y);y|x=x0=y0}●可分离变量方程。
一阶可分离变量方程:dy/dx=f(x)/g(y),可分离变量为:∫g(y)dy=∫f(x)dx,设g(y)、f(x)的原函数分别为G(y)、F(x),则可解出方程的通解:G(y)=F(x)+C。
●例:求微分方程dy/dx=2xy的通解。
解:方程是可分离变量的,分离变量后得dy/y=2xdx,两端积分∫dy/y=∫2xdx,得ln|y|=x2+C1,从而y=±e C1e x2。
因±e C1仍是任意常数,把它记作C,便得方程的通解y=Ce x2。
●齐次方程。
如果一阶微分方程dy/dx=f(x,y)中的函数f(x,y)可写成y/x的函数,即f(x,y)=φ(y/x),则称这方程为齐次方程。
数学实验基础 实验报告(1)常微分方程
实验一 常微分方程1. 分别用Euler 法和ode45解下列常微分方程并与解析解比较: (1) ,(0)1,13y x y y x '=+=<<Euler 法:function [t,y]=euler(Fun,tspan,y0,h) t=tspan(1):h:tspan(2); y(1)=y0;for i=1:length(t)-1y(i+1)=y(i)+h.*feval(Fun,t(i),y(i)); end t=t'; y=y';function f=Fun(x,y) % 常微分方程的右端函数 f=x+y;>> [x,y]=euler('Fun',[0,3],1,0.1)>> [x,y] ans =0 1.0000 0.1000 1.1000 0.2000 1.2200 0.3000 1.3620 0.4000 1.5282 0.5000 1.7210 0.6000 1.9431 0.7000 2.1974 0.8000 2.4872 0.9000 2.8159 1.0000 3.1875 1.1000 3.6062 1.2000 4.0769 1.3000 4.6045 1.4000 5.1950 1.5000 5.8545 1.6000 6.5899 1.7000 7.4089 1.8000 8.3198 1.9000 9.3318 2.0000 10.4550 2.1000 11.7005 2.2000 13.0805 2.3000 14.6086 2.4000 16.2995 2.5000 18.1694 2.6000 20.2364 2.7000 22.5200 2.8000 25.0420 2.9000 27.8262 3.0000 30.8988ode45:>> [x,y]=ode45('Fun',[0,3],1) ans =0 1.0000 0.0502 1.0528 0.1005 1.1109 0.1507 1.17460.2010 1.2442 0.2760 1.3596 0.3510 1.4899 0.4260 1.63610.5010 1.7996 0.5760 1.9817 0.6510 2.1838 0.7260 2.4074实验一 常微分方程0.8010 2.6544 0.8760 2.9264 0.9510 3.2254 1.0260 3.55351.1010 3.9131 1.1760 4.3065 1.2510 4.7364 1.3260 5.20561.4010 5.7172 1.4760 6.2744 1.5510 6.8810 1.6260 7.54061.7010 8.2574 1.7760 9.0359 1.8510 9.8808 1.9260 10.79742.0010 11.7912 2.0760 12.8683 2.1510 14.0351 2.2260 15.29862.3010 16.6664 2.3760 18.1466 2.4510 19.7478 2.5260 21.47962.6010 23.3522 2.6760 25.3764 2.7510 27.5641 2.8260 29.92812.9010 32.4820 2.9257 33.3694 2.9505 34.2796 2.9752 35.21343.0000 36.1711解析解:>> y=dsolve('Dy=x+y','y(0)=1','x') y =2*exp(x) - x - 1(2) 20.01()2sin(),(0)0,(0)1,05y y y t y y t ''''-+===<< Euler 法:实验一常微分方程function f=Fun(t,y)% 常微分方程的右端函数f=[y(2);0.01*y(2)^2-2*y(1)+sin(t)];>> [t,y]=euler('Fun',[0,5],[0,1],0.2)ode45:>> [t,y]=ode45('Fun',[0,5],[0,1])t =0 0.0001 0.0001 0.0002 0.0002 0.0005 0.0007 0.0010 0.0012 0.00250.0037 0.0050 0.0062 0.0125 0.0188 0.0251 0.0313 0.0627 0.0941 0.12550.1569 0.2819 0.4069 0.5319 0.6569 0.7819 0.9069 1.0319 1.1569 1.28191.4069 1.5319 1.6569 1.7819 1.90692.0319 2.1569 2.2819 2.4069 2.53192.6569 2.7819 2.90693.0319 3.1569 3.2819 3.4069 3.5319 3.6569 3.78193.90694.0319 4.1569 4.2819 4.4069 4.5319 4.6569 4.7427 4.8285 4.91425.0000y =0 1.0000 0.0001 1.0000 0.0001 1.0000 0.0002 1.0000 0.0002 1.00000.0005 1.0000 0.0007 1.0000 0.0010 1.0000 0.0012 1.0000 0.0025 1.00000.0037 1.0000 0.0050 1.0000 0.0062 1.0000 0.0125 1.0000 0.0188 1.00000.0251 0.9999 0.0313 0.9998 0.0627 0.9987 0.0941 0.9965 0.1253 0.99340.1564 0.9893 0.2786 0.9632 0.3966 0.9220 0.5085 0.8662 0.6126 0.79670.7072 0.7146 0.7908 0.6210 0.8620 0.5176 0.9198 0.4058 0.9632 0.28760.9915 0.1647 1.0043 0.0392 1.0013 -0.0869 0.9826 -0.2117 0.9485 -0.33310.8996 -0.4490 0.8365 -0.5578 0.7605 -0.6577 0.6725 -0.7471 0.5742 -0.8246实验一 常微分方程0.4669 -0.8889 0.3525 -0.9393 0.2327 -0.9748 0.1095 -0.9950 -0.0154 -0.9996-0.1398 -0.9887 -0.2619 -0.9624 -0.3798 -0.9212 -0.4916 -0.8657 -0.5957 -0.7970-0.6904 -0.7161 -0.7742 -0.6242 -0.8460 -0.5228 -0.9046 -0.4134 -0.9491 -0.2978-0.9789 -0.1777 -0.9934 -0.0549 -0.9945 0.0300 -0.9883 0.1146 -0.9748 0.1985-0.9543 0.28092. 求一通过原点的曲线,它在(,)x y 处的切线斜率等于22,0 1.57.x y x +<<若x 上限增为1.58,1.60会发生什么?function f=Fun(x,y) % 常微分方程的右端函数 f=2*x+y.^2;>> [x,y]=ode45('Fun',[0,1.57],0) x =0 0.0393 0.0785 0.1178 0.1570 0.1963 0.2355 0.2748 0.3140 0.3533 0.3925 0.4318 0.4710 0.5103 0.5495 0.5888 0.6280 0.6673 0.7065 0.7458 0.7850 0.8243 0.8635 0.9028 0.9420 0.9813 1.0205 1.0598 1.0990 1.1383 1.1775 1.2168 1.2560 1.2953 1.3345 1.3738 1.4130 1.4248 1.4367 1.4485 1.4604 1.4722 1.4840 1.4959 1.5077 1.5140 1.5203 1.5265 1.5328 1.5376 1.5424 1.5472 1.5519 1.5543 1.5567 1.5591 1.5614 1.5631 1.5647 1.5664 1.5681 1.5685 1.5690 1.5695 1.5700 y =实验一 常微分方程0 0.0015 0.0062 0.0139 0.0247 0.0386 0.0556 0.0758 0.09920.1259 0.1559 0.1895 0.2266 0.2675 0.3124 0.3615 0.4152 0.4738 0.5378 0.6076 0.6841 0.7679 0.8601 0.9620 1.0751 1.2014 1.3434 1.5045 1.6892 1.9037 2.1557 2.4577 2.8282 3.3003 3.9056 4.7317 5.9549 6.4431 7.0116 7.6832 8.4902 9.4821 10.7170 12.3090 14.4551 15.9220 17.7080 19.9390 22.8164 25.6450 29.2282 33.9673 40.5910 44.9434 50.3088 57.1229 66.1087 74.3108 84.7123 98.4901 117.7875 124.9206 132.9699 142.1268 152.641500.20.40.60.81 1.2 1.4 1.6若x 上限增为1.58,1.60,则超出运算的范围,发生溢出。
常微分方程实验报告_闵安游_U201310085
⑵、3 步 4 阶 Adams 方法: ①、画图程序:
w = exp(i*linspace(0,2*pi)); z = 24*(w.^3 - w.^2)./(1 - 5*w + 19*w.^2 + 9*w.^3); plot(z)
②、绝对稳定性判断程序:
9 3 19 2 5 h (1 h ) (1 h) h 0 此时的特征方程是: 24 24 24 24
fprintf('误差为: %10.10f\n',m)
测试程序:
f = @(t,y)-2*y + 2*t + 1; F = @(t)exp(-2*t) + t; Threesteps_Adams1(0,2,f,F,0.1)
第 2 题: ⑴、2 步 3 阶 Adams 方法: ①、画图程序:
w = exp(i*linspace(0,2*pi)); z = 12*w.*(w - 1)./(5*w.*w + 8*w -1); plot(z)
t = a + i*h; y2 = y1 + h/2*(3*f1 - f0); plot(t,y2,'.'); f2 = f(t,y2); y1 = y2; f0 = f1; f1 = f2; end y = linspace(a,b); plot(y,F(y),'r'); hold off m = abs(y1-F(b)); fprintf('数值解: %10.10f\n',y1) fprintf('误差为:%10.10f\n',m)
' u2 24u1 51u2 9 cos t 1/ 3sin t
实验报告1
实验报告1 函数与极限院系 班号姓名学号成绩一、实验内容函数图形的显示,极限的运算,最值的计算.二、预期目标1.熟悉Mathematica 软件的基本操作.2.掌握函数与极限的有关操作命令.3.学会利用Mathematica 软件对函数进行分析研究.三、常用命令1. 作图命令: 2. 参数作图命令: 3. 图形显示命令: 4. 求极限命令: 5. 求极值名命令:四、练习内容1.画出下列函数的图形: (1) y=cos3x作图命令:(2) f (x )=x 5+3e x+log 3(3-x ) x ∈[-2,2]作图命令:(3)⎪⎩⎪⎨⎧=+=ty t t x 2sin作图命令:(4)⎪⎩⎪⎨⎧==ty t x 33sin cos t ∈[0,2π]作图命令:2.求下列极限:(1)110002lim+∞>-n nn (2)113)2(3)2(lim ++∞>-+-+-n n n n n (3)35)3)(2)(1(limnn n n n +++∞>- (4)3522lim -+>-x x x (5)131lim +->-x x x(6)x e xx arctan lim -+∞>-(7)156182221lim +-->-x x x x (8))sin 11sin (lim x x x x x -∞>-计算结果:(1)(2)(3)(4)(5)(6)(7)(8)3.讨论函数f(x)=2x3-6x2-18x+7在点2.8附近的极值.命令:结果:五、思考与提高1.怎样对隐函数的图形进行显示?2.怎样利用软件对函数极限存在性进行判断?3.如何利用软件对函数的连续性进行判断?4.如何求函数的最大(小)值?实验报告2 微分及其应用院系 班号姓名学号成绩一、实验内容导数的运算法则,复合函数求导法及参数方程求导法等.二、预期目标1.进一步理解导数及其几何应用.2.学习Mathematica 的求导命令与求导法.三、常用命令1.求导命令: 2.求微分命令: 3.隐函数求导命令: 4. 参数方程所确定的函数求导命令:四、练习内容1.求下列函数的导数: (1)x y 2ln 1+=求导命令: 求导结果: (2)21121xx y +++=求导命令:求导结果:(3) y=cos 2(cos2x )求导命令: 求导结果:(4)y=2x/lnx求导命令: 求导结果: (5)y=ln[ln(lnt)]求导命令:求导结果: (6)xxy arccos arcsin =求导命令: 求导结果:(7)y=e arcsinx +arctane x求导命令: 求导结果:(8)xey 1sin 2-=求导命令:求导结果:2.求下列函数的二阶导数:(1) y=tanx 计算结果:y ” =(2)y=(1+x 2)arctanx 计算结果:y ” =(3)y=xtanx -cscx 计算结果:y ” =(4)y=21ln (x -1)-21ln (x+1) 计算结果:y ” = (5)⎪⎩⎪⎨⎧-==21arcsin ty t x 计算结果:y ” =(6)⎪⎩⎪⎨⎧==tb y t a x sin cos计算结果:y ” =3.求下列方程所确定的隐函数y=y (x )的导数xyd d : (1) sin (xy )+cosy=0 计算结果:xyd d =(2)arctan x y =ln 22y x + 计算结果:xyd d =(3)x y =y x计算结果:xy d d =4.验证参数方程⎪⎩⎪⎨⎧==te y t e x tt cos sin 所确定的函数y 满足关系:)d d (2)(d d 222y x y x y x xy -=+ 程序:五、思考与提高1. 如何利用函数的导数判定函数的单调性、凹凸性?2.如何求由方程F(x,y,z)=0确定的隐函数z=z(x,y)的偏导数?实验报告3 积分及其应用院系 班号姓名学号成绩一、实验内容一元函数的不定积分与定积分二、预期目标1.加深理解积分理论中分割、近似、求和、取极限的思想方法.2.学习求积分的命令Integrate 与NIntegrate .3.熟悉Mathematica 软件在积分运算的重要作用.三、常用命令1.求和命令: 2.求不定积分命令: 3.求定积分命令:四、练习内容1.求下列函数的一个原函数:(1)41x (2)212x +积分命令: 积分命令: 积分结果: 积分结果:(3))1()1(22x x x ++ (4)4211xx -+ 积分命令: 积分命令: 积分结果:积分结果:(5)x x 22sin cos 2cos (6)xxe e +1积分命令: 积分命令: 积分结果: 积分结果:(7))tan 1(cos 12x x + (8)x e x 32 积分命令: 积分命令: 积分结果:积分结果:(9)x cos 1+ (10))34cos()23sin(+⋅+x x 积分命令: 积分命令: 积分结果: 积分结果:2.计算下列定积分:(1)⎰2/6/2d cos ππx x (2)⎰+-4/02sin 12sin 1πxxdx计算结果: 计算结果:(3)⎰-2/0d cos 351πx x(4)⎰30d cot arc x x x计算结果: 计算结果:(5)⎰---222d 11x x (6)t t e td cos 2/02⎰π计算结果:计算结果:(7)⎰+12/3d 1x xx (8)⎰π222d sin x x 计算结果: 计算结果:3.计算下列积分,并求其结果关于变量x 的导数:(1)⎰+02d 1x t t (2)⎰-xt t te 0d 2积分结果: 积分结果: 关于x 的导数:关于x 的导数:(3)⎰0sin 2d )cos(x t t (4)⎰+203d 11x t t 积分结果: 积分结果: 关于x 的导数: 关于x 的导数: 4.判定广义积分⎰∞++12)1(1x x dx 及⎰--2022)2(x exdx 的敛散性,收敛时计算出积分值. ⎰∞++12)1(1x x dx ⎰--2022)2(x exdx 程序: 程序: 结果: 结果: 5.求积分⎰-102)43(x x dx 具有6位、10位有效数的近似值. 命令: 五、思考与提高1. Mathematica 系统对分段函数的积分能否进行自动处理?2.《高等数学》中所学的积分换元法在软件系统里如何应用?3.怎样用Mathematica 中动画来演示定积分的定义?实验报告4 三角级数院系 班号 姓名 学号 成绩一、实验内容级数敛散性的判定.二、预期目标1.掌握级数的展开与求和命令.2. 学习使用Mathematica 进行级数敛散性的判定.三、常用命令1.求taylor 展式命令:四、练习内容1.求下列泰勒展开式,并在同一坐标系下画出函数图形及展开式图形. (1) ln (1+x ) 在x0=0点的8阶Taylor 展开. 程序:(2) P (x )=x 4-5x 3+x 2-3x+4 在x0=4点的4阶Taylor 展开. 程序:(3) f (x )=x1在x0=-1点的n 阶Taylor 展开. 程序:2.求下列级数的和函数:(1)∑∞=--112121n n x n (2)∑∞=+1)1(1n n x n n (-1≤x ≤1) 命令:命令:结果: 结果:(3)∑∞=-+112)1(n n x n n 命令: 结果:3.判定下列级数的敛散性:(1)∑∞=12n n n(2)∑∞=++13211n nn 结论:结论:(3)∑∞=1!2n nn (4)∑∞=1)(sin n n n n n π结论: 结论:(5)∑∞=+112tann n n π(6)∑∞=12)!(n nn n 结论: 结论:4.判定下列级数是否收敛,收敛时请指出是绝对收敛,还是条件收敛? (1)∑∞=---11121)1(n n n (2)∑∞=+-122)1(n n n 结论:结论:(3)∑∞=--1ln )1(n nn n (4)∑∞=12sin n n na (a 为常数) 结论: 结论:五、思考与提高用判别法可以判别级数的敛散性,但在实际应用时,往往要使用其和,原则上可用Sum 语句求和,但许多数项级数仅仅使Sum 语句求不出其和,而另-Mathematica 命令NSum 可与判别结果一起用来求出其近似值,问:是否对任一级数均可用NSum 来求其近似值?试以∑∞=-1)1(n n 为例观察.实验报告5 空间解析几何院系 班号 姓名 学号 成绩一、实验内容空间图形的显示,简单动画的制作.二、预期目标1.能正确显示空间图形.2.能用Mathematica 制作简单的动画.三、常用命令1.三维作图命令: 2.参数方程作图命令(三维曲线): (曲面): 3.动画命令:四、练习内容1.显示下列函数图形:(1) 椭球面⎪⎩⎪⎨⎧===v z v u y v u x cos 3sin sin 5sin cos 2,),0(),2,0(ππ∈∈v u作图命令:(2) 椭球抛物面⎪⎩⎪⎨⎧===23sin 3cos 3u z v u y v u x ,其中)2,0(),2,0(π∈∈v u作图命令:(3) 双曲抛物面⎪⎩⎪⎨⎧-===3/)(22v u z v y u x ,其中)4,4(,-∈v u作图命令:(4) 圆柱螺线⎪⎩⎪⎨⎧===t z t y t x 4sin 34cos 3,其中)5,0(∈t作图命令:3. 制作平面振动动画(利用函数y x y x f 3sin 2cos ),(=,其中x,y 均属于(-1,1)).程序:五、思考与提高用参量函数与直接函数显示图形有什么区别?比较谁更容易作出图形?实验报告6 多元微分学院系 班号 姓名 学号 成绩一、 实验题目隐函数的导数,函数的偏导数,函数的极值.二、 预期目标1.求隐函数的导数.2.求函数的偏导数和全微分.3.用微分知识求函数的极值.三、常用命令1.求偏导命令: 2.求全微分命令: 3.解方程(组)命令:四、练习内容1. 设xx xy tan =,求dxdy . 命令: 结果:2. 设),(y x f z =由方程02=+--z xye z e所确定,求xz ∂∂. 程序: 结果:3. 设0ln 2=--xyz xy xz 确定函数),(y x f z =,求z 的全微分. 程序: 结果:4.求下列函数的偏导数:(1)yz x z y y x y y x z ∂∂∂∂-=,sin cos sin cos 2323,求结果:(2)yzx z v u y v u x y x z ∂∂∂∂+=-==,2,22,求,其中结果:4. 求函数22y x z +=在平面x+y=1上的最小值.程序: 结果:五、思考与提高1. 隐函数的二阶(偏)导数应如何求?2.函数的方向导数怎样求?实验报告7 多元积分学院系 班号 姓名 学号 成绩一、 实验题目空间立体体积和表面积.二、 预期目标1.用Mathematica 软件计算重积分.2.能解决空间立体体积和表面积的计算.三、常用命令1.求二重积分命令:四、练习内容1.计算下列重积分:(1)⎰⎰1D dxdy x y,其中D 1是由y=2x ,y=x ,x=4,x=2所围成的区域 . 积分命令:计算结果:(2)⎰⎰+2)(22D dxdy y x,其中D 2是由y=x ,y=x+2,y=2,y=6所围成的区域.积分命令:计算结果:(3)⎰⎰++3)1ln(22D dxdy y x ,其中D 3:0,0,122≥≥≤+y x y x . 积分命令:计算结果:(4)⎰⎰⎰Ω++3)(z y x dxdydz,其中Ω:21≤≤x ,21≤≤y ,21≤≤z . 积分命令:计算结果:(5)⎰⎰⎰Ω++222zy x dxdydz ,其中Ω是由222z y x =+及1=z 所围成的区域. 积分命令:计算结果:2.求抛物面x y x y 2,==及平面z=0,z+x=6所围成的物体(密度为1)的质量.程序: 结果:五、思考与练习1.在实验步骤1中{x,0,1}与{y,2*x,x^2+1}能不能交换次序?为什么?2.在重积分中,如果可以用换元法,也可以用Integrate直接积分时,用哪一种方法好,为什么?3.曲线积分和曲面积分如何计算?实验报告8 常微分方程院系 班号 姓名 学号 成绩一、 实验题目常微分方程(组)的精确解.二、 预期目标1.求一阶常微分方程的精确解.2.求解简单的微分方程组和高阶方程.三、常用命令1.求解微分方程命令: 2.求解微分方程组命令: 3.求微分方程数值解命令:四、练习内容1. 求x y x y tan cos '2=+的通解.命令:结果:2. 求13232=-+y xx dx dy ,且满足y(1)=0的特解. 命令:结果:3. 求⎩⎨⎧=--=++03'5'y x y e y x x t ,满足⎩⎨⎧==0)0(1)0(y x 的特解.命令:结果:五、思考与提高如果遇见无法直接用DSolve 求解的常微分方程,如22112'x y y +=+,怎么办?院系 班号 姓名 学号 成绩一、实验内容矩阵的运算(加法、数乘、乘法、转置、逆)二、预期目标熟悉Mathematica 软件中关于矩阵运算的各种命令.三、常用命令1.矩阵显示命令: 2.求矩阵转置命令: 3.求逆矩阵命令: 4.求矩阵和差命令: 5.求矩阵数乘命令: 6.求矩阵乘命令:四、思考与练习已知矩阵 ⎥⎥⎥⎥⎥⎥⎥⎦⎤⎢⎢⎢⎢⎢⎢⎢⎣⎡------------=031948118763812654286174116470561091143A⎥⎥⎥⎥⎥⎥⎥⎦⎤⎢⎢⎢⎢⎢⎢⎢⎣⎡------=503642237253619129113281510551201187851697236421B求:(1) A'; (2)A-1;(3)A*B .(1)求A'的命令: (2)求A-1的命令:A'= A-1=(3)求A*B 的命令:A*B =(请用矩阵形式表示计算结果)院系 班号 姓名 学号 成绩一、实验内容对矩阵作各种变化,初等变换.二、预期目标1.复习并掌握矩阵初等变换的方法.2.掌握Mathematic 软件中关于矩阵初等变换的相关命令.三、常用命令1.取矩阵元素命令: 2.取矩阵的子矩阵命令: 3.求矩阵维数命令:四、练习内容1.已知矩阵;302 150311101⎥⎥⎥⎦⎤⎢⎢⎢⎣⎡---=A(1)求A 的行向量组a 1,a 2,a 3, 以及列向量组b 1,b 2,b 3,b 4程序:(2)求A 的一,三,五行,二,三,四列交叉点上的元素做出子矩阵.程序:结果: 2.判断下列向量组是否线性相关(1)⎥⎥⎥⎦⎤⎢⎢⎢⎣⎡-=1211a ⎥⎥⎥⎦⎤⎢⎢⎢⎣⎡-=1302a ⎥⎥⎥⎦⎤⎢⎢⎢⎣⎡-=3123a 程序:结论:(2)⎥⎥⎥⎦⎤⎢⎢⎢⎣⎡-=1121a ⎥⎥⎥⎦⎤⎢⎢⎢⎣⎡-=1112a ⎥⎥⎥⎦⎤⎢⎢⎢⎣⎡-=1353a 程序:结论:实验报告11 行列式运算院系 班号 姓名 学号 成绩一、实验内容行列式的计算.二、预期目标1. 复习矩阵的行列式的求法,矩阵初等变换方法.2. 熟悉Mathematic 软件中关于求一个矩阵的行列式的命令把矩阵进行初等变换的命令以及与其相关的其它命令.三、常用命令1.求矩阵行列式命令:四、练习内容 1.求行列式βααββααββα+++100001000(共10阶)的值计算结果:2.利用克莱姆法则求解下列线性方程组(1)⎪⎪⎪⎩⎪⎪⎪⎨⎧-=-+--=++++=-+-+=+-+--=-+-+3322224343238243214225432154321543215422153321x x x x x x x x x x x x x x x x x x x x x x x x x程序:结果:(2) ⎪⎪⎪⎩⎪⎪⎪⎨⎧=+=++=++=++=+150650650651655454343232121x x x x x x x x x x x x x结果:2.已知⎥⎥⎥⎦⎤⎢⎢⎢⎣⎡----=876174114A ,⎥⎥⎥⎦⎤⎢⎢⎢⎣⎡---=614475914B验证:|A×B|=|A|×|B|.程序:实验报告12 求解方程组院系 班号 姓名 学号 成绩一、实验题目求AX=B 的通解.二、实验目的通过本实验,使学生认识到虽然在《线性代数》中求AX=B 的通解比较繁,但在Mathematica 软件中却是比较简单的. 三、常用命令1.矩阵化简命令: 2.解线性方程组命令: 3.求AX=0的基础解系命令:四、练习内容1.求下列矩阵的秩:(1)⎥⎥⎥⎦⎤⎢⎢⎢⎣⎡-----=253414312311112A 命令: 结果: (2)⎥⎥⎥⎦⎤⎢⎢⎢⎣⎡-----=189513411314311B 命令: 结果:2.解下列线性方程组:(1) ⎥⎥⎥⎦⎤⎢⎢⎢⎣⎡----512111211121⎥⎥⎥⎦⎤⎢⎢⎢⎣⎡-=⎥⎥⎥⎥⎦⎤⎢⎢⎢⎢⎣⎡5514321x x x x 程序:结果:(2) ⎥⎥⎥⎥⎦⎤⎢⎢⎢⎢⎣⎡----1111145212142121⎥⎥⎥⎥⎦⎤⎢⎢⎢⎢⎣⎡=⎥⎥⎥⎥⎦⎤⎢⎢⎢⎢⎣⎡3/10324321x x x x结果:实验报告13 特征值、特征向量院系 班号 姓名 学号 成绩一、实验题目计算已知矩阵的特征值和属于每一个特征值的特征向量.二、实验目的1.复习线代中的特征值与特征向量的求法.2.比较Mathematic 软件与普通方法的异同之处.三、常用命令1.求矩阵特征值命令: 2.求矩阵特征向量命令:四、练习内容求出下列矩阵的全部特征值与特征向量:1.⎥⎦⎤⎢⎣⎡-=00a a A ; 程序:结果:2.⎥⎥⎥⎦⎤⎢⎢⎢⎣⎡=001010100B ; 结果:3. ⎥⎥⎥⎥⎦⎤⎢⎢⎢⎢⎣⎡-----=1111111111111111C . 结果:实验报告14 离散型随机变量及其相关知识院系班号姓名学号成绩一、实验内容排列、组合的计算,几种离散型随机变量的产生及其相关内容.二、预期目标1.熟练掌握Mathematical软件的基本操作.2.熟悉与排列、组合、离散型随机变量有关的操作命令.3.掌握利用Mathematical软件处理简单的概率问题.三、常用命令1.(双)阶乘运算命令:2.组合数的计算命令:3.排列数的计算命令:4.服从二项分布的随机变量的生成命令:5.服从泊松分布的随机变量命令:6.将离散型随机变量的分布律拟合为函数的命令:四、练习内容1.计算下列结果(1)15!(2)15!!命令:命令:结果:结果:2.计算下列排列组合式的结果(1)P510(2)C510(3)!6!4!2!12⨯⨯命令:命令:命令:结果:结果:结果:3.生成以n=20,p=0.3为参数服从二项分布的随机变量bdist,将其分布律图形显示.程序:4.生成以p=0.4为参数服从几何分布的随机变量bdist,将其分布律图形显示.程序:5.生成以p=0.2为参数服从泊松分布的随机变量bdist,将其分布律图形显示.程序:五、思考与提高1.试分析几种离散型随机变量分布律的最值情况?2.怎样求解离散型随机变量有关的事件概率?实验报告15 连续型随机变量及其相关知识院系班号姓名学号成绩一、实验内容连续型随机变量的产生及其相关内容.二、预期目标1.熟练掌握几种连续型随机变量产生的有关操作命令.2.掌握利用软件对连续型随机变量进行分析的方法.3.掌握利用软件处理简单的概率问题.三、常用命令1.服从均匀分布的随机变量的生成命令:2.服从正态分布的随机变量的生成命令:3.服从t分布的随机变量的生成命令:4.服从χ2分布的随机变量的生成命令:5.服从F分布的随机变量的生成命令:6.求连续型随机变量的概率密度函数的命令:7.求连续型随机变量的分布函数的命令:四、练习内容1.生成以μ=10.05和σ=0.06为参数服从正态分布的连续型随机变量gdist及其概率密度函数、分布函数并图形显示;试求概率P{9.9<gdist<10.17}.程序:2.生成以a=0,b=1为参数服从柯西分布的连续型随机变量gdist及其概率密度函数、分布函数并图形显示.程序:3.生成以n1=4,n2=8为自由度服从F分布的连续型随机变量gdist及其概率密度函数、分布函数并图形显示.程序:4.生成以α=1,β=3为服从威布尔分布的连续型随机变量gdist及其概率密度函数、分布函数并图形显示.程序:五、思考与提高怎样利用软件对随机变量函数的分布进行分析,以及有关事件概率的求解?实验报告16 数字特征院系班号姓名学号成绩一、实验内容随机变量的数字特征及其相关内容.二、预期目标1.熟练掌握随机变量数字特征的有关操作命令.2.掌握利用软件对随机变量的特征函数(母函数)的求解.3.掌握利用软件处理简单的概率问题.三、常用命令1.求随机变量的期望的命令:2.求随机变量的方差的命令:3.求随机变量的标准差的命令:4.求随机变量的函数的方差的命令:5.求数据的协方差的命令:6.求数据的协方差矩阵的命令:7.求两随机变量的相关系数的命令:8.求两数据的相关系数矩阵的命令:四、练习内容1.(1)求以λ为参数服从泊松分布的随机变量的数学期望和方差.(2)求上述随机变量函数(f(x)=x2)的数学期望.(3)求服从参数λ=0.1的指数分布的随机变量的特征函数.程序:结果:2.(1)若样本data={16.5,13.8, 16.6, 15.7, 16.0, 16.4, 15.3},求样本均值、调和均值和中位数.结果:(2)若二维总体的样本data={{1612, 7627}, {1598, 6954},{1804, 8365},{1752, 9469}, {2067, 6410}, {2365, 10327},{1646, 7320}, {1579, 8196}, {1880, 9709}, {1773, 10370},{1712, 7749}, {1932, 6818}, {1820, 9307}, {1900, 6457},{1587, 8309}, {2208, 9559}, {1487, 6255}},求样本均值向量、中位数向量、方差向量和协方差矩阵.程序:结果:实验报告17 估计理论院系班号姓名学号成绩一、实验内容单个和两个总体均值、方差的估计.二、预期目标1.熟练掌握估计理论的相关操作命令.2.熟练掌握利用Mathematical软件对总体均值、方差进行估计.3.掌握利用Mathematical软件处理估计理论相关的实际问题.三、常用命令1.求总体均值的无偏估计的命令:2.求总体方差的无偏估计的命令:3.求总体方差的极大似然估计的命令:4.求单个总体均值的区间估计的命令:5.求两个总体均值之差的区间估计的命令:6.求单个总体方差的区间估计的命令:7.求两个总体方差之比的区间估计的命令:四、练习内容1.若样本data1={4506,4508,4499,4503,4504,4510,4497,4512,4514, 4505,4493,4496,4506,4502,4509,4496}来自正态总体,方差未知(置信度为0.95):求出总体均值、方差的置信区间.程序:结果:2.若样本data2={4507,4507,4497,4506,4503,4511,4498,4510,4514,4510,4493,4491,4507,4501,4510,4495}来自正态总体,设置信度为0.95:(1)若data1与data2的总体方差都未知,均值之差的置信区间;程序:结果:(2)若data1与data2的总体方差都为40,均值之差的置信区间.程序:结果:3.data1与data2的总体方差之比值的置信区间(置信度为0.95).程序:结果:实验报告18 假设检验院系班号姓名学号成绩一、实验内容对单个和两个总体均值、方差的假设检验.二、预期目标1.熟练掌握假设检验有关的操作命令.2.熟练掌握利用Mathematical软件对单个总体均值、方差的假设检验.3.掌握利用Mathematical软件对两个总体均值、方差有关的假设检验.三、常用命令1.求单个总体对均值的假设检验的命令:2.求两个总体对均值之差的假设检验的命令:3.求单个总体方差的假设检验的命令:4.求两个总体方差之比值的假设检验的命令:5.求标准正态分布有关概率的命令:6.求t分布有关概率的命令:7.求χ2分布有关概率的命令:8.求F分布有关概率的命令:四、练习内容设有甲、乙两种安眠药,比较其治疗效果.X表示服用甲药后睡眠时间延长时数,Y 表示服用乙药后睡眠时间延长时数,独立观察20个病人,其中10人服用甲药,另10人服用乙药,数据如下表:试就下列两种情况分析这两种药物的疗效有无显著性的差异.(显著性水平为0.05)(1)X与Y的方差相同;(2)X与Y的方差不同.程序:程序:结论:结论:五、思考与提高针对概率论与数理统计中左边、右边假设检验的问题,如何利用软件加以实现?31。
计算物理实验报告 常微分方程
常微分方程的边值问题和本征值问题一、问题描述利用搜索法和弦割法,得到该常微分方程的本征值,再利用打靶法计算多个本征值。
二、解决方法(一)搜索法1.先随便猜测k的一个试验值,程序中令k=12.由Numerov算法根据本题的条件,kn+1=kn=kn-1=k,s=0,得到yn+2,yn+1,yn间的迭代公式令con=(k*h)^2/12yn+2=2*(1-5*con)*yn+1/(1+con)-yn3自己给定φ的初始条件,然后利用公式得到边界值φ(1)4.然后以小的步长dk增加k值,这里令dk=1,每当φ(1)改变符号时,就将步长减半后倒退回来重复5.当步长小于所要求的容许误差时终止程序,此时的k值即为所求。
(二)弦割法1.随便猜测两个k值,这里令k0=1,k1=22.自己给定φ的初始条件,对两个k值分别利用上述公式进行迭代,得到边界值y1(1)和y2(1)。
3.比较y1(1)和y2(1)的绝对值大小。
若绝对值大,说明对应的k值距离本征值距离较远。
4.将(k0+k1)/2赋给k2,边界值绝对值小的对应的k值保持不变,边界值绝对值大的对应k值重新定位k2的值。
5.重复进行实验,当y1(1)和y(2)的差的绝对值小于容许误差时终止程序。
此时k1的值即为所求。
当搜索法和弦割法大致求出了一个本征值后,利用打靶法,调整k值再度进行搜索,得到多个本征值,绘出其中一个本征值对应的函数图像,观察其性质。
三、程序实现1.搜索法subroutine add(t,y0,y1) !利用子程序表示函数值的迭代implicit nonereal(8)::t,h,con,y0,y1,y2integer::i,nn=10000h=1.0/ncon=(t*h)**2/12do i=1,n-1y2=2*(1-5*con)*y1/(1+con)-y0 !利用Numerov算法,得到迭代公式y0=y1 !向前迭代y1=y2end doreturnend subroutine addprogram zy3implicit nonereal(8)::diffk,dk,yold,k,b0,b1integer::sb0=0.01 !取初始值,根据题目条件,令y0=y1,来保证x=0的位置导数为0b1=0.01s=1k=s !给定一个猜测的k值,此为搜索的初值dk=1 !给定步长diffk=0.0000001 !给定步长最后达到的误差范围call add(k,b0,b1)yold=b1 !通过运行子程序,得到由初始值积到x=1时的不为0的函数值do while(abs(dk)>diffk) !开始搜索k=k+dk !在k中走一步b0=0.01b1=0.01call add(k,b0,b1)if(yold*b1<0)then!若果y1变号k=k-dk !后退dk=dk/2.0 !步长减半end ifend dowrite(*,*)k !写出求得的本征值end2.弦割法subroutine add(t,b) !利用子程序表示函数值的迭代implicit nonereal(8)::t,h,con,y0,y1,y2,binteger::i,nb=0y0=0.01y1=0.01n=10000h=1.0/ncon=(t*h)**2/12do i=1,n-1y2=2*(1-5*con)*y1/(1+con)-y0 !利用Numerov算法,得到迭代公式y0=y1 !向前迭代y1=y2end dob=abs(y1) !得到x=1处函数值的绝对值,为确定k2点的位置做准备returnend subroutine addprogram zy3real(8)::a,k0,k1,k2,dk,m1,m2,dminteger::i,nk0=1 !给两个启动值k1=2k2=0dm=0.00000001 !表示k0和k1对应的函数值相等时允许的误差 m1=0 !此值表示k值取k0时,x=1处函数值的绝对值 m2=0 !此值表示k值取k1时,x=1处函数值的绝对值do while (.true.)call add(k0,m1)call add(k1,m2) !运行子程序,分别得到两个k值对应的x=1处的函数值的绝对值if(abs(m1-m2)<dm)exit!当k0k1对应的绝对值近似相等时退出循环if(m1>m2)then!如果k0对应的函数值绝对值较大k2=(k1+k0)/2.0 !k2点取在k1和k0的平均值k0=k2 !当k0对应的函数值绝对值较大时,表示其离本征值较远,而将其舍弃不用,赋k2值k1=k1 !此时k1距离本征值较近,不用变elsek2=(k1+k0)/2.0 !反之,舍去k1,k0不变k0=k0k1=k2end ifend dowrite(*,*)k2 !得到本征值kend program zy33.打靶法得多个本征值subroutine add(t,y0,y1) !利用子程序表示函数值的迭代implicit nonereal(8)::t,h,con,y0,y1,y2integer::i,nn=10000h=1.0/ncon=(t*h)**2/12do i=1,n-1y2=2*(1-5*con)*y1/(1+con)-y0 !利用Numerov算法,得到迭代公式y0=y1 !向前积分y1=y2end doreturnend subroutine addprogram zy3implicit nonereal(8)::diffk,dk,yold,k,b0,b1integer::sb0=0.01 !取初始值,根据题目条件,令y0=y1,来保证x=0的位置导数为0b1=0.01do s=1,100,2 !改变k的初值k=sdk=1 !给定步长diffk=0.0000001 !给定步长最后达到的误差范围call add(k,b0,b1)yold=b1 !通过运行子程序,得到由初始值积到x=1时的不为0的函数值do while(abs(dk)>diffk) !开始搜索k=k+dk !在k中走一步b0=0.01b1=0.01call add(k,b0,b1)if(yold*b1<0)then!若果y1变号k=k-dk !后退dk=dk/2.0 !步长减半end ifend dowrite(*,*)k !写出求得的本征值end doend4.选取一个本征值,看函数图像(k=1.5708)program zy3implicit nonereal(8)::diffk,dk,yold,k,y0,y1,y2,h,coninteger::i,nk=1.5708y0=0.001y1=0.001open(unit=10,file='p.txt')n=1000h=10.0/ncon=(k*h)**2/12do i=1,n-1write(10,*)y1,i*hy2=2*(1-5*con)*y1/(1+con)-y0 y0=y1y1=y2end doclose(10)end四、程序结果1.搜索法2.弦割法3.打靶法4.k=1.5708时的函数五、实验中出现的问题在利用打靶法找本征值的过程中,如果步长太小,会出现如下情况而当笔者将3代入k值用搜索法的时候,得到笔者认为出线该情况的可能是循环程序彼此间出现了问题,而当调整了步长之后,冲突减小,得到比较正常的特征值。
常微分方程数值解实验报告
常微分方程数值解实验报告学院:数学与信息科学专业:信息与计算科学:思义学号:201216524 课程:常微分方程数值解实验一:常微分方程的数值解法1、分别用Euler 法、改进的Euler 法(预报校正格式)和S —K 法求解初值问题。
(h=0.1)并与真解作比较。
⎩⎨⎧=++-=10(1y')y x y 1.1实验代码:%欧拉法function [x,y]=naeuler(dyfun,xspan,y0,h)%dyfun 是常微分方程,xspan 是x 的取值围,y0是初值,h 是步长 x=xspan(1):h:xspan(2); y(1)=y0; for n=1:length(x)-1y(n+1)=y(n)+h*feval(dyfun,x(n),y(n)); end%改进的欧拉法function [x,m,y]=naeuler2(dyfun,xspan,y0,h)%dyfun 是常微分方程,xspan 是x 的取值围,y0是初值,h 是步长。
%返回值x 为x 取值,m 为预报解,y 为校正解 x=xspan(1):h:xspan(2); y(1)=y0;m=zeros(length(x)-1,1); for n=1:length(x)-1 k1=feval(dyfun,x(n),y(n)); y(n+1)=y(n)+h*k1; m(n)=y(n+1);k2=feval(dyfun,x(n+1),y(n+1));y(n+1)=y(n)+h*(k1+k2)/2;end%四阶S—K法function [x,y]=rk(dyfun,xspan,y0,h)%dyfun是常微分方程,xspan是x的取值围,y0是初值,h是步长。
x=xspan(1):h:xspan(2);y(1)=y0;for n=1:length(x)-1k1=feval(dyfun,x(n),y(n));k2=feval(dyfun,x(n)+h/2,y(n)+(h*k1)/2);k3=feval(dyfun,x(n)+h/2,y(n)+(h*k2)/2);k4=feval(dyfun,x(n)+h,y(n)+h*k3);y(n+1)=y(n)+(h/6)*(k1+2*k2+2*k3+k4);end%主程序x=[0:0.1:1];y=exp(-x)+x;dyfun=inline('-y+x+1');[x1,y1]=naeuler(dyfun,[0,1],1,0.1);[x2,m,y2]=naeuler2(dyfun,[0,1],1,0.1);[x3,y3]=rk(dyfun,[0,1],1,0.1);plot(x,y,'r',x1,y1,'+',x2,y2,'*',x3,y3,'o');xlabel('x');ylabel('y');legend('y为真解','y1为欧拉解','y2为改进欧拉解','y3为S—K解','Location','NorthWest');1.2实验结果:x 真解y 欧拉解y1 预报值m 校正值y2 S—K解y30.0 1.0000 1.0000 1.0000 1.00000.1 1.0048 1.0000 1.0000 1.0050 1.00480.2 1.0187 1.0100 1.0145 1.0190 1.01870.3 1.0408 1.0290 1.0371 1.0412 1.04080.4 1.0703 1.0561 1.0671 1.0708 1.07030.5 1.1065 1.0905 1.1037 1.1071 1.10650.6 1.1488 1.1314 1.1464 1.1494 1.14880.7 1.1966 1.1783 1.1945 1.1972 1.19660.8 1.2493 1.2305 1.2475 1.2500 1.24930.9 1.3066 1.2874 1.3050 1.3072 1.30661.0 1.3679 1.3487 1.3665 1.3685 1.36792、选取一种理论上收敛但是不稳定的算法对问题1进行计算,并与真解作比较。
高等数学数学实验报告(两篇)2024
引言概述:高等数学数学实验报告(二)旨在对高等数学的相关实验进行探究与研究。
本次实验报告共分为五个大点,每个大点讨论了不同的实验内容。
在每个大点下,我们进一步细分了五到九个小点,对实验过程、数据收集、数据分析等进行了详细描述。
通过本次实验,我们可以更好地理解高等数学的概念和应用。
正文内容:一、微分方程实验1.利用欧拉法求解微分方程a.介绍欧拉法的原理和步骤b.详细阐述欧拉法在实际问题中的应用c.给出具体的实例,展示欧拉法的计算步骤2.应用微分方程建立模型求解实际问题a.介绍微分方程模型的建立方法b.给出一个具体的实际问题,使用微分方程建立模型c.详细阐述模型求解步骤和结果分析3.使用MATLAB求解微分方程a.MATLAB求解微分方程的基本语法和函数b.给出一个具体的微分方程问题,在MATLAB中进行求解c.分析结果的准确性和稳定性二、级数实验1.了解级数的概念和性质a.简要介绍级数的定义和基本概念b.阐述级数收敛和发散的判别法c.讨论级数的性质和重要定理2.使用级数展开函数a.介绍级数展开函数的原理和步骤b.给出一个函数,使用级数展开进行近似计算c.分析级数近似计算的精确度和效果3.级数的收敛性与运算a.讨论级数收敛性的判别法b.介绍级数的运算性质和求和法则c.给出具体的例题,进行级数的运算和求和三、多元函数极值与最值实验1.多元函数的极值点求解a.介绍多元函数的极值点的定义和求解方法b.给出一个多元函数的实例,详细阐述求解过程c.分析极值点对应的函数值和意义2.多元函数的条件极值与最值a.讨论多元函数的条件极值的判定法b.给出一个具体的多元函数,求解其条件极值和最值c.分析条件极值和最值对应的函数值和意义3.利用MATLAB进行多元函数极值与最值的计算a.MATLAB求解多元函数极值与最值的基本语法和函数b.给出一个多元函数的具体问题,在MATLAB中进行求解c.分析结果的准确性和可行性四、曲线积分与曲面积分实验1.曲线积分的计算方法与应用a.介绍曲线积分的定义和计算方法b.给出一个具体的曲线积分问题,详细阐述计算过程c.分析曲线积分结果的几何意义2.曲线积分的应用举例a.讨论曲线积分在实际问题中的应用b.给出一个实际问题,使用曲线积分进行求解c.分析曲线积分结果的实际意义和应用价值3.曲面积分的计算方法与应用a.介绍曲面积分的定义和计算方法b.给出一个具体的曲面积分问题,详细阐述计算过程c.分析曲面积分结果的几何意义五、空间解析几何实验1.空间曲线的参数方程表示与性质a.介绍空间曲线的参数方程表示和性质b.给出一个具体的空间曲线,转化为参数方程表示c.分析参数方程对应的几何意义和性质2.平面与空间直线的位置关系a.讨论平面与空间直线的位置关系的判定方法b.给出一个具体的平面与空间直线的问题,判定其位置关系c.分析位置关系对应的几何意义和应用实例3.空间直线与平面的夹角和距离计算a.介绍空间直线与平面的夹角和距离的计算方法b.给出一个具体的空间直线和平面,计算其夹角和距离c.分析夹角和距离计算结果的几何意义总结:通过本次高等数学数学实验报告(二),我们深入了解了微分方程、级数、多元函数极值与最值、曲线积分、曲面积分以及空间解析几何的相关概念和应用。
常微分方程的数值解法实验报告
常微分方程的数值解法实验报告实验报告:常微分方程的数值解法摘要:常微分方程(ODE)是描述动力学系统中物理量随时间变化的数学方程,广泛应用于自然科学和工程领域。
然而,对于一些复杂的非线性ODE,很难找到解析解。
因此,我们需要数值解法来求解这些方程。
本实验报告将介绍四种常见的常微分方程数值解法,分别是欧拉法、改进的欧拉法、四阶龙格-库塔法和自适应步长的龙格-库塔法,并通过数值实验比较它们的精度和效率。
1.引言在实际问题中,许多物理量的变化规律可以由常微分方程描述。
然而,对于复杂的非线性ODE,很难找到解析解。
因此,为了解决这类问题,我们需要借助数值方法来求解。
2.方法本实验采用四种常见的常微分方程数值解法:欧拉法、改进的欧拉法、四阶龙格-库塔法和自适应步长的龙格-库塔法。
(1)欧拉法是最简单的数值解法,通过将微分方程转化为差分方程,使用离散的步长来近似微分方程。
(2)改进的欧拉法在欧拉法的基础上进行了改进,使用预估-校正的方法来提高精度。
(3)四阶龙格-库塔法是一种经典的数值解法,通过利用不同步长处的斜率来近似微分方程,具有较高的精度。
(4)自适应步长的龙格-库塔法是在四阶龙格-库塔法的基础上改进而来的,根据步长的大小自适应地选择不同的步长,同时保证精度和效率。
3.实验设计为了比较这四种数值解法的精度和效率,我们设计了两个实验。
实验一是求解一阶常微分方程:dy/dx = -2x,初始条件y(0) = 1,解析解为y = 1 - x^2、实验二是求解二阶常微分方程:d^2y/dx^2 + y = 0,初始条件y(0) = 0,dy/dx(0) = 1,解析解为y = sin(x)。
4.结果与分析实验一中,比较四种数值解法在不同步长下的近似解和解析解,计算其误差。
实验结果表明,四阶龙格-库塔法和自适应步长的龙格-库塔法具有较高的精度,而欧拉法和改进的欧拉法的精度较低。
实验二中,我们比较四种数值解法在不同步长下的近似解和解析解,并计算其误差。
常微分方程--酒驾问题
东南大学数学建模实验报告实验内容:酒驾问题一实验目的(1)掌握常微分方程建模问题(2)学会使用Matlab进行常微分方程的求解二实验内容与要求国家质量监督检验检疫局 2004年5月31日发布了新的《车辆驾驶人员血液、呼气酒 精含量阈值与检验》国家标准,新标准规定,车辆驾驶人 员血液中的酒精含量大于或等于20毫克/百毫升,小于80 毫克/百毫升为饮酒驾车(原标准是小于100毫克/百毫 升),血液中的酒精含量大于或等于80毫克/百毫升为醉 酒驾车(原标准是大于或等于100毫克/百毫升 )。
在中某人午12点喝了一瓶啤酒,下午6点检查时符合 新的驾车标准,紧接着他在吃晚饭时又喝了一瓶啤酒,为 了保险起见他呆到凌晨2点才回家,又一次遭遇检查时却 被定为饮酒驾车,这让他懊恼又困惑,为什么喝了同样多 的酒,两次检查结果会不一样呢?请你参考下面的数据建立饮酒后 血液中酒精含量的数学模型,并讨论以下问题: 1、对某人碰到的情况作出解释; 2、假设酒是在很短时间内喝的,在喝了3瓶啤酒或半斤低度白酒后多长时间内驾车就会违反上述标准.3、怎样估计血液中酒精含量在什么时候最高。
4、根据你的模型论证:如果天天喝酒,是否能开车? 以下是某人喝了两瓶啤酒后血液酒精浓度(毫克/百毫升)三 假设及建模假设一:机体分为中心室和周边室,两个室的容积在过程中保持不变。
假设二:药物从一室向另一室的转移速率,及向体外的排除速率,与该室的酒精浓度成正比。
假设三:只在中心室一体外有酒精交换,即酒精从体外进入中心室,最后又从中心室排出体外,与转移和排除的数量相比,酒精的吸收可以忽略。
建模:二室模型的示意图如下图所示:饮酒()t f 0两个房室中酒精量)(),(21t x t x 满足的微分方程。
)(1t x 的变化率由一室向二室的转移112x k -,一室向体外排除113xk -,二室向一室的转移221x k 及酒精)(0t f 组成;)(2t x的变化率由一室向二室的转移112x k 及二室向一室的转移221x k -组成,于是有: )(022********t f x k x k x k dtdx ++--=2211122x k x k dtdx -= (1) )(t x i 与血液中酒精含量)(t c i 、房室容积i V 显然有关系式2,1.................................),........()(==i t c V t x i i i (2)将(2)式代入(1)式可得:2211122121022112113121)()(c k c k V V dt dc V t f c k V Vc k k dt dc -=+++-= (3)喝酒相当于在酒精进入中心室之前先有一个将酒精吸收入血液的过程,可以简化为有一个吸收室,如下图,)(0t x 为吸收室的酒精,酒精由吸收室进入中心室的转移速率系数为01k ,于是)(0t x 满足:00010)0(D x x k dt dx =-= (4)当0)0(,)0(,0)(2110===c V D c t f 时,(3)可以化为: t t Be Ae t c βα--+=)(1四 代码及结果format short g% 题中提供的某人喝了两瓶啤酒后血液酒精浓度随时间变化表t=[ 0.25; 0.5; 0.75; 1; 1.5; 2; 2.5; 3; 3.5; 4; 4.5; 5; 6; 7; 8; 9; 10; 11; 12; 13; 14; 15; 16 ];c=[ 30; 68; 75; 82; 84; 77; 70; 68; 58; 51; 50; 41; 38; 35; 28; 25; 18; 15; 12; 10; 7; 7; 4 ];% 根据此变化表拟合求解相关系数ft =fittype('A1*exp(-a*x)+B1*exp(-b*x)');options = fitoptions('Method','NonlinearLeastSquares');options.StartPoint = [0 -1000 0 0];cfit = fit(t,c,ft,options);plot( cfit, t, c, 'o' );A1=cfit.A1B1=cfit.B1a=cfit.ab=cfit.b由此解得:(数值见右图,拟合曲线见下图)A1 = 110.55B1 = -151.46a = 0.17949b = 2.8243%---1---%%问题:某人中午12点喝了一瓶啤酒,下午6点检查合格,晚饭又喝一瓶,次日凌晨2点检查未通过,请对此情况做出解释。
(完整版)Mathematica——常微分方程、拉氏变换与级数实验
§13.5常微分方程、拉氏变换与级数实验[学习目标]1. 会用Mathematica 求解微分方程(组);2. 能用Mathematica 求微分方程(组)的数值解;3. 会利用Mathematica 进行拉氏变换与逆变换;4. 能进行幂级数和傅里叶级数的展开.一、 常微分方程(组)Mathematica 能求常微分方程(组)的准确解,能求解的类型大致覆盖了人工求解的范围,功能很强。
但不如人灵活(例如在隐函数和隐方程的处理方面),输出的结果与教材上的答案可能在形式上不同。
另外,Mathematica 求数值解也很方便,且有利于作出解的图形。
在本节中,使用Laplace 变换解常微分方程(组)的例子也是十分成功的,过去敬而远之的方法如今可以轻而易举的实现了.求准确解的函数调用格式如下:DSolve [eqn ,y [x ],x ] 求方程eqn 的通解y (x ),其中自变量是x 。
DSolve[{eqn ,y[x 0]= =y 0},y[x ],x ] 求满足初始条件y (x 0)= y 0的特解y (x )。
DSolve[{eqn1,eqn2,…},{y 1[x ],y 2[x ],…},x] 求方程组的通解. DSolve[{equ1,…,y 1[x 0]= =y 10,…},{y 1[x],y 2[x ],…},x] 求方程组的特解。
说明:应当特别注意,方程及各项参数的表述方式很严格,容易出现输入错误。
微分方程的表示法只有通过例题才能说清楚.例1 解下列常微分方程(组):(1)25)1(12+++='x x yy ,(2)y x x y y )(132++=', (3) ⎩⎨⎧-='='yz z y , (4)⎩⎨⎧-='='yz zy 的通解及满足初始条件y (0)=0,z (0)=1的特解。
解:In [1]:=DSolve [y ′[x ]= =2y[x ]/(x+1)+(x+1)^(5/2), y [x],x ]Out [1]=⎭⎬⎫⎩⎨⎧⎭⎬⎫⎩⎨⎧+++→]1[)1()1(32][22/7c x x x yIn [2]:=DSolve [y ′[x]= =(1+y [x ]^2)/((x+x^3)y[x ]),y[x],x]Out [2]={{2211]1[11][x c x x y ++---→}, {2211]1[11][xc x x y ++--→}} In[3]:=DSolve [{y ′[x]= =z [x ],z ′[x]= = —y [x ]}, {y [x],z[x ]},x ]Out[3]={{y[x ]→C[1]Cos[x]+ C [2]Sin [x ], z [x ]→C[2]Cos[x ]— C [1]Sin[x]}}In [4]:=DSolve [{y ′[x]= =z [x],z ′[x ]= = —y[x],y[0]= =0,z [0]= =1}, {y[x],z [x ]},x]Out[4]={{y [x ]→Sin [x ],z[x]→Cos[x]}}提示:认真观察上例,可以从中学习输入格式,未知函数总带有自变量,等号用连续键入两个等号表示,这两点由于不习惯会出错!导数符号用键盘上的撇号,连续两撇表示二阶导数,这与习惯相同。
- 1、下载文档前请自行甄别文档内容的完整性,平台不提供额外的编辑、内容补充、找答案等附加服务。
- 2、"仅部分预览"的文档,不可在线预览部分如存在完整性等问题,可反馈申请退款(可完整预览的文档不适用该条件!)。
- 3、如文档侵犯您的权益,请联系客服反馈,我们会尽快为您处理(人工客服工作时间:9:00-18:30)。
实验一 常微分方程
1. 分别用Euler 法和ode45解下列常微分方程并与解析解比较: (1) ,(0)1,13y x y y x '=+=<<
odefun=inline (‘x+y ’,‘x ’,‘y ’); [x1,y1]=ode45(odefun,[0,3],1); [x2,y2]=euler(odefun,[0,3],1,0.01); S=dsolve(‘Dy=x+y’,’y(0)=1’,’x’); hold on
plot(x1,y1,’o’); plot(x2,y2,’r’);
fplot(‘2*exp(x )-x-1’,[0,3]);
title(‘Euler 法和ode45法与解析解’); Legend(‘ode45法’,’Euler 法’,’解析解’); ylabel(‘x 轴’); ylabel(‘y 轴’); hold off
(2) 20.01()2sin(),(0)0,(0)1,05y y y t y y t ''''-+===<< odefun=@(t,y)[y(2);0.01*y(2)^2-2*y(1)+sin(t)]; [t,y]=ode45(odefun,[0 5],[0:1]); plot(t,y(:,1));
title(‘第一题(2)图像’); xlabel(‘x 轴’); ylabel(‘y 轴’);
实验一 常微分方程
2. 求一通过原点的曲线,它在(,)x y 处的切线斜率等于2
2,0 1.57.x y x +<<若x 上限增为1.58,1.60会发生什么?
t2(1.57)
实验一 常微分方程
3. 求解刚性方程组:
11212
121100.25999.750.5,(0)1,050.
999.751000.250.5,(0)1,y y y y x y y y y '=-++=⎧<<⎨
'=-+=-⎩
4. (温度过程)夏天把开有空调的室内一支读数为20℃的温度计放到户外,10分钟后读2
5.2℃, 再过10分钟后读数28.32℃。
建立一个较合理的模型来推算户外温度。
实验一常微分方程
5. (广告效应)某公司生产一种耐用消费品,市场占有率为5%时开始做广告,一段时间的市场跟踪调查后,该公司发现:单位时间内购买人口百分比的相对增长率与当时还没有买的百分比成正比,且估得此比例系数为0.5。
(1) 建立该问题的数学模型,并求其数值解与模拟结果作以比较;
(2) 厂家问:要做多少时间广告,可使市场购买率达到80%?
实验一常微分方程
6. (肿瘤生长) 肿瘤大小V生长的速率与V的a次方成正比,其中a为形状参数,0≤a≤1;而其比例系数K随时间减小,减小速率又与当时的K值成正比,比例系数为环境参数b。
设某肿瘤参数a=1, b=0.1, K的初始值为2,V的初始值为1。
问
(1)此肿瘤生长不会超过多大?
(2)过多长时间肿瘤大小翻一倍?
(3)何时肿瘤生长速率由递增转为递减?
(4)若参数a=2/3呢?
ans=2/exp(t/10)
ans=exp(20)/exp(20/exp(t/10))
fanbei=36
实验一常微分方程Zhuanbian=2034930319768065/(2097152*exp(20/exp(t/10)*exp(t/10))
实验一常微分方程
fanbei=40
zhuanbian=(2*(20/exp(t/10)-23)^2/(9*exp(t/10))
选做题:
1.(生态系统的振荡现象)第一次世界大战中,因为战争很少捕鱼,按理战后应能捕到更多的鱼才是。
可是大战后,在地中海却捕不到鲨鱼,因而渔民大惑不解。
令x1为鱼饵的数量,x2为鲨鱼的数量,t为时间。
微分方程为
(5.20)
式中a1, a2, b1, b2都是正常数。
第一式鱼饵x1的增长速度大体上与x1成正比,即按a1x1比率增加, 而被鲨鱼吃掉的部分按b1x1x2的比率减少;第二式中鲨鱼的增长速度由于生存竞争的自然死亡或互相咬食按a2x2的比率减少,但又根据鱼饵的量的变化按b2x1x2的比率增加。
对a1=3, b1=2, a2=2.5, b2=1, x1(0)=x2(0)=1求解。
画出解曲线图和相轨线图,可以观察到鱼饵和鲨鱼数量的周期振荡现象。
实验一常微分方程
实验报告总结:
本次题中多次运用欧拉方法和ode算法我不是很懂,通过书本和同学的讨论,最后才弄懂。
常微分方程在现实中的应用很重要,用处很广。