通过在Matlab中计算ODE函数中的ODE值来检测稳态 [英] Detecting steady-state by calculating ODE values in ODE function in Matlab

问题描述

我有一个ODE方程系统，我想解决，但有一个棘手的部分，当系统达到稳态时，我想改变一个（或多个）参数的值。例如，请考虑以下内容：

I have a system of ODE equations which I want to solve, but there is a tricky part that when the system reaches steady-state, I would like to change the value of one (or more) parameters. For example, consider the following:

function dydt = diff(t,x,params)
   F = params(1);
   G = params(2);
   dydt = zeros(2,1);
   dydt(1) = F*x(1) - G*x(1)*x(2);
   dydt(2) = (F-G)*x(2);
end

我希望我的代码工作，以便在系统达到稳定时例如，状态，F的值变为10并且G的值变为2。我想通过使用，例如，

I would like my code to work such that when the system has reached steady-state, the value of F is changed to 10 and the value of G is changed to 2, for example. I was thinking of detecting the values of dydt(1) and dydt(2) by using, for example,

if norm(dydt)<1
   F = 10;
   G = 2;
end

我如何在Matlab中为ODE表达式做到这一点？如果我把这条if条件放在ODE表达式之前，那么dydt的值将始终为零。但是如果我在ODE表达式之后放置这个If条件，那么If条件将没有用来纠正ODE表达式。

How do I do that for ODE expression in Matlab? If I put this if condition before the ODE expression, I the value of dydt will always be zero. But if I put this If condition after the ODE expression, the If conditions will have no use to correct the ODE expression.

谢谢！

推荐答案

假设ODE的参数是固定的，并且不依赖于系统的状态。你要做的是模拟分段连续ODE。这类问题通常通过活动地点解决（假设存在解决方案）。您需要在关键点停止模拟，更改参数，然后使用与上一个条件相同的初始条件重新启动新模拟。

A parameter of an ODE is assumed to be fixed and does not depend on the state of the system. What you're trying to do is simulate a piecewise continuous ODE. This sort of problem is usually solved with event location (assuming a solution exists). You need to stop the simulation at the key point, change your parameters, and restart a new simulation with initial conditions the same as those at the end of the previous.

这里是ODE功能的一个小例子。我不知道您的初始条件，初始参数值或其他设置，所以这只是为了证明这个方案：

Here's a little example with your ODE function. I don't know your initial conditions, initial parameter values, or other settings, so this is just to demonstrate this scheme:

function eventsdemo
params = [-1.5 1];
tspan = [0 10];
x0 = [1;1];
opts = odeset('Events',@events);
[t1,x1] = ode45(@(t,x)f(t,x,params),tspan,x0,opts); % Simulate with events

% Change parameters, set initial conditions based on end of previous run
params = [1.5 1];
x0 = x1(end,:);
tspan = [t1(end) 10];
[t2,x2] = ode45(@(t,x)f(t,x,params),tspan,x0); % Simulate again

% Concatenate results, removing duplicate points
t = [t1;t2(2:end)];
x = [x1;x2(2:end,:)];

figure;
plot(t,x);
hold on;
plot(t2(1),x2(1,:),'k*'); % Plot event location


function dxdt = f(t,x,params) %#ok<INUSL>
F = params(1);
G = params(2);
dxdt = [F*x(1) - G*x(1)*x(2);
        (F-G)*x(2)];


function [value,isterminal,direction] = events(t,x) %#ok<INUSL>
value = norm(x)-1e-3; % Don't try to detect exact zero for asymptotics
isterminal = true;
direction = -1;

在这种情况下，解决方案渐近地接近（0,0）处的稳定定点。重要的是不要尝试完全检测（0,0）解决方案可能永远不会到达那一点。您可以使用小容差。但是，根据您的系统，更改参数后，选择此容差可能会影响行为。您也可以考虑将此问题重新表述为边界值问题。我不知道你正在尝试用这个系统做什么，所以我不能说太多其他的东西（这可能是这个网站的主题）。

In this case the solution asymptotically approaches a stable fixed point at (0,0). It's important to not try to detect (0,0) exactly as the solution may never reach that point. You can instead use a small tolerance. However, depending on your system, it's possible that choosing this tolerance could impact the behavior after you change parameters. You could also consider reformulating this problem as a boundary value problem. I don't know what you're trying to do with this system so I can't say much else (and it would probably be off-topic for this site).

这篇关于通过在Matlab中计算ODE函数中的ODE值来检测稳态的文章就介绍到这了，希望我们推荐的答案对大家有所帮助，也希望大家多多支持IT屋！