There are many simulators out there, some of which are quite good, for integrating ordinary differential equation (ODE) models and saving/displaying the results as a function of time. Some do much more than this. One limitation of these simulators is that they do not interact directly with Igor, i.e. data needs to be saved/loaded/imported/exported. I propose a simulator framework for Igor to integrate arbitrary numbers of ODEs and output the results directly to Igor waves for graph/table display. Furthermore, for rapid model prototyping, the model should be extensible and the model parameters easily editable within Igor. The simulator could use loop-based, fixed-step numeric integration, or IntegrateODE, depending on what is faster or more convenient for a given model.
I have solved pieces of this problem many times, and before I re-invent the wheel again, I would like to solicit feedback/comments/code from anyone who is doing or has done something similar. Eventually I would like this project to culminate in a standard package for release on IgorExchange.
I have developed some code for simulating chemical reaction systems using IntegrateODE. The input is through a table containing the chemical reactions, initial concentrations, and rate coefficients. Currently, up to 3 reactants and 2 products are possible per reaction, and I have run a model with around 50 reactions successfully. I've written procedures to analyze the chemical mechanism, develop the differential equations, modify the igor procedures, and run the model. Igor waves for each model component are created automatically.
Probably the worst developed part of this is that I could not find a more elegant way of updating the differential equations rather than to automatically modify igor procedures by first removing them from the experiment, dynamically editing them as text files (in igor), and reloading them into the experiment. I would be very interested to see if there are better solutions to this.
I am willing to share the code, however, I don't know if it currently is worth a project here, maybe after some testing by someone else. Please email me if you're interested in trying it out and giving me some feedback.
Harald
....I am willing to share the code, however, I don't know if it currently is worth a project here, maybe after some testing by someone else. Please email me if you're interested in trying it out and giving me some feedback.
My impression (rightly or wrongly) is, folks seem to see this forum as a place only to post "completed" projects. That is a bit of a shame - I think the user community might benefit by seeing ongoing developments too (not to mention the PR side for WaveMetrics). In this regard, a few nice things about this forum are that a) you are free to post a project with any number of "buyer beware" notices, b) you can open up your project to other being able to edit it directly (rather than through email exchanges of files) via SVN, c) you can use this forum rather than email to track and post notices (bug reports, issues, ...), and d) you have an archival place where the development on the files are stored off residence from your own hard drive (should it crash, as they are apt to do).
In other words, I would say don't be bashful ... post what you have and see what feedback comes from it.
--
J. J. Weimer
Chemistry / Chemical & Materials Engineering, UAH
I'm not familiar with other simulator packages, and would love to see a listing of features.
I take it that the user can write out the derivatives and these get transferred into the code, maybe the way Harald suggested?
A feature I think would be great would be to build in some sensitivity analysis methods. A closely related general package that would be nice to have would be a dynamical systems toolbox.
I'd love to contribute but can't really promise to have time for this these days.
function differentials(pw,xx,yw,dydx) Wave pw // parameter wave (input) Variable xx // x value at which to calculate derivatives Wave yw // wave containing y[i] (input) Wave dydx // wave to receive dy[i]/dx (output)
There is something to be said for being able to write out differential equations like you would on paper or on a blackboard. Some programs (like XPP by Bard Ermentrout) allow you to do this, and that is one of their greatest strengths.
Perhaps such a parser would be easier to write when it would be asked to parse within an Igor notebook for code defined as ...
d: x, y
i: t
v: tmax = 1000
v: dt = 0.001
b: x[0] = 1
b: y[0] = 3
e: x = cos(x) + 2*y
e: y = sin(y) - 3*x
Where d (dependent), i (independent), v (integration variables), b (boundary conditions), and e (equations) set the framework.
I can envision a special "IDE Editor" notebook/scratchpad embedded in a control panel wherein one would write such as the above, click a button, and have the code appear below (for final tweaking and then copy/pasting).
(NOTE ADDED as AFTERTHOUGHT)
In fact, as I think about this, everything above the e: lines (number of dependent parameters; the integration start, limit, and step size; and the boundary conditions on the parameters) can be set by SetVariable type controls in the panel that has the embedded notebook. A "scratchpad" area would be used to enter the equations.
Just my thoughts.
--
J. J. Weimer
Chemistry / Chemical & Materials Engineering, UAH
I've added the procedure file kinsim_ex.ipf to the project page http://www.igorexchange.com/node/1333. Feel free to test it and let me know about any issues you may have with it.
As said in the description, load the package and execute "Display_help" from the "KinSim" menu. Make sure to save the experiment in a separate folder before any further steps.
It's written in a pretty crude form, but does the job for my application, so it would probably be hard to modify.
I hope it'll be useful for some...
Harald
Your release notes include the statement, "Doesn't seem to handle well with very stiff ODE's, which could possibly be improved by using a different integrator (Gear?)."
IntegrateODE/M=3 uses the BDF (Backwards Differentiation Formula, or Gear method). It was included specifically for stiff systems. Have you tried it? It shouldn't be too hard to include a selector for the integration method.
John,
I've tried it but wasn't successful, there were some extra parameters that needed to be set and I wasn't sure what they were. So instead, I combined fast reactions into slower ones, basically treating them in the steady-state assumption, which made the model faster overall. But I'm sure the Gear method would be useful if you can't make these assumptions, I'd be happy to get some examples where it would make sense to use and get an explanation on how.
Harald
Jeff Weimer (jjweimer) and I have started a project called "ODEpad" which incorporates some of these ideas. We will let the list know when we have something worth releasing.
Following up on this old thread, as of Sep 2019: We have recently released a free Igor-based easy-to-use kinetic solver (KinSim v4) for teaching and research, which a lot of people are using now:
At least 15 research papers have been published using KinSim so far. We have verified that KinSim runs as fast as the FACSIMILE commercial solver which some people in our community use.
#include <FitODE>
I don't think they get much use, and they're kind of old. I'm sure that having you use them would improve them!
John Weeks
WaveMetrics, Inc.
support@wavemetrics.com
January 26, 2010 at 12:07 pm - Permalink
Probably the worst developed part of this is that I could not find a more elegant way of updating the differential equations rather than to automatically modify igor procedures by first removing them from the experiment, dynamically editing them as text files (in igor), and reloading them into the experiment. I would be very interested to see if there are better solutions to this.
I am willing to share the code, however, I don't know if it currently is worth a project here, maybe after some testing by someone else. Please email me if you're interested in trying it out and giving me some feedback.
Harald
January 26, 2010 at 12:24 pm - Permalink
My impression (rightly or wrongly) is, folks seem to see this forum as a place only to post "completed" projects. That is a bit of a shame - I think the user community might benefit by seeing ongoing developments too (not to mention the PR side for WaveMetrics). In this regard, a few nice things about this forum are that a) you are free to post a project with any number of "buyer beware" notices, b) you can open up your project to other being able to edit it directly (rather than through email exchanges of files) via SVN, c) you can use this forum rather than email to track and post notices (bug reports, issues, ...), and d) you have an archival place where the development on the files are stored off residence from your own hard drive (should it crash, as they are apt to do).
In other words, I would say don't be bashful ... post what you have and see what feedback comes from it.
--
J. J. Weimer
Chemistry / Chemical & Materials Engineering, UAH
January 26, 2010 at 03:51 pm - Permalink
I take it that the user can write out the derivatives and these get transferred into the code, maybe the way Harald suggested?
A feature I think would be great would be to build in some sensitivity analysis methods. A closely related general package that would be nice to have would be a dynamical systems toolbox.
I'd love to contribute but can't really promise to have time for this these days.
January 26, 2010 at 09:00 pm - Permalink
or this:
tmax=100000
dt=0.0001
x(0)=1
y(0)=3
dx/dt = cos(x) + 2y
dy/dt = sin(y) - 3x
into this (for dumb integration):
variable i,xx=1,yy=3
make /o/n=(tmax/dt) xHistory=nan,yHistory=nan
for(i=0;i<100000;i+=1)
xx += dt*(cos(xx)+2*yy)
yy += dt*(sin(yy)-3*xx)
xHistory[i]=xx
yHistory[i]=yy
endfor
or this (for integrateODE based integration):
Wave pw // parameter wave (input)
Variable xx // x value at which to calculate derivatives
Wave yw // wave containing y[i] (input)
Wave dydx // wave to receive dy[i]/dx (output)
dydx[0] = cos(yw[0])+2*yw[1])
dydx[1] = sin(yw[1])-3*yw[0]
return 0
end
There is something to be said for being able to write out differential equations like you would on paper or on a blackboard. Some programs (like XPP by Bard Ermentrout) allow you to do this, and that is one of their greatest strengths.
Rick
January 27, 2010 at 09:27 am - Permalink
Perhaps such a parser would be easier to write when it would be asked to parse within an Igor notebook for code defined as ...
d: x, y
i: t
v: tmax = 1000
v: dt = 0.001
b: x[0] = 1
b: y[0] = 3
e: x = cos(x) + 2*y
e: y = sin(y) - 3*x
Where d (dependent), i (independent), v (integration variables), b (boundary conditions), and e (equations) set the framework.
I can envision a special "IDE Editor" notebook/scratchpad embedded in a control panel wherein one would write such as the above, click a button, and have the code appear below (for final tweaking and then copy/pasting).
(NOTE ADDED as AFTERTHOUGHT)
In fact, as I think about this, everything above the e: lines (number of dependent parameters; the integration start, limit, and step size; and the boundary conditions on the parameters) can be set by SetVariable type controls in the panel that has the embedded notebook. A "scratchpad" area would be used to enter the equations.
Just my thoughts.
--
J. J. Weimer
Chemistry / Chemical & Materials Engineering, UAH
January 27, 2010 at 10:27 am - Permalink
As said in the description, load the package and execute "Display_help" from the "KinSim" menu. Make sure to save the experiment in a separate folder before any further steps.
It's written in a pretty crude form, but does the job for my application, so it would probably be hard to modify.
I hope it'll be useful for some...
Harald
January 27, 2010 at 03:56 pm - Permalink
IntegrateODE/M=3 uses the BDF (Backwards Differentiation Formula, or Gear method). It was included specifically for stiff systems. Have you tried it? It shouldn't be too hard to include a selector for the integration method.
John Weeks
WaveMetrics, Inc.
support@wavemetrics.com
January 28, 2010 at 09:36 am - Permalink
I've tried it but wasn't successful, there were some extra parameters that needed to be set and I wasn't sure what they were. So instead, I combined fast reactions into slower ones, basically treating them in the steady-state assumption, which made the model faster overall. But I'm sure the Gear method would be useful if you can't make these assumptions, I'd be happy to get some examples where it would make sense to use and get an explanation on how.
Harald
January 28, 2010 at 10:07 pm - Permalink
Rick
January 31, 2010 at 02:29 pm - Permalink
Following up on this old thread, as of Sep 2019: We have recently released a free Igor-based easy-to-use kinetic solver (KinSim v4) for teaching and research, which a lot of people are using now:
https://pubs.acs.org/doi/10.1021/acs.jchemed.9b00033 (paper describing it for teaching purposes)
http://tinyurl.com/kinsim-release (Pls always download the code from here, as there are multiple updates since the paper version)
At least 15 research papers have been published using KinSim so far. We have verified that KinSim runs as fast as the FACSIMILE commercial solver which some people in our community use.
Regards,
-Jose L Jimenez <jose.jimenez@colorado.edu>
September 24, 2019 at 01:59 pm - Permalink
Thanks for the update. I am amazed!
September 24, 2019 at 02:05 pm - Permalink