版块导航: 正在加载中...

登录注册

应《网络安全法》要求，自2017年10月1日起，未进行实名认证将不得使用互联网跟帖服务。为保障您的帐号能够正常使用，请尽快对帐号进行手机号验证，感谢您的理解与支持！

24小时热门版块排行榜

南方科技大学公共卫生及应急管理学院2026级博士研究生招生报考通知（长期有效）

返回列表

zhaoqian59

捐助贵宾 (小有名气)

应助: 1 (幼儿园)
金币: 1514
散金: 8
红花: 5
帖子: 194
在线: 57.4小时
虫号: 2656296
注册: 2013-09-15
性别: MM
专业: 爆炸与冲击动力学

[求助] 这个程序是不是哪里有问题？初始值跟实验值是一样的，但是模拟结果不对已有1人参与

f(x)
pl=1000;
phyd=110000;
P0=29.39E9;
a0=0.0341;
v=1.5;

x1=x(1);
x2=x(2);
pg=P0*(a0/x1)^(3*v);
y=[x2 (1/pl*(pg-phyd)-3/2*x2^2)/x1]';

Ts = 0.001; % 仿真步长
Tn = 1;    % 仿真终止时间
t  = 0:Ts:Tn; % 仿真时间范围
N  = length(t);
x=[0.034 0]'; % 系统初值：x=[a a的导数]
%% Variable Declaration and Initialization
a=0.034*ones(1,N);
a_dt=zeros(1,N);
%x = [a a_d]'; % intial states
%% Simulation Iteration
for n = 1:N
a(n)=x(1);
a_dt(n)=x(2);
k1 = fx(x);
k2 = fx(x+Ts/2*k1);
k3 = fx(x+Ts/2*k2);
k4 = fx( x+Ts *k3);
x = x + Ts/6 * ( k1 + 2*k2 + 2*k3 + k4 );
end
%% Plot Response
plot(t,a)
xlabel('t/(s)')
ylabel('a')
是一个二阶非线性常微分方程

回复此楼

» 收录本帖的淘帖专辑推荐

matlab编程绘图

» 猜你喜欢

三甲基碘化亚砜的氧化反应已经有4人回复
请问下大家为什么这个铃木偶联几乎不反应呢已经有5人回复
请问有评职称，把科研教学业绩算分排序的高校吗已经有5人回复
孩子确诊有中度注意力缺陷已经有12人回复
2025冷门绝学什么时候出结果已经有3人回复
天津工业大学郑柳春团队欢迎化学化工、高分子化学或有机合成方向的博士生和硕士生加入已经有4人回复
康复大学泰山学者周祺惠团队招收博士研究生已经有6人回复
AI论文写作工具：是科研加速器还是学术作弊器？已经有3人回复
论文投稿，期刊推荐已经有4人回复
硕士和导师闹得不愉快已经有13人回复

爆炸力学

1楼2015-08-11 09:29:52

已阅回复此楼关注TA 给TA发消息送TA红花 TA的回帖

月只蓝

主管区长 (职业作家)

专家经验: +1059
应助: 1712 (讲师)
贵宾: 8.888
金币: 68120.7
散金: 1938
红花: 443
沙发: 4
帖子: 4373
在线: 3291.4小时
虫号: 1122189
注册: 2010-10-14
专业: 宇宙学
管辖: 计算模拟区

【答案】应助回帖

感谢参与，应助指数 +1

楼主自己编写的四阶龙格库塔来求解常微分方程？
能给出方程具体形式么，初值，以及各常数？

赞一下

回复此楼

MATLAB、MS小问题、普通问题请发帖求助！时间精力有限，恕不接受无偿私信求助。

2楼2015-08-11 10:37:11

已阅回复此楼关注TA 给TA发消息送TA红花 TA的回帖

zhaoqian59

捐助贵宾 (小有名气)

应助: 1 (幼儿园)
金币: 1514
散金: 8
红花: 5
帖子: 194
在线: 57.4小时
虫号: 2656296
注册: 2013-09-15
性别: MM
专业: 爆炸与冲击动力学

引用回帖:

2楼: Originally posted by 月只蓝 at 2015-08-11 10:37:11
楼主自己编写的四阶龙格库塔来求解常微分方程？
能给出方程具体形式么，初值，以及各常数？

aa''+1.5a'=1/1000*(P-Pwater)
P=P0*(a0/a)^v
a是关于t的函数，Pwater=1.1*100000，P0=29E9，a0=0.0341，v=3.谢谢你啦~

赞一下

回复此楼

爆炸力学

3楼2015-08-11 15:03:26

已阅回复此楼关注TA 给TA发消息送TA红花 TA的回帖

月只蓝

主管区长 (职业作家)

专家经验: +1059
应助: 1712 (讲师)
贵宾: 8.888
金币: 68120.7
散金: 1938
红花: 443
沙发: 4
帖子: 4373
在线: 3291.4小时
虫号: 1122189
注册: 2010-10-14
专业: 宇宙学
管辖: 计算模拟区

引用回帖:

3楼: Originally posted by zhaoqian59 at 2015-08-11 15:03:26
aa''+1.5a'=1/1000*(P-Pwater)
P=P0*(a0/a)^v
a是关于t的函数，Pwater=1.1*100000，P0=29E9，a0=0.0341，v=3.谢谢你啦~...

t=0时，对应a a'的初值呢？

赞一下

回复此楼

MATLAB、MS小问题、普通问题请发帖求助！时间精力有限，恕不接受无偿私信求助。

4楼2015-08-11 15:08:12

已阅回复此楼关注TA 给TA发消息送TA红花 TA的回帖

zhaoqian59

捐助贵宾 (小有名气)

应助: 1 (幼儿园)
金币: 1514
散金: 8
红花: 5
帖子: 194
在线: 57.4小时
虫号: 2656296
注册: 2013-09-15
性别: MM
专业: 爆炸与冲击动力学

引用回帖:

4楼: Originally posted by 月只蓝 at 2015-08-11 15:08:12
t=0时，对应a a'的初值呢？...

t=0时，a(t)=0.0341,a'=0

赞一下

回复此楼

爆炸力学

5楼2015-08-11 15:14:44

已阅回复此楼关注TA 给TA发消息送TA红花 TA的回帖

月只蓝

主管区长 (职业作家)

专家经验: +1059
应助: 1712 (讲师)
贵宾: 8.888
金币: 68120.7
散金: 1938
红花: 443
沙发: 4
帖子: 4373
在线: 3291.4小时
虫号: 1122189
注册: 2010-10-14
专业: 宇宙学
管辖: 计算模拟区

【答案】应助回帖

引用回帖:

3楼: Originally posted by zhaoqian59 at 2015-08-11 15:03:26
aa''+1.5a'=1/1000*(P-Pwater)
P=P0*(a0/a)^v
a是关于t的函数，Pwater=1.1*100000，P0=29E9，a0=0.0341，v=3.谢谢你啦~...

记 a=u1，a'=u2；
原二阶方程可将阶为方程组：
u1'=u2;
u2'=1/u1* ( (P-Pwater)/1000 - 1.5*u2 );
其中，P=P0*(a0/u1)^v;

MATLAB代码如下：

CODE:

function solve_odes_two_oder

clear all;clc

u0=[0.0341 0];

tspan=linspace(0,0.1,100);

[t u]=ode45(@odefun,tspan,u0);

figure(1)

plot(t,u(:,1),'bo--',t,u(:,2),'r-*'),legend('a',' da/dt ')

function f=odefun(t,u)

Pwater=1.1*100000;

P0=29E9;

a0=0.0341;

v=3;

P=P0*(a0/u(1))^v;

f(1)=u(2);

f(2)=1/u(1)*   (  (P-Pwater)/1000 - 1.5*u(2)     );

f=f';

结果中，a和a'迅速增大，即便在0.1的时间之内，见附图1。

附图1.png

赞一下

回复此楼

MATLAB、MS小问题、普通问题请发帖求助！时间精力有限，恕不接受无偿私信求助。

6楼2015-08-11 15:47:29

已阅回复此楼关注TA 给TA发消息送TA红花 TA的回帖

zhaoqian59

捐助贵宾 (小有名气)

应助: 1 (幼儿园)
金币: 1514
散金: 8
红花: 5
帖子: 194
在线: 57.4小时
虫号: 2656296
注册: 2013-09-15
性别: MM
专业: 爆炸与冲击动力学

引用回帖:

6楼: Originally posted by 月只蓝 at 2015-08-11 15:47:29
记 a=u1，a'=u2；
原二阶方程可将阶为方程组：
u1'=u2;
u2'=1/u1* ( (P-Pwater)/1000 - 1.5*u2 );
其中，P=P0*(a0/u1)^v;

MATLAB代码如下：

function solve_odes_two_oder
clear all;clc
u0=;
...

复制过去之后运行不出来.....

赞一下

回复此楼

爆炸力学

7楼2015-08-12 09:40:23

已阅回复此楼关注TA 给TA发消息送TA红花 TA的回帖