@ The 64bit Maple2015 comes with its own compiler. Unfortunately we can use only Digits =15 when you call dsolve with compile=true. I think this needs attention from Maplesoft folks. As of now, no gain in 64bit compared to 32 bit compiler. I was hoping that 64 bit with compile will allow for higher number of digits to be used for dsolve with compile=true, may be I am missing some information/theory.

For lsode, even lsode[backfull] works for this problem. I would suggest that you use this if this problem is stiff. Backward solvers are better for stiff problems.

Also, method = mebdfi works. But won't give result at x = 0.1, but works from 0.1 to 1. 

Both of these methods are slow and won't take compile=true.

Stiff =true works if use Digits:=20. This means that you can't use compile=true option for the default solvers as well. This might have changed in recent 64 bit versions, not sure.

Your problem seems to very badly scaled. The Gamma term seems to be nearly zero till x = 0.8 and then increases. Rescaling to a different x variable might help.

@cr-ms 

Try to change the DAE solver used. Some solvers might consume too much RAM. I am not sure 64GB will help for any general case. If you keep increasing the size, you might face the same problem.

To add to this, Maple has provided lsode, gear and other methods just for educational purposes. They are not implemented in an efficient manner. (They can only solve ODEs while the original Gear or LSODE methods can do DAEs).

For stiff problems, using stiff=true or rosenbrock methods provide the same answers. For some problems, using stiff=true, implicit =true provides much more robust solver (the most robust solver in Maple based on my experience).

I would strongly suggest scaling the variables and using stiff solvers (use stiff=true). Also, try different tolerances (1e-20 or so if you are interested in very small numbers).

It is not a good idea to integrate past or near singualrity. Perhaps DAE formulation combined with events will help. If you post the mws fomrat code, I will try.

@ 

Carl Love posted similar approach. Sorry for the duplication.

Unfortunately help files are not good for this. The examples given are very easy and simple. Dealing with matrices, vectors or arrays and converting external C file to Maple is not well explained.

While Compiler Compile might work, this is not useful for future sessions and it is not trivial to use it for future purposes. 

Providing an example for converting external Fortran or C files for Newton Raphson method (not single variable, but multiple variables), Runge Kutta methods, etc will help a lot. My strong belief is that instead of adding more or more routines in Maple, Maplesoft should rather make it easy to use existing C or Fortran solvers or routines by providng more examples on wrappers, linking multiple DLL files or C files, etc. 

Sorry for the rant, of course if you were to post your Maple or C code, I can find out if I can do it. So far, I have been successful in accessing very few C files in dll form (exe format is much easier to do).

@Allan Wittkopf 

Thanks for your attention and help. Your solution though unelegant, is great as it combines fsolve in a procedural form. (In fact I was surprised that even compile = true works).  I think this is the approach you will need for DAEs as well where in you can find initial conditions for algebraic variables.

Though fsolve works for this case, it won't work for largescale problems as the ODE solution might need only 1e-6 or crude tolerances but fsolve needs more.

If I might I will add that the solvers convert the system to dy/dt = f form, not just linear (say Mdy/dt = f). If M is allowed to a function of y or at least constant, we can solve more problems.

@joha 

I would appreciate if you can provide meaningful ranges for the parameters of h. There is way to solve this without solving for dh/dt at t = 0. I would like to check my codes. If you have a reference, please mention the same.

Thanks

@Carl Love 

Thanks. Agreed. Perhaps time lapse in posts and I am always confused using the mapleprimes not sure of which post is current or latest.

@Preben Alsholm 

Even in the DAE format one can specify u(0)=alpha, and then add alpha to the parameters list. But alpha should come from the fsolve.

@Preben Alsholm 

If there are multiple solutions for the nonlinear equation, then there will be a different solution for each case. The user will have to decide what makes sense.

If all the roots are meaningul, then there is more than solution for the nonlinear ODE just like we have multiple roots for quadratic or nonlinear equations. 

There is no ODE solver in the world which can solve and give results for all possible solutions numerically (as of now in my experience and search). If the model is known to have multiple solutions, there are approaches to get those solutions and software available for the same (only the final steady state values).

@Carl Love 

The values of the parameters should be the set by the user. Non zero values also work. See code attached. But fsolve needs to be used to find the initial condition.

Download convertimplicitODEtoDAE.mws

@acer 

It works for me in Maple 17 with the commands/code I posted. I was trying to warn the user to not expect this method to work for large scale implicit ODEs.

@Dave L 

I think it will be great to have an fsolve option for a very simple Newton Raphson method, wherein we are guaranteed an ouput with Digits =15  and specified absolute and relative tolerances and a maximum number of iterations. Displaying the iterations is also a great addtion to have so that users will know if there is a progress being made.