Solucionadores de ode Matlab: mudança de estado e tempo especificado

Estou resolvendo um conjunto de ODEs (dy / dt) em t = 0, todas as condições iniciais t = 0 y_0 = (0,0,0). Posso adicionar algum número aos valores de y em momentos diferentes (por exemplo, em t = 10, y1 deve ser adicionado a esse número; em t = 20, y2 deve ser adicionado a esse número, etc.) e resolver as equações?

Ok, como Simon McKenzie diz que realmente precisamos de mais informações sobre o seu problema urgente , mas acho que posso ajudar. Pelo que você nos deu, eu vou assumir que você tem uma function myfun que você passa para algo como ode45

 y_0 = [0,0,0]; % here Tfinal is the time in seconds that you want to simulate to % and specify the tspan so that you will get a solution at each whole % number of seconds, I believe [t,y] = ode45(@myfun,[0:0.1:Tfinal],y_0); 

Em algum lugar você definiu sua function, aqui eu chamarei de myfun

 function dy = myfun(t,y) % here let's check to see if it's time to add your offsets % still, you may want to put a little fudge factor for the time, but % you'll have to experiment, I'll set it up though EPS = 1e-12; if( t < 10 + EPS || t > 10 - EPS ) y(1) = y(1) + 10; . . . . % this only makes sense if you then use y(1) in the compuatation 

Caso contrário, basta adicionar o deslocamento ao vetor da solução retornada, ou seja,

 idx10 = find( t == 10 ); % should be there since it's in the tspan y(idx10:end) = y(idx10:end) + 10; % I guess you add it from that point on? 

Inserindo grandes descontinuidades em seu ODE da maneira que você sugere (e da maneira ilustrada por @macduff) pode levar a menos precisão e tempos de computação mais longos (especialmente com ode45ode15s pode ser uma opção melhor ou pelo menos se certificar de que seu absoluto e relativo tolerâncias são adequadas). Você produziu efetivamente um sistema muito rígido . Se você quiser adicionar algum número aos ODEs iniciando em um horário específico, lembre-se de que o solucionador só avalia essas equações em pontos específicos no tempo de sua própria escolha. (Não se deixe enganar pelo fato de poder obter saídas de tamanho de passo fixas especificando tspan como mais de dois elementos – todos os solucionadores do Matlab são solucionadores de tamanho de passo variables ​​e escolhem seus passos verdadeiros com base em critérios de erro.)

Uma opção melhor é integrar o sistema por partes e append as saídas resultantes de cada execução em conjunto:

 % t = 0 to t = 10, pass parameter a = 0 to add to ODEs a = 0; tspan = [0 10]; [T,Y] = ode45(@(t,y)myfun(t,y,a),tspan,y_0); % t = 10 to t = 20, pass parameter a = 10 to add to ODEs a = 10; [t,y] = ode45(@(t,y)myfun(t,y,a),tspan+T(end),Y(end,:)); T = [T;t(2:end)]; Y = [Y;y(2:end,:)]; % t = 20 to t = 30, pass parameter a = 20 to add to ODEs a = 20; [t,y] = ode45(@(t,y)myfun(t,y,a),tspan+T(end),Y(end,:)); T = [T;t(2:end)]; Y = [Y;y(2:end,:)]; 

O editor do Matlab pode reclamar que o array T e Y não está sendo pré-alocado e / ou em crescimento, mas, neste caso, não há problema, já que eles estão crescendo em grandes partes apenas algumas vezes. Como alternativa, se você quiser tamanhos de etapas de saída fixos, poderá fazer isso:

 dt = 0.01; T = 0:dt:30; Y = zeros(length(T),length(y_0)); % t = 0 to t = 10, pass parameter a = 0 to add to ODEs a = 0; [~,Y(1:10/dt+1,:)] = ode45(@(t,y)myfun(t,y,a),T(1:10/dt+1),y_0); % t = 10 to t = 20, pass parameter a = 10 to add to ODEs a = 10; [~,Y(10/dt+1:20/dt+1,:)] = ode45(@(t,y)myfun(t,y,a),T(10/dt+1:20/dt+1),Y(10/dt+1,:)); % t = 20 to t = 30, pass parameter a = 20 to add to ODEs a = 20; [~,Y(20/dt+1:end,:)] = ode45(@(t,y)myfun(t,y,a),T(20/dt+1:end),Y(20/dt+1,:)); 

Pode-se facilmente converter os dois blocos de código acima for laços mais compactos, se desejado.

Em ambos os casos, sua function ODE myfun incorpora o parâmetro da seguinte maneira:

 function ydot = myfun(t,y,a) y(1) = ... % add a however you like ...