This is by splitting your discontinuous domain to its three chunks for which it is continuous in this case from 0-15, 15-25 and 25-30 seconds. With that said, as 'annoying' as it is (and trust be I have been in the same place), the only way to tackle the problem in MATLAB is by the way that you have done. Unfortunately matlab does not have any standard solver that handle discontinuities, and even though to the best of my knowledge this is not directly mentioned in MATLAB's documentation, it is something that is mentioned in this excellent article that outlines and compares the capabilities of MATLAB's ODE solver suits with those of other ODE suits such as Python, Julia, etc. You will notice that this input is discontinuous, and thus non-differentiable.
#MATLAB ODE45 CODE#
I was unable to run your code as it must be missing some input, but I see that your input to the differential equations is the current $I$ which you control as a user, and it looks something like this:
![matlab ode45 matlab ode45](https://www.mathworks.com/help/examples/matlab/win64/ODEWithTimeDependentTermsExample_01.png)
Gbarl=0.003 % mS/cm^2 Leakage conductanceÄydt = [((1/Cm)*(I(chunk)-(INa+IK+Il))) % Normal caseįunction = DE(varargin) Title('Voltage Change for Hodgkin-Huxley Model') Y = V(end, :) % Final position is initial value for next interval H=alphah(V)/(alphah(V)+betah(V)) % Initial h-value N=alphan(V)/(alphan(V)+betan(V)) % Initial n-value M=alpham(V)/(alpham(V)+betam(V)) % Initial m-value % Plotting purposes (set I(idx) equal to last value of I) H1=alphah(V)/(alphah(V)+betah(V)) % Initial h-value N1=alphan(V)/(alphan(V)+betan(V)) % Initial n-value M1=alpham(V)/(alpham(V)+betam(V)) % Initial m-value So is this correct? Can I keep my way of approaching the problem? This means that documentation suggests a method which is driving a numerical method outside the specified limits.
![matlab ode45 matlab ode45](https://i.stack.imgur.com/sYHou.png)
The Dormand Prince Runge Kutta integrator with step size control is designed to operator on differentiable functions. Someone has suggested that it is not correct because the example code in the documentation of ODE45 uses INTERP1 to calculate a parameter in the function to be calculated. However, I am aware that under ODE45 for MATLAB there is a section for time-dependent terms.
![matlab ode45 matlab ode45](https://i.ytimg.com/vi/fx3bl4oA_0U/hqdefault.jpg)
This code that I'm posting works: I integrate in chunks, as many as I define at the beginning, and depending on the number n of timesteps in which current I changes its value (which is known because the user defines it), I will call ode45 n times, every time using the last values of previous iteration as starting values. I have tried in different ways to see what happens to voltage V and gating conductances m, n and h when, at time step x, current I switched from 0 to 0.1, and then at time step x + n it gets back to 0.