Dear Prof. Preben,

I check your code with some other integrals, and comparing their integral plots with first 2  proposal methods. Your idea, and code work perfectly. You are a Maple genius. Thanks you so much.

May God bless you with good health, lots of happiness and love.

Dear Prof. Preben,

I check your code with some other integrals, and comparing their integral plots with first 2  proposal methods. Your idea, and code work perfectly. You are a Maple genius. Thanks you so much.

May God bless you with good health, lots of happiness and love.

Thanks Preben, your code looks like awesome.

Would you please send me your worksheet? I have the following error when running dsolve

"Error, (in dsolve/numeric/process_input) system must be entered as a set/list of expressions/equations"

Thanks Preben, your code looks like awesome.

Would you please send me your worksheet? I have the following error when running dsolve

"Error, (in dsolve/numeric/process_input) system must be entered as a set/list of expressions/equations"

@Carl Love : thanks you Carl for your effort. I'm still thinking how to calculate the integral. As you know, we only need to solve first 9 differential equations to obtain the solutions, and the last equation doesn't affect other variables. Is there any way that Maple only solves first 9 equations and then substitutes the resluts to the last equation at each instant? 

Update: In addition, without last ODE, by using "eliminate" function in order to remove lambda 1 and lambda 2 from the equations, Maple can solve first 9 ODEs without any errors or warning. 

@Carl Love : thanks you Carl for your effort. I'm still thinking how to calculate the integral. As you know, we only need to solve first 9 differential equations to obtain the solutions, and the last equation doesn't affect other variables. Is there any way that Maple only solves first 9 equations and then substitutes the resluts to the last equation at each instant? 

Update: In addition, without last ODE, by using "eliminate" function in order to remove lambda 1 and lambda 2 from the equations, Maple can solve first 9 ODEs without any errors or warning. 

[I am very sorry, STHence, I accidentally overwrote your reply when I was trying to put in my own. But the essence was that you said that you had forgotten to include the with(plots) command in your posted worksheet. Feel free to post your reply again, if it matters to you, and I will edit it back into this position. ---Carl Love]

[I am very sorry, STHence, I accidentally overwrote your reply when I was trying to put in my own. But the essence was that you said that you had forgotten to include the with(plots) command in your posted worksheet. Feel free to post your reply again, if it matters to you, and I will edit it back into this position. ---Carl Love]

Thanks Markiyan so much. Your answer is great.

In addition, this DAE solver may have the following error when we increase the simulation time:

"Warning, cannot evaluate the solution further right of 3.1304860, maxfun limit exceeded (see ?dsolve,maxfun for details)"

However, if "maxfun" increases, this error can be eliminated.

For example,

> odeplot(solution, [t, y[A](t)], 0 .. 10, numpoints = 100);
Warning, cannot evaluate the solution further right of 3.1304860, maxfun limit exceeded (see ?dsolve,maxfun for details)

>

proc(x_rkf45)  ...  end;

>odeplot(solution3, [t, y[A](t)], 0 .. 10, numpoints = 100)

Thanks Markiyan so much. Your answer is great.

In addition, this DAE solver may have the following error when we increase the simulation time:

"Warning, cannot evaluate the solution further right of 3.1304860, maxfun limit exceeded (see ?dsolve,maxfun for details)"

However, if "maxfun" increases, this error can be eliminated.

For example,

> odeplot(solution, [t, y[A](t)], 0 .. 10, numpoints = 100);
Warning, cannot evaluate the solution further right of 3.1304860, maxfun limit exceeded (see ?dsolve,maxfun for details)

>

proc(x_rkf45)  ...  end;

>odeplot(solution3, [t, y[A](t)], 0 .. 10, numpoints = 100)

Thanks you Acer and Markiyan so much. All two answers are great, and having the same result for the equation in my question.

For Markiyan's answer, it's more convenient to calculate any functions from the solutions without re-solving the differential equations. However, when applying this method to my differential equations (different from the equations in my question) which having 9 variables, it takes very long time to calculate the integral of solution.

For Preben's answer, less time comsumption is its most advantage.

Thanks you Acer and Markiyan so much. All two answers are great, and having the same result for the equation in my question.

For Markiyan's answer, it's more convenient to calculate any functions from the solutions without re-solving the differential equations. However, when applying this method to my differential equations (different from the equations in my question) which having 9 variables, it takes very long time to calculate the integral of solution.

For Preben's answer, less time comsumption is its most advantage.

Thanks you Preben for your help. It's a very good idea.

Thanks you Preben for your help. It's a very good idea.

Thanks you Markiyan, your code works very well.