一.实验题目
1.(必做题)解微分方程(组)
(1) 322232(1)(0)0,'(0)1,''(0)1,d y d y dy y dx dx dx y y y ?=---???===-?
(提示可以考虑[0,20]x ∈,以特解函数及其一阶、二阶导数曲线图形来表示)
解: ①将高阶微分方程化为一阶微分方程组,设y y y y y y ''='==321,,,则有
()???
????---='='='2122333221
1y y y y y y y y ②建立函数文件
function dy=myfun(x,y)
dy=[y(2);y(3);(y(3)-1)^2-y(2)-y(1)^2];
③主程序:
[x,y]=ode45('myfun',[0,20],[0;1;-1]);
plot(x,y(:,1),'*',x,y(:,2),'+',x,y(:,3),'o')
%legend('y','y 的一阶导数','y 的二阶导数');
④结果
注意此题得不到解析解,只能用数值解,解法可参看PPT 中数值解例题3 (2)运用数值解手段描述下面常微分方程组在初值0[0;0;110]x e ∈-下的相空间的相轨线.
11232233
1223'()8()/3()()'()10()10()'()()()28()()x t x t x t x t x t x t x t x t x t x t x t x t =-+??=-+??=-+-?
解:①建立函数文件
function dx=lorenz(t,x)
dx=[-8/3*x(1)+x(2)*x(3);-10*x(2)+10*x(3);-x(1)*x(2)+28*x(2)-x(3)];
②主程序文件
[t,x]=ode45('lorenz',[0,100],[0;0;1e-10]);
axis([0 40 -20 20 -20 20]);
plot3(x(:,1),x(:,2),x(:,3));
grid on
③结果
(3)求0.0199.99100(0)2,(0)1dy y z
dx dz z dx y z ?=--???=-??==???
的数值解,并画出图像 解:首先建立odefun1 .m 如下:
function dy=odefun1(x,y);
dy=[-0.01*y(1)-99.99*y(2);-100*y(2)];
然后建立主程序shiyan2_3.m
clc
clear
close all
[x,y]=ode15s('odefun1',[0 100],[2;1])
plot(x,y(:,1) ,'*',x,y(:,2),'r*')
结果为
(4)求下列方程的通解及特解
222'''()022,'22x y xy x n y y y πππ?++-=??????==- ? ????
???(Bessel 方程,令12n =) 解:求通解的主程序为(syms n )
diff_y='x^2*D2y+x*Dy+(x^2-(1/2)^2)*y=0';
y=dsolve(diff_y,'x')
结果为:
y=C1/x^(1/2)*sin(x)+C2/x^(1/2)*cos(x)
y =
(2^(1/2)*C12*cos(x))/(pi^(1/2)*x^(1/2)) + (2^(1/2)*C13*sin(x))/(pi^(1/2)*x^(1/2)) 求特解的主程序为
diff_y='x^2*D2y+x*Dy+(x^2-(1/2)^2)*y=0';
y=dsolve(diff_y,'y(pi/2)=2,Dy(pi/2)=-2/pi','x')
结果为:
y =2^(1/2)*pi^(1/2)/x^(1/2)*sin(x)
y =
(2*sin(x)*(pi/2)^(1/2))/x^(1/2) + (cos(x)*(2/(pi/2)^(1/2) - pi/(pi/2)^(3/2)))/(2*x^(1/2))
2.(必做题)凶杀案作案时间问题:受害者的尸体于晚上7:30被发现,法医于晚上8:20赶到凶案现场,测得尸体温度为32.6℃;一小时后,当尸体即将被抬走时,测得尸体温度为31.4℃,室温在几个小时内始终保持21.1℃。此案最大的嫌疑犯张某声称自己是无罪的,并有证人说:“下午张某一直在办公室上班,5:00时打完电话后就离开了办公室”。从张某到受害者家(凶案现场)步行需5分钟,现在的问题是,张某不在凶案现场的证言能否被采信,使他排除在嫌疑犯之外。(提示:Newton 冷却定理告诉我们“物体在介质中冷却速度同该物体温度与介质温度之差成正比”)
解:首先应确定凶案的发生时间,若死亡时间在下午5点5分之前,则张某就不是嫌疑犯,否则不能将张某排除。
设T(t)表示t 时刻尸体的温度,并记晚上8:20为t=0,则T(0)=32.6℃,T(1)=31.4℃。假设受害者死亡时体温是正常的,即T=37℃(查资料)是要确定受害者死亡的时间,也就是求T(t)=37℃的时刻,进而确定张某是否是嫌疑犯。
人体体温受大脑神经中枢调节。人死亡后体温调节的功能消失,尸体的温度受外界环境温度的影响。假设尸体温度的变化率服从牛顿冷却定律,即尸体温度的变化律与他同周围的温度差成正比。即
模型:由Newton 冷却定理可得一阶线性微分方程模型 (21.1)(0)32.6
t dT T dt T λ?=--???=?
求解:(1)首先用dsolve 求解该方程的解析解
程序:
syms lamd
sy3d11='DT+lamd*(T-21.1)=0';
T=dsolve(sy3d11,'T(0)=32.6','t')
结果:
T =211/10+23/2*exp(-lamd*t)
或
T =23/(2*exp(lamd*t)) + 211/10
(2)求解参数lamd
可以利用初始条件“1小时后,当尸体即将被抬走时,测得尸体温度为31.4℃” 由上式可以得到:31.4=21.1+11.5*exp(-1*lamd)
lamd 的值为0.11020314013361429463890984998294
(程序:lamd=solve('31.4-21.1-11.5*exp(-1*lamd)=0','lamd')
(3)求解t0
当T=37℃时,有t=-2.95 小时
(程序:t0=solve('37-21.1-11.5*exp(-0.11*t)','t'))
=-2小时57分,
8小时20分-2小时57分=5小时23分。即死亡时间大约在下午5:23,
因此张某不能被排除在嫌疑犯之外。
北京交通大学海滨学院考试试题 课程名称:数学实验2010-2011第一学期出题教师:数学组适用专业: 09机械, 物流, 土木, 自动化 班级:学号:姓名: 选做题目序号: 1.一对刚出生的幼兔经过一个月可以长成成兔, 成兔再经过一个月后可以 繁殖出一对幼兔. 如果不计算兔子的死亡数, 请用Matlab程序给出在未来24个月中每个月的兔子对数。 解: 由题意每月的成兔与幼兔的数量如下表所示: 1 2 3 4 5 6 ··· 成兔0 1 1 2 3 5··· 幼兔 1 0 1 1 2 3··· 运用Matlab程序: x=zeros(1,24); x(1)=1;x(2)=1; for i=2:24 x(i+1)=x(i)+x(i-1); end x 结果为x = 1 1 2 3 5 8 13 21 3 4 5 5 89 144 233 377 610 987 1597 2584 4181 6765 1094 6 7711 2865 7 46368 2.定积分的过程可以分为分割、求和、取极限三部分, 以1 x e dx 为例, 利用
已学过的Matlab 命令, 通过作图演示计算积分的过程, 并与使用命令int() 直接积分的结果进行比较. 解:根据求积分的过程,我们先对区间[0,1]进行n 等分, 然后针对函数x e 取和,取和的形式为10 1 i n x i e e dx n ξ=≈ ∑ ? ,其中1[ ,]i i i n n ξ-?。这里取i ξ为区间的右端点,则当10n =时,1 x e dx ?可用10 101 1.805610 i i e ==∑ 来近似计算, 当10n =0时,100 100 1 01 =1.7269100 i x i e e dx =≈ ∑?,当10n =000时,10000 10000 1 1 =1.718410000 i x i e e dx =≈ ∑ ?. 示意图如下图,Matlab 命令如下: x=linspace (0,1,21); y=exp(x); y1=y(1:20); s1=sum(y1)/20 y2=y(2:21); s2=sum(y2)/20 plot(x,y); hold on for i=1:20 fill([x(i),x(i+1),x(i+1),x(i),x(i)],[0,0,y(i),y(i),0],'b') end syms k;symsum(exp(k/10)/10,k,1,10);%n=10 symsum(exp(k/100)/100,k,1,100);%n=100 symsum(exp(k/10000)/10000,k,1,10000);%n=10000
东华大学20 ~ 20 学年第__ __学期期_末_试题A 踏实学习,弘扬正气;诚信做人,诚实考试;作弊可耻,后果自负 课程名称______高等数学实验___________使用专业____ 班级_____________姓名________________学号__________ 机号 要求:写出M 函数(如果需要的话)、MATLAB 指令和计算结果。1.设矩阵A = 6 14230215 1 0321 21----, 求A 的行列式和特征值。 2. 设 f (x ,y ) =2x cos (xy 2 ),求 21,2 x y f x y ==???。
3. 求积分? --1 2 2 1)2(x x xdx 的数值解。 4. 求解微分方程0.5e - x d y -sin x d x=0, y (0)=0, 要求写出x =2 时的y 值。 5. 求解下列方程在k=6,θ=π/3附近的解???=-=-1)sin (3 )cos 1(θθθk k
6. 取k 7. 编写一个M 函数文件,使对任意给定的精度ε, 求N 使得 επ≤-∑=612 1 2 N n n 并对ε= 0.001求解。
8. 在英国工党成员的第二代加入工党的概率为0.5,加入保守党的概率为0.4,加入自由党的概率为0.1。而保守党成员的第二代加入保守党的概率为0.7,加入工党的概率为0.2,加入自由党的概率为0.1。而自由党成员的第二代加入保守党的概率为0.2,加入工党的概率为0.4,加入自由党的概率为0.4。求自由党成员的第三代加入工党的概率是多少?假设这样的规律保持不变,在经过很多代后,英国政党大致分布如何?
大学数学实验 项目一 矩阵运算与方程组求解 实验1 行列式与矩阵 实验目的 掌握矩阵的输入方法. 掌握利用Mathematica (4.0以上版本) 对矩阵进行转置、加、减、数乘、相乘、乘方等运算, 并能求矩阵的逆矩阵和计算方阵的行列式. 基本命令 在Mathematica 中, 向量和矩阵是以表的形式给出的. 1. 表在形式上是用花括号括起来的若干表达式, 表达式之间用逗号隔开. 如输入 {2,4,8,16} {x,x+1,y,Sqrt[2]} 则输入了两个向量. 2. 表的生成函数 (1) 最简单的数值表生成函数Range, 其命令格式如下: Range[正整数n]—生成表{1,2,3,4,…,n }; Range[m, n]—生成表{m ,…,n }; Range[m, n, dx]—生成表{m ,…,n }, 步长为d x . (2) 通用表的生成函数Table. 例如,输入命令 Table[n^3,{n,1,20,2}] 则输出 {1,27,125,343,729,1331,2197,3375,4913,6859} 输入 Table[x*y,{x,3},{y,3}] 则输出 {{1,2,3},{2,4,6},{3,6,9}} 3. 表作为向量和矩阵 一层表在线性代数中表示向量, 二层表表示矩阵. 例如,矩阵 ??? ? ??5432 可以用数表{{2,3},{4,5}}表示. 输入 A={{2,3},{4,5}} 则输出 {{2,3},{4,5}} 命令MatrixForm[A]把矩阵A 显示成通常的矩阵形式. 例如, 输入命令: MatrixForm[A] 则输出 ??? ? ??5432 但要注意, 一般地, MatrixForm[A]代表的矩阵A 不能参与运算. 输入 B={1,3,5,7} 输出为 {1,3,5,7} 输入 MatrixForm[B] 输出为
一、实验内容 P206第六题 function f=wuyan2(c) y=[3.9 5.3 7.2 9.6 12.9 17.1 23.2 31.41 38.6 50.2 62.9 76.0 92.0 106.5 123.2 131.7 150.7 179.3 204.0 226.5 251.4 281.4] t=[0 10 20 30 40 50 60 70 80 90 100 110 120 130 140 150 160 170 180 190 200 210] f=y-c(1)/(1+c(1)/3.9-1)*exp^(-c(2)*t) c0=[1 1] c=lsqnonlin('wuyan2',c0) P206第七题 function f=wuyan1(c) q=[0.4518 0.4862 0.5295 0.5934 0.7171 0.8964 1.0202 1.1963 1.4928 1.6909 1.8548 2.1618 2.6638 3.4634 4.6759 5.8478 6.7885 7.4463 7.8345 8.2068 8.9468 9.7315 10.5172 11.7390 13.6876 ]; k=[0.0911 0.0961 0.1230 0.1430 0.1860 0.2543 0.3121 0.3792 0.4754 0.4410 0.4517 0.5595 0.8080 1.3072 1.7042 2.0019 2.2914 2.4941 2.8406 2.9855 3.2918 3.7214 4.3500 5.5567 7.0477]; l=[4.2361 4.3725 4.5295 4.6436 4.8179 4.9873 5.1282 5.2783 5.4334 5.5329 6.4749 6.5491 6.6152 6.6808 6.7455 6.8065 6.8950 6.9820 7.0637 7.1394 7.2085 7.3025 7.3470 7.4432 7.5200]; f=q-c(1)*k.^c(2).*l.^c(3) c0=[1 1 1] c=lsqnonlin('wuyan1',c0) c = 0.4091 0.6401 1.1446 a=0.4091 α=0.6401 β=1.1446 P239第五题 c=[-20 -30]; A=[1 2;5 4]; b=[20 70]; v1=[0 0]; [x,f,ef,out,lag]=linprog(c,A,b,[],[],v1) z=-f x = 10.0000 5.0000
清华大学数学实验报告4
————————————————————————————————作者: ————————————————————————————————日期: ?
电13 苗键强2011010645
一、实验目的 1.掌握用 MATLAB 软件求解非线性方程和方程组的基本用法, 并对结果作初步分析; 2.练习用非线性方程和方程组建立实际问题的模型并进行求解。 二、实验内容 题目1 【问题描述】 (Q1)小张夫妇以按揭方式贷款买了1套价值20万元的房子,首付了5万元,每月还款1000元,15年还清。问贷款利率是多少? (Q2)某人欲贷款50 万元购房,他咨询了两家银行,第一家银行 开出的条件是每月还4500元,15 年还清;第二家银行开出的条件是每年还45000 元,20年还清。从利率方面看,哪家银行较优惠(简单假设:年利率=月利率×12)? 【分析与解】 假设初始贷款金额为x0,贷款利率为p,每月还款金额为x,第i 个月还完当月贷款后所欠银行的金额为x i,(i=1,2,3,......,n)。由题意可知: x1=x0(1+p)?x x2=x0(1+p)2?x(1+p)?x x3=x0(1+p)3?x(1+p)2?x(1+p)?x ……
x n=x0(1+p)n?x(1+p)n?1???x(1+p)?x =x0(1+p)n?x (1+p)n?1 p =0 因而有: x0(1+p)n=x (1+p)n?1 p (1) 则可以根据上述方程描述的函数关系求解相应的变量。 (Q1) 根据公式(1),可以得到以下方程: 150p(1+p)180?(1+p)180+1=0 设 f(p)=150p(1+p)180?(1+p)180+1,通过计算机程序绘制f(p)的图像以判断解p的大致区间,在Matlab中编程如下: fori = 1:25 t = 0.0001*i; p(i) = t; f(i) =150*t*(1+t).^180-(1+t).^180+1; end; plot(p,f),hold on,grid on; 运行以上代码得到如下图像:
1.(1) [1 2 3 4;0 2 -1 1;1 -1 2 5;]+(1/2).*([2 1 4 10;0 -1 2 0;0 2 3 -2]) 2. A=[3 0 1;-1 2 1;3 4 2],B=[1 0 2;-1 1 1;2 1 1] X=(B+2*A)/2 3. A=[-4 -2 0 2 4;-3 -1 1 3 5] abs(A)>3 % 4. A=[-2 3 2 4;1 -2 3 2;3 2 3 4;0 4 -2 5] det(A),eig(A),rank(A),inv(A) 求计算机高手用matlab解决。 >> A=[-2,3,2,4;1,-2,3,2;3,2,3,4;0,4,-2,5] 求|A| >> abs(A) ans = ( 2 3 2 4 1 2 3 2 3 2 3 4 0 4 2 5 求r(A) >> rank(A) ans =
4 求A-1 《 >> A-1 ans = -3 2 1 3 0 -3 2 1 2 1 2 3 -1 3 -3 4 求特征值、特征向量 >> [V,D]=eig(A) %返回矩阵A的特征值矩阵D 与特征向量矩阵V , V = - + + - - + - + - + - + D = { + 0 0 0 0 - 0 0 0 0 + 0 0 0 0 - 将A的第2行与第3列联成一行赋给b >> b=[A(2,:),A(:,3)'] b = 《 1 - 2 3 2 2 3 3 -2