I'm trying to solve a system of ODE's with MATLAB's numerical ode solver e.g. ode45. The first thing to do is to implement a function f(t,y) returns the time derivatives of the state vector:
function dydt = myFunc(t,y)%system of ODE's. y is a vector of three statesP = y(1)*y(2);dy1 = P - y(2);dy2 = 4*y(2) +2*y(3)*y(1);dy3 = y(1)+y(3);dydt = [dy1;dy2;dy3];end
Then, the time evolution of y can be solved with the ode23 solver using:
% [0 100] is the time interval in which the solution is sought; [1; 2; 3] is the initial condiiton of states y(1), y(2), and y(3).[t,y] = ode23(@myFunc, [0 100], [1; 2; 3;]);
My problem is that in my function defined above (just an example, not my actual system), the expression of P is not simply y(1)*y(2), but y(1)*y(2)*(_some low pass filter_). This low pass filter is causing the major trouble: I can't find a way to implement it, be it in time domain or in laplace form. Using a filter in laplace form, matlab will return an error saying laplace transform cannot be performed on data type: double.
% with a filter in laplace formP = ilaplace(laplace(y(1)*y(2)) * omega / (s + omega));
Using a time-domain filter, I couldn't even retrieve the historical data of the states y(1) and y(2) in order to do that.
- Using a state-space model instead of solving system of ode's: this can not be done because state space models in matlab only support linear systems but my model is nonlinear.
- Using some vector global variable that could track the trajectory of y(1) and y(2) by appending new values every time myFunc is called: Possibly feasible, but ode23/45 uses variable step sizes and furthermore I'm not sure if myFunc will be called in the sequency of time when running ode23/45. If feasible, how can this be implemented exactly?
- Any other thoughts?
Thanks for your help.