Systems of Ordinary Differential Equations
by Dennis Darland
Mail To: pal at dennisdarland dot com
"True laws can only be expressed in differential equations." - Bertrand Russell - Human Knowledge: Its Scope and Limits, Page 166.
Table of Contents
- What omnisode is.
- Where it came from.
- Test Results.
- Known Bugs.
- Bugs Fixed.
- Plans for the Future.
- Items Accomplished.
- Taylor Series - Bibliography.
- Old version of this page. (Before extensive revision on 12/21/2012)
- It is mainly a Ruby program which generates a program to solve a set of one or more ordinary differential equations.
- It currently can generate either c, c++, ruby, Maxima or Maple.
- It uses a specifiable number of terms of the Taylor series of the equations.
- The Taylor series is iteratively computed from recursive relations determined by the equations.
- It runs under either Linux or cygwin, or windows.
- The program is now run by a ruby program tomni.rb which replaces the script files and works on windows as well as linux and cygwin.
- A ".ode" file specifies any problem to be solved.
- A list of ".ode" files are placed in a ".odes" file the name of which is given to tomni.rb.
- You run it simply with "ruby tomni.rb"
- ".odes" files - "easy.odes" and "any.odes" are provided along with many ".ode" files.
- Of course, you must obtain, c, c++, ruby, maxima or maple separately.
- All can be obtained for free except maple.
- Note: where I say analytic solution, I really mean a solution given as a closed-form expression and evaluated at a point numerically.
- It is open source software licensed under the GNU General Public License v3.0
- Note: The peculiar way some of the code is written is due to the way arrays in maxima work (they seem to be call by value). After learning how to use "modedeclare" in maxima, this
no longer seems necessary, but actually the generated code is probably more efficient in time, and not take up much more space.
Omnisode - Where it came from
- I became interested in this method from Y. F. Chang at the University of Nebraska, Lincoln in about 1977-1978.
- Unfortunately, I became ill and was unable to complete my Masters degree.
- I retained an interest in it, even though it was many years before I had a computer actually powerful enough for this method.
- I intermittently remained in contact with Prof. Chang and Prof. George Corliss who had also been at the UNL.
- Prof. Corliss is now at Marquette. I believe Prof. Chang is deceased.
- The latest version of Professors Chang and Corliss's that they gave me is at Atomft311
- The version of Professors Chang and Corliss's that I made minimal changes to for SUSE Linux 13.1 is at Atomft311 for current Linux
- I only changed makefiles and scripts and have Prof. Corliss's permission (9/6/2014) to make it public.
- The Atomft311 User Manual is at Atom311 Manual
- In the Makefile, Pof. Corliss states:
An 80+ page printed ATOMFT User Manual containing 38 examples
of ODE systems is available at no charge from Corliss.
ATOMFT is distributed at no cost, but it carries two obligations:
- If you have trouble, do not just bitch. contact Corliss!
- If you use ATOMFT in your research, send Corliss a copy of any resulting papers or technical reports.
- I suggest, if you have trouble with the makefiles or scripts for linux, you contact me, Dennis J. Darland, before contacting Prof. Corliss.
- There are paths that will need to be modified in the scripts and makefiles, depending on your setup - they generally include 'dennis'
- I am not very familiar with the FORTRAN code.
Omnisode - Documentation
- There is documentation on the ode files. (Revised 1/19/2014)
- There is documentation on the script files used to run omnisode. THESE ARE OBSOLETE. See the README below. (Revised 1/19/2014)
- I reworked the README file to make it more useful. It covers the new tomni.rb whch replaces all the script files in rev0049 (2/11/2014)
- rdoc (Ruby Documentation) for release 0045 of the ode file preprocessor.
- rdoc (Ruby Documentation) for release 0045 of omnisode. (diffeq generator)
- rdoc (Ruby Documentation) for release 0045 of the diffeq post processor.
- rdoc (Ruby Documentation) for release 0045 of diffeq.
Omnisode - Latest Test Results
- Old test results are at here. They are mainly obsolete
- September 21-23, 2014 - Tests with SUSE Linux 13.1 with new Laptop (One I had been using had some hardware issues due to heat.)
What became rev0057 of omnisode.
See Notes below for changes included.
- maple - cases involving Si or erf
Comparisons with other tools - ATOMFT, Maple dsolve and maxima runga-kutta
Comparisons for First Painleve transcendent y''(t) = 6*y(t)*y(t) + t, y(0) = 1, y'(0) = 0
- y(1.1)=87.774060162627567911Maple dsolve solution (numeric/taylorseries)
- y(1.1)=87.7740601626275705(h = 0.0011)Omnisode solutions using different increments using Maple as the target language.
- y(1.1)=87.774060162627570239(h = 0.01)
- y(1.1)=87.77405966362392 (h = 0.1)
- y(1.1)=8.77739E+01Single Precision ATOMFT input gives ATOMFT output
- y(1.1)=8.77740601591047D+01Double Precision ATOMFT input gives ATOMFT output
- y(1.1)=8.77739E+01Single Precision ATOMFT input gives ATOMFT output - equal spacing
- y(1.1)=87.77028865886018Maxima (Runge-Kutta) Input gives Maxima output
- Chang's example 3_1 - Multiple Equations
x''(t) = -5.8 * x(t) * (x(t)*x(t) + y(t) * y(t)) ** (-1.5)
y''(t) = -5.8 * y(t) * (x(t)*x(t) + y(t) * y(t)) ** (-1.5)
x(1) = -1
y(1) = 0
- Maple dsolve solution (numeric/taylorseries)
- Omnisode solutions using increment of 0.01 using maple as the target language.
- Omnisode solutions using increment of 0.0001 using gcc as the target language.
- Single Precision ATOMFT input gives ATOMFT output
According to comments should give:
Omnisode - Known Bugs
- Sometimes computed increment (glob_h) does not give desired accuracy. (this may be impossible - I need to think more about it).
- Problem with comments in the 1st (equation) section
- Seems display interval is required - no default.
- Solution not always displayed at display interval (boundary problem) - Improved - may want to rethink this.
- I think I need to convert systems of higher order equations to systems of 1st order equations - what I tried to do requires too small increment it seems.
Will probably implement this via another preprocessor.
Omnisode - Bugs Fixed
- First ode file in a batch is not copied to results directory. - Fixed.
- Constant to full series power not accurate. (expt_c_sin_new.ode) - Fixed.
- Linear function to full series power not accurate. (expt_lin_sin_new.ode) OK except h determined is too large for desired accuracy.
- Full series function to full series power not accurate. (expt_sin_sin_new.ode) OK except h determined is too large for desired accuracy.
- Fixed missing "tail.html" from release. (in rev0027)
- Fixed memory problem in c and c++
- January 4, 2014 Removed contents of temp subdirectory in scripts to avoid any confusion with previous tests.
- January 4, 2014 Made match test scripts use "if" to terminate after any problem fails - it had been easy to overlook errors.
- January 4, 2014 There were places where I suspect blanks were needed as separators of tokens in equation block. I made preodein.rb add them to temp output.
- January 4, 2014 I had introduced an error in reorganizing the structure of omnisode.rb, which is now much better.
I set the default MAX_TERMS to be 30 - it should be over written later with the value of the ode file anyway.
But accidentally in maxima, I allocated the array before updating with the value of the ode file.
Almost all tests worked - except for sin.ode. The default was 30. sin.ode set it to 40. - Array out of bounds!
I examined every subscript in the program very carefully. I put range checks in. A couple of those were in tricky places.
But none of that helped. Then I figured it out - it was an easy fix.
- January 4, 2014 I made log10 only compute log(10.0) only once - only needed for maxima
- January 4, 2014 Equations no longer must all be on one line each. There was an easy fix in preodein.rb. I did this after the 4th set of tests but this has no other effect.
- January 5, 2014 Added modedeclare for a local variable I missed -for Maxima.
- January 5, 2014 Used alias instead of functions I needed different names for. (Improved performance by about 20% for Maxima)
- January 7, 2014 Cleaned up code - removed dead code, but left methods that would be useful if reimplanting Maple
- January 7, 2014 Made it possible to have comments or blank lines in equation section.
- January 8, 2014 Made it possible to have comments after end of any line anywhere in ode file (except ! which must be on line by itself)
- September 22, 2014 Fixed "/" screw up.
Omnisode - Download
- The current omnisode (0057) was published on the web September 23, 2014.
- I reworked the README file to make it more useful.
- Both tar.gz and zip versions - contents should be identical.
- The latest published version may be obtained at The Sourceforge Download page for Omnisode
- Omnisode download stats as of 8/29/2014
Omnisode - Plans for the future
- Make omnisode, atsgen, atsinc, preodein, and preindent have comments for rdoc/
- Need better error messages on syntax errors is equations.
- Investigate use of interval arithmetic. Now (12/15/2012) think this is probably impractical.
- Ruby support using APFP Almost certainly impractical.
- Perl support - Little reason to do - except I like playing with languages.
- Python support - Little reason to do - except I like playing with languages.
- Note that the support would be for the programs generated, omnisode will remain a Ruby program.
Omnisode - Items accomplished
- May 22, 2012 - Language independent ".ode" files (they define a set of equations to solve)
- May 22, 2012 - Documentation of ode file format
- May 24, 2012 - Capability to add test data to a html table
- May 25-31, 2012 - Links to detail files added to tables generated.
- May 25-31, 2012 - Automatic generation of test scripts.
- May 25-31, 2012 - Automatic indentation of programs generated.
- May 25-31, 2012 - Work files moved to "temp" subdirectory.
- May 25-31, 2012 - Added "ode.over" file to override values in all ".ode" files for batch tests.
- June 5-12, 2012 - Improved significantly handling of systems of (multiple) equations.
- June 12, 2012 - More documentation.
- June 15-18, 2012 - Optimized factorial calculations.
- July 2-10, 2012 Treat constants and linear expressions as special cases of series in order to optimize these common cases. More efficient bothe in time and space. (Mostly Done - except difficulties in a few cases).
- July 11, 2012 Put in warnings for special cases which seem to have problems - or have trouble checking. (Also need to check these cases again).
- August 13, 2012 C++ code generation functioning (Almost c also).
- August 21, 2012 Now C as well as C++.
- September 6, 2012 Made to support negative increment (H).
- September 6, 2012 Added capability to display values calculated only at specified intervals.
- November 7, 2012 Now Ruby code can be generated.
- December 14, 2012 Attempting to control error using estimate of radius of convergence and other factors - still needs much testing.
- December 22, 2012 Revised documentation.
- January 16, 2013 Worked on special cases with difficulties - devised tests for which I had analytic solutions and fixed any problems.
- January 26, 2013 Examined how estimation of singularities is working. Radius works well. Order less accurate.
- January 27, 2013 Added link to ode.over above rest of data in output table.
- January 27, 2013 Devised way to set maximum h smaller than display_interval.
- January 27, 2013 Cleaned up code (Mainly deleting dead (commented out) code).
- January 27, 2013 Cleaned up ode files (Mainly removing obsolete glob_h)
- February 24, 2013 Added support for Windows - just for Ruby.
- April 4, 2013 Support for Si and erf functions (Maple only) released.
- May 15, 2013 Added logic for ratio test for singularities.
- May 20, 2013 Results file now contains adjacent output data of singularities - actual (if location given in ode file), ratio test, three term test, and six term test.
- December 28, 2013 The last few months the main thing I have done is rewrite the singularity code..
The calculations for 2, 3 and six term tests is separated from a function (which I think is not yet quite right) to see if the radius is reasonably approaching a limit.
- December 28, 2013 Also made corrections to printing of times - especially for c, c++ and Maxima. In c and c++ I had missed warning messages due to how I handled stderr in the scripts.
In Maxima there is no function 'trunc' - functions which do not exist in maxima just seem to do nothing - no error. Improved use of singularity calculations.
- December 30, 2013 Improved attempt to see if radius is approaching a limit. Improved use of radius calculated in adjusting increment.
- January 8, 2014 Removed the two relative error columns from output table - the number of correct digits is sufficient.
- January 8, 2014 Added four columns to table for the smallest of the given, ratio, three term and six term tests for singularities.
- January 8, 2014 Changed default values of glob_upper_ratio_limit and glob_lower_ratio_limit - will be more restrictive on what is considered a singularity.
- January 9, 2014 More Singularity Info now in Table.
- January 9-19, 2014 Singularity analysis studied thoroughly. Now make sure tests are tending toward limit.
- January 9-19, 2014 Have program selecting which analysis to use to adjust step size.
- January 9-19, 2014 Improve logic to select initial step size.
- January 9-19, 2014 Fixed - It had seemed in some cases spaces were needed to separate tokens in equations in ode file.
- January 9-19, 2014 Fixed = Individual equations in ode file (1st section) had needed to be on one line. Now permit multiple lines.
- January 9-19, 2014 Comments may now be anywhere in ode file.
- January 23, 2014 Checked & already use last 4 terms - (Use last 3 terms of series to estimate error instead of only last term).
- February 5-6, 2014 Reworked main part of omnisode.rb and support Maple again
- February 7-11, 2014 Devised tomni.rb to replace shell scripts - works also for windows.
Taylor Series - Bibliography
- D. Barton, I. M. Willers, and R. V. M. - "The automatic solution of ordinary differential equations by the method of Taylor series", in The Computer Journal, 14(3): 243-248, (1970)
- Peter Henrici - Applied and Computational Complex Analysis: Volume 1: Power Series-Integration-Conformal Mapping-Location of Zeros(1973)
- Y. F. Chang - "Automatic Solution of Differential Equations" in Constructive and Computational
Methods for Differential and Integral Equations, Springer-Verlag, (1974)
- Peter Henrici - Applied and Computational Complex Analysis: Volume 2: Special Functions-Integral Transforms-Asymptotics-Continued Fractions(1977)
- George F. Corliss - "Integrating ODE's Along Paths on Riemann Surfaces", in Proceedings of the Seventh Manitoba Conference on
Numerical Mathematics and Computing (September 29 - October 1, 1977)
- George Corliss and David Lowery - "Choosing a stepsize for Taylor series methods for solving ODE's" in Journal of Computational and
Applied Mathematics, volume 3 and 4, (1977)
- Y. F. Chang - Start of Draft of book on Taylor Series (1978)
- D. Barton - "On Taylor Series and stiff equations" - ACM Trans. Math. Softw. 6, 3 (Sept. 1980).
- Louis B. Rall - Automatic Differentiation: Techniques and Applications, Springer-Verlag, (1981)
(this was hard to get, I finally found a discarded copy from the author's university's library on a web bookstore)
- George F, Corliss - "Integrating ODE's in the Complex Plane - Pole Vaulting", in Mathematics of Computation, Volume 15, Number 152,
- Y. F. Chang and G. Corliss - "Ratio-Like and Recurrence Relation Tests for Convergence of Series", in J. Inst. Maths Applics 25, (1980)
- George Corliss and Y. F. Chang - "Solving Ordinary Differential Equations Using Taylor Series" in ACM Transactions on
Mathematical Software Volume 8 Number 2, June 1982, pp. 114-144
- Y. F. Chang - "The ATOMCC Toolbox" in BYTE, Volume 11, Number 4, pp 215-224, April, 1986
- Peter Henrici - Applied and Computational Complex Analysis: Volume 3: Discrete Fourier Analysis-Cauchy Integrals-Construction of Conformal Maps-Univalent Functions(1986)
- Andreas Griewank and George F. Corliss - Editors -
Automatic Differentiation of Algorithms: Theory, Implementation and Application, siam (1992)
- Y. F. Chang - ATOMFT User Manual, ver. 3.11 (1994)
- Martin Berz, Christian Bischof, George Corliss and Andreas Griewank - Editors -
Computational Differentiation: Techniques, Applications and Tools, siam (1996)
- Angel Jorba and Maorong Zou Taylor User's Manual, November 13, 2001
- George Corliss, Christele Faure, Andreas Griewank, Laurent Hascoet and Uwe Naumann - Editors,
Automatic Differentiation of Algorithms: From Simulation to Optimization , Springer, (2002)
- Angel Jorba and Maorong Zou A software package for the numerical integration of ODEs by means of high-order Taylor methods, (January 27, 2004)
- Roberto Barrio - Performance of the Taylor series method for ODEs/DAEs in Applied Mathematics and Computation, 163, no. 2, 525--545 (2005)
- R. Barrio, F. Blesa, M. Lara, - VSVO formulation of the Taylor method for the numerical solution of ODEs in Computers and Mathematics and Applications 50 (1-2) 93--111. (2005)
- R. Barrio - Sensitivity analysis of ODE's/DAE's using the Taylor series method in SIAM J. on Scientific Computing 27 (6) 1929--1947. (2006)
- Martin Bucker, George Corliss, Paul Hovland, Uwe Naumann and Boyana Norris - Editors Automatic Differentiation: Applications, Theory, and Implementations, Springer,(2006)
- A. Abad, R. Barrio, F. Blesa, M. Rodríguez - TIDES tutorial: Integrating ODEs by using the Taylor Series Method., Monografías de la Academia de Ciencias de la Universidad de Zaragoza. 36 pp. 1 - 116 , (2011)
- A. Abad, R. Barrio, F. Blesa, M. Rodríguez - Algorithm 924: TIDES, a Taylor series Integrator for Differential EquationS in ACM Transactions on Mathematical Software. 39, no. 1, article 5. (2012)
- A. Abad, R. Barrio, F. Blesa, M. Rodríguez - TIDES: a Taylor Integrator for Differential Equations How to get the program (2012)
- September 23, 2013
I had had a hard time believing that the program worked as well as it had with the "/" for "*" error.
I looked at the code again. I hadn't changed the code - I must have changed the generated code.
I changed the code and it did much worse.
But anyway I discovered - in many cases - I can use a larger (I previously mis-stated smaller) h.
The omnisode code was actually the same for 0055 and 0056 and will be for 0057.
But the .ode and .odes files will be current in 0057.
- September 21-22, 2014
I was reminded in comparing my results to Dr Y. F. Chang & Dr. George Corliss's ATOMFT that I should be able to use a smaller h (increment).
After thinking where such a error most likely could have occurred I studied the "jump" code generated in one case.
I discovered an "/" where there should be an "*". I have corrected this ran tests for gcc and made a release.
Many cases work fine with larger increments now.
There are two cases - mtest2.ode and mtest6.ode which are giving problems and which are insensitive to lowering h.
- September 8-9, 2014
Testing comparing to Maple 18, ATOMFT, and Maxima (Runge-Kutta)
Have noted some bugs
- September 4, 2014
Worked on getting Chang's ATOMFT working on SUSE Linux 13.1
Seemed only changes required were to makefiles and scripts.
I hope I can make useful comparisons between ATOMFT and omnisode.
After I have it working, I may make modifications to make comparisons easier.
I will also preserve the original and version with only makefiles and scripts modified od ATOMFT.
Corliss specified easy rules to comply with for using his and Chang's work.
- August 29, 2014
Uploaded test results for maple.
Released rev0055 of omnisode.
- August 27-28, 2014
Added evalf in Maple for estimation of error.
Changed "analytic" to "closed_form".
Used glob_prec instead of literals in check in error calculations.
Simplified iteration count - eliminating unneeded variables.
Removed some commented out code.
Put put test results for c,c++, and Ruby on web site.
I expect this to become rev0055 of omnisode.
Still need to test maple.
I know there are some problems with maxima which I have been unable to diagnose, but it mostly works, although very slow.
- August 25-26, 2014
Worked on occasional division by zero error.
Tracked down and found error in handling of error reporting if error was exactly zero.
I hadn't completely overlooked this - but was not handling it carefully enough.
It shows up more often in maxima - but caused in c also.
Fixed much of it - but still some problem - perhaps elsewhere.
Quitting for now.
- August 24, 2014
Undid what I had done for division by zero - didn't fix problem and caused other problems.
Made maple use local variable (float) in factorial_1. (didn't seem to help)
Made maxima return 1.0 for expt to 0 power.
Still getting division by 0.0 error with mtest7_sm_h.ode. (maxima only)
Spent several hours with no success.
- August 23-24, 2014
Made improvement to iteration count (eliminated unmercenary variable)
Avoided division by zero in calculation of needed step size in rare cases.
Made display interval such that at least 10 steps will be displayed. (convenient because of way ode.over file works - for easier testing in many cases)>
Ran tests on c, c++, and Ruby.
I still believe there is some problem with maxima - but testing it takes much longer.
- August 18-20, 2014
Was going to run complete tests on new computer & make slightly updated release (0051)
Had made change in preodein.rb and preindent.rb.
Then noticed big errors in c and c++ results (I thought I looked at all this before).
Ruby and Maple seemed OK on a quick look.
NEED to look into this.
Now I have fixed problem (it was introduced in preodein.rb fix). I am rerunning all tests.
- August 13, 2014
Continued support for python.
Fixed problem with indentation - in pre-python code also.
- August 12, 2014
Fixed error involving absolute value in preodein.rb
Ran new tests for c, c++, ruby, maxima under Linux on new computer.
Not studied results thoroughly - no time.
Also ran test with Maple under Windows.
- August 11, 2014
Started adding support for python.
Ran test for release 0050 with new laptop.
- May 17, 2014
Uploaded release 0050 with fix for Digits affecting Maple estimate of error only.
- February 24, 2014
I went to study the Barrio, et al. material and either there is a temporary server problem or it is gone.
- February 23, 2014
Next I need to study Ahab & Barrio's work (which I added a while back to my bibliography but have not studied yet.
The next goal for my own program is studying error estimation and control.
- February 15, 2014
Fixed on problem where 'Digits' was not taken into account (for Maple) in error estimate.
Otherwise just running tests.
- February 14, 2014
Some testing on cases with less accuracy.
Seems smaller h helps - but hard to judge - only covers part of domain of original problem and new h over whole domain would take a long time.
Also - I think is estimate of error haven't taken account of 'Digits' in Maple.
- February 12, 2014
Added to automatically open tables in browser if possible after calculations (not released yet)
Ran test with erf and Si in Maple.
Added option to delete target directory if it exists. (not released yet)
- February 9-11, 2014
Replaced shell scripts to run omnisode with ruby program - tomni.rb.
Works with linux, cygwin and windows//
Works with c, c++, ruby, maxima and maple.
See Test Results.
- February 7-8, 2014
Working on ruby program (starting with tr.rb) to consolidate the logic in the various scripts.
tomni.rb works with languages c, c++, ruby, maxima, and maple.
and works with os's linux, windows, and cygwin.
- February 5-6, 2014
Got Maple 17 on February 5, 2014
Have omnisde working for Maple on February 6, 2014
More testing etc. needed
- January 31, 2014
I released Rev0064 of omnisode:
And uploaded test results:
- January 30, 2014 (Later)
More tests (include better comments in table) with additional cases with smaller increment in cases where results were poor.
This seems to show the logic is correct in those cases except for selection of increment.
Also I noticed using a smaller increment sometimes resulted in singularities not being detected, although the solution itself was more accurate.
The remaining exceptions (where smaller h does not help) are both tests involving multiple equations (mtest2 and mtest6).
I have not identified the problems with these.
More tests for c with estimate of error
- January 29-30,2014
Tests with estimate of error (not relying on actual solution) based on intermediate calculations
Fixed some errors in ode files - mostly about given location of singularities.
Also adjusted estimate of error to account for complexity of problem and number of terms used.
Also adjusted calculation of increment to by used for the same things.
Tests for c with estimate of error
- January 28-29, 2014
Added estimate of error based on intermediate calculations (not solution)
A lot of the tests with poor accuracy involve sqrt.
I should look into that.
Many others involve expt
Which sqrt is a special case of.
I should look into that also.
Also any others with poor accuracy.
I am too tired now!!!
- January 18-19, 2014
Working more on automatic selection of increment and singularity detection.
After more testing and debugging released rev0045 of omnisode.
See Test Results.
- January 15-17, 2014
Working on automatic selection of increment.
Fixed a few other bugs (don't remember details)
Removed some more dead code.
- January 13-14, 2014
Some experimenting with maxima:
A problem solved using higher level maxima functions that omnisode originally would have taken 25 hours to solve but was cut to 3.5 hours after some optimizing
The output of the maxima program
I talked to Dr Corliss (on the phone) today.
He mentioned Barrio's work & I located it with google and added it to my bibliography.
I haven't had a chance to study it yet.
- January 11, 2014
Uploaded Latest Test Results.
Reworked the README file and made release 0044 available.
- January 10, 2014
I've been studying singularity data. I provided "given" data for the ode files I could. Been comparing.
I found if often helped to increase "max_terms", but that you could also go too far.
Running tests and will probably upload results tomorrow. I find I think better with more sleep!
- January 9, 2014
Tests with just Maxims (1/9/2014).
More Singularity Info now in Table.
- January 8, 2014 - Extensive tests after many important changes.
See Bugs Fixed. and other Items Accomplished.
Tests with just c.
Tests with just c++.
Tests with just Ruby.
Tests with just Maxima coming soon
Later: More tests:
Tests with just c.
Tests with just c++.
Tests with just Ruby
- January 7, 2014
Corrected another problem with automatic negative h - actually there wee 3 mistakes
Fixed mtest7.ode - impossible problem.
Corrected that good digits were reposte 2 low in table - compared with output which was correct.
I made them match.
Used some macros in c and c++ - which should speed them up.
Correct in c where I had used abs instead of fabs.
I suspect equations with less accuracy are due to nature of the problem.
Next thing to do is add 3 more entries to the table - whether any singularities were detect by the ratio, three term and six term tests.
That will make examining results easier.
Running tests before going to bed.
- January 6, 2013
In testing I discovered that I had recently introduced a bug in negative increment
It is determined automatically if the end point is less than the start point.
The example was expt_c_sin_new.ode
I have fixed it, and am starting more testing.
Later I noticed mtest7.ode taking a long time
I looked it and realized it has no unique solution with those initial conditions.
I'm waiting till tomorrow to try changing the .ode file - other cases with accuracy problems may have similar causes.
- January 5, 2014
More Testing - Then made release 0041
Later I made the following improvements at the suggestion of someone on the Maxima list
Added modedeclare for a local variable I missed -for Maxima.
Used alias instead of functions I needed different names for. (Improved performance by about 20% for Maxima)
- January 4, 2014 Fixed Several Bugs Fixed.
I had messes the parsing all up - and then ended up putting all the old code back from an old copy. I did make preodein.rb add spaces between tokens. I think that may have sometimes mattered.
I can make another simple change there so equations don't have to be all on one line.
Current Test Results.
Fixed so equations each can be on more than one line. After test results.
Now I am going to study test results - line change has no effect on them.
The problems I found in studying test results mainly appeared to be caused by bad initial conditions causing such things as sqrt of negative numbers.
I did find an error in defaults being applied after the 2nd & override blocks instead of before. I fixed that.
I am VERY TIRED - going to bed now (8:00pm)
- January 3, 2014 My friends on the maxima list helped me arrive at this simple maxima program to do part of solving a diffeq
It does a lot for its size. I intend not to become much more familiar with maxima!
- January 2, 2014
Got Ruby version to work in a preliminary test.
I expect the others will be a little easier.
Later, I have the c, c++ and Maxma versions running again. They need much more testng, and there are still many ways they can be approved. I amy putting the
results on the web, as I can study them easier there.
- January 1, 2014
I decided to dig into the organization of the code in omnisode.
It has become disorganized as it developed.
It will make it much easier to understand and to work on - hence more reliable.
I have attempted to make it so that input equation can break on a line. I have not tested that yet.
The Ruby version is functioning again.
I am also trying to improve the code in other ways.
- December 31, 2013
Ran test for maxima.
Worked on how to use modedeclare and compile
Resulting test source code: Maxima sample code for modedeclare, compile, declare_variable, and array
The test results
Later that day I I found how to use madedeclare in a slightly different way - which will be much easier to implement in omnisode:
Resulting test source code: Maxima sample code (modifed slightly) again for modedeclare, compile, declare_variable, and array
The modified test results
I started working implementing the necessaty ckanges to omnisode
It is a big job.
- December 30, 2013
Improved attempt to see if radius is approaching a limit. Improved use of radius calculated in adjusting increment. Other minor fixes.
More Test Results.
- December 27-28, 2013
Worked on singularity logic.
No longer have any computer with Maple.
Al least temporarily stopped testing with Maxima also - it's not really good for number crunching - but it sometimes detects syntax other languages miss.
Ran across warnings in c and c=+ I had missed before.
They mostly concerned time measurement (type mismatches)
I believe all that is now correct for c and c++ - should look over more for Ruby soon and Maxima someday.
Thee low lever singularity functions seem finisher. Almost completely rewritten since last release - I simplified to coding.
The radii of convergence work well. The order does not seem to I have been over it carefully.
But there is a higher order function which tries to see if they seem to be going toward a limit.
That needs more woek.
After sleep and more work uploaded December 28, 2013 - Tests with just c. Improvements to coding for singularities for ratio, three term and six term tests. Good results there for radii but not order. Also corrections to time display functions especially for c, c++ and Maxima. I am limiting big test runs to c to save time. Ruby, c+, and Maxima seem OK with limited testing. I no longer have Maple available..
More Test Results.
- December 16, 2013
Worked about 3 hrs in early morning.
Logic revisions for singularities seems basically functional.
Get good results (in preliminary testing) for radius for ratio, three term and six term tests.
There is trouble for the order (as before) for each of these tests.
I wonder if there is catastrophic subtraction error.
I may disable that code tomorrow and proceed with testing without it.
I can then best figure out how to use that information.
- December 15, 2012
Worked on revised singularity code.
Massive recoding on this part.
I saw what looked like some errors.
But just trying to finish recoding.
Too hard to figure out what the old code was doing!!!
- December 12, 2013
Have coded new structure to singularity detection.
Still most details to fill in
Much easier to understand.
- December 10, 2013
I am working on seeing if the three & six term tests are going toward a limit (like ratio test)
I am trying to simplify the code & use more common code - but it gets tricky.
I gave up for now -- I need to think more about it.
- December 2, 2013
Realized I made a mistake yesterday.
I have realized the relation  in Chang's typescript (part 4)
Can also be more directly derived from 4 in Appendix A of Chang's typescript (end)
I have used Maple to derive the radius of convergence & order
I have gone over this part of the code thoroughly and am now getting better results.
I still need to go over it more, but need a break.
I screw things all up if I try to work too long.
There are errors in my equations of May 30, 2013 in Notes below.
Apologies to Dr Corliss!!!!
No matter how I try, even though I can find no error, the equations for the order do not work well.
I have "if'd" some of the code out for now - I seem to be getting nothing useful from it.
I am going to add code to the three and six term tests to see if they seem to be reasonably going toward a limit (as I did with the ratio test).
- December 1, 2013
Managed to solve problem - but solution involves 'RootOf'
So using it will be limited to Maple and possibly Maxima (I don't know Maxima that well)
(I am using Maple 11)
- November 29, 2013
Worked on problem - unable to solve.
- November 28, 2013
Worked on solving for complex singularity using my modified three term test.
Couldn't get Maple to solve it.
- October 23., 2013
Worked on my_test.mxt
- August 22, 2013
Installed latest maxima for windows.
I believe in open source and cannot afford the cost of Maple
- August 18, 2013
Have ideas. Waiting till I can afford new Maple, or have time to learn more maxima.
Also I have free version of latex, and have learned a little.
- July 3, 2013
Minor corrections for syntax for Maple and Maxima. Then started comprehensive tests. Then added to test section.
- June 24, 2013
Added derivatives w.r.t (n-1)th and (n-2)th terms and added to guess at error in 3 term analysis.
I worked out the equations a couple weeks ago but hadn't finished coding.
- June 13, 2013
Went over June 10 and 11 work and found error in calculation of "part1" - corrected and re-ran tests.
Tried to estimate error in radius and order by multiplying by taylor term - really need error in taylor term - so would be fraction of this.
- June 11, 2013
Optimized computing derivatives of radius and order for singularities r.r.t. last Taylor term. Also decreased number of terms and increased maximum increment. See test results.
- June 10, 2013
Tried computing derivatives of radius and order for singularities r.r.t. last Taylor term. See test results.
- June 9, 2013
Ran comprehensive test, and uploaded release corresponding to May 30 version.
- May 30, 2013
After studying the three term analysis, I decided to use the formulae in Chang's appendix, which I have tested.
For f(x) = (x-a)**p
F(M) = (p - M + 2) * H * F(M-1)/Rc/(M-1)
F(M-1) = (p - M + 1) * H * F(M-2)/Rc/(M-2)
I Used Maple to solve and got:
temp1 = 2 F[m-1] * F[m-1] + m F[m-2] * F[m] - m * F[m-1] F[m-1] - F[m] * F[m-2]
temp2 = H * F[m-1] F[m-2]
temp3 = F[m-2] * F[m] + 4 * m * F[m-1] * F[m-1] + m*m * F[m-2] * F[m] - 2*m * F[m-2] * F[m] - 4 * F[m-1] * F[m-1] - m*m * F[m-1] * F[m-1]
Rc = temp2/temp1
p = temp3/temp1
It seemed to still do very well (even better that Coriss' equations for the radius and not quite so bad for the order - but still not good).
- May 29, 2013
Worked a little on determining order of singularities. Looked at Chang & Corliss' paper on it & noticed formula that goes with ratio test and added it to my program. Also tried to put three term test closer to form in paper. I had tried, because of problems with it, to solve the eq's with Maple and come up with different formulation. Nothing I have tried has worked. I suspect root cause is the problem I found on April 24, 2013. I have some ideas on how to improve results. But I probably cannot get to it for a few days.
- May 27, 2013
Uploaded results from test with larger display interval.
- May 26, 2013
Uploaded latest test results.
Tests which place ratio test and three term and six term analyses adjacent to each other.
I have uploaded corresponding release (0036) today.
I also started another general test with display interval 0.1 instead of 0.01 for comparison of speed and accuracy. (The display interval is also a minimum increment.)
- May 25, 2013
Found and fixed problem introduced when debugging memory problem in c and c+. Not quite so much faster now. Still quite a bit faster than it had been. Also fixed problem with missing comma in code generated for maxima in cases of multiple equations. Because of that am re-running general test overnight. It used to take about 16 hours - now about 6 hours.
- May 24-25, 2013
In debugging logic to optimize size of increment, glob_h, I found and corrected a lot - basically it had not been working correctly at all. The framework is now functioning correctly - the details can, I hope, be improved on. Also I tracked down a memory problem which only affected c and c+. The tests are now running much much faster, which will make progress faster as well. I am too tired to look over the latest results closely now - I am tired. I will after some sleep - and hopefully have results to upload.
- May 23-24, 2013
In going over tests I saw there is inaccuracy in c and c++ on add_c_lin, although ruby, maxima, and maple work OK on it. Also c and c++ wor OK on sin. In debugging I was lred to working on the opimal determination of the increment, glob_h. This led to a much larger glob_h being used, and thus much faster computations, with very good accuracy. But the inaccuracy in the one case of c and c++ remained. Anyway I am very fatigued. Thus I started a general test to see if I can find a pattern. The error almost has to be in the generated atomall function. The results should narrow it down. Also I eliminated tracking the 1st detection of singularities. Only the detailed results will help there.
- May 21, 2013
I ran general tests yesterday and today. It was a comedy of errors. I had changed the path to the results, from html/results to html/omniresults, in order to better specify that these file need not be backed up. I thought I had changed it everywhere. But I hadn't in the 1st line of gentest.sh, which removes the prior test. Then all the results were not saved. Then I corrected that and reran it. This time, in the paths in the t*l.sh files I had accidentally changed it in the results file name itself. I wrote a quick ruby program and renamed all the files. Then I looked at the results. C and C++ were giving inaccurate results on add_c_lin.ode, but Ruby, Maple and Maxima gave good results. I am too tired to work more on it now. I remember I haven't changed the path in the ruby script for windows either.
Later - Have new laptop set up ruby works OK - but duplicated problem with c and c++ - maxima and maple not installed yet.
- May 20, 2013
More work on preparing for testing of location and order of singularities.
The results look promising. But order has trouble in three term test. Perhaps using more terms could help. I made a release of this version at Download omnisode (0035).
- May 19-20, 2013
Modified omnisode to accept (fully implemented for time being only for generating ruby - will need to preprocess subscripts) location of closest singularity and compare actual (given) radius
and order with those computed by the ratio, three term and six term tests. These are output in the results files to which there are links in the results table
Now that these results are easily examined, I can start some serious debugging!!!
- May 19, 2013
I am getting ideas on what to do next. I am considering adding columns to the results table for the actual radius and order of singularities, and also the results of the equations derived
by Corliss and Chang. The radii given by those seemed close, even though I think I found an error in their derivation. (also retaining the radii given by the ratio test.)
- May 15-16, 2013
Worked on implementing ratio test (see Henrici Vol 1) instead of Corliss' and Chang's tests. See Test Results..
- May 14-15, 2013
Worked on debugging ratio test.
Main difficulty was I was only keeping minimum radius - instead of radius for each series.
I was making progress, but became exhausted and went to bed.
- May 13-14, 2013
Modified omnisode to use ratio test as found in Henrici, Volume 1
Debugged and started comprehensive test (on older slower computer - I need this one for other things)
NOTE: Also found I was getting bus error with emacs - switched to emacs-x11 and resolved the issue.
- May 3, 2013
Version 0033, matching latest test results, was published 4/30/2013.
This page was not updated till today due to technical difficulties.
- April 29, 2013
Someone has recommended Chebyshev polynomials.
Some links are:
I ran the problem given in the 1st link with omnisode - Test from Chebfun
- One hour tutoriial(I couldn't get to work)
- Approximation Theory and Approximation Practice
It has trouble with expt near zero, but does OK, it seems, with the products.
The link gave an exact symbolic solution but not the numerical value.
- April 24, 2013
Looked over test results and uploaded them. See Test Results.
Cases intended for testing singularities (singN.ode (where N is a digit) and nonlinearN.ode (same) gave accurate results for radius.
I did not get good results for order and it now says "Unknown" - that was easier than removing it altogether.
Some functions which I would have expected to be entire gave finite positive radii.
The equation will always give an answer and my only check is that the radius is positive.
I had been checking for positive order, but all orders except negative integers are permissible, and how can you check that with floating point?
After more study:
I kept getting wrong answers for the order, though the radii were close in my test cases. I kept studying the analysis in Chang's typescript more & more closely.
The trouble is going from equation  to equation 
d f(x)/dx = - p * (u(x)**(-p-1)) * d u(x)/dx 
F(2)*U(1) / h = - p * F(1) * U(2) / h 
I think there is an error in that one needs to evaluate the derivatives of u(x)**(-p-1) at x=0. and I think Corliss took the derivatives of u(x) at x= 0 and took them to the (-p-1) power. Though I cannot quite get his result. (I have not tried very hard as it would have to be wrong). But if one takes the derivatives of , then the terms of u(x)**(-p-1) will not disappear as they do in , unless p is a negative integer.
[Added 4/25/2013 - the equations (on which the "three term analysis" (p. 79.)) referred to are in part 4 (pp. 74-75) of Chang's typescript.
Also it should be noted that the "six term test" (p. 96) is based on the relations used in the three term test (p. 82).
A variation of the results also were published in Y. F. Chang and G. Corliss - "Ratio-Like and Recurrence Relation Tests for Convergence of Series", in J. Inst. Maths Applics 25, (1980), though not all of the derivation.]
- April 23, 2013
Checked over test results.
Reduced glob_max_h for a few cases.
Removed mtest6 cases with larger glob_max_h
Adjusted bounds for a few cases.
Started another test.
- April 22, 2013
Studied results of yesterday's test.
Hard to check results, because with different speeds of different languages, they timeout at different points.
Modified program to report data at point where singularity is first detected in order to have data in table easier to compare.
Started another test.
- April 21, 2013
Went over the derivation of equations for poles.
I believe there is an error going from eq  to eq  in Chang's typescript (see bibliography - it's in part 4 ).
There eq  is:
F(2) * U(1) / h = - p * F(1) * U(2) /h
I believe it should be:
F(2) * U(1) / h = - p * F(1) * (U(2) ** (- p - 1)) / h
The derived equation for the radius of convergence seems to work anyway.
I haven't figured how to derive results from the corrected equation .
I wasn't using the order anyway.
I started a new comprehensive test. The order is now given as "Unknown".
- April 19, 2013
Worked much trying to get radius to work for real poles.
Redid formulae and couldn't make it work any way I tried.
Worked many hours while following news on search for Boston Marathon bomber.
Did however make some fixes, though they did not fix basic problem.
Not enough for a release.
- April 18, 2013
Worked on calculations of poles.
Made some corrections.
Complex poles - both radius and order working very well.
Real poles - radius works very well - but not order.
Suspect problem may not be in calculation itself - need to check it - and use of arrays - etc.
- April 4, 2013
rev0030 released - Si and erf (Maple only)
- March 30, 2013
Added and tested Si (Maple only)
- March 29, 2013
Finished and tested erf (Maple only)
- March 28, 2013
Added logic for erf for constant and linear arguments. (Will only work for Maple - needs erf built-in function.)
Started testing it.
Still need full series case - shell added.
Also have shell for Si function added. (Also for Maple only - requires Si built-in function.)
- March 19, 2013
Studied some generated code and did some testing. (All OK)
Thinking about how to go about handling singularities.
- March 2, 2013
Decided to focus next on testing and checking logic for locating singularities.
- February 25, 2013
Created zip archive and added to sode sourceforge download section.
Contents should be identical with tar.gz file.
Added for convenience of windows users.
- February 24, 2013
Finish changes to also support windows (rev0029) - but just Ruby code generation.
I don't have the other languages under Windows.
- February 23, 2013
Modified omnisode.rb to use "File.join"
Had difficulty with generated code so "entry.html" is now generated in current directory."
Scripts will have to be modified for this.
Have ruby program tr.rb working under windows (functions as tr.sh under Linux and cygwin)
Will have to modify (and test) shell scripts to account for entry.html being moved.
Assumes editor in path named "e" - my name for emacs under windows.
- February 22, 2013
Started working on making omnisode work under windows.
I'm using windows 7 and installed the latest Ruby today.
Installed omnisode0027 on windows machine.
First just working on generating ruby. (I don't have any other languages installed on windows.)
Only change to code so far is path separator (which would also work with other languages).
Used 'SEPARATOR' to indicate '\' where linux (& cygwin) would use '/'.
Haven't decided how program should determine how to set that variable.
Made bare minimum batch file to test with.
Omnisode worked with sin.ode and ode.over.
Omnisode failed to solve Sin.ode did not work right without ode.over. (Needs investigation - probable happens under linux also.)
Need substitute for linux shell scripts (Ruby program?)
- February 20, 2013
Added rdoc (ruby documentation) for the various programs in used in omnisode.
- February 3, 2013
Test downloaded latest version of omnisode on a different machine (which uses cygwin).
Only difficulty was emacs is not set up same on this machine.
Only tested with Ruby, but "ruby diffeq.rb" worked after trying "tr.sh" which had created diffeq.rb.
I don't have any of the other languages on that machine.
- January 30, 2013
The code seems stable, and I need to attend to other things.
I will let it incubate in the back of my mind before I do more - this is how I get my best ideas.
- January 28-29, 2013
Did testing and uploaded table of results.
The radius of convergence seem to be calculated accurately (which is all I use for now) but not the order as well especially for real singularities.
I have not been able to figure out why, although I have spent quite a bit of time on it.
Also the increment for a desired accuracy is not always calculated correctly. A maximum h may be set to help with this.
Made release rev0027 of omnisode 1/29/2013.
- January 27, 2013
Added link to copy of ode.over to output table.
Made first ode file copy correctly.
Modified singularity code so it tells which equation.
Made so glob_max_h can set h less than would be otherwise computed - and which may be less than display_interval.
Cleaned up both omnisode.rb and ode files.
Started comprehensive test in all 5 languages.
- January 26, 2013
Ended up using equation from Chang and Corliss (1980)
Results added to test results section.
- January 25, 2013
Attempted to fix problem with order of pole for real singularities.
No luck, although I did get re-familiarized with the code.
Spent about three hours.
Later: Worked on this.
Studied - Y. F. Chang and G. Corliss - "Ratio-Like and Recurrence Relation Tests for Convergence of Series", in J. Inst. Maths Applics 25, (1980)
Derived related equation (using Maple). Different because they used different convention for indices. Better results but still not very good.
- January 24, 2013
Studied test results for radius of convergence and order of poles.
Both appeared correct for complex poles.
Radius of convergence but not order of pole seemed correct for real poles.
I could not locate an error in the program.
- January 17, 2013
Released rev0025 of omnisode.
- January 15-16, 2013
Devised tests for cases I didn't have analytic solutions for comparison. Five of the eight worked well. The other three need attention.
New Tests for eight special cases.
Later: Initially 3 had problems - due to problems in entering the equations in the ode files. Two required a smaller h for better accuracy.
New Tests for seven special cases.
Two cases requiring smaller increment for good accuracy.
- January 14, 2013
A Summary Table (1/14/2013) With a few bugs fixed - Extensive tests for single equations using c. c++, Ruby, Maple, and Maxima.
Single equations (not yet systems) work fairly well.
Tests with multiple equations (1/14/2013) still one problem is less accurate than others.
These tests just using Ruby - to save time.
Tests for multiple equations using Ruby.
The current omnisode (0024) was published on the web January 14, 2013.
The latest published version may be obtained at The Sourceforge Download page for Omnisode
- January 12, 2013
Had time to take a quick look and found the problem in c.
It was introduced by one of the last of the other fixes.
I corrected it and started another big test.
Later: Looked over test results carefully.
All lines in test script resulted in rows in results table.
mult_sin_lin.ode was not as accurate as others.
Found mistake and started another test.
- January 11, 2013
Test finished - some sort of problem with script (it seems) for c.
Too busy today to work on it.
- January 10, 2013
Fixed remaining 3 problems. Actually just 2 hard to find bugs.
Put back all commented out tests in script and started test on single equations.
Need to look closely at results - especially results on singularities.
- January 9, 2013
I had previously just checked for test results in the output tables mainly for low numbers of correct digits.
I had noticed some results missing - I knew these could be caused my syntax errors in the generated code - but did not think there were many.
Also there were a few (just for Maxima) which I had removed from the script which would get "stuck."
Today I went though and noted all rows for which there were some results, but some rows were missing.
I found 49 rows missing.
I am going to start making corrections tonight, but it will probably take a few days to make corrections, run tests, and make a new release.
After working through the evening, I had fixed all except 3 cases of missing rows. Required only 4 fixes I think.
There is also, though, a problem with estimation of singularities. I may have introduced it in one of the 4 fixes.
That area has been less tested anyway. Everything else has to work right before it can.
- January 6, 2013
Uploaded release 0023
- January 5, 2013
Uploaded results with fix from Jan 3 uploaded.
Will make release with fix soon.
- January 3, 2013
Found and fixed problem in emit_pre_assign and emit_assign.
- December 22, 2012
Uploaded rev0022 of omnisode.
Also extensive test results on single equations uploaded.
- December 21, 2012
Test equation is:
y''''' = -y'''
Error appears in y''''
Somehow must miss its calculation or using some uninitialized value in its calculation.
Will sleep on it!
NO SLEEP!!! - Worked ALL Night!!!
I think I fixed the problems I had detected!
8:30am Started General Test.
Will sleep (hopefully) while test runs.
Single equations seem to work fairly well.
Some problems with multiple equations.
- December 20, 2012
I studied and improved debug code.
Located and fixed errors in pre_emit_assign and emit_assign.
Can tell these are correct 1st time through.
There must be additional error in jump for higher order equations.
Began inserting debugging there.
- December 19, 2012
I decided to do testing on a problem that had higher order derivatives on both lhs and rhs, but a single equation.
y''''' = - y'''
I put in debug to give both calculated and evaluated values of analytic solutions of first several terms.
At the beginning of a run, they were very accurate, but then suddenly a term goes wildly astray.
Thought it might be hardware.
Tried on a different computer - had to install ruby and program.
Had similar results - but not exactly same. - Sudden deviation - but not same place.
Must be software bug - like memory clobber - but some of the languages I am using are supposed to check that.
Enough for now!!!
- December 17, 2012
Tests with different display interval for a problem with multiple equations.
However computed glob_h for desired accuracy was generally smaller than interval anyway.
See Special tests 0002
- December 14-15, 2012
Completed general test and one special test comparing different interval sizes for the same problem.
See General test 0001
In the following table results are given for various increment sizes.
The increment is no longer set directly - the value used is computed to try to achieve a desired accuracy.
But the increment will not be larger than an interval set to display the results at. So there is still some control.
It turns out larger increments are better - at least up to some point. I still need a lot more testing on this.
See Special tests 0001
Note: Results files larger than 10MB have been removed. Many are still larger than 1MB. Some had been over 500MB.
- December 13, 2012 (Later)
- Fixed time in c and c++.
- Discovered larger error when using larger interval seems more due to larger range tested than increased step error.
- Started work on testing specific cases for c only (ones showing larger error in general test).
- Rerunning general test with time fixed in c and c++.
- December 12-13, 2012
Looked at test results and ran a few more smaller tests.
- Because some tests reached max time, Maxima has greater observed accuracy, but only because it is slower and had fewer iterations to accumulate error.
- Times in c and c++ are too large. I can tell the tests are running much faster than measured time. Maybe what I am using as a second is really a millisecond?
- Most tests run accurately.
- Most tests that are not accurate, are accurate when smaller increment (h) is used.
- When the smaller h is used, the results seem too accurate to be mistaken.
- Need more systematic testing for cases where error seems larger (with similar h).
- Note: because results often time out, before completing range, increasing h not only increases error on individual step, but also greatly increased the total range covered - thus magnifying the error reported.
- Regarding this last note: the last fact may be overcome with more systematic settings in ".ode" and "ode.over" files.
- December 11, 2012
Looked at problem again with correct debugging.
It appeared that I had already solved problem,
and that remaining appearance of problem was in error in debug code.
Started another general test for overnight.
- December 10, 2012
Fixed code added for debugging.
Then too tired to debug problem.
Also updated project for SourceForge.
- December 9, 2012
Isolated one error, but could not figure it out - could not seem able to think clearly.
- December 4, 2012
Fixed a couple things related to problem below.
Seems I have, at least, a problem with a sign somewhere.
- November 28, 2012
Located one problem and fixed.
Emit_pre_diff needed modification for the special case handling.
Still something wrong.
Too tired to think any more.
- November 27, 2012
Narrowed bug down.
I believe it is in emit_pre_assign and emit_assign (at least).
But I could not think clearly about how it should work.
- November 22, 2012
Added debug code for 1 specific problem (diff2.ode)
Debug code seems OK, but still haven't solved problem.
Seem unable to think (clearly)
- November 21, 2012
Added Total Expected time to solve problem to output and replaced expected lime left with it in table produced.
I think it is more useful.
- November 19, 2012
Seems there definitely is a problem (see 11/18/2012)
Studied code and didn't find the problem.
May have to resurrect some debugging code.
Too tired tonight.
- November 18, 2012
I went over test results.
It seems that there are accuracy problems iff the derivative on the lhs is order greater than 1.
That is an isolated section of code to examine, but I am too tired today.
- November 16, 2012 (early morning)
Fixed problem in ode files with lines containing only blanks - introduced with logic for comments.
Improved logic for determining size of increment - but still not really correct.
Started another general test.
(later) Went over results, fixed some problems and started retest.
- November 11, 2012
Comments now OK in ode files (except for top equations section)
Increment automatically calculated to obtain desired accuracy for specified range. (Including negative increments.)
The last item is coded and runs but it seems not to always give the desired accuracy and thus needs much more work and testing.
- November 10, 2012
Continued working on the same things - some progress
- November 9, 2012
Worked mostly on automatically determining optimal increment size.
It started working, but then I got sidetracked on getting comments working correctly in preprocessor (preodein.rb).
Making progress there also.
- November 8, 2012
Working on testing.
Tweaking ode files and not using "ode.over" - thorough testing requires that for now. (Especially size of increment.)
Would be good to make program select values itself - hopefully that will come later.
- November 7, 2012
Worked on logic to adjust step size based on estimate of radius of convergence.
Started test on this including c, c++, Maple, Maxima and Ruby.
- November 6, 2012
Indenting now OK.
Started lengthy general test.
Intend to make some support to adjust increment size for singularities.
Also to take larger step size as long as sufficiently accurate.
Also to support Ruby with my APFP library.
- November 5, 2012
Have Ruby code generation working.
Also working on making output code indented properly for every supported language. (Maple, maxima, c , c++, Ruby).
- November 4, 2012
Studying articles relating to increment size and pole-vaulting around singularities.
Interesting - but seems (especially pole-vaulting) to require special analysis of individual cases.
Started adding support for generating Ruby code.
- October 23, 2012
The tests with the larger increment (below) have catastrophic errors if they pass over singularities. (I have known this,
but failed to comment on it here.) I have yet to work out how to handle this - I am researching and thinking about it.
- September 21, 2012 Tests with larger increment (0.005) - Most results OK - But some very bad - No changes in code.
Report of tests for c. c++, Maple, and Maxima (WITH LARGER H - 0.005).
Expect to mostly be studying problem of error control and experimenting with code - no more releases likely for a while.
- September 8, 2012
The current omnisode (0018) was published on the web.
- September 7, 2012
Improved September 6th work and made test below:
Tests with negative H and/or display interval.
- September 6, 2012
Made to support negative increment (H)
Added capability to display values calculated only at specified intervals.
- September 3, 2012
Improved correct digits reporting (had been incorrectly reporting value for last equation for all in systems of equations)
Better report of tests for c. c++, Maple, and Maxima.
Released this version (0017).
This version (0017) may be obtained at The Sourceforge Download page for Omnisode
- September 2, 2012
Fixed problems in my work done on August 30, 2012. Started general test.
- August 30. 2012
Fixed problem with "Correct Digits" (see below) and started general test.
- August 22, 2012
I have realized the "Correct digits" reported in the cases of multiple equations is that of the last equation - for all of the equations.
This will be fixed in the next release.
- August 21, 2012
Ran general batch test with fixes from yesterday
There were some problems, so there is more to do before I post results.
However I can give better time comparisons.
I discovered the time reported by Maple is CPU time - so I am using the difference in the times in the time stamps instead.
Maxima can report either clock or CPU time - there are 2 functions - I am using clock time.
C and C++ report clock time.
The following was obtained with mtest6.ode
Maple did 100 iterations in 24 seconds (17 seconds CPU)
C++ did 69 iterations in 60 seconds
Maxima did 8 iterations in 67 seconds
Maple 4.17 iterations/second
C++ 1.15 iterations/second
Maxima 0.12 iterations/second
NOTE: The comparison is not altogether fair. Digits was set to 32 for Maple, while only 16 digit precision was obtained for C++ and Maxima
Fixed problems with changes to support local variables in ode files for maxima.
Fixed Digits reported in c and c++ (got it right this time).
Reversed order of entries is results table.
Added script files for c.
A Summary Table (8/21/2012) Now including c and c++ (really just c code but compiled with c++) is at:
Extensive tests for c. c++, Maple, and Maxima.
Published current omnisode (0016) on the web.
The latest published version may be obtained at The Sourceforge Download page for Omnisode
- August 20, 2012
Fixed a number of problems with last release (omnisode0014)
- Made to work with c as well as c++ (using gnu gcc and g++)
- Fixed time measurement in c and c++ (Also note Maple time is CPU time - Maxima, C and C++ time are clock time)
- Fixed digits field in table for c and c++
- Added number of correct digits to output and table when analytic solution given
- Added beautify program for tcppl.sh script.
- Made preodein.rb handle local declarations for maxima (only affects analytic solutions)
- August 12-13, 2012
Testing and debugging c++. Test results at Extensive tests for Maple, Maxima, and c++
For some reason the tests "expt_c_sin", "expt_lin_sin", "expt_sin_lin", and "expt_sin_sin" got stuck in Maxima ands were removed from the test.
Need to study results thoroughly.
- August 11, 2012
Ordered the 3 volumes of Peter Henrici - Applied and Computational Complex Analysis
Finished work on generating c and c++
Started testing c and c++ - a few corrections, but going well.
- August 10, 2012
Added blocks and declared all local variables in maximas. Maxima was a little slower with the blocks added, but not very significantly. I knew where they all went, as I had them (local declarations) for Maple, and mint checks that so I should not have missed any. I have found out how to use big floats in maxima. It will be quite a bit of work, and as maxima is sort of slow, and presumably this would make it slower, I am not pursuing this now. Most of the time, I suspect input data to real problems is not known that accurately anyway. I will be looking into interval arithmetic though - probably not for maxima. Maple has something but I don't know if it will do what I need. I understand what mode_declare does now, but it doesn't support big floats, so I will avoid it for now. I don't have not had problems without it, and I think it is only effective if the code is compiled. I am not interested in compiling as then I could not use big floats, and I want to do that someday. Next, I guess, is finishing up C and C++.
- August 9, 20121
I ran a timing test comparison for c++, Maple and Maxima.
Set time limit to 1 minute.
Maximum iteration was sufficient to solve entire problem.
Problem was sin.ode for all three.
Output was redirected to file only - none to screen
c++ solved entire problem in 23 seconds (140001 iterations)
Maple performed 11542 iterations in 60 seconds
Maxima performed 372 iterations in 60 seconds
That amounts to:
c++ 6,087 iterations/second
Maple 192 iterations/second
Maxima 6.2 iterations/second
c++ is almost 32 times faster than Maple which is almost 31 times faster than Maxima.
So c++ is almost 982 times faster than Maxima.
Later (I think my dates may be off the last few days):
Removed mode_declare's from maxima
Then maxima performed 394 iterations in 60 seconds
Was able to get bloats (with trig functions) to work in maxima. I had misunderstood syntax before.
Am going to try to compile maxima code - I realize I cannot do this with if I use bfloats - which I don't now.
- August 8, 2012
c++ mostly working. So much faster will have to modify output.
Files produced are too large to view, if enough iterations are used in order measure time.
- August 7, 2012
c++ version running and giving good answers with very limited testing.
Some problem with displaying time data, but it appears to be much faster than either Maple or Maxima.
- August 3-6, 2012
Working on adding c and c++ to languages generated.
Going much easier than anticipated.
Mostly just tedious parts left - a few tricky ones also.
- July 27, 2012
Checked coding of Chang's case 5 for expt_lin_lin.ode and found nothing wrong. Tried using Taylor series of 250 terms but this did not help. The analytic solution and derivation of Chang statements are beyond my skill level currently.
I got to thinking that I remembered something on this in Rall (1981). Looked and found it. Not sure if it will
help or not. No time now.
- July 11, 2012
Made release 0013 It may be downloaded from The Sourceforge Download page for Omnisode- Main change is to handling constants and linear expressions as special cases.
See test results above
- July 1, 2012
Did sqrt, and modified array allocation to save memory as made possible with optimized logic.
- June 27-30, 2012
All cases of exponentiation coded, except sqrt.
Next - a lot of testing!
- June 25, 2012
Continued changes for constants and linear expressions - exp, and log done.
Log has truncation problem
- June 24, 2012
What I used here was the time left - I should have used the time elapsed + time left. - I did not think of that before.
|Language||Compute every time||Pre-compute all||Compute as needed and save||With Linear Opt|
|Maple ||1 min 40 sec (16 digits) ||29 sec (16 digits) ||41 sec (32 digits)||20 sec (32 digits)|
|Maxima ||40 min 33 sec (16 digits) ||31 min 58 sec (16 digits) ||31 min 25 sec (16 digits)||8 min 52 sec (16 digits)|
But I don't think it really matters - they should be proportional to each other.
By linear opt, I mean taking advantage of the fact that the operand to the sin function is linear (In this particular example).
Did more testing on arccos & arctan. Fixed 2 typos in arctan. Both did better when dividing increment by 100 even though
multiplying iterations by 100. Thus it must be truncation error. (Getting about 20 place accuracy.)
- June 23, 2012
I ran into less accuracy with arccos and arctan than I have been typically getting. (I know the domain of values matters and
that I need to do much more thorough testing. I am to mostly trying to catch blatant errors now. But I was getting significantly less
accuracy (for the domain I was using) for arccos and arctan than other functions. So I decided to look at the Fortran code generated by Chang's
Atomft. First I tried it on Linux, but was unable to find g77. So I tried cygwin, and after a little work (mainly running d2u on the files
and reworking the make files to my liking. - I had used Atomft before but only retained the original zip file.) ran it on the arc functions.
I got error messages - they are not supported by Atomft. Then I went back to Linux and, with a little web browsing, found out that g77 is now
called gfortran - and got the same results.
- June 21, 2012
Continued changes for constants and linear expressions - sinh, cosh, tanh, and arcsin done.
- June 20, 2012
Continued changes for constants and linear expressions - division, cosine and tangent done.
- June 19, 2012
Continued changes for constants and linear expressions - subtraction and multiplication done.
Also looking at some Taylor series capabilities in Maxima - should also look at Maple
- June 18, 2012
Same tests, but factorials computed only as needed, but then saved for when needed again
should be same results except faster(June 18, 2012)
Started revisions to treat constants and linear expressions as special cases of series in order to optimize these common cases.
I have modified the basic structure for this and implemented and tested for addition and sine.
- June 17, 2012
Same tests, but factorials pre-computed
should be same results except faster(June 17, 2012)
I have revised the factorial logic again. Instead of pre-computing all the factorials, I am just computing them the first time needed.
Also I noticed (Affecting maple only) I had Digits set to 16 in sin.ode, to 50 in h2sin.ode, h3sin.ode and sing4.ode, and to 100 in sing4.ode.
These are all now set to 32. (like everything else).
The timing results for sin.ode were:
I have another general test running. It will take several hours and is on another computer so my editing here does not affect it
|Language||Compute every time||Pre-compute all||Compute as needed and save|
|Maple ||1 min 40 sec (16 digits) ||29 sec (16 digits) ||41 sec (32 digits)|
|Maxima ||40 min 33 sec (16 digits) ||31 min 58 sec (16 digits) ||31 min 25 sec (16 digits)|
- June 16, 2012
Same tests, but just increment (H) increased another factor of 10
new results (increment increased by factor of 10 )(June 16, 2012)
NOTE: sometimes maxima's error is less - but because it reached the maximum time (set to 15 minutes) before it reached maximum iterations.
Pre-calculated factorials to increase efficiency. Running test.
- June 13, 2012
Same tests, but increment and number of iterations both increased by factors of 10
new results (mult H & and iterations by factors of 10 )(June 13, 2012)
Same tests, but just number of iterations increased another factor of 10
new results (iterations increased by factor of 10 )(June 13, 2012)
- June 12, 2012
After much thought and work on systems of equations I started getting
The poorest are mtest5.ode, mtest9.ode, mtest9_rev.ode and log.ode.
BTW: mtest6.ode and mtest6_rev.ode were previously called complicated.ode and complicatedrev.ode
- June 8, 2012
I added a capability to give initial conditions of order for a dependent variable higher than what appears on the LHS, to go along
with permitting them on the RHS. This, I think is needed for some systems of equations. For such terms, there values will not be updated
in the "atomall" loop. Rather, they will be projected in the "jump_series" function. I have all this coded and running, but have not
yet made any examples where it is really used. I had to make a change to the "ode" files in order to accomplish this. Next I need to
update the existing ode files (I have only done one.), and come up with some real examples.
- June 6, 2012
Made global variable in omnisode to set which of 2 ways to jump - relating to the subiter methods and multiple equations.
It runs but too tired to really test. Need better sets of equations - I have some ideas - but too tired now.
- June 5, 2012
Had to make changes to omnisode to handle the sort of problem I was thinking of - I think I need a yet harder problem!
- Table comparing tests using three methods on difficult system of equations
I realized I made a mistake on the prior table - I had the method overridden in "ode.over." I reran it and got these results:
Table comparing tests (2nd try) using three methods on difficult system of equations
I need interdependent equations at higher orders I think. I think maybe I have now made it so that one need not have the highest derivatives on the LHS!
- June 4, 2012
Uploaded omnisode0010.tar.gz to The Sourceforge Download page for Omnisode
Modified omnisode to use either one sub-iteration or total order sub-iterations or total order + max terms sub-iterations according to the variable glob_subiter_method is set to 1, 2 or 3 in the ".ode" file or "ode.over" file. This is to make testing easier. I need to invent an example that would seem to require the higher number of sub-iterations.
- June 3, 2012
Working on improving my HTML and omnisode documentation.
I made documentation on the script files used to run omnisode.
- June 2, 2012
Testing various numbers of sub-iterations with various increments:
- H = 0.00001
- H = 0.00005
- H = 0.0001
- H = 0.0002
- H = 0.0005
- H = 0.001
- June 1, 2012
I am uncertain of something having to do with the solution of multiple (simultaneous) equations. Obviously some either have no
solution or no unique solution. E.g.:
diff(y1,x,1) = diff(y2,x,1)
diff(y2,x,1) = diff(y1,x,1)
Has as solutions any differentiable equations where y1 and y2 differ by s constant.
When I first started I first tried sorting the equations - but this cannot always be done. Then I tried executing the code (atomall)
in a loop, until the solution stopped changing much, but this was inadequate. Then I tried executing the loop the number of times
there were equations, but I decided this was insufficient. I currently have the loop executing the sum of the orders of the equations
times. I am not sure if this is adequate or not. Perhaps one needs to add the number of terms of the series to this. Of course, one
still cannot arrive at a solution when there is none or a unique solution when there is none, but perhaps one can usually (except when
the function is badly behaving for some other reason) obtain a solution in cases where there is, indeed, a unique solution.
I did some testing on this while watching the Celtics
- May 31, 2012
Uploaded new test results with links to details.
Waiting till I have things better tested and documented before another release.
- May 27, 2012
Testing a nonlinear equation.
Updated information on omnisode on Sourceforge and added latest tar file to its download page.
- Msy 26, 2012
Ran automated test.
Also ran tests on systems of equations.
- May 25, 2012
I have realized a slight problem with how I track the "revision" of the code in my table.html. Fairly often, especially when I
reach a staple point, I run a script to save the three main pieces of ruby code. They are saved with numbers I specify, and the script
prevents overwriting. It also saves the revision numbers in a file, used in the table.html when a test is recorded. Well, in testing
I noticed when a test completed all the way to the end point (no timeout or maximum iterations reached), the cell containing
the expected time remaining was not filled, and everything after it is shifted a cell to the left. So I worked fixing it. But I did not
save my work after every attempt, so the "revision" codes only indicate the last time the code was saved - not the last time it
I've worked more on the html table The prior results are at
html table 0001
- May 24, 2012
I spent the last two days implementing automatic testing - constructing an
html table automatically.
I haven't made any more changes to the algorithm. It will take a little time to decide how to package the automatic table
- May 23, 2012
The new tar file is omnisode0007.tar.gz which may be obtained at the Sourceforge download page for Omnisode.
There is one ode file format for both Maple and Maxima.
The ode file format was documented yesterday.
I briefly tested all ode files with maxima, and all seemed OK except complicated.ode and log.ode. Complicated seemed stuck. Maybe
it was just very slow. Complicated.ode seemed OK with Maple. Also log.ode was not very accurate. But logs are difficult.
- May 22, 2012
Maple: Sample generated Maple program
Maple: Sample Maple program results
Maxima: Sample generated Maxima program
Maple: Sample Maxima program results
I have written and tested making the ode files language independent.
I have documented the ode files.
It will be tomorrow till I get the tar file made - less than a half hour.
- May 21, 2012
More testing of new version (tar file from May 20).
I believe Maple is limited to about 32 digits precision. That's all I seem to get, and there is a negligible time difference in
calculations assigning either 32 or 100 to "Digits".
- May 20, 2012
After considerable effort trying to get my software to use maxima big floats, I realized that maxima does not support many functions
(e.g. trig) for big floats. My software, thus, only supports 16 digits for maxima. I did, however, fix a few other things in trying.
I have had time to do only limited testing, but I know of no problems, other than that the order of poles seems a little off, although
the radius of convergence seems ok. Next I plan to do extensive testing and work on documentation.
omnisode0006.tar.gz was made.
- May 19, 2012
It works for both Maple and Maxima now.
The biggest difficulty had to do with arrays in Maxima. It seems they are call by value.
I ended up using global arrays, and having the Ruby program generate code for each separate array.
I haven't had enough time to test enough, but the code is greatly improved!
The latest results are at:
omnisode0004.tar.gz was made.
I always worry that something will happen to me, before I finish.
- May 15, 2012
It is working pretty well now.
I need to refamiliarize myself with sourceforge's tools but there is a tar file on my own web site for now:
omnitar0001.tar.gz was made. (removed as obsolete)
The testrecords.html has also been updated.
There is one sample output (and Maple echos the program) at Sample Solution
- May 14, 2012
I have fixed a couple problems including multiple equations.
The results have been added to the testrecords.html
- May 12, 2012
I am making some progress on this again the past few weeks.
I talk about this some on my blog Dennis' Daze
I am not releasing anything yet, but may start using the source version control tool soon.
A table of results of testing is at testrecords.html
- December 12, 2010
I have abandoned this project.
There are problems I am unable to solve and my mind is not working well.
I intend to devote myself to philosophy - and may do some prolog programming for that.
- These Ruby programs generate programs in Maple or Ruby to solve Systems of Ordinary Differential Equations.
A long Taylor series method, pioneered by Prof. Y.F. Chang, who taught at the University of Nebraska
in the late 1970's when I was a graduate student there, is used. The number of
Taylor series terms can be specified in the problem file, though it is usually
30. The Maple version is gone as, at least for now, I do not have Maple available to me.
I have started on a maxima version. In it I am taking time to remove some of the biggest
inefficiencies. Also I may have found an error that I am trying to resolve. The Ruby version below
would possibly also contain the error. But no one has reported any problems to me.
- Coming soon here - a maxima version
- There are problems with the Ruby and Maple versions which will have to wait until the maxima version is released.
Y. F. Chang's Draft Typescript(1978) on Taylor Series
Chang on ATS
- The Taylor series terms are used to calculate the values of the dependent
variables, and also (optionally) the radius of convergence and order of any
singularities. These can be used to reduce the size
of the next increment. These can be used to reduce the size
of the next increment. Increasing it in the middle results in large
errors because each successive term is multiplied by increasing powers of
the ratio of the new h to the old h. In the beginning I determine h by starting with a very
small h and increase it until the estimated truncation error is as large as specified as permissible.
- To visit the sourceforge project web page try Sode Project Page
- For info on RubyApfp:
RubyApfp Home Page