蹦极运动员自由落体微分方程求解(尝试)
题目描述:
一个质量为68.1kg的蹦极运动员从一个静止的热气球上滑落,计算蹦极运动员前12s自由落体过程中的速度。如果所使用的绳索无限长,确定蹦极运动员能达到的最终速度。假设阻尼系数为0.25kg/m。
数学模型:
d
v
d
t
=
g
−
c
d
m
v
2
\frac{dv}{dt} = g - \frac{c_d}{m}v^2
dtdv=g−mcdv2
一、完整代码
实质为已知初值的简单微分方程求解
%自由落体微分方程求解
%% 1、求解微分方程通解
syms v(t) g cd m %定义变量
eqn = diff(v,t) == g - cd/m*v^2; %列写微分方程
cond = v(0)==0; %微分方程初始条件
sol = dsolve(eqn,cond); %求解微分方程
%% 2、符号表达式求值
m0 = 68.1;
g0 = 9.81;
cd0 = 0.25; %常值变量条件
vt = subs(subs(subs(sol,'m',m0),'g',g0),'cd',cd0); %代入常值变量
i = 1;
%% 3、循环求解具体速度值,步长0.5,时间区间[0,12]
for ti = 0:0.5:12
vval(i) = double(subs(vt,'t',ti));%计算数值用double
i = i+1;
end
%% 4、根据得出的数据点绘图plot
tf = 0:0.5:12; %和vval数据量相同
plot(tf,vval,'--b','LineWidth',2)
title('自由落体速度随时间变化曲线'),xlabel('时间t(s)'),ylabel('速度(m/s)');
%可以直接在图上修改数据&链接到原数据表+显示数据点
grid on %显示网格线
hold on %添加新绘图时保留当前绘图
plot(11,vval(23),'ro','MarkerSize',8)%标记(11,23)对应的点
text(11,vval(23)+3,num2str(double(vval(23))),'HorizontalAlignment','center')
%num2str将数值数组转换为表示数字的字符数组
%添加文本说明,显示当前对应点具体数值
%vval需要显示为数值,前两项定位标记,为11的居中
%% 5、蹦极运动员速度求极限值limit
limval = double(limit(vt,'t',inf));
%inf表示t趋向于无穷大
%vt为代入其他常值之后的仅与t有关的表达式
%% 6、输出极限值
%显示变量或者文字 dsip,仅有一个参数,需要提前将不同类型字符串合并
outstr = strcat('蹦极运动员自由下落的极限速度是:',num2str(limval));
disp(outstr)
二、 涉及的知识点
1、微分方程的直接表示+已知初值desolve求解简单微分方程
此处变量并非微分方程的变量,而是方程未知量,防止报错
syms v(t) g cd m %定义变量
eqn = diff(v,t) == g - cd/m*v^2; %列写微分方程
cond = v(0)==0; %微分方程初始条件
sol = dsolve(eqn,cond);
2、利用subs代入已知量
由于m、g、cd被定义为未知量,因此不可直接赋值求解,需要subs代入数值
vt = subs(subs(subs(sol,'m',m0),'g',g0),'cd',cd0); %代入常值变量
3、double变计算公式为数值
double用于计算最终数值,若无则输出为代入数值后 vt 的计算式
for ti = 0:0.5:12 %可不用“;”
vval(i) = double(subs(vt,'t',ti));
i = i+1;
end
4、plot绘图、title图名、xlabel坐标说明、text添加文本说明(留坑)
可突出显示具体的点,注意title等必须在绘图plot命令之后
plot(tf,vval,'--b','LineWidth',2)%对应数据量相同
title('自由落体速度随时间变化曲线'),xlabel('时间t(s)'),ylabel('速度(m/s)');
%可以直接在图上修改数据&链接到原数据表+显示数据点
grid on %显示网格线
hold on %添加新绘图时保留当前绘图
plot(11,vval(23),'ro','MarkerSize',8)%标记(11,23)对应的点
text(11,vval(23)+3,num2str(double(vval(23))),'HorizontalAlignment','center')
三、碎碎念
1、句末加“;”,该句不输出
a = 2, b = 3, c = 8;
命令行窗口输出为 a = 2, b = 3,仅c = 8该句不输出
2、两个“%”+空格即可划分模块;单个%为注释
效果如图
3、clc,clf,clear
clc 清除命令行(之前运行的东西还在);
clf 清空图表;
clear 清空工作区;
命令行键盘方向键“↑”扔可调用之前命令;
4、grid on 添加网格线、hold on 在原图上继续添加曲线
5、limit取极限、inf无穷大
6、 num2str转换成字符表示;strcat字符拼接
outstr = strcat('蹦极运动员自由下落的极限速度是:',num2str(limval));
disp(outstr)
7、disp单变量显示
disp(outstr)
8、save\load
save(‘Bengjiathlete.mat’)保存当前文件的变量;
load(‘Bengjiathlete.mat’)重新调出当前文件的变量,或者当前文件夹目录双击可以自动执行load
9、快捷键Alt+Backspace恢复之前删除(直接Ctrl+Z即可)
10、Ctrl+→向右移动一个单词
11、格式显示format short/long
12、help(‘’)很重要
尾末留坑(那些孩子遇到未解决的疑惑)
①绘图弹出窗口“交互式探查绘图数据”
②以实时脚本打开
③工作区变量会出现显示1x1sym的形式