如何在ode45中插入条件if ?

2次查看(最近30天)

西蒙Maucci 2021年2月1日

评论道: 明星黾 2021年2月1日

答:接受明星黾

你好，我想模拟一个漏洞百出的集成火灾模型。我使用的是Ode45，它工作得很好，但我想插入一个If条件，以便当y值小于特定值(-50)时，y必须重置为-65。我不知道如何实现它。有什么建议吗?非常感谢。

                    x0 = -70;
                   
                    t0 = 0;
                   
                    dt = 0.1;
                   
                    达峰时间= 2;
                   
                    tvec = t0: dt:最高温度;
                   
                    (t, x) =数值(@test_ode tvec, x0)
                   
                    情节(t, x)
                   
                    函数[dydt] = test_ode (t, y)
                   
                    dydt = -1/10 * y + 13
                   
                    结束

0评论
显示隐藏1以上的评论

登录评论。

登录回答这个问题。

接受的答案

明星黾 2021年2月1日

0
链接

直接链接到这个答案

//www.tatmou.com/matlabcentral/answers/732928-how-to-insert-an-if-conditional-in-ode45#answer_611808

试试这个:

                        x0 = -70;
                       
                        t0 = 0;
                       
                        dt = 0.1;
                       
                        达峰时间= 20;
                       
                        % tvec = t0: dt:最高温度;
                       
                        Tvec = [t0 tmax];
                       
                        选择= odeset (“事件”, @depolEvents);
                       
                        为k = 1:10
                       
                        [t{k}，x{k}]=ode45(@test_ode, tvec, x0, opts);
                       
                        x0 = -65;
                       
                        Tvec = [t{k}(end) t{k}(end)+tmax]
                       
                        % tvec = t{k}(end): dt: t{k}(end)+tmax;
                       
                        结束
                       
                        电视=猫(1 t {:});
                       
                        十五=猫(1,x {:});
                       
                        数字
                       
                        情节(电视、十五)
                       
                        网格
                       
                        函数[dydt] = test_ode (t, y)
                       
                        dydt = -1/10 * y + 13;
                       
                        结束
                       
                        函数(位置、isterminal方向)= depolEvents (t, y)
                       
                        Position = y(1) + 50;%我们想要的值为零
                       
                        isterminal = 1;%停止集成
                       
                        方向= 1;0可以从任意方向接近
                       
                        结束

这使用了技术歌唱活动的位置停止在循环的下一个增量中重新启动的集成。

如果它不像你写的那样做(因为我可能不完全理解你想做什么)，用它来试验一下，得到你想要的结果。

4评论
显示隐藏 3年长的评论

明星黾 2021年2月1日

我的荣幸!

我发布的代码运行没有错误 R2020b ．我在发布之前测试了一下。

我不知道是什么导致了这个问题，因为“左边”是单元格数组，所以可以是任何东西，或任何大小的向量或矩阵。的组件 “t” 和 “x” 在代码中是(或应该是)列向量，在测试时也是。

当我使用发布的代码运行它时，它生成了如下图:

编辑- (2021年2月1日19:40)

我还是不明白你的单元格阵列有什么问题。

这是一种避免它们的方法:

                             电视= [];
                            
                             十五= [];
                            
                             为k = 1:10
                            
                             [t,x]=ode45(@test_ode, tvec, x0, opts);
                            
                             x0 = -65;
                            
                             Tvec = [t(end) t(end)+tmax]
                            
                             电视=[电视;t];
                            
                             十五=[十五;x];
                            
                             结束
                            
                             数字
                            
                             情节(电视、十五)
                            
                             网格

用它代替，然后删除这段完整的代码:

                             为k = 1:10
                            
                             [t{k}，x{k}]=ode45(@test_ode, tvec, x0, opts);
                            
                             x0 = -65;
                            
                             Tvec = [t{k}(end) t{k}(end)+tmax]
                            
                             结束
                            
                             电视=猫(1 t {:});
                            
                             十五=猫(1,x {:});

因为部分替换会抛出错误，所以使用单元格数组删除块，并用这个编辑过的注释中的第一个块替换它。

修订后的代码如下:

                             x0 = -70;
                            
                             t0 = 0;
                            
                             dt = 0.1;
                            
                             达峰时间= 20;
                            
                             Tvec = [t0 tmax];
                            
                             选择= odeset (“事件”, @depolEvents);
                            
                             电视= [];
                            
                             十五= [];
                            
                             为k = 1:10
                            
                             [t,x]=ode45(@test_ode, tvec, x0, opts);
                            
                             x0 = -65;
                            
                             Tvec = [t(end) t(end)+tmax]
                            
                             电视=[电视;t];
                            
                             十五=[十五;x];
                            
                             结束
                            
                             数字
                            
                             情节(电视、十五)
                            
                             网格
                            
                             函数[dydt] = test_ode (t, y)
                            
                             dydt = -1/10 * y + 13;
                            
                             结束
                            
                             函数(位置、isterminal方向)= depolEvents (t, y)
                            
                             Position = y(1) + 50;%我们想要的值为零
                            
                             isterminal = 1;%停止集成
                            
                             方向= 1;0可以从任意方向接近
                            
                             结束